73 character (len=255),
parameter:: subcategory=
"volgrid6d_class"
77 type(griddim_def) :: griddim
78 TYPE(datetime),
pointer :: time(:)
79 TYPE(vol7d_timerange),
pointer :: timerange(:)
80 TYPE(vol7d_level),
pointer :: level(:)
81 TYPE(volgrid6d_var),
pointer :: var(:)
82 TYPE(grid_id),
POINTER :: gaid(:,:,:,:)
83 REAL,
POINTER :: voldati(:,:,:,:,:,:)
84 integer :: time_definition
85 integer :: category = 0
90 MODULE PROCEDURE volgrid6d_init
96 MODULE PROCEDURE volgrid6d_delete, volgrid6dv_delete
102 MODULE PROCEDURE volgrid6d_read_from_file
103 MODULE PROCEDURE import_from_gridinfo, import_from_gridinfovv, &
104 volgrid6d_import_from_file
110 MODULE PROCEDURE volgrid6d_write_on_file
111 MODULE PROCEDURE export_to_gridinfo, export_to_gridinfov, export_to_gridinfovv,&
112 volgrid6d_export_to_file
118 MODULE PROCEDURE volgrid6d_transform_compute, volgrid6d_v7d_transform_compute,&
119 v7d_volgrid6d_transform_compute, v7d_v7d_transform_compute
125 MODULE PROCEDURE volgrid6d_transform, volgrid6dv_transform,&
126 volgrid6d_v7d_transform, volgrid6dv_v7d_transform, v7d_volgrid6d_transform, &
131 MODULE PROCEDURE vg6d_wind_rot
135 MODULE PROCEDURE vg6d_wind_unrot
141 MODULE PROCEDURE display_volgrid6d,display_volgrid6dv
156 MODULE PROCEDURE vg6d_rounding, vg6dv_rounding
162 wind_rot,wind_unrot,vg6d_c2a,
display,volgrid6d_alloc,volgrid6d_alloc_vol, &
163 volgrid_get_vol_2d, volgrid_set_vol_2d, volgrid_get_vol_3d, volgrid_set_vol_3d
173 SUBROUTINE volgrid6d_init(this, griddim, time_definition, categoryappend)
174 TYPE(volgrid6d) :: this
175 TYPE(griddim_def),
OPTIONAL :: griddim
176 INTEGER,
INTENT(IN),
OPTIONAL :: time_definition
177 CHARACTER(len=*),
INTENT(in),
OPTIONAL :: categoryappend
179 character(len=512) :: a_name
181 if (
present(categoryappend))
then
182 call l4f_launcher(a_name,a_name_append=trim(subcategory)//
"."//trim(categoryappend))
184 call l4f_launcher(a_name,a_name_append=trim(subcategory))
186 this%category=l4f_category_get(a_name)
192 call init(this%griddim)
194 if (
present(griddim))
then
195 call copy (griddim,this%griddim)
198 CALL vol7d_var_features_init()
200 if(
present(time_definition))
then
201 this%time_definition = time_definition
203 this%time_definition = 0
206 nullify (this%time,this%timerange,this%level,this%var)
207 nullify (this%gaid,this%voldati)
209 END SUBROUTINE volgrid6d_init
222 SUBROUTINE volgrid6d_alloc(this, dim, ntime, nlevel, ntimerange, nvar, ini)
223 TYPE(volgrid6d),
INTENT(inout) :: this
224 TYPE(grid_dim),
INTENT(in),
OPTIONAL :: dim
225 INTEGER,
INTENT(in),
OPTIONAL :: ntime
226 INTEGER,
INTENT(in),
OPTIONAL :: nlevel
227 INTEGER,
INTENT(in),
OPTIONAL :: ntimerange
228 INTEGER,
INTENT(in),
OPTIONAL :: nvar
229 LOGICAL,
INTENT(in),
OPTIONAL :: ini
238 IF (
PRESENT(ini))
THEN
245 if (
present(dim))
call copy (dim,this%griddim%dim)
248 IF (
PRESENT(ntime))
THEN
250 IF (
ASSOCIATED(this%time))
DEALLOCATE(this%time)
254 ALLOCATE(this%time(ntime),stat=stallo)
257 CALL raise_fatal_error()
261 this%time(i) = datetime_miss
268 IF (
PRESENT(nlevel))
THEN
269 IF (nlevel >= 0)
THEN
270 IF (
ASSOCIATED(this%level))
DEALLOCATE(this%level)
274 ALLOCATE(this%level(nlevel),stat=stallo)
277 CALL raise_fatal_error()
286 IF (
PRESENT(ntimerange))
THEN
287 IF (ntimerange >= 0)
THEN
288 IF (
ASSOCIATED(this%timerange))
DEALLOCATE(this%timerange)
292 ALLOCATE(this%timerange(ntimerange),stat=stallo)
295 CALL raise_fatal_error()
304 IF (
PRESENT(nvar))
THEN
306 IF (
ASSOCIATED(this%var))
DEALLOCATE(this%var)
310 ALLOCATE(this%var(nvar),stat=stallo)
313 CALL raise_fatal_error()
317 CALL init(this%var(i))
323 end SUBROUTINE volgrid6d_alloc
334 SUBROUTINE volgrid6d_alloc_vol(this, ini, inivol, decode)
335 TYPE(volgrid6d),
INTENT(inout) :: this
336 LOGICAL,
INTENT(in),
OPTIONAL :: ini
337 LOGICAL,
INTENT(in),
OPTIONAL :: inivol
338 LOGICAL,
INTENT(in),
OPTIONAL :: decode
347 IF (
PRESENT(inivol))
THEN
353 IF (this%griddim%dim%nx > 0 .AND. this%griddim%dim%ny > 0)
THEN
356 IF (.NOT.
ASSOCIATED(this%var))
CALL volgrid6d_alloc(this, nvar=1, ini=ini)
357 IF (.NOT.
ASSOCIATED(this%time))
CALL volgrid6d_alloc(this, ntime=1, ini=ini)
358 IF (.NOT.
ASSOCIATED(this%level))
CALL volgrid6d_alloc(this, nlevel=1, ini=ini)
359 IF (.NOT.
ASSOCIATED(this%timerange))
CALL volgrid6d_alloc(this, ntimerange=1, ini=ini)
361 IF (optio_log(decode))
THEN
362 IF (.NOT.
ASSOCIATED(this%voldati))
THEN
367 ALLOCATE(this%voldati(this%griddim%dim%nx,this%griddim%dim%ny,&
368 SIZE(this%level),
SIZE(this%time), &
369 SIZE(this%timerange),
SIZE(this%var)),stat=stallo)
372 CALL raise_fatal_error()
377 IF (linivol) this%voldati = rmiss
382 IF (.NOT.
ASSOCIATED(this%gaid))
THEN
386 ALLOCATE(this%gaid(
SIZE(this%level),
SIZE(this%time), &
387 SIZE(this%timerange),
SIZE(this%var)),stat=stallo)
390 CALL raise_fatal_error()
404 this%gaid = grid_id_new()
410 &trying to allocate a volume with an invalid or unspecified horizontal grid')
411 CALL raise_fatal_error()
414 END SUBROUTINE volgrid6d_alloc_vol
430 SUBROUTINE volgrid_get_vol_2d(this, ilevel, itime, itimerange, ivar, voldati)
432 INTEGER,
INTENT(in) :: ilevel
433 INTEGER,
INTENT(in) :: itime
434 INTEGER,
INTENT(in) :: itimerange
435 INTEGER,
INTENT(in) :: ivar
436 REAL,
POINTER :: voldati(:,:)
438 IF (
ASSOCIATED(this%voldati))
THEN
439 voldati => this%voldati(:,:,ilevel,itime,itimerange,ivar)
442 IF (.NOT.
ASSOCIATED(voldati))
THEN
443 ALLOCATE(voldati(this%griddim%dim%nx,this%griddim%dim%ny))
445 CALL grid_id_decode_data(this%gaid(ilevel,itime,itimerange,ivar), voldati)
448 END SUBROUTINE volgrid_get_vol_2d
464 SUBROUTINE volgrid_get_vol_3d(this, itime, itimerange, ivar, voldati)
466 INTEGER,
INTENT(in) :: itime
467 INTEGER,
INTENT(in) :: itimerange
468 INTEGER,
INTENT(in) :: ivar
469 REAL,
POINTER :: voldati(:,:,:)
473 IF (
ASSOCIATED(this%voldati))
THEN
474 voldati => this%voldati(:,:,:,itime,itimerange,ivar)
477 IF (.NOT.
ASSOCIATED(voldati))
THEN
478 ALLOCATE(voldati(this%griddim%dim%nx,this%griddim%dim%ny,
SIZE(this%level)))
482 DO ilevel = 1,
SIZE(this%level)
484 CALL grid_id_decode_data(this%gaid(ilevel,itime,itimerange,ivar), &
492 END SUBROUTINE volgrid_get_vol_3d
506 SUBROUTINE volgrid_set_vol_2d(this, ilevel, itime, itimerange, ivar, voldati)
508 INTEGER,
INTENT(in) :: ilevel
509 INTEGER,
INTENT(in) :: itime
510 INTEGER,
INTENT(in) :: itimerange
511 INTEGER,
INTENT(in) :: ivar
512 REAL,
INTENT(in) :: voldati(:,:)
514 IF (
ASSOCIATED(this%voldati))
THEN
517 CALL grid_id_encode_data(this%gaid(ilevel,itime,itimerange,ivar), voldati)
520 END SUBROUTINE volgrid_set_vol_2d
534 SUBROUTINE volgrid_set_vol_3d(this, itime, itimerange, ivar, voldati)
536 INTEGER,
INTENT(in) :: itime
537 INTEGER,
INTENT(in) :: itimerange
538 INTEGER,
INTENT(in) :: ivar
539 REAL,
INTENT(in) :: voldati(:,:,:)
543 IF (
ASSOCIATED(this%voldati))
THEN
548 DO ilevel = 1,
SIZE(this%level)
550 CALL grid_id_encode_data(this%gaid(ilevel,itime,itimerange,ivar), &
558 END SUBROUTINE volgrid_set_vol_3d
564 SUBROUTINE volgrid6d_delete(this)
567 INTEGER :: i, ii, iii, iiii
573 if (
associated(this%gaid))
then
575 DO i=1 ,
SIZE(this%gaid,1)
576 DO ii=1 ,
SIZE(this%gaid,2)
577 DO iii=1 ,
SIZE(this%gaid,3)
578 DO iiii=1 ,
SIZE(this%gaid,4)
579 CALL delete(this%gaid(i,ii,iii,iiii))
584 DEALLOCATE(this%gaid)
595 if (
associated( this%time ))
deallocate(this%time)
596 if (
associated( this%timerange ))
deallocate(this%timerange)
597 if (
associated( this%level ))
deallocate(this%level)
598 if (
associated( this%var ))
deallocate(this%var)
600 if (
associated(this%voldati))
deallocate(this%voldati)
604 call l4f_category_delete(this%category)
606 END SUBROUTINE volgrid6d_delete
618 subroutine volgrid6d_write_on_file (this,unit,description,filename,filename_auto)
621 integer,
optional,
intent(inout) :: unit
622 character(len=*),
intent(in),
optional :: filename
623 character(len=*),
intent(out),
optional :: filename_auto
624 character(len=*),
INTENT(IN),
optional :: description
627 character(len=254) :: ldescription,arg,lfilename
628 integer :: ntime, ntimerange, nlevel, nvar
630 logical :: opened,exist
642 call date_and_time(values=tarray)
645 if (
present(description))
then
646 ldescription=description
648 ldescription=
"Volgrid6d generated by: "//trim(arg)
651 if (.not.
present(unit))
then
662 lfilename=trim(arg)//
".vg6d"
663 if (
index(arg,
'/',back=.true.) > 0) lfilename=lfilename(
index(arg,
'/',back=.true.)+1 : )
665 if (
present(filename))
then
666 if (filename /=
"")
then
671 if (
present(filename_auto))filename_auto=lfilename
674 inquire(unit=lunit,opened=opened)
675 if (.not. opened)
then
676 inquire(file=lfilename,exist=exist)
677 if (exist)
CALL raise_error(
'file exist; cannot open new file')
678 if (.not.exist)
open (unit=lunit,file=lfilename,form=
"UNFORMATTED")
682 if (
associated(this%time)) ntime=
size(this%time)
683 if (
associated(this%timerange)) ntimerange=
size(this%timerange)
684 if (
associated(this%level)) nlevel=
size(this%level)
685 if (
associated(this%var)) nvar=
size(this%var)
688 write(unit=lunit)ldescription
689 write(unit=lunit)tarray
692 write(unit=lunit) ntime, ntimerange, nlevel, nvar
695 if (
associated(this%time))
call write_unit(this%time, lunit)
696 if (
associated(this%level))
write(unit=lunit)this%level
697 if (
associated(this%timerange))
write(unit=lunit)this%timerange
698 if (
associated(this%var))
write(unit=lunit)this%var
703 if (
associated(this%voldati))
write(unit=lunit)this%voldati
705 if (.not.
present(unit))
close(unit=lunit)
707 end subroutine volgrid6d_write_on_file
716 subroutine volgrid6d_read_from_file (this,unit,filename,description,tarray,filename_auto)
719 integer,
intent(inout),
optional :: unit
720 character(len=*),
INTENT(in),
optional :: filename
721 character(len=*),
intent(out),
optional :: filename_auto
722 character(len=*),
INTENT(out),
optional :: description
723 integer,
intent(out),
optional :: tarray(8)
725 integer :: ntime, ntimerange, nlevel, nvar
727 character(len=254) :: ldescription,lfilename,arg
728 integer :: ltarray(8),lunit
729 logical :: opened,exist
737 if (.not.
present(unit))
then
748 lfilename=trim(arg)//
".vg6d"
749 if (
index(arg,
'/',back=.true.) > 0) lfilename=lfilename(
index(arg,
'/',back=.true.)+1 : )
751 if (
present(filename))
then
752 if (filename /=
"")
then
757 if (
present(filename_auto))filename_auto=lfilename
760 inquire(unit=lunit,opened=opened)
761 if (.not. opened)
then
762 inquire(file=lfilename,exist=exist)
763 IF (.NOT. exist)
CALL raise_fatal_error(
'file '//trim(lfilename)//
' does not exist, cannot open')
764 open (unit=lunit,file=lfilename,form=
"UNFORMATTED")
767 read(unit=lunit)ldescription
768 read(unit=lunit)ltarray
770 call l4f_log(l4f_info,
"Info: reading volgrid6d from file: "//trim(lfilename))
771 call l4f_log(l4f_info,
"Info: description: "//trim(ldescription))
774 if (
present(description))description=ldescription
775 if (
present(tarray))tarray=ltarray
779 read(unit=lunit) ntime, ntimerange, nlevel, nvar
782 call volgrid6d_alloc (this, &
783 ntime=ntime, ntimerange=ntimerange, nlevel=nlevel, nvar=nvar)
785 call volgrid6d_alloc_vol (this)
787 if (
associated(this%time))
call read_unit(this%time, lunit)
788 if (
associated(this%level))
read(unit=lunit)this%level
789 if (
associated(this%timerange))
read(unit=lunit)this%timerange
790 if (
associated(this%var))
read(unit=lunit)this%var
795 if (
associated(this%voldati))
read(unit=lunit)this%voldati
797 if (.not.
present(unit))
close(unit=lunit)
799 end subroutine volgrid6d_read_from_file
821 SUBROUTINE import_from_gridinfo(this, gridinfo, force, dup_mode, clone, &
825 LOGICAL,
INTENT(in),
OPTIONAL :: force
826 INTEGER,
INTENT(in),
OPTIONAL :: dup_mode
827 LOGICAL ,
INTENT(in),
OPTIONAL :: clone
828 LOGICAL,
INTENT(IN),
OPTIONAL :: isanavar
830 CHARACTER(len=255) :: type
831 INTEGER :: itime0, itimerange0, itime1, itimerange1, itime, itimerange, &
832 ilevel, ivar, ldup_mode
836 REAL,
ALLOCATABLE :: tmpgrid(:,:)
838 IF (
PRESENT(dup_mode))
THEN
844 call get_val(this%griddim,proj_type=type)
847 call l4f_category_log(this%category,l4f_debug,
"import_from_gridinfo: "//trim(type))
850 if (.not.
c_e(type))
then
851 call copy(gridinfo%griddim, this%griddim)
855 CALL volgrid6d_alloc_vol(this, ini=.true.)
857 else if (.not. (this%griddim == gridinfo%griddim ))
then
860 "volgrid and gridinfo grid type or size are different, gridinfo rejected")
867 ilevel =
index(this%level, gridinfo%level)
868 IF (ilevel == 0 .AND. optio_log(force))
THEN
869 ilevel =
index(this%level, vol7d_level_miss)
870 IF (ilevel /= 0) this%level(ilevel) = gridinfo%level
873 IF (ilevel == 0)
THEN
875 "volgrid6d: level not valid for volume, gridinfo rejected")
880 IF (optio_log(isanavar))
THEN
882 itime1 =
SIZE(this%time)
884 itimerange1 =
SIZE(this%timerange)
886 correctedtime = gridinfo%time
887 IF (this%time_definition == 1 .OR. this%time_definition == 2) correctedtime = correctedtime + &
888 timedelta_new(sec=gridinfo%timerange%p1)
889 itime0 =
index(this%time, correctedtime)
890 IF (itime0 == 0 .AND. optio_log(force))
THEN
891 itime0 =
index(this%time, datetime_miss)
892 IF (itime0 /= 0) this%time(itime0) = correctedtime
894 IF (itime0 == 0)
THEN
896 "volgrid6d: time not valid for volume, gridinfo rejected")
902 correctedtimerange = gridinfo%timerange
903 IF (this%time_definition == 2) correctedtimerange%p1 = 0
904 itimerange0 =
index(this%timerange, correctedtimerange)
905 IF (itimerange0 == 0 .AND. optio_log(force))
THEN
906 itimerange0 =
index(this%timerange, vol7d_timerange_miss)
907 IF (itimerange0 /= 0) this%timerange(itimerange0) = gridinfo%timerange
909 IF (itimerange0 == 0)
THEN
911 "volgrid6d: timerange not valid for volume, gridinfo rejected")
915 itimerange1 = itimerange0
918 ivar =
index(this%var, gridinfo%var)
919 IF (ivar == 0 .AND. optio_log(force))
THEN
920 ivar =
index(this%var, volgrid6d_var_miss)
921 IF (ivar /= 0) this%var(ivar) = gridinfo%var
925 "volgrid6d: var not valid for volume, gridinfo rejected")
930 DO itimerange = itimerange0, itimerange1
931 DO itime = itime0, itime1
932 IF (
ASSOCIATED(this%gaid))
THEN
934 IF (
c_e(this%gaid(ilevel,itime,itimerange,ivar)))
THEN
938 IF (optio_log(
clone))
CALL delete(this%gaid(ilevel,itime,itimerange,ivar))
941 IF (optio_log(
clone))
THEN
942 CALL copy(gridinfo%gaid, this%gaid(ilevel,itime,itimerange,ivar))
947 this%gaid(ilevel,itime,itimerange,ivar) = gridinfo%gaid
950 IF (
ASSOCIATED(this%voldati))
THEN
951 IF (.NOT.dup .OR. ldup_mode == 0)
THEN
952 this%voldati(:,:,ilevel,itime,itimerange,ivar) =
decode_gridinfo(gridinfo)
953 ELSE IF (ldup_mode == 1)
THEN
956 this%voldati(:,:,ilevel,itime,itimerange,ivar) = tmpgrid(:,:)
958 ELSE IF (ldup_mode == 2)
THEN
959 WHERE(.NOT.
c_e(this%voldati(:,:,ilevel,itime,itimerange,ivar)))
960 this%voldati(:,:,ilevel,itime,itimerange,ivar) =
decode_gridinfo(gridinfo)
967 "gaid not allocated, you probably need to call volgrid6d_alloc_vol first")
975 END SUBROUTINE import_from_gridinfo
982 SUBROUTINE export_to_gridinfo(this, gridinfo, itime, itimerange, ilevel, ivar, &
983 gaid_template, clone)
987 INTEGER :: itimerange
990 TYPE(
grid_id),
INTENT(in),
OPTIONAL :: gaid_template
991 LOGICAL,
INTENT(in),
OPTIONAL :: clone
994 LOGICAL :: usetemplate
995 REAL,
POINTER :: voldati(:,:)
1002 IF (.NOT.
c_e(this%gaid(ilevel,itime,itimerange,ivar)))
THEN
1004 CALL l4f_category_log(this%category,l4f_debug,
"empty gaid found, skipping export")
1009 usetemplate = .false.
1010 IF (
PRESENT(gaid_template))
THEN
1011 CALL copy(gaid_template, gaid)
1013 CALL l4f_category_log(this%category,l4f_debug,
"template cloned to a new gaid")
1015 usetemplate =
c_e(gaid)
1018 IF (.NOT.usetemplate)
THEN
1020 CALL copy(this%gaid(ilevel,itime,itimerange,ivar), gaid)
1022 CALL l4f_category_log(this%category,l4f_debug,
"original gaid cloned to a new one")
1025 gaid = this%gaid(ilevel,itime,itimerange,ivar)
1029 IF (this%time_definition == 1 .OR. this%time_definition == 2)
THEN
1030 correctedtime = this%time(itime) - &
1031 timedelta_new(sec=this%timerange(itimerange)%p1)
1033 correctedtime = this%time(itime)
1036 CALL init(gridinfo,gaid, this%griddim, correctedtime, this%timerange(itimerange), &
1037 this%level(ilevel), this%var(ivar))
1040 CALL export(gridinfo%griddim, gridinfo%gaid)
1042 IF (
ASSOCIATED(this%voldati))
THEN
1043 CALL encode_gridinfo(gridinfo, this%voldati(:,:,ilevel,itime,itimerange,ivar))
1044 ELSE IF (usetemplate)
THEN
1046 CALL volgrid_get_vol_2d(this, ilevel, itime, itimerange, ivar, voldati)
1051 END SUBROUTINE export_to_gridinfo
1071 SUBROUTINE import_from_gridinfovv(this, gridinfov, dup_mode, clone, decode, &
1072 time_definition, anavar, categoryappend)
1075 INTEGER,
INTENT(in),
OPTIONAL :: dup_mode
1076 LOGICAL ,
INTENT(in),
OPTIONAL :: clone
1077 LOGICAL,
INTENT(in),
OPTIONAL :: decode
1078 INTEGER,
INTENT(IN),
OPTIONAL :: time_definition
1079 CHARACTER(len=*),
INTENT(IN),
OPTIONAL :: anavar(:)
1080 CHARACTER(len=*),
INTENT(in),
OPTIONAL :: categoryappend
1082 INTEGER :: i, j, stallo
1083 INTEGER :: ngrid, ntime, ntimerange, nlevel, nvar, ltime_definition
1085 CHARACTER(len=512) :: a_name
1086 TYPE(
datetime),
ALLOCATABLE :: correctedtime(:)
1087 LOGICAL,
ALLOCATABLE :: isanavar(:)
1088 TYPE(vol7d_var) :: lvar
1092 if (
present(categoryappend))
then
1093 call l4f_launcher(a_name,a_name_append=trim(subcategory)//
"."//trim(categoryappend))
1095 call l4f_launcher(a_name,a_name_append=trim(subcategory))
1097 category=l4f_category_get(a_name)
1103 IF (
PRESENT(time_definition))
THEN
1104 ltime_definition = max(min(time_definition, 2), 0)
1106 ltime_definition = 0
1109 ngrid=count_distinct(gridinfov%array(1:gridinfov%arraysize)%griddim,back=.true.)
1111 ' different grid definition(s) found in input data')
1113 ALLOCATE(this(ngrid),stat=stallo)
1114 IF (stallo /= 0)
THEN
1116 CALL raise_fatal_error()
1119 IF (
PRESENT(categoryappend))
THEN
1120 CALL init(this(i), time_definition=ltime_definition, categoryappend=trim(categoryappend)//
"-vol"//
t2c(i))
1122 CALL init(this(i), time_definition=ltime_definition, categoryappend=
"vol"//
t2c(i))
1126 this(:)%griddim=pack_distinct(gridinfov%array(1:gridinfov%arraysize)%griddim, &
1130 ALLOCATE(isanavar(gridinfov%arraysize))
1131 isanavar(:) = .false.
1132 IF (
PRESENT(anavar))
THEN
1133 DO i = 1, gridinfov%arraysize
1134 DO j = 1,
SIZE(anavar)
1135 lvar =
convert(gridinfov%array(i)%var)
1136 IF (lvar%btable == anavar(j))
THEN
1137 isanavar(i) = .true.
1143 t2c(gridinfov%arraysize)//
' constant-data messages found in input data')
1146 IF (ltime_definition == 1 .OR. ltime_definition == 2)
THEN
1147 ALLOCATE(correctedtime(gridinfov%arraysize))
1148 correctedtime(:) = gridinfov%array(1:gridinfov%arraysize)%time
1149 DO i = 1, gridinfov%arraysize
1150 correctedtime(i) = correctedtime(i) + &
1151 timedelta_new(sec=gridinfov%array(i)%timerange%p1)
1154 IF (ltime_definition == 2)
THEN
1155 ALLOCATE(correctedtimerange(gridinfov%arraysize))
1156 correctedtimerange(:) = gridinfov%array(1:gridinfov%arraysize)%timerange
1157 correctedtimerange(:)%p1 = 0
1161 IF (
PRESENT(anavar))
THEN
1162 j = count((this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1163 .AND. .NOT.isanavar(:))
1166 ' has only constant data, this is not allowed')
1168 CALL raise_fatal_error()
1171 IF (ltime_definition == 1 .OR. ltime_definition == 2)
THEN
1172 ntime = count_distinct(correctedtime, &
1173 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1174 .AND. .NOT.isanavar(:), back=.true.)
1176 ntime = count_distinct(gridinfov%array(1:gridinfov%arraysize)%time, &
1177 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1178 .AND. .NOT.isanavar(:), back=.true.)
1180 IF (ltime_definition == 2)
THEN
1181 ntimerange = count_distinct(correctedtimerange, &
1182 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1183 .AND. .NOT.isanavar(:), back=.true.)
1185 ntimerange = count_distinct(gridinfov%array(1:gridinfov%arraysize)%timerange, &
1186 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1187 .AND. .NOT.isanavar(:), back=.true.)
1189 nlevel = count_distinct(gridinfov%array(1:gridinfov%arraysize)%level, &
1190 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1192 nvar = count_distinct(gridinfov%array(1:gridinfov%arraysize)%var, &
1193 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1200 CALL volgrid6d_alloc(this(i),this(i)%griddim%dim,ntime=ntime, &
1201 ntimerange=ntimerange,nlevel=nlevel,nvar=nvar)
1203 IF (ltime_definition == 1 .OR. ltime_definition == 2)
THEN
1204 this(i)%time = pack_distinct(correctedtime, ntime, &
1205 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1206 .AND. .NOT.isanavar(:), back=.true.)
1208 this(i)%time = pack_distinct(gridinfov%array(1:gridinfov%arraysize)%time, ntime, &
1209 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1210 .AND. .NOT.isanavar(:), back=.true.)
1212 CALL sort(this(i)%time)
1214 IF (ltime_definition == 2)
THEN
1215 this(i)%timerange = pack_distinct(correctedtimerange, ntimerange, &
1216 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1217 .AND. .NOT.isanavar(:), back=.true.)
1219 this(i)%timerange = pack_distinct(gridinfov%array(1:gridinfov%arraysize)%timerange, &
1220 ntimerange, mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1221 .AND. .NOT.isanavar(:), back=.true.)
1223 CALL sort(this(i)%timerange)
1225 this(i)%level=pack_distinct(gridinfov%array(1:gridinfov%arraysize)%level, &
1226 nlevel,mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1228 CALL sort(this(i)%level)
1230 this(i)%var=pack_distinct(gridinfov%array(1:gridinfov%arraysize)%var, nvar, &
1231 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1237 CALL volgrid6d_alloc_vol(this(i), decode=decode)
1241 IF (ltime_definition == 1 .OR. ltime_definition == 2)
DEALLOCATE(correctedtime)
1242 IF (ltime_definition == 2)
DEALLOCATE(correctedtimerange)
1244 DO i = 1, gridinfov%arraysize
1248 CALL l4f_category_log(category,L4F_INFO, &
1249 "to
volgrid6d index:
"//t2c(index(this%griddim, gridinfov%array(i)%griddim)))
1252 CALL import(this(index(this%griddim, gridinfov%array(i)%griddim)), &
1253 gridinfov%array(i), dup_mode=dup_mode, clone=clone, isanavar=isanavar(i))
1257 !chiudo il logger temporaneo
1258 CALL l4f_category_delete(category)
1260 END SUBROUTINE import_from_gridinfovv
1263 !> Export a \a volgrid6d object to an \a arrayof_gridinfo object.
1264 !! The multidimensional \a volgrid6d structure is serialized into a
1265 !! one-dimensional array of gridinfo_def objects, which is allocated
1266 !! to the proper size if not already allocated, or it is extended
1267 !! keeping the old data if any.
1268 SUBROUTINE export_to_gridinfov(this, gridinfov, gaid_template, clone)
1269 TYPE(volgrid6d),INTENT(inout) :: this !< volume to be exported
1270 TYPE(arrayof_gridinfo),INTENT(inout) :: gridinfov !< output array of gridinfo_def objects
1271 TYPE(grid_id),INTENT(in),OPTIONAL :: gaid_template !< \a grid_id template to be used for output data replacing the one contained in \a this
1272 LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE., clone the grid_id included in \a this rather than making a shallow copy
1274 INTEGER :: i ,itime, itimerange, ilevel, ivar
1275 INTEGER :: ntime, ntimerange, nlevel, nvar
1276 TYPE(gridinfo_def) :: gridinfol
1279 CALL l4f_category_log(this%category,L4F_DEBUG,"start export_to_gridinfov
")
1282 ! this is necessary in order not to repeatedly and uselessly copy the
1283 ! same array of coordinates for each 2d grid during export, the
1284 ! side-effect is that the computed projection in this is lost
1285 CALL dealloc(this%griddim%dim)
1288 ntime=size(this%time)
1289 ntimerange=size(this%timerange)
1290 nlevel=size(this%level)
1294 DO itimerange=1,ntimerange
1298 CALL init(gridinfol)
1299 CALL export(this, gridinfol, itime, itimerange, ilevel, ivar, &
1300 gaid_template=gaid_template, clone=clone)
1301 IF (c_e(gridinfol%gaid)) THEN
1302 CALL insert(gridinfov, gridinfol)
1304 CALL delete(gridinfol)
1312 END SUBROUTINE export_to_gridinfov
1315 !> Export an array of \a volgrid6d objects to an \a arrayof_gridinfo object.
1316 !! The multidimensional \a volgrid6d structures are serialized into a
1317 !! one-dimensional array of gridinfo_def objects, which is allocated
1318 !! to the proper size if not already allocated, or it is extended
1319 !! keeping the old data if any.
1320 SUBROUTINE export_to_gridinfovv(this, gridinfov, gaid_template, clone)
1323 TYPE(volgrid6d),INTENT(inout) :: this(:) !< volume array to be exported
1324 TYPE(arrayof_gridinfo),INTENT(inout) :: gridinfov !< output array of gridinfo_def objects
1325 TYPE(grid_id),INTENT(in),OPTIONAL :: gaid_template !< \a grid_id template to be used for output data replacing the one contained in \a this
1326 LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE., clone the grid_id included in \a this rather than making a shallow copy
1330 DO i = 1, SIZE(this)
1332 CALL l4f_category_log(this(i)%category,L4F_DEBUG, &
1333 "export_to_gridinfovv grid
index:
"//t2c(i))
1335 CALL export(this(i), gridinfov, gaid_template=gaid_template, clone=clone)
1338 END SUBROUTINE export_to_gridinfovv
1341 Import the content of a supported file (like grib or gdal-supported
!>
1342 !! format) into an array of \a volgrid6d objects. This method imports
1343 !! a set of gridded fields from a file object into a suitable number
1344 !! of \a volgrid6d objects. The data are imported by creating a
1345 !! temporary \a gridinfo object, importing it into the \a volgrid6d
1346 !! object cloning the gaid's and then destroying the gridinfo, so it
1347 !! works similarly to volgrid6d_class::import_from_gridinfovv() method.
1348 !! For a detailed explanation of the \a anavar argument, see the
1349 !! documentation of volgrid6d_class::import_from_gridinfovv() method.
1350 SUBROUTINE volgrid6d_import_from_file(this, filename, dup_mode, decode, &
1351 time_definition, anavar, categoryappend)
1352 importTYPE(volgrid6d),POINTER :: this(:) !< object in which to
1353 importCHARACTER(len=*),INTENT(in) :: filename !< name of file from which to
1354 INTEGER,INTENT(in),OPTIONAL :: dup_mode !< determines the behavior in case of duplicate metadata: if \a dup_mode is not provided or 0, a duplicate field overwrites, if \a dup_mode is 1, duplicate fields are merged with priority to the last
1355 LOGICAL,INTENT(in),OPTIONAL :: decode !< if provided and \a .FALSE. the data volume in the elements of \a this is not allocated and successive work will be performed on grid_id's
1356 INTEGER,INTENT(IN),OPTIONAL :: time_definition !< 0=time is reference time; 1=time is validity time
1357 CHARACTER(len=*),INTENT(IN),OPTIONAL :: anavar(:) !< list of variables (B-table code equivalent) to be treated as time-independent data
1358 character(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
1360 TYPE(arrayof_gridinfo) :: gridinfo
1362 CHARACTER(len=512) :: a_name
1366 IF (PRESENT(categoryappend))THEN
1367 CALL l4f_launcher(a_name,a_name_append= &
1368 TRIM(subcategory)//".
"//TRIM(categoryappend))
1370 CALL l4f_launcher(a_name,a_name_append=TRIM(subcategory))
1372 category=l4f_category_get(a_name)
1374 CALL import(gridinfo, filename=filename, categoryappend=categoryappend)
1376 IF (gridinfo%arraysize > 0) THEN
1378 CALL import(this, gridinfo, dup_mode=dup_mode, clone=.TRUE., decode=decode, &
1379 time_definition=time_definition, anavar=anavar, &
1380 categoryappend=categoryappend)
1382 CALL l4f_category_log(category,L4F_INFO,"deleting gridinfo
")
1383 CALL delete(gridinfo)
1386 CALL l4f_category_log(category,L4F_INFO,"file does not contain gridded data
")
1390 CALL l4f_category_delete(category)
1392 END SUBROUTINE volgrid6d_import_from_file
1395 !> High level method for exporting a volume array to file.
1396 !! All the information contained into an array of \a volgrid6d
1397 !! objects, i.e. dimension descriptors and data, is exported to a file
1398 !! using the proper output driver (typically grib_api for grib
1399 !! format). If a template is provided, it will determine the
1400 !! characteristic of the output file, otherwise the \a grid_id
1401 !! descriptors contained in the volgrid6d object will be used
1402 SUBROUTINE volgrid6d_export_to_file(this, filename, gaid_template, categoryappend)
1403 TYPE(volgrid6d) :: this(:) !< volume(s) to be exported
1404 CHARACTER(len=*),INTENT(in) :: filename !< output file name
1405 TYPE(grid_id),INTENT(in),OPTIONAL :: gaid_template !< template for the output file, if provided the grid_id information stored in the volgrid6d objects will be ignored
1406 character(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
1408 TYPE(arrayof_gridinfo) :: gridinfo
1410 CHARACTER(len=512) :: a_name
1412 IF (PRESENT(categoryappend)) THEN
1413 CALL l4f_launcher(a_name,a_name_append=TRIM(subcategory)//".
"//TRIM(categoryappend))
1415 CALL l4f_launcher(a_name,a_name_append=TRIM(subcategory))
1417 category=l4f_category_get(a_name)
1420 CALL l4f_category_log(category,L4F_DEBUG,"start
export to file
")
1423 CALL l4f_category_log(category,L4F_INFO,"writing
volgrid6d to grib file:
"//TRIM(filename))
1425 !IF (ASSOCIATED(this)) THEN
1426 CALL export(this, gridinfo, gaid_template=gaid_template, clone=.TRUE.)
1427 IF (gridinfo%arraysize > 0) THEN
1428 CALL export(gridinfo, filename)
1429 CALL delete(gridinfo)
1432 ! CALL l4f_category_log(category,L4F_INFO,"volume
volgrid6d is not associated
")
1436 CALL l4f_category_delete(category)
1438 END SUBROUTINE volgrid6d_export_to_file
1441 !> Array destructor for \a volgrid6d class.
1442 !! Delete an array of \a volgrid6d objects and deallocate the array
1444 SUBROUTINE volgrid6dv_delete(this)
1445 TYPE(volgrid6d),POINTER :: this(:) !< vector of volgrid6d object
1449 IF (ASSOCIATED(this)) THEN
1450 DO i = 1, SIZE(this)
1452 CALL l4f_category_log(this(i)%category,L4F_DEBUG, &
1455 CALL delete(this(i))
1460 END SUBROUTINE volgrid6dv_delete
1463 ! Internal method for performing grid to grid computations
1464 SUBROUTINE volgrid6d_transform_compute(this, volgrid6d_in, volgrid6d_out, &
1465 lev_out, var_coord_vol, clone)
1466 TYPE(grid_transform),INTENT(in) :: this ! oggetto di trasformazione per il grigliato
1467 type(volgrid6d), INTENT(in) :: volgrid6d_in ! oggetto da trasformare
1468 type(volgrid6d), INTENT(inout) :: volgrid6d_out ! oggetto trasformato; deve essere completo (init, alloc, alloc_vol)
1469 TYPE(vol7d_level),INTENT(in),OPTIONAL :: lev_out(:) ! vol7d_level object defining target vertical grid, for vertical interpolations
1470 INTEGER,INTENT(in),OPTIONAL :: var_coord_vol ! index of variable defining vertical coordinate values in input volume
1471 LOGICAL,INTENT(in),OPTIONAL :: clone ! se fornito e \c .TRUE., clona i gaid da volgrid6d_in a volgrid6d_out
1473 INTEGER :: ntime, ntimerange, inlevel, onlevel, nvar, &
1474 itime, itimerange, ilevel, ivar, levshift, levused, lvar_coord_vol, spos
1475 REAL,POINTER :: voldatiin(:,:,:), voldatiout(:,:,:), coord_3d_in(:,:,:)
1476 TYPE(vol7d_level) :: output_levtype
1480 call l4f_category_log(volgrid6d_in%category,L4F_DEBUG,"start volgrid6d_transform_compute
")
1488 lvar_coord_vol = optio_i(var_coord_vol)
1490 if (associated(volgrid6d_in%time))then
1491 ntime=size(volgrid6d_in%time)
1492 volgrid6d_out%time=volgrid6d_in%time
1495 if (associated(volgrid6d_in%timerange))then
1496 ntimerange=size(volgrid6d_in%timerange)
1497 volgrid6d_out%timerange=volgrid6d_in%timerange
1500 IF (ASSOCIATED(volgrid6d_in%level))THEN
1501 inlevel=SIZE(volgrid6d_in%level)
1503 IF (PRESENT(lev_out)) THEN
1504 onlevel=SIZE(lev_out)
1505 volgrid6d_out%level=lev_out
1506 ELSE IF (ASSOCIATED(volgrid6d_in%level))THEN
1507 onlevel=SIZE(volgrid6d_in%level)
1508 volgrid6d_out%level=volgrid6d_in%level
1511 if (associated(volgrid6d_in%var))then
1512 nvar=size(volgrid6d_in%var)
1513 volgrid6d_out%var=volgrid6d_in%var
1515 ! allocate once for speed
1516 .NOT.
IF (ASSOCIATED(volgrid6d_in%voldati)) THEN
1517 ALLOCATE(voldatiin(volgrid6d_in%griddim%dim%nx, volgrid6d_in%griddim%dim%ny, &
1520 .NOT.
IF (ASSOCIATED(volgrid6d_out%voldati)) THEN
1521 ALLOCATE(voldatiout(volgrid6d_out%griddim%dim%nx, volgrid6d_out%griddim%dim%ny, &
1525 CALL get_val(this, levshift=levshift, levused=levused)
1527 IF (c_e(lvar_coord_vol)) THEN
1528 CALL get_val(this%trans, output_levtype=output_levtype)
1529 .OR.
IF (output_levtype%level1 == 103 output_levtype%level1 == 108) THEN
1530 spos = firsttrue(volgrid6d_in%level(:) == vol7d_level_new(1))
1532 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1533 'output level '//t2c(output_levtype%level1)// &
1534 ' requested, but height/press of surface not provided in volume')
1536 .NOT..AND..NOT.
IF (c_e(levshift) c_e(levused)) THEN
1537 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1538 'internal inconsistence, levshift and levused undefined when they should be')
1544 ! IF (c_e(var_coord_vol)) THEN
1545 ! IF (ivar == var_coord_vol) CYCLE ! skip coordinate variable in output
1547 DO itimerange=1,ntimerange
1549 ! skip empty columns where possible, improve
1550 .AND.
IF (c_e(levshift) c_e(levused)) THEN
1551 .NOT.
IF (ANY(c_e( &
1552 volgrid6d_in%gaid(levshift+1:levshift+levused,itime,itimerange,ivar) &
1555 DO ilevel = 1, MIN(inlevel,onlevel)
1556 ! if present gaid copy it
1557 .AND..NOT.
IF (c_e(volgrid6d_in%gaid(ilevel,itime,itimerange,ivar)) &
1558 c_e(volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))) THEN
1560 IF (optio_log(clone)) THEN
1561 CALL copy(volgrid6d_in%gaid(ilevel,itime,itimerange,ivar),&
1562 volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))
1564 CALL l4f_category_log(volgrid6d_in%category,L4F_DEBUG, &
1565 "cloning gaid, level
"//t2c(ilevel))
1568 volgrid6d_out%gaid(ilevel,itime,itimerange,ivar) = &
1569 volgrid6d_in%gaid(ilevel,itime,itimerange,ivar)
1573 ! if out level are more, we have to clone anyway
1574 DO ilevel = MIN(inlevel,onlevel) + 1, onlevel
1575 .AND..NOT.
IF (c_e(volgrid6d_in%gaid(inlevel,itime,itimerange,ivar)) &
1576 c_e(volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))) then
1578 CALL copy(volgrid6d_in%gaid(inlevel,itime,itimerange,ivar),&
1579 volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))
1581 CALL l4f_category_log(volgrid6d_in%category,L4F_DEBUG, &
1582 "forced cloning gaid, level
"//t2c(inlevel)//"->
"//t2c(ilevel))
1587 IF (c_e(lvar_coord_vol)) THEN
1588 NULLIFY(coord_3d_in)
1589 CALL volgrid_get_vol_3d(volgrid6d_in, itime, itimerange, lvar_coord_vol, &
1591 IF (c_e(spos)) THEN ! compute difference wrt surface coordinate
1592 IF (spos == 0) THEN ! error condition, set all to missing and goodnight
1593 coord_3d_in(:,:,levshift+1:levshift+levused) = rmiss
1595 DO ilevel = levshift+1, levshift+levused
1596 .AND.
WHERE(c_e(coord_3d_in(:,:,ilevel)) c_e(coord_3d_in(:,:,spos)))
1597 coord_3d_in(:,:,ilevel) = coord_3d_in(:,:,ilevel) - &
1598 coord_3d_in(:,:,spos)
1600 coord_3d_in(:,:,ilevel) = rmiss
1606 CALL volgrid_get_vol_3d(volgrid6d_in, itime, itimerange, ivar, &
1608 IF (ASSOCIATED(volgrid6d_out%voldati)) & ! improve!!!!
1609 CALL volgrid_get_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
1611 IF (c_e(lvar_coord_vol)) THEN
1612 CALL compute(this, voldatiin, voldatiout, convert(volgrid6d_in%var(ivar)), &
1613 coord_3d_in(:,:,levshift+1:levshift+levused)) ! subset coord_3d_in
1615 CALL compute(this, voldatiin, voldatiout, convert(volgrid6d_in%var(ivar)))
1617 CALL volgrid_set_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
1623 IF (c_e(lvar_coord_vol)) THEN
1624 DEALLOCATE(coord_3d_in)
1626 .NOT.
IF (ASSOCIATED(volgrid6d_in%voldati)) THEN
1627 DEALLOCATE(voldatiin)
1629 .NOT.
IF (ASSOCIATED(volgrid6d_out%voldati)) THEN
1630 DEALLOCATE(voldatiout)
1634 END SUBROUTINE volgrid6d_transform_compute
1637 !> Performs the specified abstract transformation on the data provided.
1638 !! The abstract transformation is specified by \a this parameter; the
1639 !! corresponding specifical transformation (\a grid_transform object)
1640 !! is created and destroyed internally. The output transformed object
1641 !! is created internally and it does not require preliminary
1643 SUBROUTINE volgrid6d_transform(this, griddim, volgrid6d_in, volgrid6d_out, &
1644 lev_out, volgrid6d_coord_in, maskgrid, maskbounds, clone, decode, categoryappend)
1645 TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
1646 TYPE(griddim_def),INTENT(in),OPTIONAL :: griddim !< griddim specifying the output grid (required by most transformation types)
1647 TYPE(volgrid6d),INTENT(inout) :: volgrid6d_in !< object to be transformed, it is not modified, despite the INTENT(inout)
1648 TYPE(volgrid6d),INTENT(out) :: volgrid6d_out !< transformed object, it does not require initialisation
1649 TYPE(vol7d_level),INTENT(in),OPTIONAL,TARGET :: lev_out(:) !< vol7d_level object defining target vertical grid, for vertical interpolations
1650 TYPE(volgrid6d),INTENT(in),OPTIONAL :: volgrid6d_coord_in !< object providing time constant input vertical coordinate for some kind of vertical interpolations
1651 REAL,INTENT(in),OPTIONAL :: maskgrid(:,:) !< 2D field to be used for defining subareas according to its values, it must have the same shape as the field to be interpolated (for transformation subtype 'maskfill')
1652 REAL,INTENT(in),OPTIONAL :: maskbounds(:) !< array of boundary values for defining a subset of valid points where the values of \a maskgrid are within the first and last value of \a maskbounds (for transformation type 'metamorphosis:maskfill')
1653 LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \a .TRUE. , clone the \a gaid's from \a volgrid6d_in to \a volgrid6d_out
1654 LOGICAL,INTENT(in),OPTIONAL :: decode !< determine whether the data in \a volgrid6d_out should be decoded or remain coded in gaid, if not provided, the decode status is taken from \a volgrid6d_in
1655 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
1657 TYPE(grid_transform) :: grid_trans
1658 TYPE(vol7d_level),POINTER :: llev_out(:)
1659 TYPE(vol7d_level) :: input_levtype, output_levtype
1660 TYPE(vol7d_var) :: vcoord_var
1661 INTEGER :: i, k, ntime, ntimerange, nlevel, nvar, var_coord_in, var_coord_vol, &
1662 cf_out, nxc, nyc, nxi, nyi, i3, i4, i5, i6, &
1663 ulstart, ulend, spos
1664 REAL,ALLOCATABLE :: coord_3d_in(:,:,:)
1665 TYPE(geo_proj) :: proj_in, proj_out
1666 CHARACTER(len=80) :: trans_type
1668 LOGICAL,ALLOCATABLE :: mask_in(:)
1671 call l4f_category_log(volgrid6d_in%category, L4F_DEBUG, "start volgrid6d_transform
")
1679 if (associated(volgrid6d_in%time)) ntime=size(volgrid6d_in%time)
1680 if (associated(volgrid6d_in%timerange)) ntimerange=size(volgrid6d_in%timerange)
1681 if (associated(volgrid6d_in%level)) nlevel=size(volgrid6d_in%level)
1682 if (associated(volgrid6d_in%var)) nvar=size(volgrid6d_in%var)
1684 .OR..OR..OR.
IF (ntime == 0 ntimerange == 0 nlevel == 0 nvar == 0) THEN
1685 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1687 ' ntimerange='//t2c(ntimerange)//' nlevel='//t2c(nlevel)//' nvar='//t2c(nvar))
1688 CALL init(volgrid6d_out) ! initialize to empty
1693 CALL get_val(this, trans_type=trans_type)
1695 ! store desired output component flag and unrotate if necessary
1697 .AND..OR.
IF (PRESENT(griddim) (trans_type == 'inter' trans_type == 'boxinter' &
1698 .OR.
trans_type == 'stencilinter')) THEN ! improve condition!!
1699 CALL get_val(volgrid6d_in%griddim, proj=proj_in)
1700 CALL get_val(griddim, component_flag=cf_out, proj=proj_out)
1701 ! if different projections wind components must be referred to geographical system
1702 IF (proj_in /= proj_out) CALL vg6d_wind_unrot(volgrid6d_in)
1703 ELSE IF (PRESENT(griddim)) THEN ! just get component_flag, the rest is rubbish
1704 CALL get_val(griddim, component_flag=cf_out)
1708 var_coord_in = imiss
1709 var_coord_vol = imiss
1710 IF (trans_type == 'vertint') THEN
1711 IF (PRESENT(lev_out)) THEN
1713 ! if volgrid6d_coord_in provided and allocated, check that it fits
1714 IF (PRESENT(volgrid6d_coord_in)) THEN
1715 IF (ASSOCIATED(volgrid6d_coord_in%voldati)) THEN
1717 ! strictly 1 time and 1 timerange
1718 .OR.
IF (SIZE(volgrid6d_coord_in%voldati,4) /= 1 &
1719 SIZE(volgrid6d_coord_in%voldati,5) /= 1) THEN
1720 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1721 'volume providing constant input vertical coordinate must have &
1722 &only 1 time and 1 timerange')
1723 CALL init(volgrid6d_out) ! initialize to empty
1728 ! search for variable providing vertical coordinate
1729 CALL get_val(this, output_levtype=output_levtype)
1730 vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
1731 .NOT.
IF (c_e(vcoord_var)) THEN
1732 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1733 'requested output level type '//t2c(output_levtype%level1)// &
1734 ' does not correspond to any known physical variable for &
1735 &providing vertical coordinate')
1736 CALL init(volgrid6d_out) ! initialize to empty
1741 DO i = 1, SIZE(volgrid6d_coord_in%var)
1742 IF (convert(volgrid6d_coord_in%var(i)) == vcoord_var) THEN
1748 .NOT.
IF (c_e(var_coord_in)) THEN
1749 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1750 'volume providing constant input vertical coordinate contains no &
1751 &variables matching output level type '//t2c(output_levtype%level1))
1752 CALL init(volgrid6d_out) ! initialize to empty
1756 CALL l4f_category_log(volgrid6d_in%category, L4F_INFO, &
1757 'Coordinate for vertint found in coord volume at position '// &
1760 ! check horizontal grid
1761 ! this is too strict (component flag and so on)
1762 ! IF (volgrid6d_coord_in%griddim /= volgrid6d_in%griddim) THEN
1763 CALL get_val(volgrid6d_coord_in%griddim, nx=nxc, ny=nyc)
1764 CALL get_val(volgrid6d_in%griddim, nx=nxi, ny=nyi)
1765 .OR.
IF (nxc /= nxi nyc /= nyi) THEN
1766 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1767 'volume providing constant input vertical coordinate must have &
1768 &the same grid as the input')
1769 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1770 'vertical coordinate: '//t2c(nxc)//'x'//t2c(nyc)// &
1771 ', input volume: '//t2c(nxi)//'x'//t2c(nyi))
1772 CALL init(volgrid6d_out) ! initialize to empty
1777 ! check vertical coordinate system
1778 CALL get_val(this, input_levtype=input_levtype)
1779 mask_in = & ! implicit allocation
1780 .AND.
(volgrid6d_coord_in%level(:)%level1 == input_levtype%level1) &
1781 (volgrid6d_coord_in%level(:)%level2 == input_levtype%level2)
1782 ulstart = firsttrue(mask_in)
1783 ulend = lasttrue(mask_in)
1784 .OR.
IF (ulstart == 0 ulend == 0) THEN
1785 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1786 'coordinate file does not contain levels of type '// &
1787 t2c(input_levtype%level1)//'/'//t2c(input_levtype%level2)// &
1788 ' specified for input data')
1789 CALL init(volgrid6d_out) ! initialize to empty
1794 coord_3d_in = volgrid6d_coord_in%voldati(:,:,ulstart:ulend,1,1,var_coord_in) ! implicit allocation
1796 .OR.
IF (output_levtype%level1 == 103 &
1797 output_levtype%level1 == 108) THEN ! surface coordinate needed
1798 spos = firsttrue(volgrid6d_coord_in%level(:) == vol7d_level_new(1))
1800 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1801 'output level '//t2c(output_levtype%level1)// &
1802 ' requested, but height/press of surface not provided in coordinate file')
1803 CALL init(volgrid6d_out) ! initialize to empty
1807 DO k = 1, SIZE(coord_3d_in,3)
1808 .AND.
WHERE(c_e(coord_3d_in(:,:,k)) &
1809 c_e(volgrid6d_coord_in%voldati(:,:,spos,1,1,var_coord_in)))
1810 coord_3d_in(:,:,k) = coord_3d_in(:,:,k) - &
1811 volgrid6d_coord_in%voldati(:,:,spos,1,1,var_coord_in)
1813 coord_3d_in(:,:,k) = rmiss
1821 .NOT.
IF (c_e(var_coord_in)) THEN ! search for coordinate within volume
1822 ! search for variable providing vertical coordinate
1823 CALL get_val(this, output_levtype=output_levtype)
1824 vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
1825 IF (c_e(vcoord_var)) THEN
1826 DO i = 1, SIZE(volgrid6d_in%var)
1827 IF (convert(volgrid6d_in%var(i)) == vcoord_var) THEN
1833 IF (c_e(var_coord_vol)) THEN
1834 CALL l4f_category_log(volgrid6d_in%category, L4F_INFO, &
1835 'Coordinate for vertint found in input volume at position '// &
1842 CALL init(volgrid6d_out, griddim=volgrid6d_in%griddim, &
1843 time_definition=volgrid6d_in%time_definition, categoryappend=categoryappend)
1844 IF (c_e(var_coord_in)) THEN
1845 CALL init(grid_trans, this, lev_in=volgrid6d_in%level, lev_out=lev_out, &
1846 coord_3d_in=coord_3d_in, categoryappend=categoryappend)
1848 CALL init(grid_trans, this, lev_in=volgrid6d_in%level, lev_out=lev_out, &
1849 categoryappend=categoryappend)
1852 CALL get_val(grid_trans, output_level_auto=llev_out) ! get levels if auto-generated
1853 .NOT.
IF (ASSOCIATED(llev_out)) llev_out => lev_out
1854 nlevel = SIZE(llev_out)
1856 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1857 'volgrid6d_transform: vertint requested but lev_out not provided')
1858 CALL init(volgrid6d_out) ! initialize to empty
1864 CALL init(volgrid6d_out, griddim=griddim, &
1865 time_definition=volgrid6d_in%time_definition, categoryappend=categoryappend)
1866 CALL init(grid_trans, this, in=volgrid6d_in%griddim, out=volgrid6d_out%griddim, &
1867 maskgrid=maskgrid, maskbounds=maskbounds, categoryappend=categoryappend)
1871 IF (c_e(grid_trans)) THEN ! transformation is valid
1873 CALL volgrid6d_alloc(volgrid6d_out, ntime=ntime, nlevel=nlevel, &
1874 ntimerange=ntimerange, nvar=nvar)
1876 IF (PRESENT(decode)) THEN ! explicitly set decode status
1878 ELSE ! take it from input
1879 ldecode = ASSOCIATED(volgrid6d_in%voldati)
1881 ! force decode if gaid is readonly
1882 decode_loop: DO i6 = 1,nvar
1883 DO i5 = 1, ntimerange
1886 IF (c_e(volgrid6d_in%gaid(i3,i4,i5,i6))) THEN
1887 .OR.
ldecode = ldecode grid_id_readonly(volgrid6d_in%gaid(i3,i4,i5,i6))
1895 IF (PRESENT(decode)) THEN
1896 .NEQV.
IF (ldecodedecode) THEN
1897 CALL l4f_category_log(volgrid6d_in%category, L4F_WARN, &
1898 'volgrid6d_transform: decode status forced to .TRUE. because driver does not allow copy')
1902 CALL volgrid6d_alloc_vol(volgrid6d_out, decode=ldecode)
1904 !ensure unproj was called
1905 !call griddim_unproj(volgrid6d_out%griddim)
1907 IF (trans_type == 'vertint') THEN
1909 CALL l4f_category_log(volgrid6d_in%category, L4F_DEBUG, &
1910 "volgrid6d_transform: vertint to
"//t2c(nlevel)//" levels
")
1912 CALL compute(grid_trans, volgrid6d_in, volgrid6d_out, lev_out=llev_out, &
1913 var_coord_vol=var_coord_vol, clone=clone)
1915 CALL compute(grid_trans, volgrid6d_in, volgrid6d_out, clone=clone)
1918 IF (cf_out == 0) THEN ! unrotated components are desired
1919 CALL wind_unrot(volgrid6d_out) ! unrotate if necessary
1920 ELSE IF (cf_out == 1) THEN ! rotated components are desired
1921 CALL wind_rot(volgrid6d_out) ! rotate if necessary
1925 ! should log with grid_trans%category, but it is private
1926 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1927 'volgrid6d_transform: transformation not valid')
1931 CALL delete (grid_trans)
1933 END SUBROUTINE volgrid6d_transform
1936 !> Performs the specified abstract transformation on the arrays of
1937 !! data provided. The abstract transformation is specified by \a this
1938 !! parameter; the corresponding specifical transformation (\a
1939 !! grid_transform object) is created and destroyed internally. The
1940 !! output transformed object is created internally and it does not
1941 !! require preliminary initialisation. According to the input data and
1942 !! to the transformation type, the output array may have of one or
1943 !! more \a volgrid6d elements on different grids.
1944 SUBROUTINE volgrid6dv_transform(this, griddim, volgrid6d_in, volgrid6d_out, &
1945 lev_out, volgrid6d_coord_in, maskgrid, maskbounds, clone, decode, categoryappend)
1946 TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
1947 TYPE(griddim_def),INTENT(in),OPTIONAL :: griddim !< griddim specifying the output grid (required by most transformation types)
1948 TYPE(volgrid6d),INTENT(inout) :: volgrid6d_in(:) !< object to be transformed, it is an array of volgrid6d objects, each of which will be transformed, it is not modified, despite the INTENT(inout)
1949 TYPE(volgrid6d),POINTER :: volgrid6d_out(:) !< transformed object, it is a non associated pointer to an array of volgrid6d objects which will be allocated by the method
1950 TYPE(vol7d_level),INTENT(in),OPTIONAL :: lev_out(:) !< vol7d_level object defining target vertical grid
1951 TYPE(volgrid6d),INTENT(in),OPTIONAL :: volgrid6d_coord_in !< object providing time constant input vertical coordinate for some kind of vertical interpolations
1952 REAL,INTENT(in),OPTIONAL :: maskgrid(:,:) !< 2D field to be used for defining subareas according to its values, it must have the same shape as the field to be interpolated (for transformation subtype 'maskfill')
1953 REAL,INTENT(in),OPTIONAL :: maskbounds(:) !< array of boundary values for defining a subset of valid points where the values of \a maskgrid are within the first and last value of \a maskbounds (for transformation type 'metamorphosis:maskfill')
1954 LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \a .TRUE. , clone the \a gaid's from \a volgrid6d_in to \a volgrid6d_out
1955 LOGICAL,INTENT(in),OPTIONAL :: decode !< if provided and \a .FALSE. the data volume is not allocated, but work is performed on grid_id's
1956 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
1958 INTEGER :: i, stallo
1961 allocate(volgrid6d_out(size(volgrid6d_in)),stat=stallo)
1962 if (stallo /= 0)then
1963 call l4f_log(L4F_FATAL,"allocating memory
")
1964 call raise_fatal_error()
1967 do i=1,size(volgrid6d_in)
1968 call transform(this, griddim, volgrid6d_in(i), volgrid6d_out(i), &
1969 lev_out=lev_out, volgrid6d_coord_in=volgrid6d_coord_in, &
1970 maskgrid=maskgrid, maskbounds=maskbounds, &
1971 clone=clone, decode=decode, categoryappend=categoryappend)
1974 END SUBROUTINE volgrid6dv_transform
1977 ! Internal method for performing grid to sparse point computations
1978 SUBROUTINE volgrid6d_v7d_transform_compute(this, volgrid6d_in, vol7d_out, &
1979 networkname, noconvert)
1980 TYPE(grid_transform),INTENT(in) :: this ! oggetto di trasformazione per grigliato
1981 type(volgrid6d), INTENT(in) :: volgrid6d_in ! oggetto da trasformare
1982 type(vol7d), INTENT(inout) :: vol7d_out ! oggetto trasformato
1983 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname ! imposta il network in vol7d_out (default='generic')
1984 LOGICAL,OPTIONAL,INTENT(in) :: noconvert !< do not try to match variable and convert values during transform
1986 INTEGER :: nntime, nana, ntime, ntimerange, nlevel, nvar, stallo
1987 INTEGER :: itime, itimerange, ivar, inetwork
1988 REAL,ALLOCATABLE :: voldatir_out(:,:,:)
1989 TYPE(conv_func),POINTER :: c_func(:)
1990 TYPE(datetime),ALLOCATABLE :: validitytime(:,:)
1991 REAL,POINTER :: voldatiin(:,:,:)
1994 call l4f_category_log(volgrid6d_in%category,L4F_DEBUG,"start volgrid6d_v7d_transform_compute
")
2003 if (present(networkname))then
2004 call init(vol7d_out%network(1),name=networkname)
2006 call init(vol7d_out%network(1),name='generic')
2009 if (associated(volgrid6d_in%timerange))then
2010 ntimerange=size(volgrid6d_in%timerange)
2011 vol7d_out%timerange=volgrid6d_in%timerange
2014 if (associated(volgrid6d_in%time))then
2015 ntime=size(volgrid6d_in%time)
2017 if (vol7d_out%time_definition == volgrid6d_in%time_definition) then
2019 ! i time sono definiti uguali: assegno
2020 vol7d_out%time=volgrid6d_in%time
2023 ! converto reference in validity
2024 allocate (validitytime(ntime,ntimerange),stat=stallo)
2026 call l4f_category_log(volgrid6d_in%category,L4F_FATAL,"allocating memory
")
2027 call raise_fatal_error()
2031 do itimerange=1,ntimerange
2032 if (vol7d_out%time_definition > volgrid6d_in%time_definition) then
2033 validitytime(itime,itimerange) = &
2034 volgrid6d_in%time(itime) + timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
2036 validitytime(itime,itimerange) = &
2037 volgrid6d_in%time(itime) - timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
2042 nntime = count_distinct(reshape(validitytime,(/ntime*ntimerange/)), back=.TRUE.)
2043 vol7d_out%time=pack_distinct(reshape(validitytime,(/ntime*ntimerange/)), nntime,back=.TRUE.)
2048 IF (ASSOCIATED(volgrid6d_in%level)) THEN
2049 nlevel = SIZE(volgrid6d_in%level)
2050 vol7d_out%level=volgrid6d_in%level
2053 IF (ASSOCIATED(volgrid6d_in%var)) THEN
2054 nvar = SIZE(volgrid6d_in%var)
2055 .NOT.
IF ( optio_log(noconvert)) THEN
2056 CALL vargrib2varbufr(volgrid6d_in%var, vol7d_out%dativar%r, c_func)
2060 nana = SIZE(vol7d_out%ana)
2062 ! allocate once for speed
2063 .NOT.
IF (ASSOCIATED(volgrid6d_in%voldati)) THEN
2064 ALLOCATE(voldatiin(volgrid6d_in%griddim%dim%nx, volgrid6d_in%griddim%dim%ny, &
2068 ALLOCATE(voldatir_out(nana,1,nlevel),stat=stallo)
2069 IF (stallo /= 0) THEN
2070 CALL l4f_category_log(volgrid6d_in%category,L4F_FATAL,"allocating memory
")
2071 CALL raise_fatal_error()
2076 do itimerange=1,ntimerange
2077 ! do ilevel=1,nlevel
2080 èè
!non chiaro se questa sezione utile o no
2081 !ossia il compute sotto sembra prevedere voldatir_out solo in out
2082 !!$ if (vol7d_out%time_definition == volgrid6d_in%time_definition) then
2083 !!$ voldatir_out=reshape(vol7d_out%voldatir(:,itime,ilevel,itimerange,ivar,inetwork),(/nana,1/))
2085 !!$ voldatir_out=reshape(vol7d_out%voldatir(:,index(vol7d_out%time,validitytime(itime,ilevel)),ilevel,itimerange,ivar,inetwork),(/nana,1/))
2088 CALL volgrid_get_vol_3d(volgrid6d_in, itime, itimerange, ivar, &
2091 CALL compute(this, voldatiin, voldatir_out, vol7d_out%dativar%r(ivar))
2093 if (vol7d_out%time_definition == volgrid6d_in%time_definition) then
2094 vol7d_out%voldatir(:,itime,:,itimerange,ivar,inetwork) = &
2097 vol7d_out%voldatir(:,index(vol7d_out%time,validitytime(itime,itimerange)),:,itimerange,ivar,inetwork)=&
2098 RESHAPE(voldatir_out,(/nana,nlevel/))
2101 ! 1 indice della dimensione "anagrafica
"
2102 ! 2 indice della dimensione "tempo
"
2103 ! 3 indice della dimensione "livello verticale
"
2104 ! 4 indice della dimensione "intervallo temporale
"
2105 ! 5 indice della dimensione "variabile
"
2106 ! 6 indice della dimensione "rete
"
2113 deallocate(voldatir_out)
2114 .NOT.
IF (ASSOCIATED(volgrid6d_in%voldati)) THEN
2115 DEALLOCATE(voldatiin)
2117 if (allocated(validitytime)) deallocate(validitytime)
2119 ! Rescale valid data according to variable conversion table
2120 IF (ASSOCIATED(c_func)) THEN
2122 CALL compute(c_func(ivar), vol7d_out%voldatir(:,:,:,:,ivar,:))
2127 end SUBROUTINE volgrid6d_v7d_transform_compute
2130 !> Performs the specified abstract transformation on the data provided.
2131 !! The abstract transformation is specified by \a this parameter; the
2132 !! corresponding specifical transformation (\a grid_transform object)
2133 !! is created and destroyed internally. The output transformed object
2134 !! is created internally and it does not require preliminary
2136 SUBROUTINE volgrid6d_v7d_transform(this, volgrid6d_in, vol7d_out, v7d, &
2137 maskgrid, maskbounds, networkname, noconvert, find_index, categoryappend)
2138 TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
2139 TYPE(volgrid6d),INTENT(inout) :: volgrid6d_in !< object to be transformed, it is not modified, despite the INTENT(inout)
2140 TYPE(vol7d),INTENT(out) :: vol7d_out !< transformed object, it does not requires initialisation
2141 TYPE(vol7d),INTENT(in),OPTIONAL :: v7d !< object containing a list of points over which transformation has to be done (required by some transformation types)
2142 REAL,INTENT(in),OPTIONAL :: maskgrid(:,:) !< 2D field to be used for defining subareas according to its values, it must have the same shape as the field to be interpolated (for transformation type 'maskinter')
2143 REAL,INTENT(in),OPTIONAL :: maskbounds(:) !< array of boundary values for defining subareas from the values of \a maskgrid, the number of subareas is SIZE(maskbounds) - 1, if not provided a default based on extreme values of \a makgrid is used
2144 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname !< set the output network name in vol7d_out (default='generic')
2145 LOGICAL,OPTIONAL,INTENT(in) :: noconvert !< do not try to match variable and convert values during transform
2146 PROCEDURE(basic_find_index),POINTER,OPTIONAL :: find_index
2147 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
2149 type(grid_transform) :: grid_trans
2150 INTEGER :: ntime, ntimerange, nlevel, nvar, nana, time_definition, nnetwork, stallo
2151 INTEGER :: itime, itimerange, inetwork
2152 TYPE(datetime),ALLOCATABLE :: validitytime(:,:)
2153 INTEGER,ALLOCATABLE :: point_index(:)
2154 TYPE(vol7d) :: v7d_locana
2157 call l4f_category_log(volgrid6d_in%category,L4F_DEBUG,"start volgrid6d_v7d_transform
")
2160 call vg6d_wind_unrot(volgrid6d_in)
2168 call get_val(this,time_definition=time_definition)
2169 .not.
if ( c_e(time_definition)) then
2170 time_definition=1 ! default to validity time
2173 IF (PRESENT(v7d)) THEN
2174 CALL vol7d_copy(v7d, v7d_locana)
2176 CALL init(v7d_locana)
2179 if (associated(volgrid6d_in%timerange)) ntimerange=size(volgrid6d_in%timerange)
2181 if (associated(volgrid6d_in%time)) then
2183 ntime=size(volgrid6d_in%time)
2185 if (time_definition /= volgrid6d_in%time_definition) then
2187 ! converto reference in validity
2188 allocate (validitytime(ntime,ntimerange),stat=stallo)
2190 call l4f_category_log(volgrid6d_in%category,L4F_FATAL,"allocating memory
")
2191 call raise_fatal_error()
2195 do itimerange=1,ntimerange
2196 if (time_definition > volgrid6d_in%time_definition) then
2197 validitytime(itime,itimerange) = &
2198 volgrid6d_in%time(itime) + timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
2200 validitytime(itime,itimerange) = &
2201 volgrid6d_in%time(itime) - timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
2206 ntime = count_distinct(reshape(validitytime,(/ntime*ntimerange/)), back=.TRUE.)
2207 deallocate (validitytime)
2213 if (associated(volgrid6d_in%level)) nlevel=size(volgrid6d_in%level)
2214 if (associated(volgrid6d_in%var)) nvar=size(volgrid6d_in%var)
2216 CALL init(grid_trans, this, volgrid6d_in%griddim, v7d_locana, &
2217 maskgrid=maskgrid, maskbounds=maskbounds, find_index=find_index, &
2218 categoryappend=categoryappend)
2219 CALL init (vol7d_out,time_definition=time_definition)
2221 IF (c_e(grid_trans)) THEN
2223 nana=SIZE(v7d_locana%ana)
2224 CALL vol7d_alloc(vol7d_out, nana=nana, ntime=ntime, nlevel=nlevel, &
2225 ntimerange=ntimerange, ndativarr=nvar, nnetwork=nnetwork)
2226 vol7d_out%ana = v7d_locana%ana
2228 CALL get_val(grid_trans, output_point_index=point_index)
2229 IF (ALLOCATED(point_index)) THEN
2230 ! check that size(point_index) == nana?
2231 CALL vol7d_alloc(vol7d_out, nanavari=1)
2232 CALL init(vol7d_out%anavar%i(1), 'B01192')
2235 CALL vol7d_alloc_vol(vol7d_out)
2237 IF (ALLOCATED(point_index)) THEN
2238 DO inetwork = 1, nnetwork
2239 vol7d_out%volanai(:,1,inetwork) = point_index(:)
2242 CALL compute(grid_trans, volgrid6d_in, vol7d_out, networkname, noconvert)
2244 CALL l4f_log(L4F_ERROR, 'vg6d_v7d_transform: transformation not valid')
2248 CALL delete(grid_trans)
2251 ! set variables to a conformal state
2252 CALL vol7d_dballe_set_var_du(vol7d_out)
2255 CALL delete(v7d_locana)
2257 END SUBROUTINE volgrid6d_v7d_transform
2260 !> Performs the specified abstract transformation on the arrays of
2261 !! data provided. The abstract transformation is specified by \a this
2262 !! parameter; the corresponding specifical transformation (\a
2263 !! grid_transform object) is created and destroyed internally. The
2264 !! output transformed object is created internally and it does not
2265 !! require preliminary initialisation. The transformation performed on
2266 !! each element of the input \a volgrid6d array object is merged into
2267 !! a single \a vol7d output object.
2268 SUBROUTINE volgrid6dv_v7d_transform(this, volgrid6d_in, vol7d_out, v7d, &
2269 maskgrid, maskbounds, networkname, noconvert, find_index, categoryappend)
2270 TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
2271 TYPE(volgrid6d),INTENT(inout) :: volgrid6d_in(:) !< object to be transformed, it is an array of volgrid6d objects, each of which will be transformed, it is not modified, despite the INTENT(inout)
2272 TYPE(vol7d),INTENT(out) :: vol7d_out !< transformed object, it does not require initialisation
2273 TYPE(vol7d),INTENT(in),OPTIONAL :: v7d !< object containing a list of points over which transformation has to be done (required by some transformation types)
2274 REAL,INTENT(in),OPTIONAL :: maskgrid(:,:) !< 2D field to be used for defining subareas according to its values, it must have the same shape as the field to be interpolated (for transformation type 'maskinter')
2275 REAL,INTENT(in),OPTIONAL :: maskbounds(:) !< array of boundary values for defining subareas from the values of \a maskgrid, the number of subareas is SIZE(maskbounds) - 1, if not provided a default based on extreme values of \a makgrid is used
2276 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname !< set the output network name in vol7d_out (default='generic')
2277 LOGICAL,OPTIONAL,INTENT(in) :: noconvert !< do not try to match variable and convert values during transform
2278 PROCEDURE(basic_find_index),POINTER,OPTIONAL :: find_index
2279 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
2282 TYPE(vol7d) :: v7dtmp
2286 CALL init(vol7d_out)
2288 DO i=1,SIZE(volgrid6d_in)
2289 CALL transform(this, volgrid6d_in(i), v7dtmp, v7d=v7d, &
2290 maskgrid=maskgrid, maskbounds=maskbounds, &
2291 networkname=networkname, noconvert=noconvert, find_index=find_index, &
2292 categoryappend=categoryappend)
2293 CALL vol7d_append(vol7d_out, v7dtmp)
2296 END SUBROUTINE volgrid6dv_v7d_transform
2299 ! Internal method for performing sparse point to grid computations
2300 SUBROUTINE v7d_volgrid6d_transform_compute(this, vol7d_in, volgrid6d_out, networkname, gaid_template)
2301 TYPE(grid_transform),INTENT(in) :: this ! object specifying the specific transformation
2302 type(vol7d), INTENT(in) :: vol7d_in ! object to be transformed
2303 type(volgrid6d), INTENT(inout) :: volgrid6d_out ! transformed object
2304 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname ! select the network to be processed from the \a vol7d input object, default the first network
2305 TYPE(grid_id),OPTIONAL,INTENT(in) :: gaid_template ! the template (typically grib_api) to be associated with output data, it also helps in improving variable conversion
2307 integer :: nana, ntime, ntimerange, nlevel, nvar
2308 INTEGER :: ilevel, itime, itimerange, ivar, inetwork
2310 REAL,POINTER :: voldatiout(:,:,:)
2311 type(vol7d_network) :: network
2312 TYPE(conv_func), pointer :: c_func(:)
2313 !TODO category sarebbe da prendere da vol7d
2315 CALL l4f_category_log(volgrid6d_out%category, L4F_DEBUG, &
2316 'start v7d_volgrid6d_transform_compute')
2324 IF (PRESENT(networkname)) THEN
2325 CALL init(network,name=networkname)
2326 inetwork = INDEX(vol7d_in%network,network)
2327 IF (inetwork <= 0) THEN
2328 CALL l4f_category_log(volgrid6d_out%category,L4F_WARN, &
2329 'network '//TRIM(networkname)//' not found, first network will be transformed')
2336 ! no time_definition conversion implemented, output will be the same as input
2337 if (associated(vol7d_in%time))then
2338 ntime=size(vol7d_in%time)
2339 volgrid6d_out%time=vol7d_in%time
2342 if (associated(vol7d_in%timerange))then
2343 ntimerange=size(vol7d_in%timerange)
2344 volgrid6d_out%timerange=vol7d_in%timerange
2347 if (associated(vol7d_in%level))then
2348 nlevel=size(vol7d_in%level)
2349 volgrid6d_out%level=vol7d_in%level
2352 if (associated(vol7d_in%dativar%r))then
2353 nvar=size(vol7d_in%dativar%r)
2354 CALL varbufr2vargrib(vol7d_in%dativar%r, volgrid6d_out%var, c_func, gaid_template)
2357 nana=SIZE(vol7d_in%voldatir, 1)
2358 ! allocate once for speed
2359 .NOT.
IF (ASSOCIATED(volgrid6d_out%voldati)) THEN
2360 ALLOCATE(voldatiout(volgrid6d_out%griddim%dim%nx, volgrid6d_out%griddim%dim%ny, &
2365 DO itimerange=1,ntimerange
2368 ! clone the gaid template where I have data
2369 IF (PRESENT(gaid_template)) THEN
2370 DO ilevel = 1, nlevel
2371 IF (ANY(c_e(vol7d_in%voldatir(:,itime,ilevel,itimerange,ivar,inetwork)))) THEN
2372 CALL copy(gaid_template, volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))
2374 volgrid6d_out%gaid(ilevel,itime,itimerange,ivar) = grid_id_new()
2380 IF (ASSOCIATED(volgrid6d_out%voldati)) & ! improve!!!!
2381 CALL volgrid_get_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
2383 ! do the interpolation
2384 CALL compute(this, &
2385 vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), voldatiout, &
2386 vol7d_in%dativar%r(ivar))
2387 ! rescale valid data according to variable conversion table
2388 IF (ASSOCIATED(c_func)) THEN
2389 CALL compute(c_func(ivar), voldatiout(:,:,:))
2392 CALL volgrid_set_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
2399 .NOT.
IF (ASSOCIATED(volgrid6d_out%voldati)) THEN
2400 DEALLOCATE(voldatiout)
2402 IF (ASSOCIATED(c_func)) THEN
2406 END SUBROUTINE v7d_volgrid6d_transform_compute
2409 !> Performs the specified abstract transformation on the data provided.
2410 !! The abstract transformation is specified by \a this parameter; the
2411 !! corresponding specifical transformation (\a grid_transform object)
2412 !! is created and destroyed internally. The output transformed object
2413 !! is created internally and it does not require preliminary
2415 SUBROUTINE v7d_volgrid6d_transform(this, griddim, vol7d_in, volgrid6d_out, &
2416 networkname, gaid_template, categoryappend)
2417 TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
2418 TYPE(griddim_def),INTENT(in),OPTIONAL :: griddim !< griddim specifying the output grid (required by most transformation types)
2419 ! TODO ripristinare intent(in) dopo le opportune modifiche in grid_class.F90
2420 TYPE(vol7d),INTENT(inout) :: vol7d_in !< object to be transformed, it is not modified, despite the INTENT(inout)
2421 TYPE(volgrid6d),INTENT(out) :: volgrid6d_out !< transformed object, it does not require initialisation
2422 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname !< select the network to be processed from the \a vol7d input object, default the first network
2423 TYPE(grid_id),OPTIONAL,INTENT(in) :: gaid_template !< the template (typically grib_api) to be associated with output data, it also helps in improving variable conversion
2424 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
2426 type(grid_transform) :: grid_trans
2427 integer :: ntime, ntimerange, nlevel, nvar
2430 !TODO la category sarebbe da prendere da vol7d
2431 !call l4f_category_log(vol7d_out%category,L4F_DEBUG,"start volgrid6d_transform
")
2433 CALL vol7d_alloc_vol(vol7d_in) ! be safe
2434 ntime=SIZE(vol7d_in%time)
2435 ntimerange=SIZE(vol7d_in%timerange)
2436 nlevel=SIZE(vol7d_in%level)
2438 if (associated(vol7d_in%dativar%r)) nvar=size(vol7d_in%dativar%r)
2440 IF (nvar <= 0) THEN ! use vol7d category once it will be implemented
2441 CALL l4f_log(L4F_ERROR, &
2442 "trying to
transform a
vol7d object incomplete or without real variables
")
2443 CALL init(volgrid6d_out) ! initialize to empty
2448 CALL init(grid_trans, this, vol7d_in, griddim, categoryappend=categoryappend)
2449 CALL init(volgrid6d_out, griddim, time_definition=vol7d_in%time_definition, &
2450 categoryappend=categoryappend)
2452 IF (c_e(grid_trans)) THEN
2454 CALL volgrid6d_alloc(volgrid6d_out, griddim%dim, ntime=ntime, nlevel=nlevel, &
2455 ntimerange=ntimerange, nvar=nvar)
2456 ! can I avoid decode=.TRUE. here, using gaid_template?
2457 CALL volgrid6d_alloc_vol(volgrid6d_out, decode=.TRUE.)
2459 CALL compute(grid_trans, vol7d_in, volgrid6d_out, networkname, gaid_template)
2461 CALL vg6d_wind_rot(volgrid6d_out)
2463 ! should log with grid_trans%category, but it is private
2464 CALL l4f_log(L4F_ERROR, 'v7d_vg6d_transform: transformation not valid')
2468 CALL delete(grid_trans)
2470 END SUBROUTINE v7d_volgrid6d_transform
2473 ! Internal method for performing sparse point to sparse point computations
2474 SUBROUTINE v7d_v7d_transform_compute(this, vol7d_in, vol7d_out, lev_out, &
2476 TYPE(grid_transform),INTENT(in) :: this ! oggetto di trasformazione per grigliato
2477 type(vol7d), INTENT(in) :: vol7d_in ! oggetto da trasformare
2478 type(vol7d), INTENT(inout) :: vol7d_out ! oggetto trasformato
2479 TYPE(vol7d_level),INTENT(in),OPTIONAL :: lev_out(:) ! vol7d_level object defining target vertical grid, for vertical interpolations
2480 INTEGER,INTENT(in),OPTIONAL :: var_coord_vol ! index of variable defining vertical coordinate values in input volume
2482 INTEGER :: itime, itimerange, ilevel, ivar, inetwork, &
2483 levshift, levused, lvar_coord_vol, spos
2484 REAL,ALLOCATABLE :: coord_3d_in(:,:,:)
2485 TYPE(vol7d_level) :: output_levtype
2487 lvar_coord_vol = optio_i(var_coord_vol)
2488 vol7d_out%time(:) = vol7d_in%time(:)
2489 vol7d_out%timerange(:) = vol7d_in%timerange(:)
2490 IF (PRESENT(lev_out)) THEN
2491 vol7d_out%level(:) = lev_out(:)
2493 vol7d_out%level(:) = vol7d_in%level(:)
2495 vol7d_out%network(:) = vol7d_in%network(:)
2496 IF (ASSOCIATED(vol7d_in%dativar%r)) THEN ! work only when real vars are available
2497 vol7d_out%dativar%r(:) = vol7d_in%dativar%r(:)
2499 CALL get_val(this, levshift=levshift, levused=levused)
2501 IF (c_e(lvar_coord_vol)) THEN
2502 CALL get_val(this%trans, output_levtype=output_levtype)
2503 .OR.
IF (output_levtype%level1 == 103 output_levtype%level1 == 108) THEN
2504 spos = firsttrue(vol7d_in%level(:) == vol7d_level_new(1))
2506 CALL l4f_log(L4F_ERROR, &
2507 'output level '//t2c(output_levtype%level1)// &
2508 ' requested, but height/press of surface not provided in volume')
2510 .NOT..AND..NOT.
IF (c_e(levshift) c_e(levused)) THEN
2511 CALL l4f_log(L4F_ERROR, &
2512 'internal inconsistence, levshift and levused undefined when they should be')
2514 ALLOCATE(coord_3d_in(SIZE(vol7d_in%ana),1,SIZE(vol7d_in%level)))
2519 DO inetwork = 1, SIZE(vol7d_in%network)
2520 DO ivar = 1, SIZE(vol7d_in%dativar%r)
2521 DO itimerange = 1, SIZE(vol7d_in%timerange)
2522 DO itime = 1, SIZE(vol7d_in%time)
2524 ! dirty trick to make voldatir look like a 2d-array of shape (nana,1)
2525 IF (c_e(lvar_coord_vol)) THEN
2526 IF (c_e(spos)) THEN ! compute difference wrt surface coordinate
2527 IF (spos == 0) THEN ! error condition, set all to missing and goodnight
2528 coord_3d_in(:,:,levshift+1:levshift+levused) = rmiss
2530 DO ilevel = levshift+1, levshift+levused
2531 .AND.
WHERE(c_e(vol7d_in%voldatir(:,itime:itime,ilevel,itimerange,lvar_coord_vol,inetwork)) &
2532 c_e(vol7d_in%voldatir(:,itime:itime,spos,itimerange,lvar_coord_vol,inetwork)))
2533 coord_3d_in(:,:,ilevel) = vol7d_in%voldatir(:,itime:itime,ilevel,itimerange,lvar_coord_vol,inetwork) - &
2534 vol7d_in%voldatir(:,itime:itime,spos,itimerange,lvar_coord_vol,inetwork)
2536 coord_3d_in(:,:,ilevel) = rmiss
2540 CALL compute(this, &
2541 vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), &
2542 vol7d_out%voldatir(:,itime:itime,:,itimerange,ivar,inetwork), &
2543 var=vol7d_in%dativar%r(ivar), &
2544 coord_3d_in=coord_3d_in)
2546 CALL compute(this, &
2547 vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), &
2548 vol7d_out%voldatir(:,itime:itime,:,itimerange,ivar,inetwork), &
2549 var=vol7d_in%dativar%r(ivar), &
2550 coord_3d_in=vol7d_in%voldatir(:,itime:itime,:,itimerange, &
2551 lvar_coord_vol,inetwork))
2554 CALL compute(this, &
2555 vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), &
2556 vol7d_out%voldatir(:,itime:itime,:,itimerange,ivar,inetwork), &
2557 var=vol7d_in%dativar%r(ivar))
2566 END SUBROUTINE v7d_v7d_transform_compute
2569 !> Performs the specified abstract transformation on the data provided.
2570 !! The abstract transformation is specified by \a this parameter; the
2571 !! corresponding specifical transformation (\a grid_transform object)
2572 !! is created and destroyed internally. The output transformed object
2573 !! is created internally and it does not require preliminary
2574 !! initialisation. The success of the transformation can be checked
2575 !! with the \a c_e method: c_e(vol7d_out).
2576 SUBROUTINE v7d_v7d_transform(this, vol7d_in, vol7d_out, v7d, maskbounds, &
2577 lev_out, vol7d_coord_in, categoryappend)
2578 TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
2579 TYPE(vol7d),INTENT(inout) :: vol7d_in !< object to be transformed, it is not modified, despite the INTENT(inout)
2580 TYPE(vol7d),INTENT(out) :: vol7d_out !< transformed object, it does not require initialisation
2581 TYPE(vol7d),INTENT(in),OPTIONAL :: v7d !< object containing a list of points over which transformation has to be done (required by some transformation types)
2582 REAL,INTENT(in),OPTIONAL :: maskbounds(:) !< array of boundary values for defining a subset of valid points where the values of \a maskgrid are within the first and last value of \a maskbounds (for transformation type 'metamorphosis:maskfill')
2583 TYPE(vol7d_level),INTENT(in),OPTIONAL,TARGET :: lev_out(:) !< vol7d_level object defining target vertical grid, for vertical interpolations
2584 TYPE(vol7d),INTENT(in),OPTIONAL :: vol7d_coord_in !< object providing time constant input vertical coordinate for some kind of vertical interpolations
2585 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
2587 INTEGER :: nvar, inetwork
2588 TYPE(grid_transform) :: grid_trans
2589 TYPE(vol7d_level),POINTER :: llev_out(:)
2590 TYPE(vol7d_level) :: input_levtype, output_levtype
2591 TYPE(vol7d_var) :: vcoord_var
2592 REAL,ALLOCATABLE :: coord_3d_in(:,:,:)
2593 INTEGER :: var_coord_in, var_coord_vol, i, k, ulstart, ulend, spos
2594 INTEGER,ALLOCATABLE :: point_index(:)
2595 TYPE(vol7d) :: v7d_locana, vol7d_tmpana
2596 CHARACTER(len=80) :: trans_type
2597 LOGICAL,ALLOCATABLE :: mask_in(:), point_mask(:)
2599 CALL vol7d_alloc_vol(vol7d_in) ! be safe
2601 IF (ASSOCIATED(vol7d_in%dativar%r)) nvar=SIZE(vol7d_in%dativar%r)
2603 CALL init(v7d_locana)
2604 IF (PRESENT(v7d)) v7d_locana = v7d
2605 CALL init(vol7d_out, time_definition=vol7d_in%time_definition)
2607 CALL get_val(this, trans_type=trans_type)
2609 var_coord_vol = imiss
2610 IF (trans_type == 'vertint') THEN
2612 IF (PRESENT(lev_out)) THEN
2614 ! if vol7d_coord_in provided and allocated, check that it fits
2616 IF (PRESENT(vol7d_coord_in)) THEN
2617 .AND.
IF (ASSOCIATED(vol7d_coord_in%voldatir) &
2618 ASSOCIATED(vol7d_coord_in%dativar%r)) THEN
2620 ! strictly 1 time, 1 timerange and 1 network
2621 .OR.
IF (SIZE(vol7d_coord_in%voldatir,2) /= 1 &
2622 .OR.
SIZE(vol7d_coord_in%voldatir,4) /= 1 &
2623 SIZE(vol7d_coord_in%voldatir,6) /= 1) THEN
2624 CALL l4f_log(L4F_ERROR, &
2625 'volume providing constant input vertical coordinate must have &
2626 &only 1 time, 1 timerange and 1 network')
2631 ! search for variable providing vertical coordinate
2632 CALL get_val(this, output_levtype=output_levtype)
2633 vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
2634 .NOT.
IF (c_e(vcoord_var)) THEN
2635 CALL l4f_log(L4F_ERROR, &
2636 'requested output level type '//t2c(output_levtype%level1)// &
2637 ' does not correspond to any known physical variable for &
2638 &providing vertical coordinate')
2643 var_coord_in = INDEX(vol7d_coord_in%dativar%r, vcoord_var)
2645 IF (var_coord_in <= 0) THEN
2646 CALL l4f_log(L4F_ERROR, &
2647 'volume providing constant input vertical coordinate contains no &
2648 &real variables matching output level type '//t2c(output_levtype%level1))
2652 CALL l4f_log(L4F_INFO, &
2653 'Coordinate for vertint found in coord volume at position '// &
2656 ! check vertical coordinate system
2657 CALL get_val(this, input_levtype=input_levtype)
2658 mask_in = & ! implicit allocation
2659 .AND.
(vol7d_coord_in%level(:)%level1 == input_levtype%level1) &
2660 (vol7d_coord_in%level(:)%level2 == input_levtype%level2)
2661 ulstart = firsttrue(mask_in)
2662 ulend = lasttrue(mask_in)
2663 .OR.
IF (ulstart == 0 ulend == 0) THEN
2664 CALL l4f_log(L4F_ERROR, &
2665 'coordinate file does not contain levels of type '// &
2666 t2c(input_levtype%level1)//'/'//t2c(input_levtype%level2)// &
2667 ' specified for input data')
2672 coord_3d_in = vol7d_coord_in%voldatir(:,1:1,ulstart:ulend,1,var_coord_in,1) ! dirty 1:1, implicit allocation
2674 IF (output_levtype%level1 == 103 &
2675 .OR.
output_levtype%level1 == 108) THEN ! surface coordinate needed
2676 spos = firsttrue(vol7d_coord_in%level(:) == vol7d_level_new(1))
2678 CALL l4f_log(L4F_ERROR, &
2679 'output level '//t2c(output_levtype%level1)// &
2680 ' requested, but height/press of surface not provided in coordinate file')
2684 DO k = 1, SIZE(coord_3d_in,3)
2685 .AND.
WHERE(c_e(coord_3d_in(:,:,k)) &
2686 c_e(vol7d_coord_in%voldatir(:,1:1,spos,1,var_coord_in,1)))
2687 coord_3d_in(:,:,k) = coord_3d_in(:,:,k) - &
2688 vol7d_coord_in%voldatir(:,1:1,spos,1,var_coord_in,1)
2690 coord_3d_in(:,:,k) = rmiss
2698 IF (var_coord_in <= 0) THEN ! search for coordinate within volume
2699 ! search for variable providing vertical coordinate
2700 CALL get_val(this, output_levtype=output_levtype)
2701 vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
2702 IF (c_e(vcoord_var)) THEN
2703 DO i = 1, SIZE(vol7d_in%dativar%r)
2704 IF (vol7d_in%dativar%r(i) == vcoord_var) THEN
2710 IF (c_e(var_coord_vol)) THEN
2711 CALL l4f_log(L4F_INFO, &
2712 'Coordinate for vertint found in input volume at position '// &
2719 IF (var_coord_in > 0) THEN
2720 CALL init(grid_trans, this, lev_in=vol7d_in%level, lev_out=lev_out, &
2721 coord_3d_in=coord_3d_in, categoryappend=categoryappend)
2723 CALL init(grid_trans, this, lev_in=vol7d_in%level, lev_out=lev_out, &
2724 categoryappend=categoryappend)
2727 CALL get_val(grid_trans, output_level_auto=llev_out) ! get levels if auto-generated
2728 .NOT.
IF (associated(llev_out)) llev_out => lev_out
2730 .AND.
IF (c_e(grid_trans)) THEN! nvar > 0) THEN
2732 CALL vol7d_alloc(vol7d_out, nana=SIZE(vol7d_in%ana), &
2733 ntime=SIZE(vol7d_in%time), ntimerange=SIZE(vol7d_in%timerange), &
2734 nlevel=SIZE(llev_out), nnetwork=SIZE(vol7d_in%network), ndativarr=nvar)
2735 vol7d_out%ana(:) = vol7d_in%ana(:)
2737 CALL vol7d_alloc_vol(vol7d_out)
2739 ! no need to check c_e(var_coord_vol) here since the presence of
2740 ! this%coord_3d_in (external) has precedence over coord_3d_in internal
2741 ! in grid_transform_compute
2742 CALL compute(grid_trans, vol7d_in, vol7d_out, llev_out, &
2743 var_coord_vol=var_coord_vol)
2745 CALL l4f_log(L4F_ERROR, 'v7d_v7d_transform: transformation not valid')
2749 CALL l4f_log(L4F_ERROR, &
2750 'v7d_v7d_transform: vertint requested but lev_out not provided')
2756 CALL init(grid_trans, this, vol7d_in, v7d_locana, maskbounds=maskbounds, &
2757 categoryappend=categoryappend)
2758 ! if this init is successful, I am sure that v7d_locana%ana is associated
2760 .AND.
IF (c_e(grid_trans)) THEN! nvar > 0) THEN
2762 CALL vol7d_alloc(vol7d_out, nana=SIZE(v7d_locana%ana), &
2763 ntime=SIZE(vol7d_in%time), ntimerange=SIZE(vol7d_in%timerange), &
2764 nlevel=SIZE(vol7d_in%level), nnetwork=SIZE(vol7d_in%network), ndativarr=nvar)
2765 vol7d_out%ana = v7d_locana%ana
2767 CALL get_val(grid_trans, point_mask=point_mask, output_point_index=point_index)
2769 IF (ALLOCATED(point_index)) THEN
2770 CALL vol7d_alloc(vol7d_out, nanavari=1)
2771 CALL init(vol7d_out%anavar%i(1), 'B01192')
2774 CALL vol7d_alloc_vol(vol7d_out)
2776 IF (ALLOCATED(point_index)) THEN
2777 DO inetwork = 1, SIZE(vol7d_in%network)
2778 vol7d_out%volanai(:,1,inetwork) = point_index(:)
2781 CALL compute(grid_trans, vol7d_in, vol7d_out)
2783 IF (ALLOCATED(point_mask)) THEN ! keep full ana
2784 IF (SIZE(point_mask) /= SIZE(vol7d_in%ana)) THEN
2785 CALL l4f_log(L4F_WARN, &
2786 'v7d_v7d_transform: inconsistency in point size: '//t2c(SIZE(point_mask)) &
2787 //':'//t2c(SIZE(vol7d_in%ana)))
2790 CALL l4f_log(L4F_DEBUG, 'v7d_v7d_transform: merging ana from in to out')
2792 CALL vol7d_copy(vol7d_in, vol7d_tmpana, &
2793 lana=point_mask, lnetwork=(/.TRUE./), &
2794 ltime=(/.FALSE./), ltimerange=(/.FALSE./), llevel=(/.FALSE./))
2795 CALL vol7d_append(vol7d_out, vol7d_tmpana)
2800 CALL l4f_log(L4F_ERROR, 'v7d_v7d_transform: transformation not valid')
2806 CALL delete (grid_trans)
2807 .NOT.
IF ( PRESENT(v7d)) CALL delete(v7d_locana)
2809 END SUBROUTINE v7d_v7d_transform
2812 !> Unrotate the wind components.
2813 !! It converts u and v components of vector quantities relative to the
2814 !! defined grid in the direction of increasing x and y coordinates to
2815 !! u and v components relative to easterly and notherly direction. The
2816 !! original fields are overwritten.
2817 !! \todo Check and correct wind component flag (to be moved in
2819 subroutine vg6d_wind_unrot(this)
2820 type(volgrid6d) :: this !< object containing wind to be unrotated
2822 integer :: component_flag
2824 call get_val(this%griddim,component_flag=component_flag)
2826 if (component_flag == 1) then
2827 call l4f_category_log(this%category,L4F_INFO, &
2828 "unrotating vector components
")
2829 call vg6d_wind__un_rot(this,.false.)
2830 call set_val(this%griddim,component_flag=0)
2832 call l4f_category_log(this%category,L4F_INFO, &
2833 "no need to unrotate vector components
")
2836 end subroutine vg6d_wind_unrot
2839 !> Rotate the wind components.
2840 !! It converts u and v components of vector quantities
2841 !! relative to easterly and notherly direction to
2842 !! defined grid in the direction of increasing x and y coordinates.
2843 !! The original fields are overwritten.
2844 subroutine vg6d_wind_rot(this)
2845 type(volgrid6d) :: this !< object containing wind to be rotated
2847 integer :: component_flag
2849 call get_val(this%griddim,component_flag=component_flag)
2851 if (component_flag == 0) then
2852 call l4f_category_log(this%category,L4F_INFO, &
2853 "rotating vector components
")
2854 call vg6d_wind__un_rot(this,.true.)
2855 call set_val(this%griddim,component_flag=1)
2857 call l4f_category_log(this%category,L4F_INFO, &
2858 "no need to rotate vector components
")
2861 end subroutine vg6d_wind_rot
2864 ! Generic UnRotate the wind components.
2865 SUBROUTINE vg6d_wind__un_rot(this,rot)
2866 TYPE(volgrid6d) :: this ! object containing wind to be (un)rotated
2867 LOGICAL :: rot ! if .true. rotate else unrotate
2869 INTEGER :: i, j, k, l, a11, a12, a21, a22, stallo
2870 double precision,pointer :: rot_mat(:,:,:)
2871 real,allocatable :: tmp_arr(:,:)
2872 REAL,POINTER :: voldatiu(:,:), voldativ(:,:)
2873 INTEGER,POINTER :: iu(:), iv(:)
2875 .NOT.
IF (ASSOCIATED(this%var)) THEN
2876 CALL l4f_category_log(this%category, L4F_ERROR, &
2877 "trying to unrotate an incomplete
volgrid6d object
")
2878 CALL raise_fatal_error()
2882 CALL volgrid6d_var_hor_comp_index(this%var, iu, iv)
2883 .NOT.
IF (ASSOCIATED(iu)) THEN
2884 CALL l4f_category_log(this%category,L4F_ERROR, &
2885 "unrotation impossible
")
2886 CALL raise_fatal_error()
2890 ! Temporary workspace
2891 ALLOCATE(tmp_arr(this%griddim%dim%nx, this%griddim%dim%ny),stat=stallo)
2892 IF (stallo /= 0) THEN
2893 CALL l4f_category_log(this%category, L4F_FATAL, "allocating memory
")
2894 CALL raise_fatal_error()
2896 ! allocate once for speed
2897 .NOT.
IF (ASSOCIATED(this%voldati)) THEN
2898 ALLOCATE(voldatiu(this%griddim%dim%nx, this%griddim%dim%ny), &
2899 voldativ(this%griddim%dim%nx, this%griddim%dim%ny))
2902 CALL griddim_unproj(this%griddim)
2903 CALL wind_unrot(this%griddim, rot_mat)
2916 DO k = 1, SIZE(this%timerange)
2917 DO j = 1, SIZE(this%time)
2918 DO i = 1, SIZE(this%level)
2920 CALL volgrid_get_vol_2d(this, i, j, k, iu(l), voldatiu)
2921 CALL volgrid_get_vol_2d(this, i, j, k, iv(l), voldativ)
2922 ! convert units forward
2923 ! CALL compute(conv_fwd(iu(l)), voldatiu)
2924 ! CALL compute(conv_fwd(iv(l)), voldativ)
2926 ! multiply wind components by rotation matrix
2927 .AND.
WHERE(voldatiu /= rmiss voldativ /= rmiss)
2928 tmp_arr(:,:) = real(voldatiu(:,:)*rot_mat(:,:,a11) + &
2929 voldativ(:,:)*rot_mat(:,:,a12))
2930 voldativ(:,:) = real(voldatiu(:,:)*rot_mat(:,:,a21) + &
2931 voldativ(:,:)*rot_mat(:,:,a22))
2932 voldatiu(:,:) = tmp_arr(:,:)
2934 ! convert units backward
2935 ! CALL uncompute(conv_fwd(iu(l)), voldatiu)
2936 ! CALL uncompute(conv_fwd(iv(l)), voldativ)
2938 CALL volgrid_set_vol_2d(this, i, j, k, iu(l), voldatiu)
2939 CALL volgrid_set_vol_2d(this, i, j, k, iv(l), voldativ)
2945 .NOT.
IF (ASSOCIATED(this%voldati)) THEN
2946 DEALLOCATE(voldatiu, voldativ)
2948 DEALLOCATE(rot_mat, tmp_arr, iu, iv)
2950 END SUBROUTINE vg6d_wind__un_rot
2953 !!$ try to understand the problem:
2957 !!$ 1) we have only one volume: we have to provide the direction of shift
2958 !!$ compute H and traslate on it
2959 !!$ 2) we have two volumes:
2960 !!$ 1) volume U and volume V: compute H and traslate on it
2961 !!$ 2) volume U/V and volume H : translate U/V on H
2962 !!$ 3) we have tree volumes: translate U and V on H
2965 !!$ 1) do not have U in volume U
2966 !!$ 2) do not have V in volume V
2967 !!$ 3) we have others variables more than U and V in volumes U e V
2970 !!$ so the steps are:
2971 !!$ 1) find the volumes
2972 !!$ 2) define or compute H grid
2973 !!$ 3) trasform the volumes in H
2976 !!$ case 1) for only one vg6d (U or V) is not managed, but
2977 !!$ the not pubblic subroutines will work but you have to know what you want to do
2980 !> Convert grids type C to type A.
2981 !! Use this to interpolate data from grid type C to grid type A
2982 !! Grids type are defined by Arakawa.
2984 !! We need to find these types of area in a vg6d array
2985 !! T area of points with temterature etc.
2986 !! U area of points with u components of winds
2987 !! V area of points with v components of winds
2989 !! this method works if it finds
2991 !! 1) volume U and volume V: compute H and traslate on it
2992 !! 2) volume U/V and volume H : translate U/V on H
2993 !! three volumes: translate U and V on H
2995 !! try to work well on more datasets at once
2996 subroutine vg6d_c2a (this)
2998 TYPE(volgrid6d),INTENT(inout) :: this(:) !< vettor of volumes volgrid6d to elaborate
3000 integer :: ngrid,igrid,jgrid,ugrid,vgrid,tgrid
3001 doubleprecision :: xmin, xmax, ymin, ymax
3002 doubleprecision :: xmin_t, xmax_t, ymin_t, ymax_t
3003 doubleprecision :: step_lon_t,step_lat_t
3004 character(len=80) :: type_t,type
3005 TYPE(griddim_def):: griddim_t
3011 call init(griddim_t)
3013 call get_val(this(igrid)%griddim,xmin=xmin_t, xmax=xmax_t, ymin=ymin_t, ymax=ymax_t,proj_type=type_t)
3014 step_lon_t=(xmax_t-xmin_t)/dble(this(igrid)%griddim%dim%nx-1)
3015 step_lat_t=(ymax_t-ymin_t)/dble(this(igrid)%griddim%dim%ny-1)
3024 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: search u/v/t points:
"//to_char(igrid)//to_char(jgrid))
3025 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: test delta:
"//to_char(step_lon_t)//to_char(step_lat_t))
3028 if (this(igrid)%griddim == this(jgrid)%griddim ) cycle
3030 .and.
if (this(igrid)%griddim%dim%nx == this(jgrid)%griddim%dim%nx &
3031 this(igrid)%griddim%dim%ny == this(jgrid)%griddim%dim%ny ) then
3033 call get_val(this(jgrid)%griddim,xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,proj_type=type)
3035 if (type_t /= type )cycle
3038 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: test u
"//&
3039 to_char(xmin)//to_char(xmax)//to_char(ymin)//to_char(ymax))
3041 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"diff coordinate lon
"//&
3042 to_char(abs(xmin - (xmin_t+step_lon_t/2.d0)))//&
3043 to_char(abs(xmax - (xmax_t+step_lon_t/2.d0))))
3044 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"diff coordinate lat
"//&
3045 to_char(abs(ymin - (ymin_t+step_lat_t/2.d0)))//&
3046 to_char(abs(ymax - (ymax_t+step_lat_t/2.d0))))
3049 .and.
if ( abs(xmin - (xmin_t+step_lon_t/2.d0)) < 1.d-3 abs(xmax - (xmax_t+step_lon_t/2.d0)) < 1.d-3 ) then
3050 .and.
if ( abs(ymin - ymin_t) < 1.d-3 abs(ymax - ymax_t) < 1.d-3 ) then
3053 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: found u
")
3062 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: test v
"//&
3063 to_char(xmin)//to_char(xmax)//to_char(ymin)//to_char(ymax))
3066 .and.
if ( abs(ymin - (ymin_t+step_lat_t/2.d0)) < 1.d-3 abs(ymax - (ymax_t+step_lat_t/2.d0)) < 1.d-3 ) then
3067 .and.
if ( abs(xmin - xmin_t) < 1.d-3 abs(xmax - xmax_t) < 1.d-3 ) then
3070 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: found v
")
3079 ! test if we have U and V
3082 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: test u and v
"//&
3083 to_char(xmin_t)//to_char(xmax_t)//to_char(ymin_t)//to_char(ymax_t)//&
3084 to_char(xmin)//to_char(xmax)//to_char(ymin)//to_char(ymax))
3086 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"uv diff coordinate lon
"//&
3087 to_char(abs(xmin_t - xmin)-step_lon_t/2.d0)//&
3088 to_char(abs(xmax_t - xmax)-step_lon_t/2.d0))
3089 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"uv diff coordinate lat
"//&
3090 to_char(abs(ymin_t - ymin) -step_lat_t/2.d0)//&
3091 to_char(abs(ymax_t - ymax)-step_lat_t/2.d0))
3094 .and.
if ( abs(ymin - (ymin_t+step_lat_t/2.d0)) < 2.d-3 abs(ymax - (ymax_t+step_lat_t/2.d0)) < 2.d-3 ) then
3095 .and.
if ( abs(xmin_t - (xmin+step_lon_t/2.d0)) < 2.d-3 abs(xmax_t - (xmax+step_lon_t/2.d0)) < 2.d-3 ) then
3098 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: found u and v case up and right
")
3104 call init(griddim_t,xmin=xmin, xmax=xmax, ymin=ymin_t, ymax=ymax_t)
3110 ! abbiamo almeno due volumi: riportiamo U e/o V su T
3111 if (c_e(ugrid)) then
3113 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid u points: interpolate u on t
"//to_char(tgrid)//to_char(ugrid))
3116 call vg6d_c2a_grid(this(ugrid),this(tgrid)%griddim,cgrid=1)
3118 call vg6d_c2a_grid(this(ugrid),griddim_t,cgrid=1)
3120 call vg6d_c2a_mat(this(ugrid),cgrid=1)
3123 if (c_e(vgrid)) then
3125 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid v points: interpolate v on t
"//to_char(tgrid)//to_char(vgrid))
3128 call vg6d_c2a_grid(this(vgrid),this(tgrid)%griddim,cgrid=2)
3130 call vg6d_c2a_grid(this(vgrid),griddim_t,cgrid=2)
3132 call vg6d_c2a_mat(this(vgrid),cgrid=2)
3137 call delete(griddim_t)
3142 end subroutine vg6d_c2a
3145 ! Convert C grid to A grid
3146 subroutine vg6d_c2a_grid(this,griddim_t,cgrid)
3148 type(volgrid6d),intent(inout) :: this ! object containing fields to be translated (U or V points)
3149 type(griddim_def),intent(in),optional :: griddim_t ! object containing grid of T points
3150 integer,intent(in) :: cgrid ! in C grid (Arakawa) we have 0=T,1=U,2=V points
3152 doubleprecision :: xmin, xmax, ymin, ymax
3153 doubleprecision :: step_lon,step_lat
3156 if (present(griddim_t)) then
3158 call get_val(griddim_t,xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax)
3159 call set_val(this%griddim,xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax)
3160 ! improve grid definition precision
3161 CALL griddim_setsteps(this%griddim)
3170 call l4f_category_log(this%category,L4F_DEBUG,"c grid: t points, nothing to do
")
3177 call l4f_category_log(this%category,L4F_DEBUG,"c grid: u points, we need interpolation
")
3180 call get_val(this%griddim, xmin=xmin, xmax=xmax)
3181 step_lon=(xmax-xmin)/dble(this%griddim%dim%nx-1)
3182 xmin=xmin-step_lon/2.d0
3183 xmax=xmax-step_lon/2.d0
3184 call set_val(this%griddim, xmin=xmin, xmax=xmax)
3185 ! improve grid definition precision
3186 CALL griddim_setsteps(this%griddim)
3191 call l4f_category_log(this%category,L4F_DEBUG,"c grid: v points, we need interpolation
")
3194 call get_val(this%griddim, ymin=ymin, ymax=ymax)
3195 step_lat=(ymax-ymin)/dble(this%griddim%dim%ny-1)
3196 ymin=ymin-step_lat/2.d0
3197 ymax=ymax-step_lat/2.d0
3198 call set_val(this%griddim, ymin=ymin, ymax=ymax)
3199 ! improve grid definition precision
3200 CALL griddim_setsteps(this%griddim)
3204 call l4f_category_log(this%category,L4F_FATAL,"c grid type not known
")
3205 call raise_fatal_error ()
3212 call griddim_unproj(this%griddim)
3215 end subroutine vg6d_c2a_grid
3217 ! Convert C grid to A grid
3218 subroutine vg6d_c2a_mat(this,cgrid)
3220 type(volgrid6d),intent(inout) :: this ! object containing fields to be translated
3221 integer,intent(in) :: cgrid ! in C grid (Arakawa) we have 0=T,1=U,2=V points
3223 INTEGER :: i, j, k, iv, stallo
3224 REAL,ALLOCATABLE :: tmp_arr(:,:)
3225 REAL,POINTER :: voldatiuv(:,:)
3228 IF (cgrid == 0) RETURN ! nothing to do
3229 .AND.
IF (cgrid /= 1 cgrid /= 2) THEN ! wrong cgrid
3230 CALL l4f_category_log(this%category,L4F_FATAL,"c grid type
"// &
3231 TRIM(to_char(cgrid))//" not known
")
3236 ! Temporary workspace
3237 ALLOCATE(tmp_arr(this%griddim%dim%nx, this%griddim%dim%ny),stat=stallo)
3239 call l4f_log(L4F_FATAL,"allocating memory
")
3240 call raise_fatal_error()
3243 ! allocate once for speed
3244 .NOT.
IF (ASSOCIATED(this%voldati)) THEN
3245 ALLOCATE(voldatiuv(this%griddim%dim%nx, this%griddim%dim%ny), stat=stallo)
3246 IF (stallo /= 0) THEN
3247 CALL l4f_log(L4F_FATAL,"allocating memory
")
3248 CALL raise_fatal_error()
3252 IF (cgrid == 1) THEN ! U points to H points
3253 DO iv = 1, SIZE(this%var)
3254 DO k = 1, SIZE(this%timerange)
3255 DO j = 1, SIZE(this%time)
3256 DO i = 1, SIZE(this%level)
3257 tmp_arr(:,:) = rmiss
3258 CALL volgrid_get_vol_2d(this, i, j, k, iv, voldatiuv)
3260 ! West boundary extrapolation
3261 .AND.
WHERE(voldatiuv(1,:) /= rmiss voldatiuv(2,:) /= rmiss)
3262 tmp_arr(1,:) = voldatiuv(1,:) - (voldatiuv(2,:) - voldatiuv(1,:)) / 2.
3266 .AND.
WHERE(voldatiuv(1:this%griddim%dim%nx-1,:) /= rmiss &
3267 voldatiuv(2:this%griddim%dim%nx,:) /= rmiss)
3268 tmp_arr(2:this%griddim%dim%nx,:) = &
3269 (voldatiuv(1:this%griddim%dim%nx-1,:) + &
3270 voldatiuv(2:this%griddim%dim%nx,:)) / 2.
3273 voldatiuv(:,:) = tmp_arr
3274 CALL volgrid_set_vol_2d(this, i, j, k, iv, voldatiuv)
3280 ELSE IF (cgrid == 2) THEN ! V points to H points
3281 DO iv = 1, SIZE(this%var)
3282 DO k = 1, SIZE(this%timerange)
3283 DO j = 1, SIZE(this%time)
3284 DO i = 1, SIZE(this%level)
3285 tmp_arr(:,:) = rmiss
3286 CALL volgrid_get_vol_2d(this, i, j, k, iv, voldatiuv)
3288 ! South boundary extrapolation
3289 .AND.
WHERE(voldatiuv(:,1) /= rmiss voldatiuv(:,2) /= rmiss)
3290 tmp_arr(:,1) = voldatiuv(:,1) - (voldatiuv(:,2) - voldatiuv(:,1)) / 2.
3294 .AND.
WHERE(voldatiuv(:,1:this%griddim%dim%ny-1) /= rmiss &
3295 voldatiuv(:,2:this%griddim%dim%ny) /= rmiss)
3296 tmp_arr(:,2:this%griddim%dim%ny) = &
3297 (voldatiuv(:,1:this%griddim%dim%ny-1) + &
3298 voldatiuv(:,2:this%griddim%dim%ny)) / 2.
3301 voldatiuv(:,:) = tmp_arr
3302 CALL volgrid_set_vol_2d(this, i, j, k, iv, voldatiuv)
3307 ENDIF ! tertium non datur
3309 .NOT.
IF (ASSOCIATED(this%voldati)) THEN
3310 DEALLOCATE(voldatiuv)
3312 DEALLOCATE (tmp_arr)
3314 end subroutine vg6d_c2a_mat
3317 !> \brief Display object on screen
3319 !! Show brief content on screen.
3320 subroutine display_volgrid6d (this)
3322 TYPE(volgrid6d),intent(in) :: this !< object to display
3326 call l4f_category_log(this%category,L4F_DEBUG,"ora mostro gridinfo
" )
3330 call display(this%griddim)
3332 IF (ASSOCIATED(this%time))then
3333 print*,"---- time vector ----
"
3334 print*,"elements=
",size(this%time)
3335 do i=1, size(this%time)
3336 call display(this%time(i))
3340 IF (ASSOCIATED(this%timerange))then
3341 print*,"---- timerange vector ----
"
3342 print*,"elements=
",size(this%timerange)
3343 do i=1, size(this%timerange)
3344 call display(this%timerange(i))
3348 IF (ASSOCIATED(this%level))then
3349 print*,"---- level vector ----
"
3350 print*,"elements=
",size(this%level)
3351 do i=1, size(this%level)
3352 call display(this%level(i))
3356 IF (ASSOCIATED(this%var))then
3357 print*,"---- var vector ----
"
3358 print*,"elements=
",size(this%var)
3359 do i=1, size(this%var)
3360 call display(this%var(i))
3364 IF (ASSOCIATED(this%gaid))then
3365 print*,"---- gaid vector (present mask only) ----
"
3366 print*,"elements=
",shape(this%gaid)
3367 PRINT*,c_e(RESHAPE(this%gaid,(/SIZE(this%gaid)/)))
3370 print*,"--------------------------------------------------------------
"
3373 end subroutine display_volgrid6d
3376 !> \brief Display vector of object on screen
3378 !! Show brief content on screen.
3379 subroutine display_volgrid6dv (this)
3381 TYPE(volgrid6d),intent(in) :: this(:) !< vector of object to display
3384 print*,"-----------------------
volgrid6d vector ---------------------
"
3386 print*,"elements=
",size(this)
3391 call l4f_category_log(this(i)%category,L4F_DEBUG,"ora mostro il vettore
volgrid6d" )
3394 call display(this(i))
3397 print*,"--------------------------------------------------------------
"
3399 end subroutine display_volgrid6dv
3402 !> Reduce some dimensions (level and timerage) for semplification (rounding).
3403 !! Vector version of vg6dv_rounding
3404 subroutine vg6dv_rounding(vg6din,vg6dout,level,timerange,nostatproc,merge)
3405 type(volgrid6d),intent(in) :: vg6din(:) !< input volume
3406 type(volgrid6d),intent(out),pointer :: vg6dout(:) !> output volume
3407 type(vol7d_level),intent(in),optional :: level(:) !< almost equal level list
3408 type(vol7d_timerange),intent(in),optional :: timerange(:) !< almost equal timerange list
3409 logical,intent(in),optional :: merge !< if there are data on more then one almost equal levels or timeranges will be merged POINT BY POINT with priority for the fird data found ordered by icreasing var index
3410 logical,intent(in),optional :: nostatproc !< do not take in account statistical processing code in timerange and P2
3414 allocate(vg6dout(size(vg6din)))
3416 do i = 1, size(vg6din)
3417 call vg6d_rounding(vg6din(i),vg6dout(i),level,timerange,nostatproc,merge)
3420 end subroutine vg6dv_rounding
3422 !> Reduce some dimensions (level and timerage) for semplification (rounding).
3423 !! You can use this for simplify and use variables in computation like alchimia
3424 !! where fields have to be on the same coordinate
3426 !! means in time for short periods and istantaneous values
3427 !! 2 meter and surface levels
3428 !! If there are data on more then one almost equal levels or timeranges, the first var present (at least one point)
3429 !! will be taken (order is by icreasing var index).
3430 !! You can use predefined values for classic semplification
3431 !! almost_equal_levels and almost_equal_timeranges
3432 !! The level or timerange in output will be defined by the first element of level and timerange list
3433 subroutine vg6d_rounding(vg6din,vg6dout,level,timerange,nostatproc,merge)
3434 type(volgrid6d),intent(in) :: vg6din !< input volume
3435 type(volgrid6d),intent(out) :: vg6dout !> output volume
3436 type(vol7d_level),intent(in),optional :: level(:) !< almost equal level list
3437 type(vol7d_timerange),intent(in),optional :: timerange(:) !< almost equal timerange list
3438 logical,intent(in),optional :: merge !< if there are data on more then one almost equal levels or timeranges
3439 !! will be merged POINT BY POINT with priority for the fird data found ordered by icreasing var index
3440 logical,intent(in),optional :: nostatproc !< do not take in account statistical processing code in timerange and P2
3442 integer :: ilevel,itimerange
3443 type(vol7d_level) :: roundlevel(size(vg6din%level))
3444 type(vol7d_timerange) :: roundtimerange(size(vg6din%timerange))
3446 roundlevel=vg6din%level
3448 if (present(level))then
3449 do ilevel = 1, size(vg6din%level)
3450 if ((any(vg6din%level(ilevel) .almosteq. level))) then
3451 roundlevel(ilevel)=level(1)
3456 roundtimerange=vg6din%timerange
3458 if (present(timerange))then
3459 do itimerange = 1, size(vg6din%timerange)
3460 if ((any(vg6din%timerange(itimerange) .almosteq. timerange))) then
3461 roundtimerange(itimerange)=timerange(1)
3466 !set istantaneous values everywere
3467 !preserve p1 for forecast time
3468 if (optio_log(nostatproc)) then
3469 roundtimerange(:)%timerange=254
3470 roundtimerange(:)%p2=0
3474 call vg6d_reduce(vg6din,vg6dout,roundlevel,roundtimerange,merge)
3476 end subroutine vg6d_rounding
3478 !> Reduce some dimensions (level and timerage).
3479 !! You can pass a volume and specify duplicated levels and timeranges in roundlevel and roundtimerange;
3480 !! you get unique levels and timeranges in output.
3481 !! If there are data on equal levels or timeranges, the first var present (at least one point)
3482 !! will be taken (order is by icreasing var index).
3483 !! you can specify merge and if there are data on equal levels or timeranges will be merged POINT BY POINT
3484 !! with priority for the first data found ordered by icreasing var index (require to decode all the data)
3485 !! Data are decoded only if needed so the output should be with or without voldata allocated
3486 subroutine vg6d_reduce(vg6din,vg6dout,roundlevel,roundtimerange,merge)
3487 type(volgrid6d),intent(in) :: vg6din !< input volume
3488 type(volgrid6d),intent(out) :: vg6dout !< output volume
3489 type(vol7d_level),intent(in) :: roundlevel(:) !< new level list to use for rounding
3490 type(vol7d_timerange),intent(in) :: roundtimerange(:) !< new timerange list to use for rounding
3491 logical,intent(in),optional :: merge !< if there are data on equal levels or timeranges will be merged POINT BY POINT with priority for the first data found ordered by icreasing var index (require to decode all the data)
3493 integer :: nlevel,ntime,ntimerange,nvar,ilevel,itimerange,ivar,indl,indt,itime,nx,ny
3494 real,allocatable :: vol2d(:,:)
3496 nx=vg6din%griddim%dim%nx
3497 ny=vg6din%griddim%dim%ny
3498 nlevel=count_distinct(roundlevel,back=.true.)
3499 ntime=size(vg6din%time)
3500 ntimerange=count_distinct(roundtimerange,back=.true.)
3501 nvar=size(vg6din%var)
3503 call init(vg6dout, vg6din%griddim, vg6din%time_definition, categoryappend="generated by vg6d_reduce
")
3504 call volgrid6d_alloc(vg6dout, vg6din%griddim%dim, ntime, nlevel, ntimerange, nvar)
3506 .or.
if ( ASSOCIATED(vg6din%voldati) optio_log(merge)) then
3507 call volgrid6d_alloc_vol(vg6dout,inivol=.true.,decode=.true.)
3508 allocate(vol2d(nx,ny))
3510 call volgrid6d_alloc_vol(vg6dout,inivol=.true.,decode=.false.)
3513 vg6dout%time=vg6din%time
3514 vg6dout%var=vg6din%var
3515 vg6dout%timerange=pack_distinct(roundtimerange,ntimerange,back=.true.)
3516 vg6dout%level=pack_distinct(roundlevel,nlevel,back=.true.)
3517 ! sort modified dimensions
3518 CALL sort(vg6dout%timerange)
3519 CALL sort(vg6dout%level)
3521 do ilevel=1,size(vg6din%level)
3522 indl=index(vg6dout%level,roundlevel(ilevel))
3523 do itimerange=1,size(vg6din%timerange)
3524 indt=index(vg6dout%timerange,roundtimerange(itimerange))
3528 if ( ASSOCIATED(vg6din%voldati)) then
3529 vol2d=vg6din%voldati(:,:,ilevel,itime,itimerange,ivar)
3532 if (optio_log(merge)) then
3534 .not.
if ( ASSOCIATED(vg6din%voldati)) then
3535 CALL grid_id_decode_data(vg6din%gaid(ilevel,itime,itimerange,ivar), vol2d)
3538 !! merge present data point by point
3539 .not.
where ( c_e(vg6dout%voldati(:,:,indl,itime,indt,ivar)))
3541 vg6dout%voldati(:,:,indl,itime,indt,ivar)=vol2d
3544 else if ( ASSOCIATED(vg6din%voldati)) then
3545 .not.
if ( any(c_e(vg6dout%voldati(:,:,indl,itime,indt,ivar))))then
3546 vg6dout%voldati(:,:,indl,itime,indt,ivar)=vol2d
3550 .and..not.
if (c_e(vg6din%gaid(ilevel,itime,itimerange,ivar)) c_e(vg6dout%gaid(indl,itime,indt,ivar)))then
3551 call copy (vg6din%gaid(ilevel,itime,itimerange,ivar), vg6dout%gaid(indl,itime,indt,ivar))
3558 .or.
if ( ASSOCIATED(vg6din%voldati) optio_log(merge)) then
3562 end subroutine vg6d_reduce
3565 end module volgrid6d_class
3569 !>\example example_vg6d_3.f90
3570 !!\brief Programma esempio semplice per gridinfo e volgrid6d.
3572 !! Programma che importa da file un vettore di gridinfo poi lo importa in volgrid6d. Da volgrid6d viene di nuovo creato un vettore di gridinfo per poi exportare su file.
3574 !>\example example_vg6d_5.f90
3575 !!\brief Programma trasformazione da volgrid6d a volgrid6d
3577 !! Legge grib da un file e li organizza in un vettore di strutture
3578 !! volgrid6d mettendoli a disposizione per eventuali elaborazioni;
3579 !! vengono poi riesportati a un file grib
3581 !>\example example_vg6d_8.f90
3582 !! \brief Programma scrittura su file vettore di anagrafica
3584 !>\example example_vg6d_6.f90
3585 !! \brief Programma trasformazione da volgrid6d a vol7d
3587 !>\example example_vg6d_7.f90
3588 !! \brief Programma trasformazione da vol7d a volgrid6d
3590 !>\example example_vg6d_9.f90
3591 !! \brief Example to create a grib editionNumber = 2 file from data generated in memory using a grib_api template
Functions that return a trimmed CHARACTER representation of the input variable.
Read the object from a formatted or unformatted file.
Write the object on a formatted or unformatted file.
Clone the object, creating a new independent instance of the object exactly equal to the starting one...
Decode and return the data array from a grid_id object associated to a gridinfo object.
Encode a data array into a grid_id object associated to a gridinfo object.
Emit log message for a category with specific priority.
Represent level object in a pretty string.
Destructor, it releases every information and memory buffer associated with the object.
Display on standard output a description of the volgrid6d object provided.
Export an object dirctly to a native file, to a gridinfo object or to a supported file format through...
Import an object dirctly from a native file, from a gridinfo object or from a supported file format t...
Constructor, it creates a new instance of the object.
Reduce some dimensions (level and timerage) for semplification (rounding).
Apply the conversion function this to values.
This module defines usefull general purpose function and subroutine.
Classi per la gestione delle coordinate temporali.
Module for describing geographically referenced regular grids.
Module for defining the extension and coordinates of a rectangular georeferenced grid.
This module defines an abstract interface to different drivers for access to files containing gridded...
Class for managing information about a single gridded georeferenced field, typically imported from an...
classe per la gestione del logging
Module for quickly interpreting the OPTIONAL parameters passed to a subprogram.
Classe per la gestione di un volume completo di dati osservati.
classe per import ed export di volumi da e in DB-All.e
Classe per la gestione dei livelli verticali in osservazioni meteo e affini.
Classe per la gestione degli intervalli temporali di osservazioni meteo e affini.
This module defines objects and methods for managing data volumes on rectangular georeferenced grids.
Class for managing physical variables in a grib 1/2 fashion.
Class for expressing an absolute time value.
Derived type associated to a block/message/record/band of gridded data coming from a file-like object...
Derived type defining a dynamically extensible array of TYPE(gridinfo_def) elements.
Object describing a single gridded message/band.
Definisce un oggetto contenente i volumi anagrafica e dati e tutti i descrittori delle loro dimension...
Definisce l'intervallo temporale di un'osservazione meteo.
Object describing a rectangular, homogeneous gridded dataset.