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()
281 CALL init(this%level(i))
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()
299 CALL init(this%timerange(i))
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)))
480 DO ilevel = 1,
SIZE(this%level)
481 CALL grid_id_decode_data(this%gaid(ilevel,itime,itimerange,ivar), &
486 END SUBROUTINE volgrid_get_vol_3d
500 SUBROUTINE volgrid_set_vol_2d(this, ilevel, itime, itimerange, ivar, voldati)
502 INTEGER,
INTENT(in) :: ilevel
503 INTEGER,
INTENT(in) :: itime
504 INTEGER,
INTENT(in) :: itimerange
505 INTEGER,
INTENT(in) :: ivar
506 REAL,
INTENT(in) :: voldati(:,:)
508 IF (
ASSOCIATED(this%voldati))
THEN
511 CALL grid_id_encode_data(this%gaid(ilevel,itime,itimerange,ivar), voldati)
514 END SUBROUTINE volgrid_set_vol_2d
528 SUBROUTINE volgrid_set_vol_3d(this, itime, itimerange, ivar, voldati)
530 INTEGER,
INTENT(in) :: itime
531 INTEGER,
INTENT(in) :: itimerange
532 INTEGER,
INTENT(in) :: ivar
533 REAL,
INTENT(in) :: voldati(:,:,:)
537 IF (
ASSOCIATED(this%voldati))
THEN
540 DO ilevel = 1,
SIZE(this%level)
541 CALL grid_id_encode_data(this%gaid(ilevel,itime,itimerange,ivar), &
546 END SUBROUTINE volgrid_set_vol_3d
552 SUBROUTINE volgrid6d_delete(this)
555 INTEGER :: i, ii, iii, iiii
561 if (
associated(this%gaid))
then
563 DO i=1 ,
SIZE(this%gaid,1)
564 DO ii=1 ,
SIZE(this%gaid,2)
565 DO iii=1 ,
SIZE(this%gaid,3)
566 DO iiii=1 ,
SIZE(this%gaid,4)
567 CALL delete(this%gaid(i,ii,iii,iiii))
572 DEALLOCATE(this%gaid)
583 if (
associated( this%time ))
deallocate(this%time)
584 if (
associated( this%timerange ))
deallocate(this%timerange)
585 if (
associated( this%level ))
deallocate(this%level)
586 if (
associated( this%var ))
deallocate(this%var)
588 if (
associated(this%voldati))
deallocate(this%voldati)
592 call l4f_category_delete(this%category)
594 END SUBROUTINE volgrid6d_delete
606 subroutine volgrid6d_write_on_file (this,unit,description,filename,filename_auto)
609 integer,
optional,
intent(inout) :: unit
610 character(len=*),
intent(in),
optional :: filename
611 character(len=*),
intent(out),
optional :: filename_auto
612 character(len=*),
INTENT(IN),
optional :: description
615 character(len=254) :: ldescription,arg,lfilename
616 integer :: ntime, ntimerange, nlevel, nvar
618 logical :: opened,exist
630 call date_and_time(values=tarray)
633 if (
present(description))
then
634 ldescription=description
636 ldescription=
"Volgrid6d generated by: "//trim(arg)
639 if (.not.
present(unit))
then
650 lfilename=trim(arg)//
".vg6d"
651 if (
index(arg,
'/',back=.true.) > 0) lfilename=lfilename(
index(arg,
'/',back=.true.)+1 : )
653 if (
present(filename))
then
654 if (filename /=
"")
then
659 if (
present(filename_auto))filename_auto=lfilename
662 inquire(unit=lunit,opened=opened)
663 if (.not. opened)
then
664 inquire(file=lfilename,exist=exist)
665 if (exist)
CALL raise_error(
'file exist; cannot open new file')
666 if (.not.exist)
open (unit=lunit,file=lfilename,form=
"UNFORMATTED")
670 if (
associated(this%time)) ntime=
size(this%time)
671 if (
associated(this%timerange)) ntimerange=
size(this%timerange)
672 if (
associated(this%level)) nlevel=
size(this%level)
673 if (
associated(this%var)) nvar=
size(this%var)
676 write(unit=lunit)ldescription
677 write(unit=lunit)tarray
680 write(unit=lunit) ntime, ntimerange, nlevel, nvar
683 if (
associated(this%time))
call write_unit(this%time, lunit)
684 if (
associated(this%level))
write(unit=lunit)this%level
685 if (
associated(this%timerange))
write(unit=lunit)this%timerange
686 if (
associated(this%var))
write(unit=lunit)this%var
691 if (
associated(this%voldati))
write(unit=lunit)this%voldati
693 if (.not.
present(unit))
close(unit=lunit)
695 end subroutine volgrid6d_write_on_file
704 subroutine volgrid6d_read_from_file (this,unit,filename,description,tarray,filename_auto)
707 integer,
intent(inout),
optional :: unit
708 character(len=*),
INTENT(in),
optional :: filename
709 character(len=*),
intent(out),
optional :: filename_auto
710 character(len=*),
INTENT(out),
optional :: description
711 integer,
intent(out),
optional :: tarray(8)
713 integer :: ntime, ntimerange, nlevel, nvar
715 character(len=254) :: ldescription,lfilename,arg
716 integer :: ltarray(8),lunit
717 logical :: opened,exist
725 if (.not.
present(unit))
then
736 lfilename=trim(arg)//
".vg6d"
737 if (
index(arg,
'/',back=.true.) > 0) lfilename=lfilename(
index(arg,
'/',back=.true.)+1 : )
739 if (
present(filename))
then
740 if (filename /=
"")
then
745 if (
present(filename_auto))filename_auto=lfilename
748 inquire(unit=lunit,opened=opened)
749 if (.not. opened)
then
750 inquire(file=lfilename,exist=exist)
751 IF (.NOT. exist)
CALL raise_fatal_error(
'file '//trim(lfilename)//
' does not exist, cannot open')
752 open (unit=lunit,file=lfilename,form=
"UNFORMATTED")
755 read(unit=lunit)ldescription
756 read(unit=lunit)ltarray
758 call l4f_log(l4f_info,
"Info: reading volgrid6d from file: "//trim(lfilename))
759 call l4f_log(l4f_info,
"Info: description: "//trim(ldescription))
762 if (
present(description))description=ldescription
763 if (
present(tarray))tarray=ltarray
767 read(unit=lunit) ntime, ntimerange, nlevel, nvar
770 call volgrid6d_alloc (this, &
771 ntime=ntime, ntimerange=ntimerange, nlevel=nlevel, nvar=nvar)
773 call volgrid6d_alloc_vol (this)
775 if (
associated(this%time))
call read_unit(this%time, lunit)
776 if (
associated(this%level))
read(unit=lunit)this%level
777 if (
associated(this%timerange))
read(unit=lunit)this%timerange
778 if (
associated(this%var))
read(unit=lunit)this%var
783 if (
associated(this%voldati))
read(unit=lunit)this%voldati
785 if (.not.
present(unit))
close(unit=lunit)
787 end subroutine volgrid6d_read_from_file
809 SUBROUTINE import_from_gridinfo(this, gridinfo, force, dup_mode, clone, &
813 LOGICAL,
INTENT(in),
OPTIONAL :: force
814 INTEGER,
INTENT(in),
OPTIONAL :: dup_mode
815 LOGICAL ,
INTENT(in),
OPTIONAL :: clone
816 LOGICAL,
INTENT(IN),
OPTIONAL :: isanavar
818 CHARACTER(len=255) :: type
819 INTEGER :: itime0, itimerange0, itime1, itimerange1, itime, itimerange, &
820 ilevel, ivar, ldup_mode
823 REAL,
ALLOCATABLE :: tmpgrid(:,:)
825 IF (
PRESENT(dup_mode))
THEN
831 call get_val(this%griddim,proj_type=type)
834 call l4f_category_log(this%category,l4f_debug,
"import_from_gridinfo: "//trim(type))
837 if (.not.
c_e(type))
then
838 call copy(gridinfo%griddim, this%griddim)
842 CALL volgrid6d_alloc_vol(this, ini=.true.)
844 else if (.not. (this%griddim == gridinfo%griddim ))
then
847 "volgrid and gridinfo grid type or size are different, gridinfo rejected")
854 ilevel =
index(this%level, gridinfo%level)
855 IF (ilevel == 0 .AND. optio_log(force))
THEN
856 ilevel =
index(this%level, vol7d_level_miss)
857 IF (ilevel /= 0) this%level(ilevel) = gridinfo%level
860 IF (ilevel == 0)
THEN
862 "volgrid6d: level not valid for volume, gridinfo rejected")
867 IF (optio_log(isanavar))
THEN
869 itime1 =
SIZE(this%time)
871 itimerange1 =
SIZE(this%timerange)
873 correctedtime = gridinfo%time
874 IF (this%time_definition == 1) correctedtime = correctedtime + &
875 timedelta_new(sec=gridinfo%timerange%p1)
876 itime0 =
index(this%time, correctedtime)
877 IF (itime0 == 0 .AND. optio_log(force))
THEN
878 itime0 =
index(this%time, datetime_miss)
879 IF (itime0 /= 0) this%time(itime0) = correctedtime
881 IF (itime0 == 0)
THEN
883 "volgrid6d: time not valid for volume, gridinfo rejected")
889 itimerange0 =
index(this%timerange,gridinfo%timerange)
890 IF (itimerange0 == 0 .AND. optio_log(force))
THEN
891 itimerange0 =
index(this%timerange, vol7d_timerange_miss)
892 IF (itimerange0 /= 0) this%timerange(itimerange0) = gridinfo%timerange
894 IF (itimerange0 == 0)
THEN
896 "volgrid6d: timerange not valid for volume, gridinfo rejected")
900 itimerange1 = itimerange0
903 ivar =
index(this%var, gridinfo%var)
904 IF (ivar == 0 .AND. optio_log(force))
THEN
905 ivar =
index(this%var, volgrid6d_var_miss)
906 IF (ivar /= 0) this%var(ivar) = gridinfo%var
910 "volgrid6d: var not valid for volume, gridinfo rejected")
915 DO itimerange = itimerange0, itimerange1
916 DO itime = itime0, itime1
917 IF (
ASSOCIATED(this%gaid))
THEN
919 IF (
c_e(this%gaid(ilevel,itime,itimerange,ivar)))
THEN
923 IF (optio_log(
clone))
CALL delete(this%gaid(ilevel,itime,itimerange,ivar))
926 IF (optio_log(
clone))
THEN
927 CALL copy(gridinfo%gaid, this%gaid(ilevel,itime,itimerange,ivar))
932 this%gaid(ilevel,itime,itimerange,ivar) = gridinfo%gaid
935 IF (
ASSOCIATED(this%voldati))
THEN
936 IF (.NOT.dup .OR. ldup_mode == 0)
THEN
937 this%voldati(:,:,ilevel,itime,itimerange,ivar) =
decode_gridinfo(gridinfo)
938 ELSE IF (ldup_mode == 1)
THEN
941 this%voldati(:,:,ilevel,itime,itimerange,ivar) = tmpgrid(:,:)
943 ELSE IF (ldup_mode == 2)
THEN
944 WHERE(.NOT.
c_e(this%voldati(:,:,ilevel,itime,itimerange,ivar)))
945 this%voldati(:,:,ilevel,itime,itimerange,ivar) =
decode_gridinfo(gridinfo)
952 "gaid not allocated, you probably need to call volgrid6d_alloc_vol first")
960 END SUBROUTINE import_from_gridinfo
967 SUBROUTINE export_to_gridinfo(this, gridinfo, itime, itimerange, ilevel, ivar, &
968 gaid_template, clone)
972 INTEGER :: itimerange
975 TYPE(
grid_id),
INTENT(in),
OPTIONAL :: gaid_template
976 LOGICAL,
INTENT(in),
OPTIONAL :: clone
979 LOGICAL :: usetemplate
980 REAL,
POINTER :: voldati(:,:)
987 IF (.NOT.
c_e(this%gaid(ilevel,itime,itimerange,ivar)))
THEN
989 CALL l4f_category_log(this%category,l4f_debug,
"empty gaid found, skipping export")
994 usetemplate = .false.
995 IF (
PRESENT(gaid_template))
THEN
996 CALL copy(gaid_template, gaid)
998 CALL l4f_category_log(this%category,l4f_debug,
"template cloned to a new gaid")
1000 usetemplate =
c_e(gaid)
1003 IF (.NOT.usetemplate)
THEN
1004 IF (optio_log(
clone))
THEN
1005 CALL copy(this%gaid(ilevel,itime,itimerange,ivar), gaid)
1007 CALL l4f_category_log(this%category,l4f_debug,
"original gaid cloned to a new one")
1010 gaid = this%gaid(ilevel,itime,itimerange,ivar)
1014 IF (this%time_definition == 1)
THEN
1015 correctedtime = this%time(itime) - &
1016 timedelta_new(sec=this%timerange(itimerange)%p1)
1018 correctedtime = this%time(itime)
1021 CALL init(gridinfo,gaid, this%griddim, correctedtime, this%timerange(itimerange), &
1022 this%level(ilevel), this%var(ivar))
1025 CALL export(gridinfo%griddim, gridinfo%gaid)
1027 IF (
ASSOCIATED(this%voldati))
THEN
1028 CALL encode_gridinfo(gridinfo, this%voldati(:,:,ilevel,itime,itimerange,ivar))
1029 ELSE IF (usetemplate)
THEN
1031 CALL volgrid_get_vol_2d(this, ilevel, itime, itimerange, ivar, voldati)
1036 END SUBROUTINE export_to_gridinfo
1056 SUBROUTINE import_from_gridinfovv(this, gridinfov, dup_mode, clone, decode, &
1057 time_definition, anavar, categoryappend)
1060 INTEGER,
INTENT(in),
OPTIONAL :: dup_mode
1061 LOGICAL ,
INTENT(in),
OPTIONAL :: clone
1062 LOGICAL,
INTENT(in),
OPTIONAL :: decode
1063 INTEGER,
INTENT(IN),
OPTIONAL :: time_definition
1064 CHARACTER(len=*),
INTENT(IN),
OPTIONAL :: anavar(:)
1065 CHARACTER(len=*),
INTENT(in),
OPTIONAL :: categoryappend
1067 INTEGER :: i, j, stallo
1068 INTEGER :: ngrid, ntime, ntimerange, nlevel, nvar
1070 CHARACTER(len=512) :: a_name
1071 TYPE(
datetime),
ALLOCATABLE :: correctedtime(:)
1072 LOGICAL,
ALLOCATABLE :: isanavar(:)
1073 TYPE(vol7d_var) :: lvar
1076 if (
present(categoryappend))
then
1077 call l4f_launcher(a_name,a_name_append=trim(subcategory)//
"."//trim(categoryappend))
1079 call l4f_launcher(a_name,a_name_append=trim(subcategory))
1081 category=l4f_category_get(a_name)
1087 ngrid=count_distinct(gridinfov%array(1:gridinfov%arraysize)%griddim,back=.true.)
1089 ' different grid definition(s) found in input data')
1091 ALLOCATE(this(ngrid),stat=stallo)
1092 IF (stallo /= 0)
THEN
1094 CALL raise_fatal_error()
1097 IF (
PRESENT(categoryappend))
THEN
1098 CALL init(this(i), time_definition=time_definition, categoryappend=trim(categoryappend)//
"-vol"//
t2c(i))
1100 CALL init(this(i), time_definition=time_definition, categoryappend=
"vol"//
t2c(i))
1104 this(:)%griddim=pack_distinct(gridinfov%array(1:gridinfov%arraysize)%griddim, &
1108 ALLOCATE(isanavar(gridinfov%arraysize))
1109 isanavar(:) = .false.
1110 IF (
PRESENT(anavar))
THEN
1111 DO i = 1, gridinfov%arraysize
1112 DO j = 1,
SIZE(anavar)
1113 lvar =
convert(gridinfov%array(i)%var)
1114 IF (lvar%btable == anavar(j))
THEN
1115 isanavar(i) = .true.
1121 t2c(gridinfov%arraysize)//
' constant-data messages found in input data')
1125 ALLOCATE(correctedtime(gridinfov%arraysize))
1126 correctedtime(:) = gridinfov%array(1:gridinfov%arraysize)%time
1127 IF (
PRESENT(time_definition))
THEN
1128 IF (time_definition == 1)
THEN
1129 DO i = 1, gridinfov%arraysize
1130 correctedtime(i) = correctedtime(i) + &
1131 timedelta_new(sec=gridinfov%array(i)%timerange%p1)
1137 IF (
PRESENT(anavar))
THEN
1138 j = count((this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1139 .AND. .NOT.isanavar(:))
1142 ' has only constant data, this is not allowed')
1144 CALL raise_fatal_error()
1147 ntime = count_distinct(correctedtime, &
1148 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1149 .AND. .NOT.isanavar(:), back=.true.)
1150 ntimerange = count_distinct(gridinfov%array(1:gridinfov%arraysize)%timerange, &
1151 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1152 .AND. .NOT.isanavar(:), back=.true.)
1153 nlevel = count_distinct(gridinfov%array(1:gridinfov%arraysize)%level, &
1154 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1156 nvar = count_distinct(gridinfov%array(1:gridinfov%arraysize)%var, &
1157 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1164 CALL volgrid6d_alloc(this(i),this(i)%griddim%dim,ntime=ntime, &
1165 ntimerange=ntimerange,nlevel=nlevel,nvar=nvar)
1167 this(i)%time = pack_distinct(correctedtime, ntime, &
1168 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1169 .AND. .NOT.isanavar(:), back=.true.)
1170 CALL sort(this(i)%time)
1172 this(i)%timerange = pack_distinct(gridinfov%array( &
1173 1:gridinfov%arraysize)%timerange, ntimerange, &
1174 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1175 .AND. .NOT.isanavar(:), back=.true.)
1176 CALL sort(this(i)%timerange)
1178 this(i)%level=pack_distinct(gridinfov%array(1:gridinfov%arraysize)%level, &
1179 nlevel,mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1181 CALL sort(this(i)%level)
1183 this(i)%var=pack_distinct(gridinfov%array(1:gridinfov%arraysize)%var, nvar, &
1184 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1190 CALL volgrid6d_alloc_vol(this(i), decode=decode)
1194 DEALLOCATE(correctedtime)
1196 DO i = 1, gridinfov%arraysize
1200 CALL l4f_category_log(category,L4F_INFO, &
1201 "to
volgrid6d index:
"//t2c(index(this%griddim, gridinfov%array(i)%griddim)))
1204 CALL import(this(index(this%griddim, gridinfov%array(i)%griddim)), &
1205 gridinfov%array(i), dup_mode=dup_mode, clone=clone, isanavar=isanavar(i))
1209 !chiudo il logger temporaneo
1210 CALL l4f_category_delete(category)
1212 END SUBROUTINE import_from_gridinfovv
1215 !> Export a \a volgrid6d object to an \a arrayof_gridinfo object.
1216 !! The multidimensional \a volgrid6d structure is serialized into a
1217 !! one-dimensional array of gridinfo_def objects, which is allocated
1218 !! to the proper size if not already allocated, or it is extended
1219 !! keeping the old data if any.
1220 SUBROUTINE export_to_gridinfov(this, gridinfov, gaid_template, clone)
1221 TYPE(volgrid6d),INTENT(inout) :: this !< volume to be exported
1222 TYPE(arrayof_gridinfo),INTENT(inout) :: gridinfov !< output array of gridinfo_def objects
1223 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
1224 LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE., clone the grid_id included in \a this rather than making a shallow copy
1226 INTEGER :: i ,itime, itimerange, ilevel, ivar
1227 INTEGER :: ntime, ntimerange, nlevel, nvar
1228 TYPE(gridinfo_def) :: gridinfol
1231 CALL l4f_category_log(this%category,L4F_DEBUG,"start export_to_gridinfov
")
1234 ! this is necessary in order not to repeatedly and uselessly copy the
1235 ! same array of coordinates for each 2d grid during export, the
1236 ! side-effect is that the computed projection in this is lost
1237 CALL dealloc(this%griddim%dim)
1240 ntime=size(this%time)
1241 ntimerange=size(this%timerange)
1242 nlevel=size(this%level)
1246 DO itimerange=1,ntimerange
1250 CALL init(gridinfol)
1251 CALL export(this, gridinfol, itime, itimerange, ilevel, ivar, &
1252 gaid_template=gaid_template, clone=clone)
1253 IF (c_e(gridinfol%gaid)) THEN
1254 CALL insert(gridinfov, gridinfol)
1256 CALL delete(gridinfol)
1264 END SUBROUTINE export_to_gridinfov
1267 !> Export an array of \a volgrid6d objects to an \a arrayof_gridinfo object.
1268 !! The multidimensional \a volgrid6d structures are serialized into a
1269 !! one-dimensional array of gridinfo_def objects, which is allocated
1270 !! to the proper size if not already allocated, or it is extended
1271 !! keeping the old data if any.
1272 SUBROUTINE export_to_gridinfovv(this, gridinfov, gaid_template, clone)
1275 TYPE(volgrid6d),INTENT(inout) :: this(:) !< volume array to be exported
1276 TYPE(arrayof_gridinfo),INTENT(inout) :: gridinfov !< output array of gridinfo_def objects
1277 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
1278 LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE., clone the grid_id included in \a this rather than making a shallow copy
1282 DO i = 1, SIZE(this)
1284 CALL l4f_category_log(this(i)%category,L4F_DEBUG, &
1285 "export_to_gridinfovv grid
index:
"//t2c(i))
1287 CALL export(this(i), gridinfov, gaid_template=gaid_template, clone=clone)
1290 END SUBROUTINE export_to_gridinfovv
1293 Import the content of a supported file (like grib or gdal-supported
!>
1294 !! format) into an array of \a volgrid6d objects. This method imports
1295 !! a set of gridded fields from a file object into a suitable number
1296 !! of \a volgrid6d objects. The data are imported by creating a
1297 !! temporary \a gridinfo object, importing it into the \a volgrid6d
1298 !! object cloning the gaid's and then destroying the gridinfo, so it
1299 !! works similarly to volgrid6d_class::import_from_gridinfovv() method.
1300 !! For a detailed explanation of the \a anavar argument, see the
1301 !! documentation of volgrid6d_class::import_from_gridinfovv() method.
1302 SUBROUTINE volgrid6d_import_from_file(this, filename, dup_mode, decode, &
1303 time_definition, anavar, categoryappend)
1304 importTYPE(volgrid6d),POINTER :: this(:) !< object in which to
1305 importCHARACTER(len=*),INTENT(in) :: filename !< name of file from which to
1306 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
1307 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
1308 INTEGER,INTENT(IN),OPTIONAL :: time_definition !< 0=time is reference time; 1=time is validity time
1309 CHARACTER(len=*),INTENT(IN),OPTIONAL :: anavar(:) !< list of variables (B-table code equivalent) to be treated as time-independent data
1310 character(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
1312 TYPE(arrayof_gridinfo) :: gridinfo
1314 CHARACTER(len=512) :: a_name
1318 IF (PRESENT(categoryappend))THEN
1319 CALL l4f_launcher(a_name,a_name_append= &
1320 TRIM(subcategory)//".
"//TRIM(categoryappend))
1322 CALL l4f_launcher(a_name,a_name_append=TRIM(subcategory))
1324 category=l4f_category_get(a_name)
1326 CALL import(gridinfo, filename=filename, categoryappend=categoryappend)
1328 IF (gridinfo%arraysize > 0) THEN
1330 CALL import(this, gridinfo, dup_mode=dup_mode, clone=.TRUE., decode=decode, &
1331 time_definition=time_definition, anavar=anavar, &
1332 categoryappend=categoryappend)
1334 CALL l4f_category_log(category,L4F_INFO,"deleting gridinfo
")
1335 CALL delete(gridinfo)
1338 CALL l4f_category_log(category,L4F_INFO,"file does not contain gridded data
")
1342 CALL l4f_category_delete(category)
1344 END SUBROUTINE volgrid6d_import_from_file
1347 !> High level method for exporting a volume array to file.
1348 !! All the information contained into an array of \a volgrid6d
1349 !! objects, i.e. dimension descriptors and data, is exported to a file
1350 !! using the proper output driver (typically grib_api for grib
1351 !! format). If a template is provided, it will determine the
1352 !! characteristic of the output file, otherwise the \a grid_id
1353 !! descriptors contained in the volgrid6d object will be used
1354 SUBROUTINE volgrid6d_export_to_file(this, filename, gaid_template, categoryappend)
1355 TYPE(volgrid6d) :: this(:) !< volume(s) to be exported
1356 CHARACTER(len=*),INTENT(in) :: filename !< output file name
1357 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
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
1364 IF (PRESENT(categoryappend)) THEN
1365 CALL l4f_launcher(a_name,a_name_append=TRIM(subcategory)//".
"//TRIM(categoryappend))
1367 CALL l4f_launcher(a_name,a_name_append=TRIM(subcategory))
1369 category=l4f_category_get(a_name)
1372 CALL l4f_category_log(category,L4F_DEBUG,"start
export to file
")
1375 CALL l4f_category_log(category,L4F_INFO,"writing
volgrid6d to grib file:
"//TRIM(filename))
1377 !IF (ASSOCIATED(this)) THEN
1378 CALL export(this, gridinfo, gaid_template=gaid_template, clone=.TRUE.)
1379 IF (gridinfo%arraysize > 0) THEN
1380 CALL export(gridinfo, filename)
1381 CALL delete(gridinfo)
1384 ! CALL l4f_category_log(category,L4F_INFO,"volume
volgrid6d is not associated
")
1388 CALL l4f_category_delete(category)
1390 END SUBROUTINE volgrid6d_export_to_file
1393 !> Array destructor for \a volgrid6d class.
1394 !! Delete an array of \a volgrid6d objects and deallocate the array
1396 SUBROUTINE volgrid6dv_delete(this)
1397 TYPE(volgrid6d),POINTER :: this(:) !< vector of volgrid6d object
1401 IF (ASSOCIATED(this)) THEN
1402 DO i = 1, SIZE(this)
1404 CALL l4f_category_log(this(i)%category,L4F_DEBUG, &
1407 CALL delete(this(i))
1412 END SUBROUTINE volgrid6dv_delete
1415 ! Internal method for performing grid to grid computations
1416 SUBROUTINE volgrid6d_transform_compute(this, volgrid6d_in, volgrid6d_out, &
1417 lev_out, var_coord_vol, clone)
1418 TYPE(grid_transform),INTENT(in) :: this ! oggetto di trasformazione per il grigliato
1419 type(volgrid6d), INTENT(in) :: volgrid6d_in ! oggetto da trasformare
1420 type(volgrid6d), INTENT(inout) :: volgrid6d_out ! oggetto trasformato; deve essere completo (init, alloc, alloc_vol)
1421 TYPE(vol7d_level),INTENT(in),OPTIONAL :: lev_out(:) ! vol7d_level object defining target vertical grid, for vertical interpolations
1422 INTEGER,INTENT(in),OPTIONAL :: var_coord_vol ! index of variable defining vertical coordinate values in input volume
1423 LOGICAL,INTENT(in),OPTIONAL :: clone ! se fornito e \c .TRUE., clona i gaid da volgrid6d_in a volgrid6d_out
1425 INTEGER :: ntime, ntimerange, inlevel, onlevel, nvar, &
1426 itime, itimerange, ilevel, ivar, levshift, levused, lvar_coord_vol, spos
1427 REAL,POINTER :: voldatiin(:,:,:), voldatiout(:,:,:), coord_3d_in(:,:,:)
1428 TYPE(vol7d_level) :: output_levtype
1432 call l4f_category_log(volgrid6d_in%category,L4F_DEBUG,"start volgrid6d_transform_compute
")
1440 lvar_coord_vol = optio_i(var_coord_vol)
1442 if (associated(volgrid6d_in%time))then
1443 ntime=size(volgrid6d_in%time)
1444 volgrid6d_out%time=volgrid6d_in%time
1447 if (associated(volgrid6d_in%timerange))then
1448 ntimerange=size(volgrid6d_in%timerange)
1449 volgrid6d_out%timerange=volgrid6d_in%timerange
1452 IF (ASSOCIATED(volgrid6d_in%level))THEN
1453 inlevel=SIZE(volgrid6d_in%level)
1455 IF (PRESENT(lev_out)) THEN
1456 onlevel=SIZE(lev_out)
1457 volgrid6d_out%level=lev_out
1458 ELSE IF (ASSOCIATED(volgrid6d_in%level))THEN
1459 onlevel=SIZE(volgrid6d_in%level)
1460 volgrid6d_out%level=volgrid6d_in%level
1463 if (associated(volgrid6d_in%var))then
1464 nvar=size(volgrid6d_in%var)
1465 volgrid6d_out%var=volgrid6d_in%var
1467 ! allocate once for speed
1468 .NOT.
IF (ASSOCIATED(volgrid6d_in%voldati)) THEN
1469 ALLOCATE(voldatiin(volgrid6d_in%griddim%dim%nx, volgrid6d_in%griddim%dim%ny, &
1472 .NOT.
IF (ASSOCIATED(volgrid6d_out%voldati)) THEN
1473 ALLOCATE(voldatiout(volgrid6d_out%griddim%dim%nx, volgrid6d_out%griddim%dim%ny, &
1477 CALL get_val(this, levshift=levshift, levused=levused)
1479 IF (c_e(lvar_coord_vol)) THEN
1480 CALL get_val(this%trans, output_levtype=output_levtype)
1481 .OR.
IF (output_levtype%level1 == 103 output_levtype%level1 == 108) THEN
1482 spos = firsttrue(volgrid6d_in%level(:) == vol7d_level_new(1))
1484 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1485 'output level '//t2c(output_levtype%level1)// &
1486 ' requested, but height/press of surface not provided in volume')
1488 .NOT..AND..NOT.
IF (c_e(levshift) c_e(levused)) THEN
1489 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1490 'internal inconsistence, levshift and levused undefined when they should be')
1496 ! IF (c_e(var_coord_vol)) THEN
1497 ! IF (ivar == var_coord_vol) CYCLE ! skip coordinate variable in output
1499 DO itimerange=1,ntimerange
1501 ! skip empty columns where possible, improve
1502 .AND.
IF (c_e(levshift) c_e(levused)) THEN
1503 .NOT.
IF (ANY(c_e( &
1504 volgrid6d_in%gaid(levshift+1:levshift+levused,itime,itimerange,ivar) &
1507 DO ilevel = 1, MIN(inlevel,onlevel)
1508 ! if present gaid copy it
1509 .AND..NOT.
IF (c_e(volgrid6d_in%gaid(ilevel,itime,itimerange,ivar)) &
1510 c_e(volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))) THEN
1512 IF (optio_log(clone)) THEN
1513 CALL copy(volgrid6d_in%gaid(ilevel,itime,itimerange,ivar),&
1514 volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))
1516 CALL l4f_category_log(volgrid6d_in%category,L4F_DEBUG, &
1517 "cloning gaid, level
"//t2c(ilevel))
1520 volgrid6d_out%gaid(ilevel,itime,itimerange,ivar) = &
1521 volgrid6d_in%gaid(ilevel,itime,itimerange,ivar)
1525 ! if out level are more, we have to clone anyway
1526 DO ilevel = MIN(inlevel,onlevel) + 1, onlevel
1527 .AND..NOT.
IF (c_e(volgrid6d_in%gaid(inlevel,itime,itimerange,ivar)) &
1528 c_e(volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))) then
1530 CALL copy(volgrid6d_in%gaid(inlevel,itime,itimerange,ivar),&
1531 volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))
1533 CALL l4f_category_log(volgrid6d_in%category,L4F_DEBUG, &
1534 "forced cloning gaid, level
"//t2c(inlevel)//"->
"//t2c(ilevel))
1539 IF (c_e(lvar_coord_vol)) THEN
1540 NULLIFY(coord_3d_in)
1541 CALL volgrid_get_vol_3d(volgrid6d_in, itime, itimerange, lvar_coord_vol, &
1543 IF (c_e(spos)) THEN ! compute difference wrt surface coordinate
1544 IF (spos == 0) THEN ! error condition, set all to missing and goodnight
1545 coord_3d_in(:,:,levshift+1:levshift+levused) = rmiss
1547 DO ilevel = levshift+1, levshift+levused
1548 .AND.
WHERE(c_e(coord_3d_in(:,:,ilevel)) c_e(coord_3d_in(:,:,spos)))
1549 coord_3d_in(:,:,ilevel) = coord_3d_in(:,:,ilevel) - &
1550 coord_3d_in(:,:,spos)
1552 coord_3d_in(:,:,ilevel) = rmiss
1558 CALL volgrid_get_vol_3d(volgrid6d_in, itime, itimerange, ivar, &
1560 IF (ASSOCIATED(volgrid6d_out%voldati)) & ! improve!!!!
1561 CALL volgrid_get_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
1563 IF (c_e(lvar_coord_vol)) THEN
1564 CALL compute(this, voldatiin, voldatiout, convert(volgrid6d_in%var(ivar)), &
1565 coord_3d_in(:,:,levshift+1:levshift+levused)) ! subset coord_3d_in
1567 CALL compute(this, voldatiin, voldatiout, convert(volgrid6d_in%var(ivar)))
1569 CALL volgrid_set_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
1575 IF (c_e(lvar_coord_vol)) THEN
1576 DEALLOCATE(coord_3d_in)
1578 .NOT.
IF (ASSOCIATED(volgrid6d_in%voldati)) THEN
1579 DEALLOCATE(voldatiin)
1581 .NOT.
IF (ASSOCIATED(volgrid6d_out%voldati)) THEN
1582 DEALLOCATE(voldatiout)
1586 END SUBROUTINE volgrid6d_transform_compute
1589 !> Performs the specified abstract transformation on the data provided.
1590 !! The abstract transformation is specified by \a this parameter; the
1591 !! corresponding specifical transformation (\a grid_transform object)
1592 !! is created and destroyed internally. The output transformed object
1593 !! is created internally and it does not require preliminary
1595 SUBROUTINE volgrid6d_transform(this, griddim, volgrid6d_in, volgrid6d_out, &
1596 lev_out, volgrid6d_coord_in, maskgrid, maskbounds, clone, decode, categoryappend)
1597 TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
1598 TYPE(griddim_def),INTENT(in),OPTIONAL :: griddim !< griddim specifying the output grid (required by most transformation types)
1599 TYPE(volgrid6d),INTENT(inout) :: volgrid6d_in !< object to be transformed, it is not modified, despite the INTENT(inout)
1600 TYPE(volgrid6d),INTENT(out) :: volgrid6d_out !< transformed object, it does not require initialisation
1601 TYPE(vol7d_level),INTENT(in),OPTIONAL,TARGET :: lev_out(:) !< vol7d_level object defining target vertical grid, for vertical interpolations
1602 TYPE(volgrid6d),INTENT(in),OPTIONAL :: volgrid6d_coord_in !< object providing time constant input vertical coordinate for some kind of vertical interpolations
1603 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')
1604 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')
1605 LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \a .TRUE. , clone the \a gaid's from \a volgrid6d_in to \a volgrid6d_out
1606 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
1607 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
1609 TYPE(grid_transform) :: grid_trans
1610 TYPE(vol7d_level),POINTER :: llev_out(:)
1611 TYPE(vol7d_level) :: input_levtype, output_levtype
1612 TYPE(vol7d_var) :: vcoord_var
1613 INTEGER :: i, k, ntime, ntimerange, nlevel, nvar, var_coord_in, var_coord_vol, &
1614 cf_out, nxc, nyc, nxi, nyi, i3, i4, i5, i6, &
1615 ulstart, ulend, spos
1616 REAL,ALLOCATABLE :: coord_3d_in(:,:,:)
1617 TYPE(geo_proj) :: proj_in, proj_out
1618 CHARACTER(len=80) :: trans_type
1620 LOGICAL,ALLOCATABLE :: mask_in(:)
1623 call l4f_category_log(volgrid6d_in%category, L4F_DEBUG, "start volgrid6d_transform
")
1631 if (associated(volgrid6d_in%time)) ntime=size(volgrid6d_in%time)
1632 if (associated(volgrid6d_in%timerange)) ntimerange=size(volgrid6d_in%timerange)
1633 if (associated(volgrid6d_in%level)) nlevel=size(volgrid6d_in%level)
1634 if (associated(volgrid6d_in%var)) nvar=size(volgrid6d_in%var)
1636 .OR..OR..OR.
IF (ntime == 0 ntimerange == 0 nlevel == 0 nvar == 0) THEN
1637 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1639 ' ntimerange='//t2c(ntimerange)//' nlevel='//t2c(nlevel)//' nvar='//t2c(nvar))
1640 CALL init(volgrid6d_out) ! initialize to empty
1645 CALL get_val(this, trans_type=trans_type)
1647 ! store desired output component flag and unrotate if necessary
1649 .AND..OR.
IF (PRESENT(griddim) (trans_type == 'inter' trans_type == 'boxinter' &
1650 .OR.
trans_type == 'stencilinter')) THEN ! improve condition!!
1651 CALL get_val(volgrid6d_in%griddim, proj=proj_in)
1652 CALL get_val(griddim, component_flag=cf_out, proj=proj_out)
1653 ! if different projections wind components must be referred to geographical system
1654 IF (proj_in /= proj_out) CALL vg6d_wind_unrot(volgrid6d_in)
1655 ELSE IF (PRESENT(griddim)) THEN ! just get component_flag, the rest is rubbish
1656 CALL get_val(griddim, component_flag=cf_out)
1660 var_coord_in = imiss
1661 var_coord_vol = imiss
1662 IF (trans_type == 'vertint') THEN
1663 IF (PRESENT(lev_out)) THEN
1665 ! if volgrid6d_coord_in provided and allocated, check that it fits
1666 IF (PRESENT(volgrid6d_coord_in)) THEN
1667 IF (ASSOCIATED(volgrid6d_coord_in%voldati)) THEN
1669 ! strictly 1 time and 1 timerange
1670 .OR.
IF (SIZE(volgrid6d_coord_in%voldati,4) /= 1 &
1671 SIZE(volgrid6d_coord_in%voldati,5) /= 1) THEN
1672 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1673 'volume providing constant input vertical coordinate must have &
1674 &only 1 time and 1 timerange')
1675 CALL init(volgrid6d_out) ! initialize to empty
1680 ! search for variable providing vertical coordinate
1681 CALL get_val(this, output_levtype=output_levtype)
1682 vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
1683 .NOT.
IF (c_e(vcoord_var)) THEN
1684 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1685 'requested output level type '//t2c(output_levtype%level1)// &
1686 ' does not correspond to any known physical variable for &
1687 &providing vertical coordinate')
1688 CALL init(volgrid6d_out) ! initialize to empty
1693 DO i = 1, SIZE(volgrid6d_coord_in%var)
1694 IF (convert(volgrid6d_coord_in%var(i)) == vcoord_var) THEN
1700 .NOT.
IF (c_e(var_coord_in)) THEN
1701 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1702 'volume providing constant input vertical coordinate contains no &
1703 &variables matching output level type '//t2c(output_levtype%level1))
1704 CALL init(volgrid6d_out) ! initialize to empty
1708 CALL l4f_category_log(volgrid6d_in%category, L4F_INFO, &
1709 'Coordinate for vertint found in coord volume at position '// &
1712 ! check horizontal grid
1713 ! this is too strict (component flag and so on)
1714 ! IF (volgrid6d_coord_in%griddim /= volgrid6d_in%griddim) THEN
1715 CALL get_val(volgrid6d_coord_in%griddim, nx=nxc, ny=nyc)
1716 CALL get_val(volgrid6d_in%griddim, nx=nxi, ny=nyi)
1717 .OR.
IF (nxc /= nxi nyc /= nyi) THEN
1718 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1719 'volume providing constant input vertical coordinate must have &
1720 &the same grid as the input')
1721 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1722 'vertical coordinate: '//t2c(nxc)//'x'//t2c(nyc)// &
1723 ', input volume: '//t2c(nxi)//'x'//t2c(nyi))
1724 CALL init(volgrid6d_out) ! initialize to empty
1729 ! check vertical coordinate system
1730 CALL get_val(this, input_levtype=input_levtype)
1731 mask_in = & ! implicit allocation
1732 .AND.
(volgrid6d_coord_in%level(:)%level1 == input_levtype%level1) &
1733 (volgrid6d_coord_in%level(:)%level2 == input_levtype%level2)
1734 ulstart = firsttrue(mask_in)
1735 ulend = lasttrue(mask_in)
1736 .OR.
IF (ulstart == 0 ulend == 0) THEN
1737 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1738 'coordinate file does not contain levels of type '// &
1739 t2c(input_levtype%level1)//'/'//t2c(input_levtype%level2)// &
1740 ' specified for input data')
1741 CALL init(volgrid6d_out) ! initialize to empty
1746 coord_3d_in = volgrid6d_coord_in%voldati(:,:,ulstart:ulend,1,1,var_coord_in) ! implicit allocation
1748 .OR.
IF (output_levtype%level1 == 103 &
1749 output_levtype%level1 == 108) THEN ! surface coordinate needed
1750 spos = firsttrue(volgrid6d_coord_in%level(:) == vol7d_level_new(1))
1752 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1753 'output level '//t2c(output_levtype%level1)// &
1754 ' requested, but height/press of surface not provided in coordinate file')
1755 CALL init(volgrid6d_out) ! initialize to empty
1759 DO k = 1, SIZE(coord_3d_in,3)
1760 .AND.
WHERE(c_e(coord_3d_in(:,:,k)) &
1761 c_e(volgrid6d_coord_in%voldati(:,:,spos,1,1,var_coord_in)))
1762 coord_3d_in(:,:,k) = coord_3d_in(:,:,k) - &
1763 volgrid6d_coord_in%voldati(:,:,spos,1,1,var_coord_in)
1765 coord_3d_in(:,:,k) = rmiss
1773 .NOT.
IF (c_e(var_coord_in)) THEN ! search for coordinate within volume
1774 ! search for variable providing vertical coordinate
1775 CALL get_val(this, output_levtype=output_levtype)
1776 vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
1777 IF (c_e(vcoord_var)) THEN
1778 DO i = 1, SIZE(volgrid6d_in%var)
1779 IF (convert(volgrid6d_in%var(i)) == vcoord_var) THEN
1785 IF (c_e(var_coord_vol)) THEN
1786 CALL l4f_category_log(volgrid6d_in%category, L4F_INFO, &
1787 'Coordinate for vertint found in input volume at position '// &
1794 CALL init(volgrid6d_out, griddim=volgrid6d_in%griddim, &
1795 time_definition=volgrid6d_in%time_definition, categoryappend=categoryappend)
1796 IF (c_e(var_coord_in)) THEN
1797 CALL init(grid_trans, this, lev_in=volgrid6d_in%level, lev_out=lev_out, &
1798 coord_3d_in=coord_3d_in, categoryappend=categoryappend)
1800 CALL init(grid_trans, this, lev_in=volgrid6d_in%level, lev_out=lev_out, &
1801 categoryappend=categoryappend)
1804 CALL get_val(grid_trans, output_level_auto=llev_out) ! get levels if auto-generated
1805 .NOT.
IF (ASSOCIATED(llev_out)) llev_out => lev_out
1806 nlevel = SIZE(llev_out)
1808 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1809 'volgrid6d_transform: vertint requested but lev_out not provided')
1810 CALL init(volgrid6d_out) ! initialize to empty
1816 CALL init(volgrid6d_out, griddim=griddim, &
1817 time_definition=volgrid6d_in%time_definition, categoryappend=categoryappend)
1818 CALL init(grid_trans, this, in=volgrid6d_in%griddim, out=volgrid6d_out%griddim, &
1819 maskgrid=maskgrid, maskbounds=maskbounds, categoryappend=categoryappend)
1823 IF (c_e(grid_trans)) THEN ! transformation is valid
1825 CALL volgrid6d_alloc(volgrid6d_out, ntime=ntime, nlevel=nlevel, &
1826 ntimerange=ntimerange, nvar=nvar)
1828 IF (PRESENT(decode)) THEN ! explicitly set decode status
1830 ELSE ! take it from input
1831 ldecode = ASSOCIATED(volgrid6d_in%voldati)
1833 ! force decode if gaid is readonly
1834 decode_loop: DO i6 = 1,nvar
1835 DO i5 = 1, ntimerange
1838 IF (c_e(volgrid6d_in%gaid(i3,i4,i5,i6))) THEN
1839 .OR.
ldecode = ldecode grid_id_readonly(volgrid6d_in%gaid(i3,i4,i5,i6))
1847 IF (PRESENT(decode)) THEN
1848 .NEQV.
IF (ldecodedecode) THEN
1849 CALL l4f_category_log(volgrid6d_in%category, L4F_WARN, &
1850 'volgrid6d_transform: decode status forced to .TRUE. because driver does not allow copy')
1854 CALL volgrid6d_alloc_vol(volgrid6d_out, decode=ldecode)
1856 !ensure unproj was called
1857 !call griddim_unproj(volgrid6d_out%griddim)
1859 IF (trans_type == 'vertint') THEN
1861 CALL l4f_category_log(volgrid6d_in%category, L4F_DEBUG, &
1862 "volgrid6d_transform: vertint to
"//t2c(nlevel)//" levels
")
1864 CALL compute(grid_trans, volgrid6d_in, volgrid6d_out, lev_out=llev_out, &
1865 var_coord_vol=var_coord_vol, clone=clone)
1867 CALL compute(grid_trans, volgrid6d_in, volgrid6d_out, clone=clone)
1870 IF (cf_out == 0) THEN ! unrotated components are desired
1871 CALL wind_unrot(volgrid6d_out) ! unrotate if necessary
1872 ELSE IF (cf_out == 1) THEN ! rotated components are desired
1873 CALL wind_rot(volgrid6d_out) ! rotate if necessary
1877 ! should log with grid_trans%category, but it is private
1878 CALL l4f_category_log(volgrid6d_in%category, L4F_ERROR, &
1879 'volgrid6d_transform: transformation not valid')
1883 CALL delete (grid_trans)
1885 END SUBROUTINE volgrid6d_transform
1888 !> Performs the specified abstract transformation on the arrays of
1889 !! data provided. The abstract transformation is specified by \a this
1890 !! parameter; the corresponding specifical transformation (\a
1891 !! grid_transform object) is created and destroyed internally. The
1892 !! output transformed object is created internally and it does not
1893 !! require preliminary initialisation. According to the input data and
1894 !! to the transformation type, the output array may have of one or
1895 !! more \a volgrid6d elements on different grids.
1896 SUBROUTINE volgrid6dv_transform(this, griddim, volgrid6d_in, volgrid6d_out, &
1897 lev_out, volgrid6d_coord_in, maskgrid, maskbounds, clone, decode, categoryappend)
1898 TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
1899 TYPE(griddim_def),INTENT(in),OPTIONAL :: griddim !< griddim specifying the output grid (required by most transformation types)
1900 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)
1901 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
1902 TYPE(vol7d_level),INTENT(in),OPTIONAL :: lev_out(:) !< vol7d_level object defining target vertical grid
1903 TYPE(volgrid6d),INTENT(in),OPTIONAL :: volgrid6d_coord_in !< object providing time constant input vertical coordinate for some kind of vertical interpolations
1904 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')
1905 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')
1906 LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \a .TRUE. , clone the \a gaid's from \a volgrid6d_in to \a volgrid6d_out
1907 LOGICAL,INTENT(in),OPTIONAL :: decode !< if provided and \a .FALSE. the data volume is not allocated, but work is performed on grid_id's
1908 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
1910 INTEGER :: i, stallo
1913 allocate(volgrid6d_out(size(volgrid6d_in)),stat=stallo)
1914 if (stallo /= 0)then
1915 call l4f_log(L4F_FATAL,"allocating memory
")
1916 call raise_fatal_error()
1919 do i=1,size(volgrid6d_in)
1920 call transform(this, griddim, volgrid6d_in(i), volgrid6d_out(i), &
1921 lev_out=lev_out, volgrid6d_coord_in=volgrid6d_coord_in, &
1922 maskgrid=maskgrid, maskbounds=maskbounds, &
1923 clone=clone, decode=decode, categoryappend=categoryappend)
1926 END SUBROUTINE volgrid6dv_transform
1929 ! Internal method for performing grid to sparse point computations
1930 SUBROUTINE volgrid6d_v7d_transform_compute(this, volgrid6d_in, vol7d_out, &
1931 networkname, noconvert)
1932 TYPE(grid_transform),INTENT(in) :: this ! oggetto di trasformazione per grigliato
1933 type(volgrid6d), INTENT(in) :: volgrid6d_in ! oggetto da trasformare
1934 type(vol7d), INTENT(inout) :: vol7d_out ! oggetto trasformato
1935 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname ! imposta il network in vol7d_out (default='generic')
1936 LOGICAL,OPTIONAL,INTENT(in) :: noconvert !< do not try to match variable and convert values during transform
1938 INTEGER :: nntime, nana, ntime, ntimerange, nlevel, nvar, stallo
1939 INTEGER :: itime, itimerange, ivar, inetwork
1940 REAL,ALLOCATABLE :: voldatir_out(:,:,:)
1941 TYPE(conv_func),POINTER :: c_func(:)
1942 TYPE(datetime),ALLOCATABLE :: validitytime(:,:)
1943 REAL,POINTER :: voldatiin(:,:,:)
1946 call l4f_category_log(volgrid6d_in%category,L4F_DEBUG,"start volgrid6d_v7d_transform_compute
")
1955 if (present(networkname))then
1956 call init(vol7d_out%network(1),name=networkname)
1958 call init(vol7d_out%network(1),name='generic')
1961 if (associated(volgrid6d_in%timerange))then
1962 ntimerange=size(volgrid6d_in%timerange)
1963 vol7d_out%timerange=volgrid6d_in%timerange
1966 if (associated(volgrid6d_in%time))then
1967 ntime=size(volgrid6d_in%time)
1969 if (vol7d_out%time_definition == volgrid6d_in%time_definition) then
1971 ! i time sono definiti uguali: assegno
1972 vol7d_out%time=volgrid6d_in%time
1975 ! converto reference in validity
1976 allocate (validitytime(ntime,ntimerange),stat=stallo)
1978 call l4f_category_log(volgrid6d_in%category,L4F_FATAL,"allocating memory
")
1979 call raise_fatal_error()
1983 do itimerange=1,ntimerange
1984 if (vol7d_out%time_definition > volgrid6d_in%time_definition) then
1985 validitytime(itime,itimerange) = &
1986 volgrid6d_in%time(itime) + timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
1988 validitytime(itime,itimerange) = &
1989 volgrid6d_in%time(itime) - timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
1994 nntime = count_distinct(reshape(validitytime,(/ntime*ntimerange/)), back=.TRUE.)
1995 vol7d_out%time=pack_distinct(reshape(validitytime,(/ntime*ntimerange/)), nntime,back=.TRUE.)
2000 IF (ASSOCIATED(volgrid6d_in%level)) THEN
2001 nlevel = SIZE(volgrid6d_in%level)
2002 vol7d_out%level=volgrid6d_in%level
2005 IF (ASSOCIATED(volgrid6d_in%var)) THEN
2006 nvar = SIZE(volgrid6d_in%var)
2007 .NOT.
IF ( optio_log(noconvert)) THEN
2008 CALL vargrib2varbufr(volgrid6d_in%var, vol7d_out%dativar%r, c_func)
2012 nana = SIZE(vol7d_out%ana)
2014 ! allocate once for speed
2015 .NOT.
IF (ASSOCIATED(volgrid6d_in%voldati)) THEN
2016 ALLOCATE(voldatiin(volgrid6d_in%griddim%dim%nx, volgrid6d_in%griddim%dim%ny, &
2020 ALLOCATE(voldatir_out(nana,1,nlevel),stat=stallo)
2021 IF (stallo /= 0) THEN
2022 CALL l4f_category_log(volgrid6d_in%category,L4F_FATAL,"allocating memory
")
2023 CALL raise_fatal_error()
2028 do itimerange=1,ntimerange
2029 ! do ilevel=1,nlevel
2032 èè
!non chiaro se questa sezione utile o no
2033 !ossia il compute sotto sembra prevedere voldatir_out solo in out
2034 !!$ if (vol7d_out%time_definition == volgrid6d_in%time_definition) then
2035 !!$ voldatir_out=reshape(vol7d_out%voldatir(:,itime,ilevel,itimerange,ivar,inetwork),(/nana,1/))
2037 !!$ voldatir_out=reshape(vol7d_out%voldatir(:,index(vol7d_out%time,validitytime(itime,ilevel)),ilevel,itimerange,ivar,inetwork),(/nana,1/))
2040 CALL volgrid_get_vol_3d(volgrid6d_in, itime, itimerange, ivar, &
2043 CALL compute(this, voldatiin, voldatir_out, vol7d_out%dativar%r(ivar))
2045 if (vol7d_out%time_definition == volgrid6d_in%time_definition) then
2046 vol7d_out%voldatir(:,itime,:,itimerange,ivar,inetwork) = &
2049 vol7d_out%voldatir(:,index(vol7d_out%time,validitytime(itime,itimerange)),:,itimerange,ivar,inetwork)=&
2050 RESHAPE(voldatir_out,(/nana,nlevel/))
2053 ! 1 indice della dimensione "anagrafica
"
2054 ! 2 indice della dimensione "tempo
"
2055 ! 3 indice della dimensione "livello verticale
"
2056 ! 4 indice della dimensione "intervallo temporale
"
2057 ! 5 indice della dimensione "variabile
"
2058 ! 6 indice della dimensione "rete
"
2065 deallocate(voldatir_out)
2066 .NOT.
IF (ASSOCIATED(volgrid6d_in%voldati)) THEN
2067 DEALLOCATE(voldatiin)
2069 if (allocated(validitytime)) deallocate(validitytime)
2071 ! Rescale valid data according to variable conversion table
2072 IF (ASSOCIATED(c_func)) THEN
2074 CALL compute(c_func(ivar), vol7d_out%voldatir(:,:,:,:,ivar,:))
2079 end SUBROUTINE volgrid6d_v7d_transform_compute
2082 !> Performs the specified abstract transformation on the data provided.
2083 !! The abstract transformation is specified by \a this parameter; the
2084 !! corresponding specifical transformation (\a grid_transform object)
2085 !! is created and destroyed internally. The output transformed object
2086 !! is created internally and it does not require preliminary
2088 SUBROUTINE volgrid6d_v7d_transform(this, volgrid6d_in, vol7d_out, v7d, &
2089 maskgrid, maskbounds, networkname, noconvert, find_index, categoryappend)
2090 TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
2091 TYPE(volgrid6d),INTENT(inout) :: volgrid6d_in !< object to be transformed, it is not modified, despite the INTENT(inout)
2092 TYPE(vol7d),INTENT(out) :: vol7d_out !< transformed object, it does not requires initialisation
2093 TYPE(vol7d),INTENT(in),OPTIONAL :: v7d !< object containing a list of points over which transformation has to be done (required by some transformation types)
2094 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')
2095 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
2096 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname !< set the output network name in vol7d_out (default='generic')
2097 LOGICAL,OPTIONAL,INTENT(in) :: noconvert !< do not try to match variable and convert values during transform
2098 PROCEDURE(basic_find_index),POINTER,OPTIONAL :: find_index
2099 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
2101 type(grid_transform) :: grid_trans
2102 INTEGER :: ntime, ntimerange, nlevel, nvar, nana, time_definition, nnetwork, stallo
2103 INTEGER :: itime, itimerange, inetwork
2104 TYPE(datetime),ALLOCATABLE :: validitytime(:,:)
2105 INTEGER,ALLOCATABLE :: point_index(:)
2106 TYPE(vol7d) :: v7d_locana
2109 call l4f_category_log(volgrid6d_in%category,L4F_DEBUG,"start volgrid6d_v7d_transform
")
2112 call vg6d_wind_unrot(volgrid6d_in)
2120 call get_val(this,time_definition=time_definition)
2121 .not.
if ( c_e(time_definition)) then
2122 time_definition=1 ! default to validity time
2125 IF (PRESENT(v7d)) THEN
2126 CALL vol7d_copy(v7d, v7d_locana)
2128 CALL init(v7d_locana)
2131 if (associated(volgrid6d_in%timerange)) ntimerange=size(volgrid6d_in%timerange)
2133 if (associated(volgrid6d_in%time)) then
2135 ntime=size(volgrid6d_in%time)
2137 if (time_definition /= volgrid6d_in%time_definition) then
2139 ! converto reference in validity
2140 allocate (validitytime(ntime,ntimerange),stat=stallo)
2142 call l4f_category_log(volgrid6d_in%category,L4F_FATAL,"allocating memory
")
2143 call raise_fatal_error()
2147 do itimerange=1,ntimerange
2148 if (time_definition > volgrid6d_in%time_definition) then
2149 validitytime(itime,itimerange) = &
2150 volgrid6d_in%time(itime) + timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
2152 validitytime(itime,itimerange) = &
2153 volgrid6d_in%time(itime) - timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
2158 ntime = count_distinct(reshape(validitytime,(/ntime*ntimerange/)), back=.TRUE.)
2159 deallocate (validitytime)
2165 if (associated(volgrid6d_in%level)) nlevel=size(volgrid6d_in%level)
2166 if (associated(volgrid6d_in%var)) nvar=size(volgrid6d_in%var)
2168 CALL init(grid_trans, this, volgrid6d_in%griddim, v7d_locana, &
2169 maskgrid=maskgrid, maskbounds=maskbounds, find_index=find_index, &
2170 categoryappend=categoryappend)
2171 CALL init (vol7d_out,time_definition=time_definition)
2173 IF (c_e(grid_trans)) THEN
2175 nana=SIZE(v7d_locana%ana)
2176 CALL vol7d_alloc(vol7d_out, nana=nana, ntime=ntime, nlevel=nlevel, &
2177 ntimerange=ntimerange, ndativarr=nvar, nnetwork=nnetwork)
2178 vol7d_out%ana = v7d_locana%ana
2180 CALL get_val(grid_trans, output_point_index=point_index)
2181 IF (ALLOCATED(point_index)) THEN
2182 ! check that size(point_index) == nana?
2183 CALL vol7d_alloc(vol7d_out, nanavari=1)
2184 CALL init(vol7d_out%anavar%i(1), 'B01192')
2187 CALL vol7d_alloc_vol(vol7d_out)
2189 IF (ALLOCATED(point_index)) THEN
2190 DO inetwork = 1, nnetwork
2191 vol7d_out%volanai(:,1,inetwork) = point_index(:)
2194 CALL compute(grid_trans, volgrid6d_in, vol7d_out, networkname, noconvert)
2196 CALL l4f_log(L4F_ERROR, 'vg6d_v7d_transform: transformation not valid')
2200 CALL delete(grid_trans)
2203 ! set variables to a conformal state
2204 CALL vol7d_dballe_set_var_du(vol7d_out)
2207 CALL delete(v7d_locana)
2209 END SUBROUTINE volgrid6d_v7d_transform
2212 !> Performs the specified abstract transformation on the arrays of
2213 !! data provided. The abstract transformation is specified by \a this
2214 !! parameter; the corresponding specifical transformation (\a
2215 !! grid_transform object) is created and destroyed internally. The
2216 !! output transformed object is created internally and it does not
2217 !! require preliminary initialisation. The transformation performed on
2218 !! each element of the input \a volgrid6d array object is merged into
2219 !! a single \a vol7d output object.
2220 SUBROUTINE volgrid6dv_v7d_transform(this, volgrid6d_in, vol7d_out, v7d, &
2221 maskgrid, maskbounds, networkname, noconvert, find_index, categoryappend)
2222 TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
2223 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)
2224 TYPE(vol7d),INTENT(out) :: vol7d_out !< transformed object, it does not require initialisation
2225 TYPE(vol7d),INTENT(in),OPTIONAL :: v7d !< object containing a list of points over which transformation has to be done (required by some transformation types)
2226 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')
2227 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
2228 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname !< set the output network name in vol7d_out (default='generic')
2229 LOGICAL,OPTIONAL,INTENT(in) :: noconvert !< do not try to match variable and convert values during transform
2230 PROCEDURE(basic_find_index),POINTER,OPTIONAL :: find_index
2231 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
2234 TYPE(vol7d) :: v7dtmp
2238 CALL init(vol7d_out)
2240 DO i=1,SIZE(volgrid6d_in)
2241 CALL transform(this, volgrid6d_in(i), v7dtmp, v7d=v7d, &
2242 maskgrid=maskgrid, maskbounds=maskbounds, &
2243 networkname=networkname, noconvert=noconvert, find_index=find_index, &
2244 categoryappend=categoryappend)
2245 CALL vol7d_append(vol7d_out, v7dtmp)
2248 END SUBROUTINE volgrid6dv_v7d_transform
2251 ! Internal method for performing sparse point to grid computations
2252 SUBROUTINE v7d_volgrid6d_transform_compute(this, vol7d_in, volgrid6d_out, networkname, gaid_template)
2253 TYPE(grid_transform),INTENT(in) :: this ! object specifying the specific transformation
2254 type(vol7d), INTENT(in) :: vol7d_in ! object to be transformed
2255 type(volgrid6d), INTENT(inout) :: volgrid6d_out ! transformed object
2256 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname ! select the network to be processed from the \a vol7d input object, default the first network
2257 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
2259 integer :: nana, ntime, ntimerange, nlevel, nvar
2260 INTEGER :: ilevel, itime, itimerange, ivar, inetwork
2262 REAL,POINTER :: voldatiout(:,:,:)
2263 type(vol7d_network) :: network
2264 TYPE(conv_func), pointer :: c_func(:)
2265 !TODO category sarebbe da prendere da vol7d
2267 CALL l4f_category_log(volgrid6d_out%category, L4F_DEBUG, &
2268 'start v7d_volgrid6d_transform_compute')
2276 IF (PRESENT(networkname)) THEN
2277 CALL init(network,name=networkname)
2278 inetwork = INDEX(vol7d_in%network,network)
2279 IF (inetwork <= 0) THEN
2280 CALL l4f_category_log(volgrid6d_out%category,L4F_WARN, &
2281 'network '//TRIM(networkname)//' not found, first network will be transformed')
2288 ! no time_definition conversion implemented, output will be the same as input
2289 if (associated(vol7d_in%time))then
2290 ntime=size(vol7d_in%time)
2291 volgrid6d_out%time=vol7d_in%time
2294 if (associated(vol7d_in%timerange))then
2295 ntimerange=size(vol7d_in%timerange)
2296 volgrid6d_out%timerange=vol7d_in%timerange
2299 if (associated(vol7d_in%level))then
2300 nlevel=size(vol7d_in%level)
2301 volgrid6d_out%level=vol7d_in%level
2304 if (associated(vol7d_in%dativar%r))then
2305 nvar=size(vol7d_in%dativar%r)
2306 CALL varbufr2vargrib(vol7d_in%dativar%r, volgrid6d_out%var, c_func, gaid_template)
2309 nana=SIZE(vol7d_in%voldatir, 1)
2310 ! allocate once for speed
2311 .NOT.
IF (ASSOCIATED(volgrid6d_out%voldati)) THEN
2312 ALLOCATE(voldatiout(volgrid6d_out%griddim%dim%nx, volgrid6d_out%griddim%dim%ny, &
2317 DO itimerange=1,ntimerange
2320 ! clone the gaid template where I have data
2321 IF (PRESENT(gaid_template)) THEN
2322 DO ilevel = 1, nlevel
2323 IF (ANY(c_e(vol7d_in%voldatir(:,itime,ilevel,itimerange,ivar,inetwork)))) THEN
2324 CALL copy(gaid_template, volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))
2326 volgrid6d_out%gaid(ilevel,itime,itimerange,ivar) = grid_id_new()
2332 IF (ASSOCIATED(volgrid6d_out%voldati)) & ! improve!!!!
2333 CALL volgrid_get_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
2335 ! do the interpolation
2336 CALL compute(this, &
2337 vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), voldatiout, &
2338 vol7d_in%dativar%r(ivar))
2339 ! rescale valid data according to variable conversion table
2340 IF (ASSOCIATED(c_func)) THEN
2341 CALL compute(c_func(ivar), voldatiout(:,:,:))
2344 CALL volgrid_set_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
2351 .NOT.
IF (ASSOCIATED(volgrid6d_out%voldati)) THEN
2352 DEALLOCATE(voldatiout)
2354 IF (ASSOCIATED(c_func)) THEN
2358 END SUBROUTINE v7d_volgrid6d_transform_compute
2361 !> Performs the specified abstract transformation on the data provided.
2362 !! The abstract transformation is specified by \a this parameter; the
2363 !! corresponding specifical transformation (\a grid_transform object)
2364 !! is created and destroyed internally. The output transformed object
2365 !! is created internally and it does not require preliminary
2367 SUBROUTINE v7d_volgrid6d_transform(this, griddim, vol7d_in, volgrid6d_out, &
2368 networkname, gaid_template, categoryappend)
2369 TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
2370 TYPE(griddim_def),INTENT(in),OPTIONAL :: griddim !< griddim specifying the output grid (required by most transformation types)
2371 ! TODO ripristinare intent(in) dopo le opportune modifiche in grid_class.F90
2372 TYPE(vol7d),INTENT(inout) :: vol7d_in !< object to be transformed, it is not modified, despite the INTENT(inout)
2373 TYPE(volgrid6d),INTENT(out) :: volgrid6d_out !< transformed object, it does not require initialisation
2374 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname !< select the network to be processed from the \a vol7d input object, default the first network
2375 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
2376 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
2378 type(grid_transform) :: grid_trans
2379 integer :: ntime, ntimerange, nlevel, nvar
2382 !TODO la category sarebbe da prendere da vol7d
2383 !call l4f_category_log(vol7d_out%category,L4F_DEBUG,"start volgrid6d_transform
")
2385 CALL vol7d_alloc_vol(vol7d_in) ! be safe
2386 ntime=SIZE(vol7d_in%time)
2387 ntimerange=SIZE(vol7d_in%timerange)
2388 nlevel=SIZE(vol7d_in%level)
2390 if (associated(vol7d_in%dativar%r)) nvar=size(vol7d_in%dativar%r)
2392 IF (nvar <= 0) THEN ! use vol7d category once it will be implemented
2393 CALL l4f_log(L4F_ERROR, &
2394 "trying to
transform a
vol7d object incomplete or without real variables
")
2395 CALL init(volgrid6d_out) ! initialize to empty
2400 CALL init(grid_trans, this, vol7d_in, griddim, categoryappend=categoryappend)
2401 CALL init(volgrid6d_out, griddim, time_definition=vol7d_in%time_definition, &
2402 categoryappend=categoryappend)
2404 IF (c_e(grid_trans)) THEN
2406 CALL volgrid6d_alloc(volgrid6d_out, griddim%dim, ntime=ntime, nlevel=nlevel, &
2407 ntimerange=ntimerange, nvar=nvar)
2408 ! can I avoid decode=.TRUE. here, using gaid_template?
2409 CALL volgrid6d_alloc_vol(volgrid6d_out, decode=.TRUE.)
2411 CALL compute(grid_trans, vol7d_in, volgrid6d_out, networkname, gaid_template)
2413 CALL vg6d_wind_rot(volgrid6d_out)
2415 ! should log with grid_trans%category, but it is private
2416 CALL l4f_log(L4F_ERROR, 'v7d_vg6d_transform: transformation not valid')
2420 CALL delete(grid_trans)
2422 END SUBROUTINE v7d_volgrid6d_transform
2425 ! Internal method for performing sparse point to sparse point computations
2426 SUBROUTINE v7d_v7d_transform_compute(this, vol7d_in, vol7d_out, lev_out, &
2428 TYPE(grid_transform),INTENT(in) :: this ! oggetto di trasformazione per grigliato
2429 type(vol7d), INTENT(in) :: vol7d_in ! oggetto da trasformare
2430 type(vol7d), INTENT(inout) :: vol7d_out ! oggetto trasformato
2431 TYPE(vol7d_level),INTENT(in),OPTIONAL :: lev_out(:) ! vol7d_level object defining target vertical grid, for vertical interpolations
2432 INTEGER,INTENT(in),OPTIONAL :: var_coord_vol ! index of variable defining vertical coordinate values in input volume
2434 INTEGER :: itime, itimerange, ilevel, ivar, inetwork, &
2435 levshift, levused, lvar_coord_vol, spos
2436 REAL,ALLOCATABLE :: coord_3d_in(:,:,:)
2437 TYPE(vol7d_level) :: output_levtype
2439 lvar_coord_vol = optio_i(var_coord_vol)
2440 vol7d_out%time(:) = vol7d_in%time(:)
2441 vol7d_out%timerange(:) = vol7d_in%timerange(:)
2442 IF (PRESENT(lev_out)) THEN
2443 vol7d_out%level(:) = lev_out(:)
2445 vol7d_out%level(:) = vol7d_in%level(:)
2447 vol7d_out%network(:) = vol7d_in%network(:)
2448 IF (ASSOCIATED(vol7d_in%dativar%r)) THEN ! work only when real vars are available
2449 vol7d_out%dativar%r(:) = vol7d_in%dativar%r(:)
2451 CALL get_val(this, levshift=levshift, levused=levused)
2453 IF (c_e(lvar_coord_vol)) THEN
2454 CALL get_val(this%trans, output_levtype=output_levtype)
2455 .OR.
IF (output_levtype%level1 == 103 output_levtype%level1 == 108) THEN
2456 spos = firsttrue(vol7d_in%level(:) == vol7d_level_new(1))
2458 CALL l4f_log(L4F_ERROR, &
2459 'output level '//t2c(output_levtype%level1)// &
2460 ' requested, but height/press of surface not provided in volume')
2462 .NOT..AND..NOT.
IF (c_e(levshift) c_e(levused)) THEN
2463 CALL l4f_log(L4F_ERROR, &
2464 'internal inconsistence, levshift and levused undefined when they should be')
2466 ALLOCATE(coord_3d_in(SIZE(vol7d_in%ana),1,SIZE(vol7d_in%level)))
2471 DO inetwork = 1, SIZE(vol7d_in%network)
2472 DO ivar = 1, SIZE(vol7d_in%dativar%r)
2473 DO itimerange = 1, SIZE(vol7d_in%timerange)
2474 DO itime = 1, SIZE(vol7d_in%time)
2476 ! dirty trick to make voldatir look like a 2d-array of shape (nana,1)
2477 IF (c_e(lvar_coord_vol)) THEN
2478 IF (c_e(spos)) THEN ! compute difference wrt surface coordinate
2479 IF (spos == 0) THEN ! error condition, set all to missing and goodnight
2480 coord_3d_in(:,:,levshift+1:levshift+levused) = rmiss
2482 DO ilevel = levshift+1, levshift+levused
2483 .AND.
WHERE(c_e(vol7d_in%voldatir(:,itime:itime,ilevel,itimerange,lvar_coord_vol,inetwork)) &
2484 c_e(vol7d_in%voldatir(:,itime:itime,spos,itimerange,lvar_coord_vol,inetwork)))
2485 coord_3d_in(:,:,ilevel) = vol7d_in%voldatir(:,itime:itime,ilevel,itimerange,lvar_coord_vol,inetwork) - &
2486 vol7d_in%voldatir(:,itime:itime,spos,itimerange,lvar_coord_vol,inetwork)
2488 coord_3d_in(:,:,ilevel) = rmiss
2492 CALL compute(this, &
2493 vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), &
2494 vol7d_out%voldatir(:,itime:itime,:,itimerange,ivar,inetwork), &
2495 var=vol7d_in%dativar%r(ivar), &
2496 coord_3d_in=coord_3d_in)
2498 CALL compute(this, &
2499 vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), &
2500 vol7d_out%voldatir(:,itime:itime,:,itimerange,ivar,inetwork), &
2501 var=vol7d_in%dativar%r(ivar), &
2502 coord_3d_in=vol7d_in%voldatir(:,itime:itime,:,itimerange, &
2503 lvar_coord_vol,inetwork))
2506 CALL compute(this, &
2507 vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), &
2508 vol7d_out%voldatir(:,itime:itime,:,itimerange,ivar,inetwork), &
2509 var=vol7d_in%dativar%r(ivar))
2518 END SUBROUTINE v7d_v7d_transform_compute
2521 !> Performs the specified abstract transformation on the data provided.
2522 !! The abstract transformation is specified by \a this parameter; the
2523 !! corresponding specifical transformation (\a grid_transform object)
2524 !! is created and destroyed internally. The output transformed object
2525 !! is created internally and it does not require preliminary
2526 !! initialisation. The success of the transformation can be checked
2527 !! with the \a c_e method: c_e(vol7d_out).
2528 SUBROUTINE v7d_v7d_transform(this, vol7d_in, vol7d_out, v7d, maskbounds, &
2529 lev_out, vol7d_coord_in, categoryappend)
2530 TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
2531 TYPE(vol7d),INTENT(inout) :: vol7d_in !< object to be transformed, it is not modified, despite the INTENT(inout)
2532 TYPE(vol7d),INTENT(out) :: vol7d_out !< transformed object, it does not require initialisation
2533 TYPE(vol7d),INTENT(in),OPTIONAL :: v7d !< object containing a list of points over which transformation has to be done (required by some transformation types)
2534 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')
2535 TYPE(vol7d_level),INTENT(in),OPTIONAL,TARGET :: lev_out(:) !< vol7d_level object defining target vertical grid, for vertical interpolations
2536 TYPE(vol7d),INTENT(in),OPTIONAL :: vol7d_coord_in !< object providing time constant input vertical coordinate for some kind of vertical interpolations
2537 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
2539 INTEGER :: nvar, inetwork
2540 TYPE(grid_transform) :: grid_trans
2541 TYPE(vol7d_level),POINTER :: llev_out(:)
2542 TYPE(vol7d_level) :: input_levtype, output_levtype
2543 TYPE(vol7d_var) :: vcoord_var
2544 REAL,ALLOCATABLE :: coord_3d_in(:,:,:)
2545 INTEGER :: var_coord_in, var_coord_vol, i, k, ulstart, ulend, spos
2546 INTEGER,ALLOCATABLE :: point_index(:)
2547 TYPE(vol7d) :: v7d_locana, vol7d_tmpana
2548 CHARACTER(len=80) :: trans_type
2549 LOGICAL,ALLOCATABLE :: mask_in(:), point_mask(:)
2551 CALL vol7d_alloc_vol(vol7d_in) ! be safe
2553 IF (ASSOCIATED(vol7d_in%dativar%r)) nvar=SIZE(vol7d_in%dativar%r)
2555 CALL init(v7d_locana)
2556 IF (PRESENT(v7d)) v7d_locana = v7d
2557 CALL init(vol7d_out, time_definition=vol7d_in%time_definition)
2559 CALL get_val(this, trans_type=trans_type)
2561 var_coord_vol = imiss
2562 IF (trans_type == 'vertint') THEN
2564 IF (PRESENT(lev_out)) THEN
2566 ! if vol7d_coord_in provided and allocated, check that it fits
2568 IF (PRESENT(vol7d_coord_in)) THEN
2569 .AND.
IF (ASSOCIATED(vol7d_coord_in%voldatir) &
2570 ASSOCIATED(vol7d_coord_in%dativar%r)) THEN
2572 ! strictly 1 time, 1 timerange and 1 network
2573 .OR.
IF (SIZE(vol7d_coord_in%voldatir,2) /= 1 &
2574 .OR.
SIZE(vol7d_coord_in%voldatir,4) /= 1 &
2575 SIZE(vol7d_coord_in%voldatir,6) /= 1) THEN
2576 CALL l4f_log(L4F_ERROR, &
2577 'volume providing constant input vertical coordinate must have &
2578 &only 1 time, 1 timerange and 1 network')
2583 ! search for variable providing vertical coordinate
2584 CALL get_val(this, output_levtype=output_levtype)
2585 vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
2586 .NOT.
IF (c_e(vcoord_var)) THEN
2587 CALL l4f_log(L4F_ERROR, &
2588 'requested output level type '//t2c(output_levtype%level1)// &
2589 ' does not correspond to any known physical variable for &
2590 &providing vertical coordinate')
2595 var_coord_in = INDEX(vol7d_coord_in%dativar%r, vcoord_var)
2597 IF (var_coord_in <= 0) THEN
2598 CALL l4f_log(L4F_ERROR, &
2599 'volume providing constant input vertical coordinate contains no &
2600 &real variables matching output level type '//t2c(output_levtype%level1))
2604 CALL l4f_log(L4F_INFO, &
2605 'Coordinate for vertint found in coord volume at position '// &
2608 ! check vertical coordinate system
2609 CALL get_val(this, input_levtype=input_levtype)
2610 mask_in = & ! implicit allocation
2611 .AND.
(vol7d_coord_in%level(:)%level1 == input_levtype%level1) &
2612 (vol7d_coord_in%level(:)%level2 == input_levtype%level2)
2613 ulstart = firsttrue(mask_in)
2614 ulend = lasttrue(mask_in)
2615 .OR.
IF (ulstart == 0 ulend == 0) THEN
2616 CALL l4f_log(L4F_ERROR, &
2617 'coordinate file does not contain levels of type '// &
2618 t2c(input_levtype%level1)//'/'//t2c(input_levtype%level2)// &
2619 ' specified for input data')
2624 coord_3d_in = vol7d_coord_in%voldatir(:,1:1,ulstart:ulend,1,var_coord_in,1) ! dirty 1:1, implicit allocation
2626 IF (output_levtype%level1 == 103 &
2627 .OR.
output_levtype%level1 == 108) THEN ! surface coordinate needed
2628 spos = firsttrue(vol7d_coord_in%level(:) == vol7d_level_new(1))
2630 CALL l4f_log(L4F_ERROR, &
2631 'output level '//t2c(output_levtype%level1)// &
2632 ' requested, but height/press of surface not provided in coordinate file')
2636 DO k = 1, SIZE(coord_3d_in,3)
2637 .AND.
WHERE(c_e(coord_3d_in(:,:,k)) &
2638 c_e(vol7d_coord_in%voldatir(:,1:1,spos,1,var_coord_in,1)))
2639 coord_3d_in(:,:,k) = coord_3d_in(:,:,k) - &
2640 vol7d_coord_in%voldatir(:,1:1,spos,1,var_coord_in,1)
2642 coord_3d_in(:,:,k) = rmiss
2650 IF (var_coord_in <= 0) THEN ! search for coordinate within volume
2651 ! search for variable providing vertical coordinate
2652 CALL get_val(this, output_levtype=output_levtype)
2653 vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
2654 IF (c_e(vcoord_var)) THEN
2655 DO i = 1, SIZE(vol7d_in%dativar%r)
2656 IF (vol7d_in%dativar%r(i) == vcoord_var) THEN
2662 IF (c_e(var_coord_vol)) THEN
2663 CALL l4f_log(L4F_INFO, &
2664 'Coordinate for vertint found in input volume at position '// &
2671 IF (var_coord_in > 0) THEN
2672 CALL init(grid_trans, this, lev_in=vol7d_in%level, lev_out=lev_out, &
2673 coord_3d_in=coord_3d_in, categoryappend=categoryappend)
2675 CALL init(grid_trans, this, lev_in=vol7d_in%level, lev_out=lev_out, &
2676 categoryappend=categoryappend)
2679 CALL get_val(grid_trans, output_level_auto=llev_out) ! get levels if auto-generated
2680 .NOT.
IF (associated(llev_out)) llev_out => lev_out
2682 .AND.
IF (c_e(grid_trans)) THEN! nvar > 0) THEN
2684 CALL vol7d_alloc(vol7d_out, nana=SIZE(vol7d_in%ana), &
2685 ntime=SIZE(vol7d_in%time), ntimerange=SIZE(vol7d_in%timerange), &
2686 nlevel=SIZE(llev_out), nnetwork=SIZE(vol7d_in%network), ndativarr=nvar)
2687 vol7d_out%ana(:) = vol7d_in%ana(:)
2689 CALL vol7d_alloc_vol(vol7d_out)
2691 ! no need to check c_e(var_coord_vol) here since the presence of
2692 ! this%coord_3d_in (external) has precedence over coord_3d_in internal
2693 ! in grid_transform_compute
2694 CALL compute(grid_trans, vol7d_in, vol7d_out, llev_out, &
2695 var_coord_vol=var_coord_vol)
2697 CALL l4f_log(L4F_ERROR, 'v7d_v7d_transform: transformation not valid')
2701 CALL l4f_log(L4F_ERROR, &
2702 'v7d_v7d_transform: vertint requested but lev_out not provided')
2708 CALL init(grid_trans, this, vol7d_in, v7d_locana, maskbounds=maskbounds, &
2709 categoryappend=categoryappend)
2710 ! if this init is successful, I am sure that v7d_locana%ana is associated
2712 .AND.
IF (c_e(grid_trans)) THEN! nvar > 0) THEN
2714 CALL vol7d_alloc(vol7d_out, nana=SIZE(v7d_locana%ana), &
2715 ntime=SIZE(vol7d_in%time), ntimerange=SIZE(vol7d_in%timerange), &
2716 nlevel=SIZE(vol7d_in%level), nnetwork=SIZE(vol7d_in%network), ndativarr=nvar)
2717 vol7d_out%ana = v7d_locana%ana
2719 CALL get_val(grid_trans, point_mask=point_mask, output_point_index=point_index)
2721 IF (ALLOCATED(point_index)) THEN
2722 CALL vol7d_alloc(vol7d_out, nanavari=1)
2723 CALL init(vol7d_out%anavar%i(1), 'B01192')
2726 CALL vol7d_alloc_vol(vol7d_out)
2728 IF (ALLOCATED(point_index)) THEN
2729 DO inetwork = 1, SIZE(vol7d_in%network)
2730 vol7d_out%volanai(:,1,inetwork) = point_index(:)
2733 CALL compute(grid_trans, vol7d_in, vol7d_out)
2735 IF (ALLOCATED(point_mask)) THEN ! keep full ana
2736 IF (SIZE(point_mask) /= SIZE(vol7d_in%ana)) THEN
2737 CALL l4f_log(L4F_WARN, &
2738 'v7d_v7d_transform: inconsistency in point size: '//t2c(SIZE(point_mask)) &
2739 //':'//t2c(SIZE(vol7d_in%ana)))
2742 CALL l4f_log(L4F_DEBUG, 'v7d_v7d_transform: merging ana from in to out')
2744 CALL vol7d_copy(vol7d_in, vol7d_tmpana, &
2745 lana=point_mask, lnetwork=(/.TRUE./), &
2746 ltime=(/.FALSE./), ltimerange=(/.FALSE./), llevel=(/.FALSE./))
2747 CALL vol7d_append(vol7d_out, vol7d_tmpana)
2752 CALL l4f_log(L4F_ERROR, 'v7d_v7d_transform: transformation not valid')
2758 CALL delete (grid_trans)
2759 .NOT.
IF ( PRESENT(v7d)) CALL delete(v7d_locana)
2761 END SUBROUTINE v7d_v7d_transform
2764 !> Unrotate the wind components.
2765 !! It converts u and v components of vector quantities relative to the
2766 !! defined grid in the direction of increasing x and y coordinates to
2767 !! u and v components relative to easterly and notherly direction. The
2768 !! original fields are overwritten.
2769 !! \todo Check and correct wind component flag (to be moved in
2771 subroutine vg6d_wind_unrot(this)
2772 type(volgrid6d) :: this !< object containing wind to be unrotated
2774 integer :: component_flag
2776 call get_val(this%griddim,component_flag=component_flag)
2778 if (component_flag == 1) then
2779 call l4f_category_log(this%category,L4F_INFO, &
2780 "unrotating vector components
")
2781 call vg6d_wind__un_rot(this,.false.)
2782 call set_val(this%griddim,component_flag=0)
2784 call l4f_category_log(this%category,L4F_INFO, &
2785 "no need to unrotate vector components
")
2788 end subroutine vg6d_wind_unrot
2791 !> Rotate the wind components.
2792 !! It converts u and v components of vector quantities
2793 !! relative to easterly and notherly direction to
2794 !! defined grid in the direction of increasing x and y coordinates.
2795 !! The original fields are overwritten.
2796 subroutine vg6d_wind_rot(this)
2797 type(volgrid6d) :: this !< object containing wind to be rotated
2799 integer :: component_flag
2801 call get_val(this%griddim,component_flag=component_flag)
2803 if (component_flag == 0) then
2804 call l4f_category_log(this%category,L4F_INFO, &
2805 "rotating vector components
")
2806 call vg6d_wind__un_rot(this,.true.)
2807 call set_val(this%griddim,component_flag=1)
2809 call l4f_category_log(this%category,L4F_INFO, &
2810 "no need to rotate vector components
")
2813 end subroutine vg6d_wind_rot
2816 ! Generic UnRotate the wind components.
2817 SUBROUTINE vg6d_wind__un_rot(this,rot)
2818 TYPE(volgrid6d) :: this ! object containing wind to be (un)rotated
2819 LOGICAL :: rot ! if .true. rotate else unrotate
2821 INTEGER :: i, j, k, l, a11, a12, a21, a22, stallo
2822 double precision,pointer :: rot_mat(:,:,:)
2823 real,allocatable :: tmp_arr(:,:)
2824 REAL,POINTER :: voldatiu(:,:), voldativ(:,:)
2825 INTEGER,POINTER :: iu(:), iv(:)
2827 .NOT.
IF (ASSOCIATED(this%var)) THEN
2828 CALL l4f_category_log(this%category, L4F_ERROR, &
2829 "trying to unrotate an incomplete
volgrid6d object
")
2830 CALL raise_fatal_error()
2834 CALL volgrid6d_var_hor_comp_index(this%var, iu, iv)
2835 .NOT.
IF (ASSOCIATED(iu)) THEN
2836 CALL l4f_category_log(this%category,L4F_ERROR, &
2837 "unrotation impossible
")
2838 CALL raise_fatal_error()
2842 ! Temporary workspace
2843 ALLOCATE(tmp_arr(this%griddim%dim%nx, this%griddim%dim%ny),stat=stallo)
2844 IF (stallo /= 0) THEN
2845 CALL l4f_category_log(this%category, L4F_FATAL, "allocating memory
")
2846 CALL raise_fatal_error()
2848 ! allocate once for speed
2849 .NOT.
IF (ASSOCIATED(this%voldati)) THEN
2850 ALLOCATE(voldatiu(this%griddim%dim%nx, this%griddim%dim%ny), &
2851 voldativ(this%griddim%dim%nx, this%griddim%dim%ny))
2854 CALL griddim_unproj(this%griddim)
2855 CALL wind_unrot(this%griddim, rot_mat)
2868 DO k = 1, SIZE(this%timerange)
2869 DO j = 1, SIZE(this%time)
2870 DO i = 1, SIZE(this%level)
2872 CALL volgrid_get_vol_2d(this, i, j, k, iu(l), voldatiu)
2873 CALL volgrid_get_vol_2d(this, i, j, k, iv(l), voldativ)
2874 ! convert units forward
2875 ! CALL compute(conv_fwd(iu(l)), voldatiu)
2876 ! CALL compute(conv_fwd(iv(l)), voldativ)
2878 ! multiply wind components by rotation matrix
2879 .AND.
WHERE(voldatiu /= rmiss voldativ /= rmiss)
2880 tmp_arr(:,:) = real(voldatiu(:,:)*rot_mat(:,:,a11) + &
2881 voldativ(:,:)*rot_mat(:,:,a12))
2882 voldativ(:,:) = real(voldatiu(:,:)*rot_mat(:,:,a21) + &
2883 voldativ(:,:)*rot_mat(:,:,a22))
2884 voldatiu(:,:) = tmp_arr(:,:)
2886 ! convert units backward
2887 ! CALL uncompute(conv_fwd(iu(l)), voldatiu)
2888 ! CALL uncompute(conv_fwd(iv(l)), voldativ)
2890 CALL volgrid_set_vol_2d(this, i, j, k, iu(l), voldatiu)
2891 CALL volgrid_set_vol_2d(this, i, j, k, iv(l), voldativ)
2897 .NOT.
IF (ASSOCIATED(this%voldati)) THEN
2898 DEALLOCATE(voldatiu, voldativ)
2900 DEALLOCATE(rot_mat, tmp_arr, iu, iv)
2902 END SUBROUTINE vg6d_wind__un_rot
2905 !!$ try to understand the problem:
2909 !!$ 1) we have only one volume: we have to provide the direction of shift
2910 !!$ compute H and traslate on it
2911 !!$ 2) we have two volumes:
2912 !!$ 1) volume U and volume V: compute H and traslate on it
2913 !!$ 2) volume U/V and volume H : translate U/V on H
2914 !!$ 3) we have tree volumes: translate U and V on H
2917 !!$ 1) do not have U in volume U
2918 !!$ 2) do not have V in volume V
2919 !!$ 3) we have others variables more than U and V in volumes U e V
2922 !!$ so the steps are:
2923 !!$ 1) find the volumes
2924 !!$ 2) define or compute H grid
2925 !!$ 3) trasform the volumes in H
2928 !!$ case 1) for only one vg6d (U or V) is not managed, but
2929 !!$ the not pubblic subroutines will work but you have to know what you want to do
2932 !> Convert grids type C to type A.
2933 !! Use this to interpolate data from grid type C to grid type A
2934 !! Grids type are defined by Arakawa.
2936 !! We need to find these types of area in a vg6d array
2937 !! T area of points with temterature etc.
2938 !! U area of points with u components of winds
2939 !! V area of points with v components of winds
2941 !! this method works if it finds
2943 !! 1) volume U and volume V: compute H and traslate on it
2944 !! 2) volume U/V and volume H : translate U/V on H
2945 !! three volumes: translate U and V on H
2947 !! try to work well on more datasets at once
2948 subroutine vg6d_c2a (this)
2950 TYPE(volgrid6d),INTENT(inout) :: this(:) !< vettor of volumes volgrid6d to elaborate
2952 integer :: ngrid,igrid,jgrid,ugrid,vgrid,tgrid
2953 doubleprecision :: xmin, xmax, ymin, ymax
2954 doubleprecision :: xmin_t, xmax_t, ymin_t, ymax_t
2955 doubleprecision :: step_lon_t,step_lat_t
2956 character(len=80) :: type_t,type
2957 TYPE(griddim_def):: griddim_t
2963 call init(griddim_t)
2965 call get_val(this(igrid)%griddim,xmin=xmin_t, xmax=xmax_t, ymin=ymin_t, ymax=ymax_t,proj_type=type_t)
2966 step_lon_t=(xmax_t-xmin_t)/dble(this(igrid)%griddim%dim%nx-1)
2967 step_lat_t=(ymax_t-ymin_t)/dble(this(igrid)%griddim%dim%ny-1)
2976 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: search u/v/t points:
"//to_char(igrid)//to_char(jgrid))
2977 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: test delta:
"//to_char(step_lon_t)//to_char(step_lat_t))
2980 if (this(igrid)%griddim == this(jgrid)%griddim ) cycle
2982 .and.
if (this(igrid)%griddim%dim%nx == this(jgrid)%griddim%dim%nx &
2983 this(igrid)%griddim%dim%ny == this(jgrid)%griddim%dim%ny ) then
2985 call get_val(this(jgrid)%griddim,xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,proj_type=type)
2987 if (type_t /= type )cycle
2990 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: test u
"//&
2991 to_char(xmin)//to_char(xmax)//to_char(ymin)//to_char(ymax))
2993 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"diff coordinate lon
"//&
2994 to_char(abs(xmin - (xmin_t+step_lon_t/2.d0)))//&
2995 to_char(abs(xmax - (xmax_t+step_lon_t/2.d0))))
2996 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"diff coordinate lat
"//&
2997 to_char(abs(ymin - (ymin_t+step_lat_t/2.d0)))//&
2998 to_char(abs(ymax - (ymax_t+step_lat_t/2.d0))))
3001 .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
3002 .and.
if ( abs(ymin - ymin_t) < 1.d-3 abs(ymax - ymax_t) < 1.d-3 ) then
3005 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: found u
")
3014 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: test v
"//&
3015 to_char(xmin)//to_char(xmax)//to_char(ymin)//to_char(ymax))
3018 .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
3019 .and.
if ( abs(xmin - xmin_t) < 1.d-3 abs(xmax - xmax_t) < 1.d-3 ) then
3022 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: found v
")
3031 ! test if we have U and V
3034 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: test u and v
"//&
3035 to_char(xmin_t)//to_char(xmax_t)//to_char(ymin_t)//to_char(ymax_t)//&
3036 to_char(xmin)//to_char(xmax)//to_char(ymin)//to_char(ymax))
3038 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"uv diff coordinate lon
"//&
3039 to_char(abs(xmin_t - xmin)-step_lon_t/2.d0)//&
3040 to_char(abs(xmax_t - xmax)-step_lon_t/2.d0))
3041 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"uv diff coordinate lat
"//&
3042 to_char(abs(ymin_t - ymin) -step_lat_t/2.d0)//&
3043 to_char(abs(ymax_t - ymax)-step_lat_t/2.d0))
3046 .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
3047 .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
3050 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid: found u and v case up and right
")
3056 call init(griddim_t,xmin=xmin, xmax=xmax, ymin=ymin_t, ymax=ymax_t)
3062 ! abbiamo almeno due volumi: riportiamo U e/o V su T
3063 if (c_e(ugrid)) then
3065 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid u points: interpolate u on t
"//to_char(tgrid)//to_char(ugrid))
3068 call vg6d_c2a_grid(this(ugrid),this(tgrid)%griddim,cgrid=1)
3070 call vg6d_c2a_grid(this(ugrid),griddim_t,cgrid=1)
3072 call vg6d_c2a_mat(this(ugrid),cgrid=1)
3075 if (c_e(vgrid)) then
3077 call l4f_category_log(this(igrid)%category,L4F_DEBUG,"c grid v points: interpolate v on t
"//to_char(tgrid)//to_char(vgrid))
3080 call vg6d_c2a_grid(this(vgrid),this(tgrid)%griddim,cgrid=2)
3082 call vg6d_c2a_grid(this(vgrid),griddim_t,cgrid=2)
3084 call vg6d_c2a_mat(this(vgrid),cgrid=2)
3089 call delete(griddim_t)
3094 end subroutine vg6d_c2a
3097 ! Convert C grid to A grid
3098 subroutine vg6d_c2a_grid(this,griddim_t,cgrid)
3100 type(volgrid6d),intent(inout) :: this ! object containing fields to be translated (U or V points)
3101 type(griddim_def),intent(in),optional :: griddim_t ! object containing grid of T points
3102 integer,intent(in) :: cgrid ! in C grid (Arakawa) we have 0=T,1=U,2=V points
3104 doubleprecision :: xmin, xmax, ymin, ymax
3105 doubleprecision :: step_lon,step_lat
3108 if (present(griddim_t)) then
3110 call get_val(griddim_t,xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax)
3111 call set_val(this%griddim,xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax)
3112 ! improve grid definition precision
3113 CALL griddim_setsteps(this%griddim)
3122 call l4f_category_log(this%category,L4F_DEBUG,"c grid: t points, nothing to do
")
3129 call l4f_category_log(this%category,L4F_DEBUG,"c grid: u points, we need interpolation
")
3132 call get_val(this%griddim, xmin=xmin, xmax=xmax)
3133 step_lon=(xmax-xmin)/dble(this%griddim%dim%nx-1)
3134 xmin=xmin-step_lon/2.d0
3135 xmax=xmax-step_lon/2.d0
3136 call set_val(this%griddim, xmin=xmin, xmax=xmax)
3137 ! improve grid definition precision
3138 CALL griddim_setsteps(this%griddim)
3143 call l4f_category_log(this%category,L4F_DEBUG,"c grid: v points, we need interpolation
")
3146 call get_val(this%griddim, ymin=ymin, ymax=ymax)
3147 step_lat=(ymax-ymin)/dble(this%griddim%dim%ny-1)
3148 ymin=ymin-step_lat/2.d0
3149 ymax=ymax-step_lat/2.d0
3150 call set_val(this%griddim, ymin=ymin, ymax=ymax)
3151 ! improve grid definition precision
3152 CALL griddim_setsteps(this%griddim)
3156 call l4f_category_log(this%category,L4F_FATAL,"c grid type not known
")
3157 call raise_fatal_error ()
3164 call griddim_unproj(this%griddim)
3167 end subroutine vg6d_c2a_grid
3169 ! Convert C grid to A grid
3170 subroutine vg6d_c2a_mat(this,cgrid)
3172 type(volgrid6d),intent(inout) :: this ! object containing fields to be translated
3173 integer,intent(in) :: cgrid ! in C grid (Arakawa) we have 0=T,1=U,2=V points
3175 INTEGER :: i, j, k, iv, stallo
3176 REAL,ALLOCATABLE :: tmp_arr(:,:)
3177 REAL,POINTER :: voldatiuv(:,:)
3180 IF (cgrid == 0) RETURN ! nothing to do
3181 .AND.
IF (cgrid /= 1 cgrid /= 2) THEN ! wrong cgrid
3182 CALL l4f_category_log(this%category,L4F_FATAL,"c grid type
"// &
3183 TRIM(to_char(cgrid))//" not known
")
3188 ! Temporary workspace
3189 ALLOCATE(tmp_arr(this%griddim%dim%nx, this%griddim%dim%ny),stat=stallo)
3191 call l4f_log(L4F_FATAL,"allocating memory
")
3192 call raise_fatal_error()
3195 ! allocate once for speed
3196 .NOT.
IF (ASSOCIATED(this%voldati)) THEN
3197 ALLOCATE(voldatiuv(this%griddim%dim%nx, this%griddim%dim%ny), stat=stallo)
3198 IF (stallo /= 0) THEN
3199 CALL l4f_log(L4F_FATAL,"allocating memory
")
3200 CALL raise_fatal_error()
3204 IF (cgrid == 1) THEN ! U points to H points
3205 DO iv = 1, SIZE(this%var)
3206 DO k = 1, SIZE(this%timerange)
3207 DO j = 1, SIZE(this%time)
3208 DO i = 1, SIZE(this%level)
3209 tmp_arr(:,:) = rmiss
3210 CALL volgrid_get_vol_2d(this, i, j, k, iv, voldatiuv)
3212 ! West boundary extrapolation
3213 .AND.
WHERE(voldatiuv(1,:) /= rmiss voldatiuv(2,:) /= rmiss)
3214 tmp_arr(1,:) = voldatiuv(1,:) - (voldatiuv(2,:) - voldatiuv(1,:)) / 2.
3218 .AND.
WHERE(voldatiuv(1:this%griddim%dim%nx-1,:) /= rmiss &
3219 voldatiuv(2:this%griddim%dim%nx,:) /= rmiss)
3220 tmp_arr(2:this%griddim%dim%nx,:) = &
3221 (voldatiuv(1:this%griddim%dim%nx-1,:) + &
3222 voldatiuv(2:this%griddim%dim%nx,:)) / 2.
3225 voldatiuv(:,:) = tmp_arr
3226 CALL volgrid_set_vol_2d(this, i, j, k, iv, voldatiuv)
3232 ELSE IF (cgrid == 2) THEN ! V points to H points
3233 DO iv = 1, SIZE(this%var)
3234 DO k = 1, SIZE(this%timerange)
3235 DO j = 1, SIZE(this%time)
3236 DO i = 1, SIZE(this%level)
3237 tmp_arr(:,:) = rmiss
3238 CALL volgrid_get_vol_2d(this, i, j, k, iv, voldatiuv)
3240 ! South boundary extrapolation
3241 .AND.
WHERE(voldatiuv(:,1) /= rmiss voldatiuv(:,2) /= rmiss)
3242 tmp_arr(:,1) = voldatiuv(:,1) - (voldatiuv(:,2) - voldatiuv(:,1)) / 2.
3246 .AND.
WHERE(voldatiuv(:,1:this%griddim%dim%ny-1) /= rmiss &
3247 voldatiuv(:,2:this%griddim%dim%ny) /= rmiss)
3248 tmp_arr(:,2:this%griddim%dim%ny) = &
3249 (voldatiuv(:,1:this%griddim%dim%ny-1) + &
3250 voldatiuv(:,2:this%griddim%dim%ny)) / 2.
3253 voldatiuv(:,:) = tmp_arr
3254 CALL volgrid_set_vol_2d(this, i, j, k, iv, voldatiuv)
3259 ENDIF ! tertium non datur
3261 .NOT.
IF (ASSOCIATED(this%voldati)) THEN
3262 DEALLOCATE(voldatiuv)
3264 DEALLOCATE (tmp_arr)
3266 end subroutine vg6d_c2a_mat
3269 !> \brief Display object on screen
3271 !! Show brief content on screen.
3272 subroutine display_volgrid6d (this)
3274 TYPE(volgrid6d),intent(in) :: this !< object to display
3278 call l4f_category_log(this%category,L4F_DEBUG,"ora mostro gridinfo
" )
3282 call display(this%griddim)
3284 IF (ASSOCIATED(this%time))then
3285 print*,"---- time vector ----
"
3286 print*,"elements=
",size(this%time)
3287 do i=1, size(this%time)
3288 call display(this%time(i))
3292 IF (ASSOCIATED(this%timerange))then
3293 print*,"---- timerange vector ----
"
3294 print*,"elements=
",size(this%timerange)
3295 do i=1, size(this%timerange)
3296 call display(this%timerange(i))
3300 IF (ASSOCIATED(this%level))then
3301 print*,"---- level vector ----
"
3302 print*,"elements=
",size(this%level)
3303 do i=1, size(this%level)
3304 call display(this%level(i))
3308 IF (ASSOCIATED(this%var))then
3309 print*,"---- var vector ----
"
3310 print*,"elements=
",size(this%var)
3311 do i=1, size(this%var)
3312 call display(this%var(i))
3316 IF (ASSOCIATED(this%gaid))then
3317 print*,"---- gaid vector (present mask only) ----
"
3318 print*,"elements=
",shape(this%gaid)
3319 PRINT*,c_e(RESHAPE(this%gaid,(/SIZE(this%gaid)/)))
3322 print*,"--------------------------------------------------------------
"
3325 end subroutine display_volgrid6d
3328 !> \brief Display vector of object on screen
3330 !! Show brief content on screen.
3331 subroutine display_volgrid6dv (this)
3333 TYPE(volgrid6d),intent(in) :: this(:) !< vector of object to display
3336 print*,"-----------------------
volgrid6d vector ---------------------
"
3338 print*,"elements=
",size(this)
3343 call l4f_category_log(this(i)%category,L4F_DEBUG,"ora mostro il vettore
volgrid6d" )
3346 call display(this(i))
3349 print*,"--------------------------------------------------------------
"
3351 end subroutine display_volgrid6dv
3354 !> Reduce some dimensions (level and timerage) for semplification (rounding).
3355 !! Vector version of vg6dv_rounding
3356 subroutine vg6dv_rounding(vg6din,vg6dout,level,timerange,nostatproc,merge)
3357 type(volgrid6d),intent(in) :: vg6din(:) !< input volume
3358 type(volgrid6d),intent(out),pointer :: vg6dout(:) !> output volume
3359 type(vol7d_level),intent(in),optional :: level(:) !< almost equal level list
3360 type(vol7d_timerange),intent(in),optional :: timerange(:) !< almost equal timerange list
3361 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
3362 logical,intent(in),optional :: nostatproc !< do not take in account statistical processing code in timerange and P2
3366 allocate(vg6dout(size(vg6din)))
3368 do i = 1, size(vg6din)
3369 call vg6d_rounding(vg6din(i),vg6dout(i),level,timerange,nostatproc,merge)
3372 end subroutine vg6dv_rounding
3374 !> Reduce some dimensions (level and timerage) for semplification (rounding).
3375 !! You can use this for simplify and use variables in computation like alchimia
3376 !! where fields have to be on the same coordinate
3378 !! means in time for short periods and istantaneous values
3379 !! 2 meter and surface levels
3380 !! If there are data on more then one almost equal levels or timeranges, the first var present (at least one point)
3381 !! will be taken (order is by icreasing var index).
3382 !! You can use predefined values for classic semplification
3383 !! almost_equal_levels and almost_equal_timeranges
3384 !! The level or timerange in output will be defined by the first element of level and timerange list
3385 subroutine vg6d_rounding(vg6din,vg6dout,level,timerange,nostatproc,merge)
3386 type(volgrid6d),intent(in) :: vg6din !< input volume
3387 type(volgrid6d),intent(out) :: vg6dout !> output volume
3388 type(vol7d_level),intent(in),optional :: level(:) !< almost equal level list
3389 type(vol7d_timerange),intent(in),optional :: timerange(:) !< almost equal timerange list
3390 logical,intent(in),optional :: merge !< if there are data on more then one almost equal levels or timeranges
3391 !! will be merged POINT BY POINT with priority for the fird data found ordered by icreasing var index
3392 logical,intent(in),optional :: nostatproc !< do not take in account statistical processing code in timerange and P2
3394 integer :: ilevel,itimerange
3395 type(vol7d_level) :: roundlevel(size(vg6din%level))
3396 type(vol7d_timerange) :: roundtimerange(size(vg6din%timerange))
3398 roundlevel=vg6din%level
3400 if (present(level))then
3401 do ilevel = 1, size(vg6din%level)
3402 if ((any(vg6din%level(ilevel) .almosteq. level))) then
3403 roundlevel(ilevel)=level(1)
3408 roundtimerange=vg6din%timerange
3410 if (present(timerange))then
3411 do itimerange = 1, size(vg6din%timerange)
3412 if ((any(vg6din%timerange(itimerange) .almosteq. timerange))) then
3413 roundtimerange(itimerange)=timerange(1)
3418 !set istantaneous values everywere
3419 !preserve p1 for forecast time
3420 if (optio_log(nostatproc)) then
3421 roundtimerange(:)%timerange=254
3422 roundtimerange(:)%p2=imiss
3426 call vg6d_reduce(vg6din,vg6dout,roundlevel,roundtimerange,merge)
3428 end subroutine vg6d_rounding
3430 !> Reduce some dimensions (level and timerage).
3431 !! You can pass a volume and specify duplicated levels and timeranges in roundlevel and roundtimerange;
3432 !! you get unique levels and timeranges in output.
3433 !! If there are data on equal levels or timeranges, the first var present (at least one point)
3434 !! will be taken (order is by icreasing var index).
3435 !! you can specify merge and if there are data on equal levels or timeranges will be merged POINT BY POINT
3436 !! with priority for the first data found ordered by icreasing var index (require to decode all the data)
3437 !! Data are decoded only if needed so the output should be with or without voldata allocated
3438 subroutine vg6d_reduce(vg6din,vg6dout,roundlevel,roundtimerange,merge)
3439 type(volgrid6d),intent(in) :: vg6din !< input volume
3440 type(volgrid6d),intent(out) :: vg6dout !< output volume
3441 type(vol7d_level),intent(in) :: roundlevel(:) !< new level list to use for rounding
3442 type(vol7d_timerange),intent(in) :: roundtimerange(:) !< new timerange list to use for rounding
3443 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)
3445 integer :: nlevel,ntime,ntimerange,nvar,ilevel,itimerange,ivar,indl,indt,itime,nx,ny
3446 real,allocatable :: vol2d(:,:)
3448 nx=vg6din%griddim%dim%nx
3449 ny=vg6din%griddim%dim%ny
3450 nlevel=count_distinct(roundlevel,back=.true.)
3451 ntime=size(vg6din%time)
3452 ntimerange=count_distinct(roundtimerange,back=.true.)
3453 nvar=size(vg6din%var)
3455 call init(vg6dout, vg6din%griddim, vg6din%time_definition, categoryappend="generated by vg6d_reduce
")
3456 call volgrid6d_alloc(vg6dout, vg6din%griddim%dim, ntime, nlevel, ntimerange, nvar)
3458 .or.
if ( ASSOCIATED(vg6din%voldati) optio_log(merge)) then
3459 call volgrid6d_alloc_vol(vg6dout,inivol=.true.,decode=.true.)
3460 allocate(vol2d(nx,ny))
3462 call volgrid6d_alloc_vol(vg6dout,inivol=.true.,decode=.false.)
3465 vg6dout%time=vg6din%time
3466 vg6dout%var=vg6din%var
3467 vg6dout%timerange=pack_distinct(roundtimerange,ntimerange,back=.true.)
3468 vg6dout%level=pack_distinct(roundlevel,nlevel,back=.true.)
3469 ! sort modified dimensions
3470 CALL sort(vg6dout%timerange)
3471 CALL sort(vg6dout%level)
3473 do ilevel=1,size(vg6din%level)
3474 indl=index(vg6dout%level,roundlevel(ilevel))
3475 do itimerange=1,size(vg6din%timerange)
3476 indt=index(vg6dout%timerange,roundtimerange(itimerange))
3480 if ( ASSOCIATED(vg6din%voldati)) then
3481 vol2d=vg6din%voldati(:,:,ilevel,itime,itimerange,ivar)
3484 if (optio_log(merge)) then
3486 .not.
if ( ASSOCIATED(vg6din%voldati)) then
3487 CALL grid_id_decode_data(vg6din%gaid(ilevel,itime,itimerange,ivar), vol2d)
3490 !! merge present data point by point
3491 .not.
where ( c_e(vg6dout%voldati(:,:,indl,itime,indt,ivar)))
3493 vg6dout%voldati(:,:,indl,itime,indt,ivar)=vol2d
3496 else if ( ASSOCIATED(vg6din%voldati)) then
3497 .not.
if ( any(c_e(vg6dout%voldati(:,:,indl,itime,indt,ivar))))then
3498 vg6dout%voldati(:,:,indl,itime,indt,ivar)=vol2d
3502 .and..not.
if (c_e(vg6din%gaid(ilevel,itime,itimerange,ivar)) c_e(vg6dout%gaid(indl,itime,indt,ivar)))then
3503 call copy (vg6din%gaid(ilevel,itime,itimerange,ivar), vg6dout%gaid(indl,itime,indt,ivar))
3510 .or.
if ( ASSOCIATED(vg6din%voldati) optio_log(merge)) then
3514 end subroutine vg6d_reduce
3517 end module volgrid6d_class
3521 !>\example example_vg6d_3.f90
3522 !!\brief Programma esempio semplice per gridinfo e volgrid6d.
3524 !! 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.
3526 !>\example example_vg6d_5.f90
3527 !!\brief Programma trasformazione da volgrid6d a volgrid6d
3529 !! Legge grib da un file e li organizza in un vettore di strutture
3530 !! volgrid6d mettendoli a disposizione per eventuali elaborazioni;
3531 !! vengono poi riesportati a un file grib
3533 !>\example example_vg6d_8.f90
3534 !! \brief Programma scrittura su file vettore di anagrafica
3536 !>\example example_vg6d_6.f90
3537 !! \brief Programma trasformazione da volgrid6d a vol7d
3539 !>\example example_vg6d_7.f90
3540 !! \brief Programma trasformazione da vol7d a volgrid6d
3542 !>\example example_vg6d_9.f90
3543 !! \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...
Object describing a rectangular, homogeneous gridded dataset.