libsim  Versione6.3.0
volgrid6d_class.F90
1 ! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
2 ! authors:
3 ! Davide Cesari <dcesari@arpa.emr.it>
4 ! Paolo Patruno <ppatruno@arpa.emr.it>
5 
6 ! This program is free software; you can redistribute it and/or
7 ! modify it under the terms of the GNU General Public License as
8 ! published by the Free Software Foundation; either version 2 of
9 ! the License, or (at your option) any later version.
10 
11 ! This program is distributed in the hope that it will be useful,
12 ! but WITHOUT ANY WARRANTY; without even the implied warranty of
13 ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 ! GNU General Public License for more details.
15 
16 ! You should have received a copy of the GNU General Public License
17 ! along with this program. If not, see <http://www.gnu.org/licenses/>.
18 #include "config.h"
19 
31 
52 MODULE volgrid6d_class
55 USE vol7d_class
57 USE geo_proj_class
58 USE grid_class
62 USE log4fortran
64 USE grid_id_class
67 !USE file_utilities
68 #ifdef HAVE_DBALLE
70 #endif
71 IMPLICIT NONE
72 
73 character (len=255),parameter:: subcategory="volgrid6d_class"
74 
76 type volgrid6d
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
86 end type volgrid6d
87 
89 INTERFACE init
90  MODULE PROCEDURE volgrid6d_init
91 END INTERFACE
92 
95 INTERFACE delete
96  MODULE PROCEDURE volgrid6d_delete, volgrid6dv_delete
97 END INTERFACE
98 
101 INTERFACE import
102  MODULE PROCEDURE volgrid6d_read_from_file
103  MODULE PROCEDURE import_from_gridinfo, import_from_gridinfovv, &
104  volgrid6d_import_from_file
105 END INTERFACE
106 
109 INTERFACE export
110  MODULE PROCEDURE volgrid6d_write_on_file
111  MODULE PROCEDURE export_to_gridinfo, export_to_gridinfov, export_to_gridinfovv,&
112  volgrid6d_export_to_file
113 END INTERFACE
114 
115 ! methods for computing transformations through an initialised
116 ! grid_transform object, probably too low level to be interfaced
117 INTERFACE compute
118  MODULE PROCEDURE volgrid6d_transform_compute, volgrid6d_v7d_transform_compute,&
119  v7d_volgrid6d_transform_compute, v7d_v7d_transform_compute
120 END INTERFACE
121 
124 INTERFACE transform
125  MODULE PROCEDURE volgrid6d_transform, volgrid6dv_transform,&
126  volgrid6d_v7d_transform, volgrid6dv_v7d_transform, v7d_volgrid6d_transform, &
127  v7d_v7d_transform
128 END INTERFACE
129 
130 INTERFACE wind_rot
131  MODULE PROCEDURE vg6d_wind_rot
132 END INTERFACE
133 
134 INTERFACE wind_unrot
135  MODULE PROCEDURE vg6d_wind_unrot
136 END INTERFACE
137 
140 INTERFACE display
141  MODULE PROCEDURE display_volgrid6d,display_volgrid6dv
142 END INTERFACE
143 
155 INTERFACE rounding
156  MODULE PROCEDURE vg6d_rounding, vg6dv_rounding
157 END INTERFACE
158 
159 private
160 
161 PUBLIC volgrid6d,init,delete,export,import,compute,transform, &
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
164 PUBLIC rounding, vg6d_reduce
165 
166 CONTAINS
167 
168 
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
178 
179 character(len=512) :: a_name
180 
181 if (present(categoryappend))then
182  call l4f_launcher(a_name,a_name_append=trim(subcategory)//"."//trim(categoryappend))
183 else
184  call l4f_launcher(a_name,a_name_append=trim(subcategory))
185 endif
186 this%category=l4f_category_get(a_name)
187 
188 #ifdef DEBUG
189 call l4f_category_log(this%category,l4f_debug,"init")
190 #endif
191 
192 call init(this%griddim)
193 
194 if (present(griddim))then
195  call copy(griddim,this%griddim)
196 end if
197 
198  ! call init(this%time)
199  ! call init(this%timerange)
200  ! call init(this%level)
201  ! call init(this%var)
202 
203 if(present(time_definition)) then
204  this%time_definition = time_definition
205 else
206  this%time_definition = 0 !default to reference time
207 end if
208 
209 nullify(this%time,this%timerange,this%level,this%var)
210 nullify(this%gaid,this%voldati)
211 
212 END SUBROUTINE volgrid6d_init
213 
214 
225 SUBROUTINE volgrid6d_alloc(this, dim, ntime, nlevel, ntimerange, nvar, ini)
226 TYPE(volgrid6d),INTENT(inout) :: this
227 type(grid_dim),INTENT(in),OPTIONAL :: dim
228 INTEGER,INTENT(in),OPTIONAL :: ntime
229 INTEGER,INTENT(in),OPTIONAL :: nlevel
230 INTEGER,INTENT(in),OPTIONAL :: ntimerange
231 INTEGER,INTENT(in),OPTIONAL :: nvar
232 LOGICAL,INTENT(in),OPTIONAL :: ini
233 
234 INTEGER :: i, stallo
235 LOGICAL :: linit
236 
237 #ifdef DEBUG
238 call l4f_category_log(this%category,l4f_debug,"alloc")
239 #endif
240 
241 IF (present(ini)) THEN
242  linit = ini
243 ELSE
244  linit = .false.
245 ENDIF
247 
248 if (present(dim)) call copy(dim,this%griddim%dim)
249 
250 
251 IF (present(ntime)) THEN
252  IF (ntime >= 0) THEN
253  IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
254 #ifdef DEBUG
255  call l4f_category_log(this%category,l4f_debug,"alloc ntime "//to_char(ntime))
256 #endif
257  ALLOCATE(this%time(ntime),stat=stallo)
258  if (stallo /=0)then
259  call l4f_category_log(this%category,l4f_fatal,"allocating memory")
260  CALL raise_fatal_error()
261  end if
262  IF (linit) THEN
263  DO i = 1, ntime
264  this%time(i) = datetime_miss
265 ! CALL init(this%time(i)) ! senza argomento inizializza a zero non missing
266  ! baco di datetime?
267  ENDDO
268  ENDIF
269  ENDIF
270 ENDIF
271 IF (present(nlevel)) THEN
272  IF (nlevel >= 0) THEN
273  IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
274 #ifdef DEBUG
275  call l4f_category_log(this%category,l4f_debug,"alloc nlevel "//to_char(nlevel))
276 #endif
277  ALLOCATE(this%level(nlevel),stat=stallo)
278  if (stallo /=0)then
279  call l4f_category_log(this%category,l4f_fatal,"allocating memory")
280  CALL raise_fatal_error()
281  end if
282  IF (linit) THEN
283  DO i = 1, nlevel
284  CALL init(this%level(i))
285  ENDDO
286  ENDIF
287  ENDIF
288 ENDIF
289 IF (present(ntimerange)) THEN
290  IF (ntimerange >= 0) THEN
291  IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
292 #ifdef DEBUG
293  call l4f_category_log(this%category,l4f_debug,"alloc ntimerange "//to_char(ntimerange))
294 #endif
295  ALLOCATE(this%timerange(ntimerange),stat=stallo)
296  if (stallo /=0)then
297  call l4f_category_log(this%category,l4f_fatal,"allocating memory")
298  CALL raise_fatal_error()
299  end if
300  IF (linit) THEN
301  DO i = 1, ntimerange
302  CALL init(this%timerange(i))
303  ENDDO
304  ENDIF
305  ENDIF
306 ENDIF
307 IF (present(nvar)) THEN
308  IF (nvar >= 0) THEN
309  IF (ASSOCIATED(this%var)) DEALLOCATE(this%var)
310 #ifdef DEBUG
311  call l4f_category_log(this%category,l4f_debug,"alloc nvar "//to_char(nvar))
312 #endif
313  ALLOCATE(this%var(nvar),stat=stallo)
314  if (stallo /=0)then
315  call l4f_category_log(this%category,l4f_fatal,"allocating memory")
316  CALL raise_fatal_error()
317  end if
318  IF (linit) THEN
319  DO i = 1, nvar
320  CALL init(this%var(i))
321  ENDDO
322  ENDIF
323  ENDIF
324 ENDIF
325 
326 end SUBROUTINE volgrid6d_alloc
327 
328 
337 SUBROUTINE volgrid6d_alloc_vol(this, ini, inivol, decode)
338 TYPE(volgrid6d),INTENT(inout) :: this
339 LOGICAL,INTENT(in),OPTIONAL :: ini
340 LOGICAL,INTENT(in),OPTIONAL :: inivol
341 LOGICAL,INTENT(in),OPTIONAL :: decode
342 
343 INTEGER :: stallo
344 LOGICAL :: linivol
345 
346 #ifdef DEBUG
347 call l4f_category_log(this%category,l4f_debug,"start alloc_vol")
348 #endif
350 IF (present(inivol)) THEN ! opposite condition, cannot use optio_log
351  linivol = inivol
352 ELSE
353  linivol = .true.
354 ENDIF
355 
356 IF (this%griddim%dim%nx > 0 .AND. this%griddim%dim%ny > 0) THEN
357 ! allocate dimension descriptors to a minimal size for having a
358 ! non-null volume
359  IF (.NOT. ASSOCIATED(this%var)) CALL volgrid6d_alloc(this, nvar=1, ini=ini)
360  IF (.NOT. ASSOCIATED(this%time)) CALL volgrid6d_alloc(this, ntime=1, ini=ini)
361  IF (.NOT. ASSOCIATED(this%level)) CALL volgrid6d_alloc(this, nlevel=1, ini=ini)
362  IF (.NOT. ASSOCIATED(this%timerange)) CALL volgrid6d_alloc(this, ntimerange=1, ini=ini)
363 
364  IF (optio_log(decode)) THEN
365  IF (.NOT.ASSOCIATED(this%voldati)) THEN
366 #ifdef DEBUG
367  CALL l4f_category_log(this%category,l4f_debug,"alloc voldati volume")
368 #endif
369 
370  ALLOCATE(this%voldati(this%griddim%dim%nx,this%griddim%dim%ny,&
371  SIZE(this%level), SIZE(this%time), &
372  SIZE(this%timerange), SIZE(this%var)),stat=stallo)
373  IF (stallo /= 0)THEN
374  CALL l4f_category_log(this%category,l4f_fatal,"allocating memory")
375  CALL raise_fatal_error()
376  ENDIF
377 
378 ! this is not really needed if we can check other flags for a full
379 ! field missing values
380  IF (linivol) this%voldati = rmiss
381  this%voldati = rmiss
382  ENDIF
383  ENDIF
384 
385  IF (.NOT.ASSOCIATED(this%gaid)) THEN
386 #ifdef DEBUG
387  CALL l4f_category_log(this%category,l4f_debug,"alloc gaid volume")
388 #endif
389  ALLOCATE(this%gaid(SIZE(this%level), SIZE(this%time), &
390  SIZE(this%timerange), SIZE(this%var)),stat=stallo)
391  IF (stallo /= 0)THEN
392  CALL l4f_category_log(this%category,l4f_fatal,"allocating memory")
393  CALL raise_fatal_error()
394  ENDIF
395 
396  IF (linivol) THEN
397 !!$ DO i=1 ,SIZE(this%gaid,1)
398 !!$ DO ii=1 ,SIZE(this%gaid,2)
399 !!$ DO iii=1 ,SIZE(this%gaid,3)
400 !!$ DO iiii=1 ,SIZE(this%gaid,4)
401 !!$ this%gaid(i,ii,iii,iiii) = grid_id_new() ! optimize?
402 !!$ ENDDO
403 !!$ ENDDO
404 !!$ ENDDO
405 !!$ ENDDO
406 
407  this%gaid = grid_id_new()
408  ENDIF
409  ENDIF
410 
411 ELSE
412  CALL l4f_category_log(this%category,l4f_fatal,'volgrid6d_alloc_vol: &
413  &trying to allocate a volume with an invalid or unspecified horizontal grid')
414  CALL raise_fatal_error()
415 ENDIF
416 
417 END SUBROUTINE volgrid6d_alloc_vol
418 
433 SUBROUTINE volgrid_get_vol_2d(this, ilevel, itime, itimerange, ivar, voldati)
434 TYPE(volgrid6d),INTENT(in) :: this
435 INTEGER,INTENT(in) :: ilevel
436 INTEGER,INTENT(in) :: itime
437 INTEGER,INTENT(in) :: itimerange
438 INTEGER,INTENT(in) :: ivar
439 REAL,POINTER :: voldati(:,:)
440 
441 IF (ASSOCIATED(this%voldati)) THEN
442  voldati => this%voldati(:,:,ilevel,itime,itimerange,ivar)
443  RETURN
444 ELSE
445  IF (.NOT.ASSOCIATED(voldati)) THEN
446  ALLOCATE(voldati(this%griddim%dim%nx,this%griddim%dim%ny))
447  ENDIF
448  CALL grid_id_decode_data(this%gaid(ilevel,itime,itimerange,ivar), voldati)
449 ENDIF
450 
451 END SUBROUTINE volgrid_get_vol_2d
452 
453 
467 SUBROUTINE volgrid_get_vol_3d(this, itime, itimerange, ivar, voldati)
468 TYPE(volgrid6d),INTENT(in) :: this
469 INTEGER,INTENT(in) :: itime
470 INTEGER,INTENT(in) :: itimerange
471 INTEGER,INTENT(in) :: ivar
472 REAL,POINTER :: voldati(:,:,:)
473 
474 INTEGER :: ilevel
475 
476 IF (ASSOCIATED(this%voldati)) THEN
477  voldati => this%voldati(:,:,:,itime,itimerange,ivar)
478  RETURN
479 ELSE
480  IF (.NOT.ASSOCIATED(voldati)) THEN
481  ALLOCATE(voldati(this%griddim%dim%nx,this%griddim%dim%ny,SIZE(this%level)))
482  ENDIF
483  DO ilevel = 1, SIZE(this%level)
484  CALL grid_id_decode_data(this%gaid(ilevel,itime,itimerange,ivar), &
485  voldati(:,:,ilevel))
486  ENDDO
487 ENDIF
488 
489 END SUBROUTINE volgrid_get_vol_3d
490 
491 
503 SUBROUTINE volgrid_set_vol_2d(this, ilevel, itime, itimerange, ivar, voldati)
504 TYPE(volgrid6d),INTENT(inout) :: this
505 INTEGER,INTENT(in) :: ilevel
506 INTEGER,INTENT(in) :: itime
507 INTEGER,INTENT(in) :: itimerange
508 INTEGER,INTENT(in) :: ivar
509 REAL,INTENT(in) :: voldati(:,:)
510 
511 IF (ASSOCIATED(this%voldati)) THEN
512  RETURN
513 ELSE
514  CALL grid_id_encode_data(this%gaid(ilevel,itime,itimerange,ivar), voldati)
515 ENDIF
516 
517 END SUBROUTINE volgrid_set_vol_2d
518 
519 
531 SUBROUTINE volgrid_set_vol_3d(this, itime, itimerange, ivar, voldati)
532 TYPE(volgrid6d),INTENT(inout) :: this
533 INTEGER,INTENT(in) :: itime
534 INTEGER,INTENT(in) :: itimerange
535 INTEGER,INTENT(in) :: ivar
536 REAL,INTENT(in) :: voldati(:,:,:)
537 
538 INTEGER :: ilevel
539 
540 IF (ASSOCIATED(this%voldati)) THEN
541  RETURN
542 ELSE
543  DO ilevel = 1, SIZE(this%level)
544  CALL grid_id_encode_data(this%gaid(ilevel,itime,itimerange,ivar), &
545  voldati(:,:,ilevel))
546  ENDDO
547 ENDIF
548 
549 END SUBROUTINE volgrid_set_vol_3d
550 
551 
555 SUBROUTINE volgrid6d_delete(this)
556 TYPE(volgrid6d),INTENT(inout) :: this
557 
558 INTEGER :: i, ii, iii, iiii
559 
560 #ifdef DEBUG
561 call l4f_category_log(this%category,l4f_debug,"delete")
562 #endif
563 
564 if (associated(this%gaid))then
565 
566  DO i=1 ,SIZE(this%gaid,1)
567  DO ii=1 ,SIZE(this%gaid,2)
568  DO iii=1 ,SIZE(this%gaid,3)
569  DO iiii=1 ,SIZE(this%gaid,4)
570  CALL delete(this%gaid(i,ii,iii,iiii))
571  ENDDO
572  ENDDO
573  ENDDO
574  ENDDO
575  DEALLOCATE(this%gaid)
576 
577 end if
578 
579 call delete(this%griddim)
580 
581 ! call delete(this%time)
582 ! call delete(this%timerange)
583 ! call delete(this%level)
584 ! call delete(this%var)
585 
586 if (associated( this%time )) deallocate(this%time)
587 if (associated( this%timerange )) deallocate(this%timerange)
588 if (associated( this%level )) deallocate(this%level)
589 if (associated( this%var )) deallocate(this%var)
590 
591 if (associated(this%voldati))deallocate(this%voldati)
592 
593 
594  !chiudo il logger
595 call l4f_category_delete(this%category)
596 
597 END SUBROUTINE volgrid6d_delete
598 
599 
609 subroutine volgrid6d_write_on_file (this,unit,description,filename,filename_auto)
610 
611 TYPE(volgrid6d),INTENT(IN) :: this
612 integer,optional,intent(inout) :: unit
613 character(len=*),intent(in),optional :: filename
614 character(len=*),intent(out),optional :: filename_auto
615 character(len=*),INTENT(IN),optional :: description
616 
617 integer :: lunit
618 character(len=254) :: ldescription,arg,lfilename
619 integer :: ntime, ntimerange, nlevel, nvar
620 integer :: tarray(8)
621 logical :: opened,exist
622 
623 #ifdef DEBUG
624 call l4f_category_log(this%category,l4f_debug,"write on file")
625 #endif
626 
627 ntime=0
628 ntimerange=0
629 nlevel=0
630 nvar=0
631 
632 !call idate(im,id,iy)
633 call date_and_time(values=tarray)
634 call getarg(0,arg)
635 
636 if (present(description))then
637  ldescription=description
638 else
639  ldescription="Volgrid6d generated by: "//trim(arg)
640 end if
641 
642 if (.not. present(unit))then
643  lunit=getunit()
644 else
645  if (unit==0)then
646  lunit=getunit()
647  unit=lunit
648  else
649  lunit=unit
650  end if
651 end if
652 
653 lfilename=trim(arg)//".vg6d"
654 if (index(arg,'/',back=.true.) > 0) lfilename=lfilename(index(arg,'/',back=.true.)+1 : )
655 
656 if (present(filename))then
657  if (filename /= "")then
658  lfilename=filename
659  end if
660 end if
662 if (present(filename_auto))filename_auto=lfilename
663 
664 
665 inquire(unit=lunit,opened=opened)
666 if (.not. opened) then
667  inquire(file=lfilename,exist=exist)
668  if (exist) CALL raise_error('file exist; cannot open new file')
669  if (.not.exist) open (unit=lunit,file=lfilename,form="UNFORMATTED")
670  !print *, "opened: ",lfilename
671 end if
672 
673 if (associated(this%time)) ntime=size(this%time)
674 if (associated(this%timerange)) ntimerange=size(this%timerange)
675 if (associated(this%level)) nlevel=size(this%level)
676 if (associated(this%var)) nvar=size(this%var)
677 
678 
679 write(unit=lunit)ldescription
680 write(unit=lunit)tarray
681 
682 call write_unit( this%griddim,lunit)
683 write(unit=lunit) ntime, ntimerange, nlevel, nvar
684 
685 !! prime 4 dimensioni
686 if (associated(this%time)) call write_unit(this%time, lunit)
687 if (associated(this%level)) write(unit=lunit)this%level
688 if (associated(this%timerange)) write(unit=lunit)this%timerange
689 if (associated(this%var)) write(unit=lunit)this%var
690 
691 
692 !! Volumi di valori dati
693 
694 if (associated(this%voldati)) write(unit=lunit)this%voldati
695 
696 if (.not. present(unit)) close(unit=lunit)
698 end subroutine volgrid6d_write_on_file
699 
700 
707 subroutine volgrid6d_read_from_file (this,unit,filename,description,tarray,filename_auto)
708 
709 TYPE(volgrid6d),INTENT(OUT) :: this
710 integer,intent(inout),optional :: unit
711 character(len=*),INTENT(in),optional :: filename
712 character(len=*),intent(out),optional :: filename_auto
713 character(len=*),INTENT(out),optional :: description
714 integer,intent(out),optional :: tarray(8)
715 
716 integer :: ntime, ntimerange, nlevel, nvar
717 
718 character(len=254) :: ldescription,lfilename,arg
719 integer :: ltarray(8),lunit
720 logical :: opened,exist
721 
722 #ifdef DEBUG
723 call l4f_category_log(this%category,l4f_debug,"read from file")
724 #endif
726 call getarg(0,arg)
727 
728 if (.not. present(unit))then
729  lunit=getunit()
730 else
731  if (unit==0)then
732  lunit=getunit()
733  unit=lunit
734  else
735  lunit=unit
736  end if
737 end if
738 
739 lfilename=trim(arg)//".vg6d"
740 if (index(arg,'/',back=.true.) > 0) lfilename=lfilename(index(arg,'/',back=.true.)+1 : )
741 
742 if (present(filename))then
743  if (filename /= "")then
744  lfilename=filename
745  end if
746 end if
747 
748 if (present(filename_auto))filename_auto=lfilename
750 
751 inquire(unit=lunit,opened=opened)
752 if (.not. opened) then
753  inquire(file=lfilename,exist=exist)
754  IF (.NOT. exist) CALL raise_fatal_error('file '//trim(lfilename)//' does not exist, cannot open')
755  open (unit=lunit,file=lfilename,form="UNFORMATTED")
756 end if
757 
758 read(unit=lunit)ldescription
759 read(unit=lunit)ltarray
760 
761 call l4f_log(l4f_info,"Info: reading volgrid6d from file: "//trim(lfilename))
762 call l4f_log(l4f_info,"Info: description: "//trim(ldescription))
763 !call l4f_log("Info: written on ",ltarray)
764 
765 if (present(description))description=ldescription
766 if (present(tarray))tarray=ltarray
767 
768 
769 call read_unit( this%griddim,lunit)
770 read(unit=lunit) ntime, ntimerange, nlevel, nvar
771 
772 
773 call volgrid6d_alloc(this, &
774  ntime=ntime, ntimerange=ntimerange, nlevel=nlevel, nvar=nvar)
775 
776 call volgrid6d_alloc_vol(this)
777 
778 if (associated(this%time)) call read_unit(this%time, lunit)
779 if (associated(this%level)) read(unit=lunit)this%level
780 if (associated(this%timerange)) read(unit=lunit)this%timerange
781 if (associated(this%var)) read(unit=lunit)this%var
782 
783 
784 !! Volumi di valori
785 
786 if (associated(this%voldati)) read(unit=lunit)this%voldati
787 
788 if (.not. present(unit)) close(unit=lunit)
789 
790 end subroutine volgrid6d_read_from_file
791 
792 
812 SUBROUTINE import_from_gridinfo(this, gridinfo, force, dup_mode, clone, &
813  isanavar)
814 TYPE(volgrid6d),INTENT(inout) :: this
815 type(gridinfo_def),INTENT(in) :: gridinfo
816 LOGICAL,INTENT(in),OPTIONAL :: force
817 INTEGER,INTENT(in),OPTIONAL :: dup_mode
818 LOGICAL , INTENT(in),OPTIONAL :: clone
819 LOGICAL,INTENT(IN),OPTIONAL :: isanavar
820 
821 CHARACTER(len=255) :: type
822 INTEGER :: itime0, itimerange0, itime1, itimerange1, itime, itimerange, &
823  ilevel, ivar, ldup_mode
824 LOGICAL :: dup
825 TYPE(datetime) :: correctedtime
826 REAL,ALLOCATABLE :: tmpgrid(:,:)
827 
828 IF (present(dup_mode)) THEN
829  ldup_mode = dup_mode
830 ELSE
831  ldup_mode = 0
832 ENDIF
833 
834 call get_val(this%griddim,proj_type=type)
835 
836 #ifdef DEBUG
837 call l4f_category_log(this%category,l4f_debug,"import_from_gridinfo: "//trim(type))
838 #endif
839 
840 if (.not. c_e(type))then
841  call copy(gridinfo%griddim, this%griddim)
842 ! ho gia` fatto init, altrimenti non potrei fare la get_val(this%griddim)
843 ! per cui meglio non ripetere
844 ! call init(this,gridinfo%griddim,categoryappend)
845  CALL volgrid6d_alloc_vol(this, ini=.true.) ! decode?
846 
847 else if (.not. (this%griddim == gridinfo%griddim ))then
848 
849  CALL l4f_category_log(this%category,l4f_error, &
850  "volgrid and gridinfo grid type or size are different, gridinfo rejected")
851  CALL raise_error()
852  RETURN
853 
854 end if
855 
856 ! Cerco gli indici del campo da inserire, se non trovo metto nel primo missing
857 ilevel = index(this%level, gridinfo%level)
858 IF (ilevel == 0 .AND. optio_log(force)) THEN
859  ilevel = index(this%level, vol7d_level_miss)
860  IF (ilevel /= 0) this%level(ilevel) = gridinfo%level
861 ENDIF
862 
863 IF (ilevel == 0) THEN
864  CALL l4f_category_log(this%category,l4f_error, &
865  "volgrid6d: level not valid for volume, gridinfo rejected")
866  CALL raise_error()
867  RETURN
868 ENDIF
869 
870 IF (optio_log(isanavar)) THEN ! assign to all times and timeranges
871  itime0 = 1
872  itime1 = SIZE(this%time)
873  itimerange0 = 1
874  itimerange1 = SIZE(this%timerange)
875 ELSE ! usual case
876  correctedtime = gridinfo%time
877  IF (this%time_definition == 1) correctedtime = correctedtime + &
878  timedelta_new(sec=gridinfo%timerange%p1)
879  itime0 = index(this%time, correctedtime)
880  IF (itime0 == 0 .AND. optio_log(force)) THEN
881  itime0 = index(this%time, datetime_miss)
882  IF (itime0 /= 0) this%time(itime0) = correctedtime
883  ENDIF
884  IF (itime0 == 0) THEN
885  CALL l4f_category_log(this%category,l4f_error, &
886  "volgrid6d: time not valid for volume, gridinfo rejected")
887  CALL raise_error()
888  RETURN
889  ENDIF
890  itime1 = itime0
891 
892  itimerange0 = index(this%timerange,gridinfo%timerange)
893  IF (itimerange0 == 0 .AND. optio_log(force)) THEN
894  itimerange0 = index(this%timerange, vol7d_timerange_miss)
895  IF (itimerange0 /= 0) this%timerange(itimerange0) = gridinfo%timerange
896  ENDIF
897  IF (itimerange0 == 0) THEN
898  CALL l4f_category_log(this%category,l4f_error, &
899  "volgrid6d: timerange not valid for volume, gridinfo rejected")
900  CALL raise_error()
901  RETURN
902  ENDIF
903  itimerange1 = itimerange0
904 ENDIF
905 
906 ivar = index(this%var, gridinfo%var)
907 IF (ivar == 0 .AND. optio_log(force)) THEN
908  ivar = index(this%var, volgrid6d_var_miss)
909  IF (ivar /= 0) this%var(ivar) = gridinfo%var
910 ENDIF
911 IF (ivar == 0) THEN
912  CALL l4f_category_log(this%category,l4f_error, &
913  "volgrid6d: var not valid for volume, gridinfo rejected")
914  CALL raise_error()
915  RETURN
916 ENDIF
917 
918 DO itimerange = itimerange0, itimerange1
919  DO itime = itime0, itime1
920  IF (ASSOCIATED(this%gaid)) THEN
921  dup = .false.
922  IF (c_e(this%gaid(ilevel,itime,itimerange,ivar))) THEN
923  dup = .true.
924  CALL l4f_category_log(this%category,l4f_warn,"gaid exist: grib duplicated")
925 ! avoid memory leaks
926  IF (optio_log(clone)) CALL delete(this%gaid(ilevel,itime,itimerange,ivar))
927  ENDIF
928 
929  IF (optio_log(clone)) THEN
930  CALL copy(gridinfo%gaid, this%gaid(ilevel,itime,itimerange,ivar))
931 #ifdef DEBUG
932  CALL l4f_category_log(this%category,l4f_debug,"cloning to a new gaid")
933 #endif
934  ELSE
935  this%gaid(ilevel,itime,itimerange,ivar) = gridinfo%gaid
936  ENDIF
937 
938  IF (ASSOCIATED(this%voldati))THEN
939  IF (.NOT.dup .OR. ldup_mode == 0) THEN
940  this%voldati(:,:,ilevel,itime,itimerange,ivar) = decode_gridinfo(gridinfo)
941  ELSE IF (ldup_mode == 1) THEN
942  tmpgrid = decode_gridinfo(gridinfo) ! f2003 automatic allocation
943  WHERE(c_e(tmpgrid))
944  this%voldati(:,:,ilevel,itime,itimerange,ivar) = tmpgrid(:,:)
945  END WHERE
946  ELSE IF (ldup_mode == 2) THEN
947  WHERE(.NOT.c_e(this%voldati(:,:,ilevel,itime,itimerange,ivar)))
948  this%voldati(:,:,ilevel,itime,itimerange,ivar) = decode_gridinfo(gridinfo)
949  END WHERE
950  ENDIF
951  ENDIF
952 
953  ELSE
954  CALL l4f_category_log(this%category,l4f_error, &
955  "gaid not allocated, you probably need to call volgrid6d_alloc_vol first")
956  CALL raise_error()
957  RETURN
958  ENDIF
959  ENDDO
960 ENDDO
961 
962 
963 END SUBROUTINE import_from_gridinfo
964 
965 
970 SUBROUTINE export_to_gridinfo(this, gridinfo, itime, itimerange, ilevel, ivar, &
971  gaid_template, clone)
972 TYPE(volgrid6d),INTENT(in) :: this
973 type(gridinfo_def),INTENT(inout) :: gridinfo
974 INTEGER :: itime
975 INTEGER :: itimerange
976 INTEGER :: ilevel
977 INTEGER :: ivar
978 type(grid_id),INTENT(in),OPTIONAL :: gaid_template
979 LOGICAL,INTENT(in),OPTIONAL :: clone
980 
981 TYPE(grid_id) :: gaid
982 LOGICAL :: usetemplate
983 REAL,POINTER :: voldati(:,:)
984 TYPE(datetime) :: correctedtime
985 
986 #ifdef DEBUG
987 CALL l4f_category_log(this%category,l4f_debug,"export_to_gridinfo")
988 #endif
989 
990 IF (.NOT.c_e(this%gaid(ilevel,itime,itimerange,ivar))) THEN
991 #ifdef DEBUG
992  CALL l4f_category_log(this%category,l4f_debug,"empty gaid found, skipping export")
993 #endif
994  RETURN
995 ENDIF
996 
997 usetemplate = .false.
998 IF (present(gaid_template)) THEN
999  CALL copy(gaid_template, gaid)
1000 #ifdef DEBUG
1001  CALL l4f_category_log(this%category,l4f_debug,"template cloned to a new gaid")
1002 #endif
1003  usetemplate = c_e(gaid)
1004 ENDIF
1005 
1006 IF (.NOT.usetemplate) THEN
1007  IF (optio_log(clone)) THEN
1008  CALL copy(this%gaid(ilevel,itime,itimerange,ivar), gaid)
1009 #ifdef DEBUG
1010  CALL l4f_category_log(this%category,l4f_debug,"original gaid cloned to a new one")
1011 #endif
1012  ELSE
1013  gaid = this%gaid(ilevel,itime,itimerange,ivar)
1014  ENDIF
1015 ENDIF
1016 
1017 IF (this%time_definition == 1) THEN
1018  correctedtime = this%time(itime) - &
1019  timedelta_new(sec=this%timerange(itimerange)%p1)
1020 ELSE
1021  correctedtime = this%time(itime)
1022 ENDIF
1023 
1024 CALL init(gridinfo,gaid, this%griddim, correctedtime, this%timerange(itimerange), &
1025  this%level(ilevel), this%var(ivar))
1026 
1027 ! reset the gridinfo, bad but necessary at this point for encoding the field
1028 CALL export(gridinfo%griddim, gridinfo%gaid)
1029 ! encode the field
1030 IF (ASSOCIATED(this%voldati)) THEN
1031  CALL encode_gridinfo(gridinfo, this%voldati(:,:,ilevel,itime,itimerange,ivar))
1032 ELSE IF (usetemplate) THEN ! field must be forced into template in this case
1033  nullify(voldati)
1034  CALL volgrid_get_vol_2d(this, ilevel, itime, itimerange, ivar, voldati)
1035  CALL encode_gridinfo(gridinfo, voldati)
1036  DEALLOCATE(voldati)
1037 ENDIF
1038 
1039 END SUBROUTINE export_to_gridinfo
1040 
1041 
1059 SUBROUTINE import_from_gridinfovv(this, gridinfov, dup_mode, clone, decode, &
1060  time_definition, anavar, categoryappend)
1061 TYPE(volgrid6d),POINTER :: this(:)
1062 type(arrayof_gridinfo),INTENT(in) :: gridinfov
1063 INTEGER,INTENT(in),OPTIONAL :: dup_mode
1064 LOGICAL , INTENT(in),OPTIONAL :: clone
1065 LOGICAL,INTENT(in),OPTIONAL :: decode
1066 INTEGER,INTENT(IN),OPTIONAL :: time_definition
1067 CHARACTER(len=*),INTENT(IN),OPTIONAL :: anavar(:)
1068 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend
1069 
1070 INTEGER :: i, j, stallo
1071 INTEGER :: ngrid, ntime, ntimerange, nlevel, nvar
1072 INTEGER :: category
1073 CHARACTER(len=512) :: a_name
1074 TYPE(datetime),ALLOCATABLE :: correctedtime(:)
1075 LOGICAL,ALLOCATABLE :: isanavar(:)
1076 TYPE(vol7d_var) :: lvar
1077 
1078 ! category temporanea (altrimenti non possiamo loggare)
1079 if (present(categoryappend))then
1080  call l4f_launcher(a_name,a_name_append=trim(subcategory)//"."//trim(categoryappend))
1081 else
1082  call l4f_launcher(a_name,a_name_append=trim(subcategory))
1083 endif
1084 category=l4f_category_get(a_name)
1085 
1086 #ifdef DEBUG
1087 call l4f_category_log(category,l4f_debug,"start import_from_gridinfovv")
1088 #endif
1089 
1090 ngrid=count_distinct(gridinfov%array(1:gridinfov%arraysize)%griddim,back=.true.)
1091 CALL l4f_category_log(category,l4f_info, t2c(ngrid)// &
1092  ' different grid definition(s) found in input data')
1093 
1094 ALLOCATE(this(ngrid),stat=stallo)
1095 IF (stallo /= 0)THEN
1096  CALL l4f_category_log(category,l4f_fatal,"allocating memory")
1097  CALL raise_fatal_error()
1098 ENDIF
1099 DO i = 1, ngrid
1100  IF (present(categoryappend))THEN
1101  CALL init(this(i), time_definition=time_definition, categoryappend=trim(categoryappend)//"-vol"//t2c(i))
1102  ELSE
1103  CALL init(this(i), time_definition=time_definition, categoryappend="vol"//t2c(i))
1104  ENDIF
1105 ENDDO
1106 
1107 this(:)%griddim=pack_distinct(gridinfov%array(1:gridinfov%arraysize)%griddim, &
1108  ngrid, back=.true.)
1109 
1110 ! mark elements as ana variables (time-independent)
1111 ALLOCATE(isanavar(gridinfov%arraysize))
1112 isanavar(:) = .false.
1113 IF (present(anavar)) THEN
1114  DO i = 1, gridinfov%arraysize
1115  DO j = 1, SIZE(anavar)
1116  lvar = convert(gridinfov%array(i)%var)
1117  IF (lvar%btable == anavar(j)) THEN
1118  isanavar(i) = .true.
1119  EXIT
1120  ENDIF
1121  ENDDO
1122  ENDDO
1123  CALL l4f_category_log(category,l4f_info,t2c(count(isanavar))//'/'// &
1124  t2c(gridinfov%arraysize)//' constant-data messages found in input data')
1125 ENDIF
1126 
1127 ! create time corrected for time_definition
1128 ALLOCATE(correctedtime(gridinfov%arraysize))
1129 correctedtime(:) = gridinfov%array(1:gridinfov%arraysize)%time
1130 IF (present(time_definition)) THEN
1131  IF (time_definition == 1) THEN
1132  DO i = 1, gridinfov%arraysize
1133  correctedtime(i) = correctedtime(i) + &
1134  timedelta_new(sec=gridinfov%array(i)%timerange%p1)
1135  ENDDO
1136  ENDIF
1137 ENDIF
1138 
1139 DO i = 1, ngrid
1140  IF (present(anavar)) THEN
1141  j = count((this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1142  .AND. .NOT.isanavar(:))
1143  IF (j <= 0) THEN
1144  CALL l4f_category_log(category, l4f_fatal, 'grid n.'//t2c(i)// &
1145  ' has only constant data, this is not allowed')
1146  CALL l4f_category_log(category, l4f_fatal, 'please check anavar argument')
1147  CALL raise_fatal_error()
1148  ENDIF
1149  ENDIF
1150  ntime = count_distinct(correctedtime, &
1151  mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1152  .AND. .NOT.isanavar(:), back=.true.)
1153  ntimerange = count_distinct(gridinfov%array(1:gridinfov%arraysize)%timerange, &
1154  mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1155  .AND. .NOT.isanavar(:), back=.true.)
1156  nlevel = count_distinct(gridinfov%array(1:gridinfov%arraysize)%level, &
1157  mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1158  back=.true.)
1159  nvar = count_distinct(gridinfov%array(1:gridinfov%arraysize)%var, &
1160  mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1161  back=.true.)
1162 
1163 #ifdef DEBUG
1164  CALL l4f_category_log(this(i)%category,l4f_debug,"alloc volgrid6d index: "//t2c(i))
1165 #endif
1166 
1167  CALL volgrid6d_alloc(this(i),this(i)%griddim%dim,ntime=ntime, &
1168  ntimerange=ntimerange,nlevel=nlevel,nvar=nvar)
1169 
1170  this(i)%time = pack_distinct(correctedtime, ntime, &
1171  mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1172  .AND. .NOT.isanavar(:), back=.true.)
1173  CALL sort(this(i)%time)
1174 
1175  this(i)%timerange = pack_distinct(gridinfov%array( &
1176  1:gridinfov%arraysize)%timerange, ntimerange, &
1177  mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1178  .AND. .NOT.isanavar(:), back=.true.)
1179  CALL sort(this(i)%timerange)
1180 
1181  this(i)%level=pack_distinct(gridinfov%array(1:gridinfov%arraysize)%level, &
1182  nlevel,mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1183  back=.true.)
1184  CALL sort(this(i)%level)
1185 
1186  this(i)%var=pack_distinct(gridinfov%array(1:gridinfov%arraysize)%var, nvar, &
1187  mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1188  back=.true.)
1189 
1190 #ifdef DEBUG
1191  CALL l4f_category_log(this(i)%category,l4f_debug,"alloc_vol volgrid6d index: "//t2c(i))
1192 #endif
1193  CALL volgrid6d_alloc_vol(this(i), decode=decode)
1194 
1195 ENDDO
1196 
1197 DEALLOCATE(correctedtime)
1198 
1199 DO i = 1, gridinfov%arraysize
1200 
1201 #ifdef DEBUG
1202  CALL l4f_category_log(category,l4f_debug,"import from gridinfov index: "//t2c(i))
1203  CALL l4f_category_log(category,l4f_info, &
1204  "to volgrid6d index: "//t2c(index(this%griddim, gridinfov%array(i)%griddim)))
1205 #endif
1206 
1207  CALL import(this(index(this%griddim, gridinfov%array(i)%griddim)), &
1208  gridinfov%array(i), dup_mode=dup_mode, clone=clone, isanavar=isanavar(i))
1209 
1210 ENDDO
1211 
1212 !chiudo il logger temporaneo
1213 CALL l4f_category_delete(category)
1214 
1215 END SUBROUTINE import_from_gridinfovv
1216 
1217 
1223 SUBROUTINE export_to_gridinfov(this, gridinfov, gaid_template, clone)
1224 TYPE(volgrid6d),INTENT(inout) :: this
1225 type(arrayof_gridinfo),INTENT(inout) :: gridinfov
1226 type(grid_id),INTENT(in),OPTIONAL :: gaid_template
1227 LOGICAL,INTENT(in),OPTIONAL :: clone
1228 
1229 INTEGER :: i ,itime, itimerange, ilevel, ivar
1230 INTEGER :: ntime, ntimerange, nlevel, nvar
1231 TYPE(gridinfo_def) :: gridinfol
1232 
1233 #ifdef DEBUG
1234 CALL l4f_category_log(this%category,l4f_debug,"start export_to_gridinfov")
1235 #endif
1236 
1237 ! this is necessary in order not to repeatedly and uselessly copy the
1238 ! same array of coordinates for each 2d grid during export, the
1239 ! side-effect is that the computed projection in this is lost
1240 CALL dealloc(this%griddim%dim)
1241 
1242 i=0
1243 ntime=size(this%time)
1244 ntimerange=size(this%timerange)
1245 nlevel=size(this%level)
1246 nvar=size(this%var)
1247 
1248 DO itime=1,ntime
1249  DO itimerange=1,ntimerange
1250  DO ilevel=1,nlevel
1251  DO ivar=1,nvar
1252 
1253  CALL init(gridinfol)
1254  CALL export(this, gridinfol, itime, itimerange, ilevel, ivar, &
1255  gaid_template=gaid_template, clone=clone)
1256  IF (c_e(gridinfol%gaid)) THEN
1257  CALL insert(gridinfov, gridinfol)
1258  ELSE
1259  CALL delete(gridinfol)
1260  ENDIF
1261 
1262  ENDDO
1263  ENDDO
1264  ENDDO
1265 ENDDO
1266 
1267 END SUBROUTINE export_to_gridinfov
1268 
1269 
1275 SUBROUTINE export_to_gridinfovv(this, gridinfov, gaid_template, clone)
1276 !, &
1277 ! categoryappend)
1278 TYPE(volgrid6d),INTENT(inout) :: this(:)
1279 type(arrayof_gridinfo),INTENT(inout) :: gridinfov
1280 type(grid_id),INTENT(in),OPTIONAL :: gaid_template
1281 LOGICAL,INTENT(in),OPTIONAL :: clone
1282 
1283 INTEGER :: i
1284 
1285 DO i = 1, SIZE(this)
1286 #ifdef DEBUG
1287  CALL l4f_category_log(this(i)%category,l4f_debug, &
1288  "export_to_gridinfovv grid index: "//t2c(i))
1289 #endif
1290  CALL export(this(i), gridinfov, gaid_template=gaid_template, clone=clone)
1291 ENDDO
1292 
1293 END SUBROUTINE export_to_gridinfovv
1294 
1295 
1305 SUBROUTINE volgrid6d_import_from_file(this, filename, dup_mode, decode, &
1306  time_definition, anavar, categoryappend)
1307 TYPE(volgrid6d),POINTER :: this(:)
1308 CHARACTER(len=*),INTENT(in) :: filename
1309 INTEGER,INTENT(in),OPTIONAL :: dup_mode
1310 LOGICAL,INTENT(in),OPTIONAL :: decode
1311 INTEGER,INTENT(IN),OPTIONAL :: time_definition
1312 CHARACTER(len=*),INTENT(IN),OPTIONAL :: anavar(:)
1313 character(len=*),INTENT(in),OPTIONAL :: categoryappend
1314 
1315 TYPE(arrayof_gridinfo) :: gridinfo
1316 INTEGER :: category
1317 CHARACTER(len=512) :: a_name
1318 
1319 nullify(this)
1320 
1321 IF (present(categoryappend))THEN
1322  CALL l4f_launcher(a_name,a_name_append= &
1323  trim(subcategory)//"."//trim(categoryappend))
1324 ELSE
1325  CALL l4f_launcher(a_name,a_name_append=trim(subcategory))
1326 ENDIF
1327 category=l4f_category_get(a_name)
1328 
1329 CALL import(gridinfo, filename=filename, categoryappend=categoryappend)
1330 
1331 IF (gridinfo%arraysize > 0) THEN
1332 
1333  CALL import(this, gridinfo, dup_mode=dup_mode, clone=.true., decode=decode, &
1334  time_definition=time_definition, anavar=anavar, &
1335  categoryappend=categoryappend)
1336 
1337  CALL l4f_category_log(category,l4f_info,"deleting gridinfo")
1338  CALL delete(gridinfo)
1339 
1340 ELSE
1341  CALL l4f_category_log(category,l4f_info,"file does not contain gridded data")
1342 ENDIF
1343 
1344 ! close logger
1345 CALL l4f_category_delete(category)
1346 
1347 END SUBROUTINE volgrid6d_import_from_file
1348 
1349 
1357 SUBROUTINE volgrid6d_export_to_file(this, filename, gaid_template, categoryappend)
1358 TYPE(volgrid6d) :: this(:)
1359 CHARACTER(len=*),INTENT(in) :: filename
1360 type(grid_id),INTENT(in),OPTIONAL :: gaid_template
1361 character(len=*),INTENT(in),OPTIONAL :: categoryappend
1362 
1363 TYPE(arrayof_gridinfo) :: gridinfo
1364 INTEGER :: category
1365 CHARACTER(len=512) :: a_name
1366 
1367 IF (present(categoryappend)) THEN
1368  CALL l4f_launcher(a_name,a_name_append=trim(subcategory)//"."//trim(categoryappend))
1369 ELSE
1370  CALL l4f_launcher(a_name,a_name_append=trim(subcategory))
1371 ENDIF
1372 category=l4f_category_get(a_name)
1373 
1374 #ifdef DEBUG
1375 CALL l4f_category_log(category,l4f_debug,"start export to file")
1376 #endif
1377 
1378 CALL l4f_category_log(category,l4f_info,"writing volgrid6d to grib file: "//trim(filename))
1379 
1380 !IF (ASSOCIATED(this)) THEN
1381  CALL export(this, gridinfo, gaid_template=gaid_template, clone=.true.)
1382  IF (gridinfo%arraysize > 0) THEN
1383  CALL export(gridinfo, filename)
1384  CALL delete(gridinfo)
1385  ENDIF
1386 !ELSE
1387 ! CALL l4f_category_log(category,L4F_INFO,"volume volgrid6d is not associated")
1388 !ENDIF
1389 
1390 ! close logger
1391 CALL l4f_category_delete(category)
1392 
1393 END SUBROUTINE volgrid6d_export_to_file
1394 
1395 
1399 SUBROUTINE volgrid6dv_delete(this)
1400 TYPE(volgrid6d),POINTER :: this(:)
1401 
1402 INTEGER :: i
1403 
1404 IF (ASSOCIATED(this)) THEN
1405  DO i = 1, SIZE(this)
1406 #ifdef DEBUG
1407  CALL l4f_category_log(this(i)%category,l4f_debug, &
1408  "delete volgrid6d vector index: "//trim(to_char(i)))
1409 #endif
1410  CALL delete(this(i))
1411  ENDDO
1412  DEALLOCATE(this)
1413 ENDIF
1414 
1415 END SUBROUTINE volgrid6dv_delete
1416 
1418 ! Internal method for performing grid to grid computations
1419 SUBROUTINE volgrid6d_transform_compute(this, volgrid6d_in, volgrid6d_out, &
1420  lev_out, var_coord_vol, clone)
1421 TYPE(grid_transform),INTENT(in) :: this ! oggetto di trasformazione per il grigliato
1422 type(volgrid6d), INTENT(in) :: volgrid6d_in ! oggetto da trasformare
1423 type(volgrid6d), INTENT(inout) :: volgrid6d_out ! oggetto trasformato; deve essere completo (init, alloc, alloc_vol)
1424 TYPE(vol7d_level),INTENT(in),OPTIONAL :: lev_out(:) ! vol7d_level object defining target vertical grid, for vertical interpolations
1425 INTEGER,INTENT(in),OPTIONAL :: var_coord_vol ! index of variable defining vertical coordinate values in input volume
1426 LOGICAL,INTENT(in),OPTIONAL :: clone ! se fornito e \c .TRUE., clona i gaid da volgrid6d_in a volgrid6d_out
1427 
1428 INTEGER :: ntime, ntimerange, inlevel, onlevel, nvar, &
1429  itime, itimerange, ilevel, ivar, levshift, levused, lvar_coord_vol, spos
1430 REAL,POINTER :: voldatiin(:,:,:), voldatiout(:,:,:), coord_3d_in(:,:,:)
1431 TYPE(vol7d_level) :: output_levtype
1432 
1433 
1434 #ifdef DEBUG
1435 call l4f_category_log(volgrid6d_in%category,l4f_debug,"start volgrid6d_transform_compute")
1436 #endif
1437 
1438 ntime=0
1439 ntimerange=0
1440 inlevel=0
1441 onlevel=0
1442 nvar=0
1443 lvar_coord_vol = optio_i(var_coord_vol)
1444 
1445 if (associated(volgrid6d_in%time))then
1446  ntime=size(volgrid6d_in%time)
1447  volgrid6d_out%time=volgrid6d_in%time
1448 end if
1449 
1450 if (associated(volgrid6d_in%timerange))then
1451  ntimerange=size(volgrid6d_in%timerange)
1452  volgrid6d_out%timerange=volgrid6d_in%timerange
1453 end if
1454 
1455 IF (ASSOCIATED(volgrid6d_in%level))THEN
1456  inlevel=SIZE(volgrid6d_in%level)
1457 ENDIF
1458 IF (present(lev_out)) THEN
1459  onlevel=SIZE(lev_out)
1460  volgrid6d_out%level=lev_out
1461 ELSE IF (ASSOCIATED(volgrid6d_in%level))THEN
1462  onlevel=SIZE(volgrid6d_in%level)
1463  volgrid6d_out%level=volgrid6d_in%level
1464 ENDIF
1465 
1466 if (associated(volgrid6d_in%var))then
1467  nvar=size(volgrid6d_in%var)
1468  volgrid6d_out%var=volgrid6d_in%var
1469 end if
1470 ! allocate once for speed
1471 IF (.NOT.ASSOCIATED(volgrid6d_in%voldati)) THEN
1472  ALLOCATE(voldatiin(volgrid6d_in%griddim%dim%nx, volgrid6d_in%griddim%dim%ny, &
1473  inlevel))
1474 ENDIF
1475 IF (.NOT.ASSOCIATED(volgrid6d_out%voldati)) THEN
1476  ALLOCATE(voldatiout(volgrid6d_out%griddim%dim%nx, volgrid6d_out%griddim%dim%ny, &
1477  onlevel))
1478 ENDIF
1479 
1480 CALL get_val(this, levshift=levshift, levused=levused)
1481 spos = imiss
1482 IF (c_e(lvar_coord_vol)) THEN
1483  CALL get_val(this%trans, output_levtype=output_levtype)
1484  IF (output_levtype%level1 == 103 .OR. output_levtype%level1 == 108) THEN
1485  spos = firsttrue(volgrid6d_in%level(:) == vol7d_level_new(1))
1486  IF (spos == 0) THEN
1487  CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1488  'output level '//t2c(output_levtype%level1)// &
1489  ' requested, but height/press of surface not provided in volume')
1490  ENDIF
1491  IF (.NOT.c_e(levshift) .AND. .NOT.c_e(levused)) THEN
1492  CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1493  'internal inconsistence, levshift and levused undefined when they should be')
1494  ENDIF
1495  ENDIF
1496 ENDIF
1497 
1498 DO ivar=1,nvar
1499 ! IF (c_e(var_coord_vol)) THEN
1500 ! IF (ivar == var_coord_vol) CYCLE ! skip coordinate variable in output
1501 ! ENDIF
1502  DO itimerange=1,ntimerange
1503  DO itime=1,ntime
1504 ! skip empty columns where possible, improve
1505  IF (c_e(levshift) .AND. c_e(levused)) THEN
1506  IF (.NOT.any(c_e( &
1507  volgrid6d_in%gaid(levshift+1:levshift+levused,itime,itimerange,ivar) &
1508  ))) cycle
1509  ENDIF
1510  DO ilevel = 1, min(inlevel,onlevel)
1511 ! if present gaid copy it
1512  IF (c_e(volgrid6d_in%gaid(ilevel,itime,itimerange,ivar)) .AND. .NOT. &
1513  c_e(volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))) THEN
1514 
1515  IF (optio_log(clone)) THEN
1516  CALL copy(volgrid6d_in%gaid(ilevel,itime,itimerange,ivar),&
1517  volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))
1518 #ifdef DEBUG
1519  CALL l4f_category_log(volgrid6d_in%category,l4f_debug, &
1520  "cloning gaid, level "//t2c(ilevel))
1521 #endif
1522  ELSE
1523  volgrid6d_out%gaid(ilevel,itime,itimerange,ivar) = &
1524  volgrid6d_in%gaid(ilevel,itime,itimerange,ivar)
1525  ENDIF
1526  ENDIF
1527  ENDDO
1528 ! if out level are more, we have to clone anyway
1529  DO ilevel = min(inlevel,onlevel) + 1, onlevel
1530  IF (c_e(volgrid6d_in%gaid(inlevel,itime,itimerange,ivar)) .AND. .NOT. &
1531  c_e(volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))) then
1532 
1533  CALL copy(volgrid6d_in%gaid(inlevel,itime,itimerange,ivar),&
1534  volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))
1535 #ifdef DEBUG
1536  CALL l4f_category_log(volgrid6d_in%category,l4f_debug, &
1537  "forced cloning gaid, level "//t2c(inlevel)//"->"//t2c(ilevel))
1538 #endif
1539  ENDIF
1540  ENDDO
1541 
1542  IF (c_e(lvar_coord_vol)) THEN
1543  nullify(coord_3d_in)
1544  CALL volgrid_get_vol_3d(volgrid6d_in, itime, itimerange, lvar_coord_vol, &
1545  coord_3d_in)
1546  IF (c_e(spos)) THEN ! compute difference wrt surface coordinate
1547  IF (spos == 0) THEN ! error condition, set all to missing and goodnight
1548  coord_3d_in(:,:,levshift+1:levshift+levused) = rmiss
1549  ELSE
1550  DO ilevel = levshift+1, levshift+levused
1551  WHERE(c_e(coord_3d_in(:,:,ilevel)) .AND. c_e(coord_3d_in(:,:,spos)))
1552  coord_3d_in(:,:,ilevel) = coord_3d_in(:,:,ilevel) - &
1553  coord_3d_in(:,:,spos)
1554  ELSEWHERE
1555  coord_3d_in(:,:,ilevel) = rmiss
1556  END WHERE
1557  ENDDO
1558  ENDIF
1559  ENDIF
1560  ENDIF
1561  CALL volgrid_get_vol_3d(volgrid6d_in, itime, itimerange, ivar, &
1562  voldatiin)
1563  IF (ASSOCIATED(volgrid6d_out%voldati)) & ! improve!!!!
1564  CALL volgrid_get_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
1565  voldatiout)
1566  IF (c_e(lvar_coord_vol)) THEN
1567  CALL compute(this, voldatiin, voldatiout, convert(volgrid6d_in%var(ivar)), &
1568  coord_3d_in(:,:,levshift+1:levshift+levused)) ! subset coord_3d_in
1569  ELSE
1570  CALL compute(this, voldatiin, voldatiout, convert(volgrid6d_in%var(ivar)))
1571  ENDIF
1572  CALL volgrid_set_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
1573  voldatiout)
1574  ENDDO
1575  ENDDO
1576 ENDDO
1577 
1578 IF (c_e(lvar_coord_vol)) THEN
1579  DEALLOCATE(coord_3d_in)
1580 ENDIF
1581 IF (.NOT.ASSOCIATED(volgrid6d_in%voldati)) THEN
1582  DEALLOCATE(voldatiin)
1583 ENDIF
1584 IF (.NOT.ASSOCIATED(volgrid6d_out%voldati)) THEN
1585  DEALLOCATE(voldatiout)
1586 ENDIF
1587 
1588 
1589 END SUBROUTINE volgrid6d_transform_compute
1590 
1591 
1598 SUBROUTINE volgrid6d_transform(this, griddim, volgrid6d_in, volgrid6d_out, &
1599  lev_out, volgrid6d_coord_in, maskgrid, maskbounds, clone, decode, categoryappend)
1600 TYPE(transform_def),INTENT(in) :: this
1601 type(griddim_def),INTENT(in),OPTIONAL :: griddim
1602 ! TODO ripristinare intent(in) dopo le opportune modifiche in grid_class.F90
1603 TYPE(volgrid6d),INTENT(inout) :: volgrid6d_in
1604 type(volgrid6d),INTENT(out) :: volgrid6d_out
1605 type(vol7d_level),INTENT(in),OPTIONAL,TARGET :: lev_out(:)
1606 type(volgrid6d),INTENT(in),OPTIONAL :: volgrid6d_coord_in
1607 REAL,INTENT(in),OPTIONAL :: maskgrid(:,:)
1608 REAL,INTENT(in),OPTIONAL :: maskbounds(:)
1609 LOGICAL,INTENT(in),OPTIONAL :: clone
1610 LOGICAL,INTENT(in),OPTIONAL :: decode
1611 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend
1612 
1613 TYPE(grid_transform) :: grid_trans
1614 TYPE(vol7d_level),POINTER :: llev_out(:)
1615 TYPE(vol7d_level) :: input_levtype, output_levtype
1616 TYPE(vol7d_var) :: vcoord_var
1617 INTEGER :: i, k, ntime, ntimerange, nlevel, nvar, var_coord_in, var_coord_vol, &
1618  cf_out, nxc, nyc, nxi, nyi, i3, i4, i5, i6, &
1619  ulstart, ulend, spos
1620 REAL,ALLOCATABLE :: coord_3d_in(:,:,:)
1621 TYPE(geo_proj) :: proj_in, proj_out
1622 CHARACTER(len=80) :: trans_type
1623 LOGICAL :: ldecode
1624 LOGICAL,ALLOCATABLE :: mask_in(:)
1625 
1626 #ifdef DEBUG
1627 call l4f_category_log(volgrid6d_in%category, l4f_debug, "start volgrid6d_transform")
1628 #endif
1629 
1630 ntime=0
1631 ntimerange=0
1632 nlevel=0
1633 nvar=0
1634 
1635 if (associated(volgrid6d_in%time)) ntime=size(volgrid6d_in%time)
1636 if (associated(volgrid6d_in%timerange)) ntimerange=size(volgrid6d_in%timerange)
1637 if (associated(volgrid6d_in%level)) nlevel=size(volgrid6d_in%level)
1638 if (associated(volgrid6d_in%var)) nvar=size(volgrid6d_in%var)
1639 
1640 IF (ntime == 0 .OR. ntimerange == 0 .OR. nlevel == 0 .OR. nvar == 0) THEN
1641  CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1642  "trying to transform an incomplete volgrid6d object, ntime="//t2c(ntime)// &
1643  ' ntimerange='//t2c(ntimerange)//' nlevel='//t2c(nlevel)//' nvar='//t2c(nvar))
1644  CALL init(volgrid6d_out) ! initialize to empty
1645  CALL raise_error()
1646  RETURN
1647 ENDIF
1648 
1649 CALL get_val(this, trans_type=trans_type)
1650 
1651 ! store desired output component flag and unrotate if necessary
1652 cf_out = imiss
1653 IF (present(griddim) .AND. (trans_type == 'inter' .OR. trans_type == 'boxinter' &
1654  .OR. trans_type == 'stencilinter')) THEN ! improve condition!!
1655  CALL get_val(volgrid6d_in%griddim, proj=proj_in)
1656  CALL get_val(griddim, component_flag=cf_out, proj=proj_out)
1657 ! if different projections wind components must be referred to geographical system
1658  IF (proj_in /= proj_out) CALL vg6d_wind_unrot(volgrid6d_in)
1659 ELSE IF (present(griddim)) THEN ! just get component_flag, the rest is rubbish
1660  CALL get_val(griddim, component_flag=cf_out)
1661 ENDIF
1662 
1663 
1664 var_coord_in = imiss
1665 var_coord_vol = imiss
1666 IF (trans_type == 'vertint') THEN
1667  IF (present(lev_out)) THEN
1668 
1669 ! if volgrid6d_coord_in provided and allocated, check that it fits
1670  IF (present(volgrid6d_coord_in)) THEN
1671  IF (ASSOCIATED(volgrid6d_coord_in%voldati)) THEN
1672 
1673 ! strictly 1 time and 1 timerange
1674  IF (SIZE(volgrid6d_coord_in%voldati,4) /= 1 .OR. &
1675  SIZE(volgrid6d_coord_in%voldati,5) /= 1) THEN
1676  CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1677  'volume providing constant input vertical coordinate must have &
1678  &only 1 time and 1 timerange')
1679  CALL init(volgrid6d_out) ! initialize to empty
1680  CALL raise_error()
1681  RETURN
1682  ENDIF
1683 
1684 ! search for variable providing vertical coordinate
1685  CALL get_val(this, output_levtype=output_levtype)
1686  vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
1687  IF (.NOT.c_e(vcoord_var)) THEN
1688  CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1689  'requested output level type '//t2c(output_levtype%level1)// &
1690  ' does not correspond to any known physical variable for &
1691  &providing vertical coordinate')
1692  CALL init(volgrid6d_out) ! initialize to empty
1693  CALL raise_error()
1694  RETURN
1695  ENDIF
1696 
1697  DO i = 1, SIZE(volgrid6d_coord_in%var)
1698  IF (convert(volgrid6d_coord_in%var(i)) == vcoord_var) THEN
1699  var_coord_in = i
1700  EXIT
1701  ENDIF
1702  ENDDO
1703 
1704  IF (.NOT.c_e(var_coord_in)) THEN
1705  CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1706  'volume providing constant input vertical coordinate contains no &
1707  &variables matching output level type '//t2c(output_levtype%level1))
1708  CALL init(volgrid6d_out) ! initialize to empty
1709  CALL raise_error()
1710  RETURN
1711  ENDIF
1712  CALL l4f_category_log(volgrid6d_in%category, l4f_info, &
1713  'Coordinate for vertint found in coord volume at position '// &
1714  t2c(var_coord_in))
1715 
1716 ! check horizontal grid
1717 ! this is too strict (component flag and so on)
1718 ! IF (volgrid6d_coord_in%griddim /= volgrid6d_in%griddim) THEN
1719  CALL get_val(volgrid6d_coord_in%griddim, nx=nxc, ny=nyc)
1720  CALL get_val(volgrid6d_in%griddim, nx=nxi, ny=nyi)
1721  IF (nxc /= nxi .OR. nyc /= nyi) THEN
1722  CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1723  'volume providing constant input vertical coordinate must have &
1724  &the same grid as the input')
1725  CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1726  'vertical coordinate: '//t2c(nxc)//'x'//t2c(nyc)// &
1727  ', input volume: '//t2c(nxi)//'x'//t2c(nyi))
1728  CALL init(volgrid6d_out) ! initialize to empty
1729  CALL raise_error()
1730  RETURN
1731  ENDIF
1732 
1733 ! check vertical coordinate system
1734  CALL get_val(this, input_levtype=input_levtype)
1735  mask_in = & ! implicit allocation
1736  (volgrid6d_coord_in%level(:)%level1 == input_levtype%level1) .AND. &
1737  (volgrid6d_coord_in%level(:)%level2 == input_levtype%level2)
1738  ulstart = firsttrue(mask_in)
1739  ulend = lasttrue(mask_in)
1740  IF (ulstart == 0 .OR. ulend == 0) THEN
1741  CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1742  'coordinate file does not contain levels of type '// &
1743  t2c(input_levtype%level1)//'/'//t2c(input_levtype%level2)// &
1744  ' specified for input data')
1745  CALL init(volgrid6d_out) ! initialize to empty
1746  CALL raise_error()
1747  RETURN
1748  ENDIF
1749 
1750  coord_3d_in = volgrid6d_coord_in%voldati(:,:,ulstart:ulend,1,1,var_coord_in) ! implicit allocation
1751 ! special case
1752  IF (output_levtype%level1 == 103 .OR. &
1753  output_levtype%level1 == 108) THEN ! surface coordinate needed
1754  spos = firsttrue(volgrid6d_coord_in%level(:) == vol7d_level_new(1))
1755  IF (spos == 0) THEN
1756  CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1757  'output level '//t2c(output_levtype%level1)// &
1758  ' requested, but height/press of surface not provided in coordinate file')
1759  CALL init(volgrid6d_out) ! initialize to empty
1760  CALL raise_error()
1761  RETURN
1762  ENDIF
1763  DO k = 1, SIZE(coord_3d_in,3)
1764  WHERE(c_e(coord_3d_in(:,:,k)) .AND. &
1765  c_e(volgrid6d_coord_in%voldati(:,:,spos,1,1,var_coord_in)))
1766  coord_3d_in(:,:,k) = coord_3d_in(:,:,k) - &
1767  volgrid6d_coord_in%voldati(:,:,spos,1,1,var_coord_in)
1768  ELSEWHERE
1769  coord_3d_in(:,:,k) = rmiss
1770  END WHERE
1771  ENDDO
1772  ENDIF
1773 
1774  ENDIF
1775  ENDIF
1776 
1777  IF (.NOT.c_e(var_coord_in)) THEN ! search for coordinate within volume
1778 ! search for variable providing vertical coordinate
1779  CALL get_val(this, output_levtype=output_levtype)
1780  vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
1781  IF (c_e(vcoord_var)) THEN
1782  DO i = 1, SIZE(volgrid6d_in%var)
1783  IF (convert(volgrid6d_in%var(i)) == vcoord_var) THEN
1784  var_coord_vol = i
1785  EXIT
1786  ENDIF
1787  ENDDO
1788 
1789  IF (c_e(var_coord_vol)) THEN
1790  CALL l4f_category_log(volgrid6d_in%category, l4f_info, &
1791  'Coordinate for vertint found in input volume at position '// &
1792  t2c(var_coord_vol))
1793  ENDIF
1794 
1795  ENDIF
1796  ENDIF
1797 
1798  CALL init(volgrid6d_out, griddim=volgrid6d_in%griddim, &
1799  time_definition=volgrid6d_in%time_definition, categoryappend=categoryappend)
1800  IF (c_e(var_coord_in)) THEN
1801  CALL init(grid_trans, this, lev_in=volgrid6d_in%level, lev_out=lev_out, &
1802  coord_3d_in=coord_3d_in, categoryappend=categoryappend)
1803  ELSE
1804  CALL init(grid_trans, this, lev_in=volgrid6d_in%level, lev_out=lev_out, &
1805  categoryappend=categoryappend)
1806  ENDIF
1807 
1808  CALL get_val(grid_trans, output_level_auto=llev_out) ! get levels if auto-generated
1809  IF (.NOT.ASSOCIATED(llev_out)) llev_out => lev_out
1810  nlevel = SIZE(llev_out)
1811  ELSE
1812  CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1813  'volgrid6d_transform: vertint requested but lev_out not provided')
1814  CALL init(volgrid6d_out) ! initialize to empty
1815  CALL raise_error()
1816  RETURN
1817  ENDIF
1818 
1819 ELSE
1820  CALL init(volgrid6d_out, griddim=griddim, &
1821  time_definition=volgrid6d_in%time_definition, categoryappend=categoryappend)
1822  CALL init(grid_trans, this, in=volgrid6d_in%griddim, out=volgrid6d_out%griddim, &
1823  maskgrid=maskgrid, maskbounds=maskbounds, categoryappend=categoryappend)
1824 ENDIF
1825 
1826 
1827 IF (c_e(grid_trans)) THEN ! transformation is valid
1828 
1829  CALL volgrid6d_alloc(volgrid6d_out, ntime=ntime, nlevel=nlevel, &
1830  ntimerange=ntimerange, nvar=nvar)
1831 
1832  IF (present(decode)) THEN ! explicitly set decode status
1833  ldecode = decode
1834  ELSE ! take it from input
1835  ldecode = ASSOCIATED(volgrid6d_in%voldati)
1836  ENDIF
1837 ! force decode if gaid is readonly
1838  decode_loop: DO i6 = 1,nvar
1839  DO i5 = 1, ntimerange
1840  DO i4 = 1, ntime
1841  DO i3 = 1, nlevel
1842  IF (c_e(volgrid6d_in%gaid(i3,i4,i5,i6))) THEN
1843  ldecode = ldecode .OR. grid_id_readonly(volgrid6d_in%gaid(i3,i4,i5,i6))
1844  EXIT decode_loop
1845  ENDIF
1846  ENDDO
1847  ENDDO
1848  ENDDO
1849  ENDDO decode_loop
1850 
1851  IF (present(decode)) THEN
1852  IF (ldecode.NEQV.decode) THEN
1853  CALL l4f_category_log(volgrid6d_in%category, l4f_warn, &
1854  'volgrid6d_transform: decode status forced to .TRUE. because driver does not allow copy')
1855  ENDIF
1856  ENDIF
1857 
1858  CALL volgrid6d_alloc_vol(volgrid6d_out, decode=ldecode)
1859 
1860 !ensure unproj was called
1861 !call griddim_unproj(volgrid6d_out%griddim)
1862 
1863  IF (trans_type == 'vertint') THEN
1864 #ifdef DEBUG
1865  CALL l4f_category_log(volgrid6d_in%category, l4f_debug, &
1866  "volgrid6d_transform: vertint to "//t2c(nlevel)//" levels")
1867 #endif
1868  CALL compute(grid_trans, volgrid6d_in, volgrid6d_out, lev_out=llev_out, &
1869  var_coord_vol=var_coord_vol, clone=clone)
1870  ELSE
1871  CALL compute(grid_trans, volgrid6d_in, volgrid6d_out, clone=clone)
1872  ENDIF
1873 
1874  IF (cf_out == 0) THEN ! unrotated components are desired
1875  CALL wind_unrot(volgrid6d_out) ! unrotate if necessary
1876  ELSE IF (cf_out == 1) THEN ! rotated components are desired
1877  CALL wind_rot(volgrid6d_out) ! rotate if necessary
1878  ENDIF
1879 
1880 ELSE
1881 ! should log with grid_trans%category, but it is private
1882  CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1883  'volgrid6d_transform: transformation not valid')
1884  CALL raise_error()
1885 ENDIF
1886 
1887 CALL delete(grid_trans)
1888 
1889 END SUBROUTINE volgrid6d_transform
1890 
1891 
1900 SUBROUTINE volgrid6dv_transform(this, griddim, volgrid6d_in, volgrid6d_out, &
1901  lev_out, volgrid6d_coord_in, maskgrid, maskbounds, clone, decode, categoryappend)
1902 TYPE(transform_def),INTENT(in) :: this
1903 type(griddim_def),INTENT(in),OPTIONAL :: griddim
1904 ! TODO ripristinare intent(in) dopo le opportune modifiche in grid_class.F90
1905 TYPE(volgrid6d),INTENT(inout) :: volgrid6d_in(:)
1906 type(volgrid6d),POINTER :: volgrid6d_out(:)
1907 type(vol7d_level),INTENT(in),OPTIONAL :: lev_out(:)
1908 type(volgrid6d),INTENT(in),OPTIONAL :: volgrid6d_coord_in
1909 REAL,INTENT(in),OPTIONAL :: maskgrid(:,:)
1910 REAL,INTENT(in),OPTIONAL :: maskbounds(:)
1911 LOGICAL,INTENT(in),OPTIONAL :: clone
1912 LOGICAL,INTENT(in),OPTIONAL :: decode
1913 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend
1914 
1915 INTEGER :: i, stallo
1916 
1917 
1918 allocate(volgrid6d_out(size(volgrid6d_in)),stat=stallo)
1919 if (stallo /= 0)then
1920  call l4f_log(l4f_fatal,"allocating memory")
1921  call raise_fatal_error()
1922 end if
1923 
1924 do i=1,size(volgrid6d_in)
1925  call transform(this, griddim, volgrid6d_in(i), volgrid6d_out(i), &
1926  lev_out=lev_out, volgrid6d_coord_in=volgrid6d_coord_in, &
1927  maskgrid=maskgrid, maskbounds=maskbounds, &
1928  clone=clone, decode=decode, categoryappend=categoryappend)
1929 end do
1930 
1931 END SUBROUTINE volgrid6dv_transform
1932 
1933 
1934 ! Internal method for performing grid to sparse point computations
1935 SUBROUTINE volgrid6d_v7d_transform_compute(this, volgrid6d_in, vol7d_out, &
1936  networkname, noconvert)
1937 TYPE(grid_transform),INTENT(in) :: this ! oggetto di trasformazione per grigliato
1938 type(volgrid6d), INTENT(in) :: volgrid6d_in ! oggetto da trasformare
1939 type(vol7d), INTENT(inout) :: vol7d_out ! oggetto trasformato
1940 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname ! imposta il network in vol7d_out (default='generic')
1941 LOGICAL,OPTIONAL,INTENT(in) :: noconvert
1942 
1943 INTEGER :: nntime, nana, ntime, ntimerange, nlevel, nvar, stallo
1944 INTEGER :: itime, itimerange, ivar, inetwork
1945 REAL,ALLOCATABLE :: voldatir_out(:,:,:)
1946 TYPE(conv_func),POINTER :: c_func(:)
1947 TYPE(datetime),ALLOCATABLE :: validitytime(:,:)
1948 REAL,POINTER :: voldatiin(:,:,:)
1949 
1950 #ifdef DEBUG
1951 call l4f_category_log(volgrid6d_in%category,l4f_debug,"start volgrid6d_v7d_transform_compute")
1952 #endif
1953 
1954 ntime=0
1955 ntimerange=0
1956 nlevel=0
1957 nvar=0
1958 nullify(c_func)
1959 
1960 if (present(networkname))then
1961  call init(vol7d_out%network(1),name=networkname)
1962 else
1963  call init(vol7d_out%network(1),name='generic')
1964 end if
1965 
1966 if (associated(volgrid6d_in%timerange))then
1967  ntimerange=size(volgrid6d_in%timerange)
1968  vol7d_out%timerange=volgrid6d_in%timerange
1969 end if
1970 
1971 if (associated(volgrid6d_in%time))then
1972  ntime=size(volgrid6d_in%time)
1973 
1974  if (vol7d_out%time_definition == volgrid6d_in%time_definition) then
1975 
1976  ! i time sono definiti uguali: assegno
1977  vol7d_out%time=volgrid6d_in%time
1978 
1979  else
1980  ! converto reference in validity
1981  allocate (validitytime(ntime,ntimerange),stat=stallo)
1982  if (stallo /=0)then
1983  call l4f_category_log(volgrid6d_in%category,l4f_fatal,"allocating memory")
1984  call raise_fatal_error()
1985  end if
1986 
1987  do itime=1,ntime
1988  do itimerange=1,ntimerange
1989  if (vol7d_out%time_definition > volgrid6d_in%time_definition) then
1990  validitytime(itime,itimerange) = &
1991  volgrid6d_in%time(itime) + timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
1992  else
1993  validitytime(itime,itimerange) = &
1994  volgrid6d_in%time(itime) - timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
1995  end if
1996  end do
1997  end do
1998 
1999  nntime = count_distinct(reshape(validitytime,(/ntime*ntimerange/)), back=.true.)
2000  vol7d_out%time=pack_distinct(reshape(validitytime,(/ntime*ntimerange/)), nntime,back=.true.)
2001 
2002  end if
2003 end if
2004 
2005 IF (ASSOCIATED(volgrid6d_in%level)) THEN
2006  nlevel = SIZE(volgrid6d_in%level)
2007  vol7d_out%level=volgrid6d_in%level
2008 ENDIF
2009 
2010 IF (ASSOCIATED(volgrid6d_in%var)) THEN
2011  nvar = SIZE(volgrid6d_in%var)
2012  IF (.NOT. optio_log(noconvert)) THEN
2013  CALL vargrib2varbufr(volgrid6d_in%var, vol7d_out%dativar%r, c_func)
2014  ENDIF
2015 ENDIF
2016 
2017 nana = SIZE(vol7d_out%ana)
2018 
2019 ! allocate once for speed
2020 IF (.NOT.ASSOCIATED(volgrid6d_in%voldati)) THEN
2021  ALLOCATE(voldatiin(volgrid6d_in%griddim%dim%nx, volgrid6d_in%griddim%dim%ny, &
2022  nlevel))
2023 ENDIF
2024 
2025 ALLOCATE(voldatir_out(nana,1,nlevel),stat=stallo)
2026 IF (stallo /= 0) THEN
2027  CALL l4f_category_log(volgrid6d_in%category,l4f_fatal,"allocating memory")
2028  CALL raise_fatal_error()
2029 ENDIF
2030 
2031 inetwork=1
2032 do itime=1,ntime
2033  do itimerange=1,ntimerange
2034 ! do ilevel=1,nlevel
2035  do ivar=1,nvar
2036 
2037  !non è chiaro se questa sezione è utile o no
2038  !ossia il compute sotto sembra prevedere voldatir_out solo in out
2039 !!$ if (vol7d_out%time_definition == volgrid6d_in%time_definition) then
2040 !!$ voldatir_out=reshape(vol7d_out%voldatir(:,itime,ilevel,itimerange,ivar,inetwork),(/nana,1/))
2041 !!$ else
2042 !!$ voldatir_out=reshape(vol7d_out%voldatir(:,index(vol7d_out%time,validitytime(itime,ilevel)),ilevel,itimerange,ivar,inetwork),(/nana,1/))
2043 !!$ end if
2044 
2045  CALL volgrid_get_vol_3d(volgrid6d_in, itime, itimerange, ivar, &
2046  voldatiin)
2047 
2048  CALL compute(this, voldatiin, voldatir_out, vol7d_out%dativar%r(ivar))
2049 
2050  if (vol7d_out%time_definition == volgrid6d_in%time_definition) then
2051  vol7d_out%voldatir(:,itime,:,itimerange,ivar,inetwork) = &
2052  voldatir_out(:,1,:)
2053  else
2054  vol7d_out%voldatir(:,index(vol7d_out%time,validitytime(itime,itimerange)),:,itimerange,ivar,inetwork)=&
2055  reshape(voldatir_out,(/nana,nlevel/))
2056  end if
2057 
2058 ! 1 indice della dimensione "anagrafica"
2059 ! 2 indice della dimensione "tempo"
2060 ! 3 indice della dimensione "livello verticale"
2061 ! 4 indice della dimensione "intervallo temporale"
2062 ! 5 indice della dimensione "variabile"
2063 ! 6 indice della dimensione "rete"
2064 
2065  end do
2066 ! end do
2067  end do
2068 end do
2069 
2070 deallocate(voldatir_out)
2071 IF (.NOT.ASSOCIATED(volgrid6d_in%voldati)) THEN
2072  DEALLOCATE(voldatiin)
2073 ENDIF
2074 if (allocated(validitytime)) deallocate(validitytime)
2075 
2076 ! Rescale valid data according to variable conversion table
2077 IF (ASSOCIATED(c_func)) THEN
2078  DO ivar = 1, nvar
2079  CALL compute(c_func(ivar), vol7d_out%voldatir(:,:,:,:,ivar,:))
2080  ENDDO
2081  DEALLOCATE(c_func)
2082 ENDIF
2083 
2084 end SUBROUTINE volgrid6d_v7d_transform_compute
2085 
2086 
2093 SUBROUTINE volgrid6d_v7d_transform(this, volgrid6d_in, vol7d_out, v7d, &
2094  maskgrid, maskbounds, networkname, noconvert, categoryappend)
2095 TYPE(transform_def),INTENT(in) :: this
2096 type(volgrid6d),INTENT(inout) :: volgrid6d_in
2097 type(vol7d),INTENT(out) :: vol7d_out
2098 type(vol7d),INTENT(in),OPTIONAL :: v7d
2099 REAL,INTENT(in),OPTIONAL :: maskgrid(:,:)
2100 REAL,INTENT(in),OPTIONAL :: maskbounds(:)
2101 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname
2102 LOGICAL,OPTIONAL,INTENT(in) :: noconvert
2103 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend
2104 
2105 type(grid_transform) :: grid_trans
2106 INTEGER :: ntime, ntimerange, nlevel, nvar, nana, time_definition, nnetwork, stallo
2107 INTEGER :: itime, itimerange, inetwork
2108 TYPE(datetime),ALLOCATABLE :: validitytime(:,:)
2109 INTEGER,ALLOCATABLE :: point_index(:)
2110 TYPE(vol7d) :: v7d_locana
2111 
2112 #ifdef DEBUG
2113 call l4f_category_log(volgrid6d_in%category,l4f_debug,"start volgrid6d_v7d_transform")
2114 #endif
2115 
2116 call vg6d_wind_unrot(volgrid6d_in)
2117 
2118 ntime=0
2119 ntimerange=0
2120 nlevel=0
2121 nvar=0
2122 nnetwork=1
2123 
2124 call get_val(this,time_definition=time_definition)
2125 if (.not. c_e(time_definition)) then
2126  time_definition=1 ! default to validity time
2127 endif
2128 
2129 IF (present(v7d)) THEN
2130  CALL vol7d_copy(v7d, v7d_locana)
2131 ELSE
2132  CALL init(v7d_locana)
2133 ENDIF
2134 
2135 if (associated(volgrid6d_in%timerange)) ntimerange=size(volgrid6d_in%timerange)
2136 
2137 if (associated(volgrid6d_in%time)) then
2138 
2139  ntime=size(volgrid6d_in%time)
2140 
2141  if (time_definition /= volgrid6d_in%time_definition) then
2142 
2143  ! converto reference in validity
2144  allocate (validitytime(ntime,ntimerange),stat=stallo)
2145  if (stallo /=0)then
2146  call l4f_category_log(volgrid6d_in%category,l4f_fatal,"allocating memory")
2147  call raise_fatal_error()
2148  end if
2149 
2150  do itime=1,ntime
2151  do itimerange=1,ntimerange
2152  if (time_definition > volgrid6d_in%time_definition) then
2153  validitytime(itime,itimerange) = &
2154  volgrid6d_in%time(itime) + timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
2155  else
2156  validitytime(itime,itimerange) = &
2157  volgrid6d_in%time(itime) - timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
2158  end if
2159  end do
2160  end do
2161 
2162  ntime = count_distinct(reshape(validitytime,(/ntime*ntimerange/)), back=.true.)
2163  deallocate (validitytime)
2164 
2165  end if
2166 end if
2167 
2168 
2169 if (associated(volgrid6d_in%level)) nlevel=size(volgrid6d_in%level)
2170 if (associated(volgrid6d_in%var)) nvar=size(volgrid6d_in%var)
2171 
2172 CALL init(grid_trans, this, volgrid6d_in%griddim, v7d_locana, &
2173  maskgrid=maskgrid, maskbounds=maskbounds, categoryappend=categoryappend)
2174 CALL init(vol7d_out,time_definition=time_definition)
2175 
2176 IF (c_e(grid_trans)) THEN
2177 
2178  nana=SIZE(v7d_locana%ana)
2179  CALL vol7d_alloc(vol7d_out, nana=nana, ntime=ntime, nlevel=nlevel, &
2180  ntimerange=ntimerange, ndativarr=nvar, nnetwork=nnetwork)
2181  vol7d_out%ana = v7d_locana%ana
2182 
2183  CALL get_val(grid_trans, output_point_index=point_index)
2184  IF (ALLOCATED(point_index)) THEN
2185 ! check that size(point_index) == nana?
2186  CALL vol7d_alloc(vol7d_out, nanavari=1)
2187  CALL init(vol7d_out%anavar%i(1), 'B01192')
2188  ENDIF
2189 
2190  CALL vol7d_alloc_vol(vol7d_out)
2191 
2192  IF (ALLOCATED(point_index)) THEN
2193  DO inetwork = 1, nnetwork
2194  vol7d_out%volanai(:,1,inetwork) = point_index(:)
2195  ENDDO
2196  ENDIF
2197  CALL compute(grid_trans, volgrid6d_in, vol7d_out, networkname, noconvert)
2198 ELSE
2199  CALL l4f_log(l4f_error, 'vg6d_v7d_transform: transformation not valid')
2200  CALL raise_error()
2201 ENDIF
2202 
2203 CALL delete(grid_trans)
2204 
2205 #ifdef HAVE_DBALLE
2206 ! set variables to a conformal state
2207 CALL vol7d_dballe_set_var_du(vol7d_out)
2208 #endif
2209 
2210 CALL delete(v7d_locana)
2211 
2212 END SUBROUTINE volgrid6d_v7d_transform
2213 
2214 
2223 SUBROUTINE volgrid6dv_v7d_transform(this, volgrid6d_in, vol7d_out, v7d, &
2224  maskgrid, maskbounds, networkname, noconvert, categoryappend)
2225 TYPE(transform_def),INTENT(in) :: this
2226 type(volgrid6d),INTENT(inout) :: volgrid6d_in(:)
2227 type(vol7d),INTENT(out) :: vol7d_out
2228 type(vol7d),INTENT(in),OPTIONAL :: v7d
2229 REAL,INTENT(in),OPTIONAL :: maskgrid(:,:)
2230 REAL,INTENT(in),OPTIONAL :: maskbounds(:)
2231 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname
2232 LOGICAL,OPTIONAL,INTENT(in) :: noconvert
2233 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend
2234 
2235 integer :: i
2236 TYPE(vol7d) :: v7dtmp
2237 
2238 
2239 CALL init(v7dtmp)
2240 CALL init(vol7d_out)
2241 
2242 DO i=1,SIZE(volgrid6d_in)
2243  CALL transform(this, volgrid6d_in(i), v7dtmp, v7d=v7d, &
2244  maskgrid=maskgrid, maskbounds=maskbounds, &
2245  networkname=networkname, noconvert=noconvert, categoryappend=categoryappend)
2246  CALL vol7d_append(vol7d_out, v7dtmp)
2247 ENDDO
2248 
2249 END SUBROUTINE volgrid6dv_v7d_transform
2250 
2251 
2252 ! Internal method for performing sparse point to grid computations
2253 SUBROUTINE v7d_volgrid6d_transform_compute(this, vol7d_in, volgrid6d_out, networkname, gaid_template)
2254 TYPE(grid_transform),INTENT(in) :: this ! object specifying the specific transformation
2255 type(vol7d), INTENT(in) :: vol7d_in ! object to be transformed
2256 type(volgrid6d), INTENT(inout) :: volgrid6d_out ! transformed object
2257 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname ! select the network to be processed from the \a vol7d input object, default the first network
2258 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 
2260 integer :: nana, ntime, ntimerange, nlevel, nvar
2261 INTEGER :: ilevel, itime, itimerange, ivar, inetwork
2262 
2263 REAL,POINTER :: voldatiout(:,:,:)
2264 type(vol7d_network) :: network
2265 TYPE(conv_func), pointer :: c_func(:)
2266 !TODO category sarebbe da prendere da vol7d
2267 #ifdef DEBUG
2268 call l4f_category_log(volgrid6d_out%category,l4f_debug,"start v7d_volgrid6d_transform_compute")
2269 #endif
2270 
2271 ntime=0
2272 ntimerange=0
2273 nlevel=0
2274 nvar=0
2275 
2276 if (present(networkname))then
2277  call init(network,name=networkname)
2278  inetwork= index(vol7d_in%network,network)
2279 else
2280  inetwork=1
2281 end if
2282 
2283 ! no time_definition conversion implemented, output will be the same as input
2284 if (associated(vol7d_in%time))then
2285  ntime=size(vol7d_in%time)
2286  volgrid6d_out%time=vol7d_in%time
2287 end if
2288 
2289 if (associated(vol7d_in%timerange))then
2290  ntimerange=size(vol7d_in%timerange)
2291  volgrid6d_out%timerange=vol7d_in%timerange
2292 end if
2293 
2294 if (associated(vol7d_in%level))then
2295  nlevel=size(vol7d_in%level)
2296  volgrid6d_out%level=vol7d_in%level
2297 end if
2298 
2299 if (associated(vol7d_in%dativar%r))then
2300  nvar=size(vol7d_in%dativar%r)
2301  CALL varbufr2vargrib(vol7d_in%dativar%r, volgrid6d_out%var, c_func, gaid_template)
2302 end if
2303 
2304 nana=SIZE(vol7d_in%voldatir, 1)
2305 ! allocate once for speed
2306 IF (.NOT.ASSOCIATED(volgrid6d_out%voldati)) THEN
2307  ALLOCATE(voldatiout(volgrid6d_out%griddim%dim%nx, volgrid6d_out%griddim%dim%ny, &
2308  nlevel))
2309 ENDIF
2310 
2311 DO ivar=1,nvar
2312  DO itimerange=1,ntimerange
2313  DO itime=1,ntime
2314 
2315 ! clone the gaid template where I have data
2316  IF (present(gaid_template)) THEN
2317  DO ilevel = 1, nlevel
2318  IF (any(c_e(vol7d_in%voldatir(:,itime,ilevel,itimerange,ivar,inetwork)))) THEN
2319  CALL copy(gaid_template, volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))
2320  ELSE
2321  volgrid6d_out%gaid(ilevel,itime,itimerange,ivar) = grid_id_new()
2322  ENDIF
2323  ENDDO
2324  ENDIF
2325 
2326 ! get data
2327  IF (ASSOCIATED(volgrid6d_out%voldati)) & ! improve!!!!
2328  CALL volgrid_get_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
2329  voldatiout)
2330 ! do the interpolation
2331  CALL compute(this, &
2332  vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), voldatiout, &
2333  vol7d_in%dativar%r(ivar))
2334 ! rescale valid data according to variable conversion table
2335  IF (ASSOCIATED(c_func)) THEN
2336  CALL compute(c_func(ivar), voldatiout(:,:,:))
2337  ENDIF
2338 ! put data
2339  CALL volgrid_set_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
2340  voldatiout)
2341 
2342  ENDDO
2343  ENDDO
2344 ENDDO
2345 
2346 IF (.NOT.ASSOCIATED(volgrid6d_out%voldati)) THEN
2347  DEALLOCATE(voldatiout)
2348 ENDIF
2349 IF (ASSOCIATED(c_func)) THEN
2350  DEALLOCATE(c_func)
2351 ENDIF
2352 
2353 END SUBROUTINE v7d_volgrid6d_transform_compute
2354 
2355 
2362 SUBROUTINE v7d_volgrid6d_transform(this, griddim, vol7d_in, volgrid6d_out, &
2363  networkname, gaid_template, categoryappend)
2364 TYPE(transform_def),INTENT(in) :: this
2365 type(griddim_def),INTENT(in),OPTIONAL :: griddim
2366 ! TODO ripristinare intent(in) dopo le opportune modifiche in grid_class.F90
2367 TYPE(vol7d),INTENT(inout) :: vol7d_in
2368 type(volgrid6d),INTENT(out) :: volgrid6d_out
2369 CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname
2370 type(grid_id),OPTIONAL,INTENT(in) :: gaid_template
2371 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend
2372 
2373 type(grid_transform) :: grid_trans
2374 integer :: ntime, ntimerange, nlevel, nvar
2375 
2376 
2377 !TODO la category sarebbe da prendere da vol7d
2378 !call l4f_category_log(vol7d_out%category,L4F_DEBUG,"start volgrid6d_transform")
2379 
2380 CALL vol7d_alloc_vol(vol7d_in) ! be safe
2381 ntime=SIZE(vol7d_in%time)
2382 ntimerange=SIZE(vol7d_in%timerange)
2383 nlevel=SIZE(vol7d_in%level)
2384 nvar=0
2385 if (associated(vol7d_in%dativar%r)) nvar=size(vol7d_in%dativar%r)
2386 
2387 IF (nvar <= 0) THEN ! use vol7d category once it will be implemented
2388  CALL l4f_log(l4f_error, &
2389  "trying to transform a vol7d object incomplete or without real variables")
2390  CALL init(volgrid6d_out) ! initialize to empty
2391  CALL raise_error()
2392  RETURN
2393 ENDIF
2394 
2395 CALL init(grid_trans, this, vol7d_in, griddim, categoryappend=categoryappend)
2396 CALL init(volgrid6d_out, griddim, time_definition=vol7d_in%time_definition, &
2397  categoryappend=categoryappend)
2398 
2399 IF (c_e(grid_trans)) THEN
2400 
2401  CALL volgrid6d_alloc(volgrid6d_out, griddim%dim, ntime=ntime, nlevel=nlevel, &
2402  ntimerange=ntimerange, nvar=nvar)
2403 ! can I avoid decode=.TRUE. here, using gaid_template?
2404  CALL volgrid6d_alloc_vol(volgrid6d_out, decode=.true.)
2405 
2406  CALL compute(grid_trans, vol7d_in, volgrid6d_out, networkname, gaid_template)
2407 
2408  CALL vg6d_wind_rot(volgrid6d_out)
2409 ELSE
2410 ! should log with grid_trans%category, but it is private
2411  CALL l4f_log(l4f_error, 'v7d_vg6d_transform: transformation not valid')
2412  CALL raise_error()
2413 ENDIF
2414 
2415 CALL delete(grid_trans)
2416 
2417 END SUBROUTINE v7d_volgrid6d_transform
2418 
2419 
2420 ! Internal method for performing sparse point to sparse point computations
2421 SUBROUTINE v7d_v7d_transform_compute(this, vol7d_in, vol7d_out, lev_out, &
2422  var_coord_vol)
2423 TYPE(grid_transform),INTENT(in) :: this ! oggetto di trasformazione per grigliato
2424 type(vol7d), INTENT(in) :: vol7d_in ! oggetto da trasformare
2425 type(vol7d), INTENT(inout) :: vol7d_out ! oggetto trasformato
2426 TYPE(vol7d_level),INTENT(in),OPTIONAL :: lev_out(:) ! vol7d_level object defining target vertical grid, for vertical interpolations
2427 INTEGER,INTENT(in),OPTIONAL :: var_coord_vol ! index of variable defining vertical coordinate values in input volume
2428 
2429 INTEGER :: itime, itimerange, ilevel, ivar, inetwork, &
2430  levshift, levused, lvar_coord_vol, spos
2431 REAL,ALLOCATABLE :: coord_3d_in(:,:,:)
2432 TYPE(vol7d_level) :: output_levtype
2433 
2434 lvar_coord_vol = optio_i(var_coord_vol)
2435 vol7d_out%time(:) = vol7d_in%time(:)
2436 vol7d_out%timerange(:) = vol7d_in%timerange(:)
2437 IF (present(lev_out)) THEN
2438  vol7d_out%level(:) = lev_out(:)
2439 ELSE
2440  vol7d_out%level(:) = vol7d_in%level(:)
2441 ENDIF
2442 vol7d_out%network(:) = vol7d_in%network(:)
2443 IF (ASSOCIATED(vol7d_in%dativar%r)) THEN ! work only when real vars are available
2444  vol7d_out%dativar%r(:) = vol7d_in%dativar%r(:)
2445 
2446  CALL get_val(this, levshift=levshift, levused=levused)
2447  spos = imiss
2448  IF (c_e(lvar_coord_vol)) THEN
2449  CALL get_val(this%trans, output_levtype=output_levtype)
2450  IF (output_levtype%level1 == 103 .OR. output_levtype%level1 == 108) THEN
2451  spos = firsttrue(vol7d_in%level(:) == vol7d_level_new(1))
2452  IF (spos == 0) THEN
2453  CALL l4f_log(l4f_error, &
2454  'output level '//t2c(output_levtype%level1)// &
2455  ' requested, but height/press of surface not provided in volume')
2456  ENDIF
2457  IF (.NOT.c_e(levshift) .AND. .NOT.c_e(levused)) THEN
2458  CALL l4f_log(l4f_error, &
2459  'internal inconsistence, levshift and levused undefined when they should be')
2460  ENDIF
2461  ALLOCATE(coord_3d_in(SIZE(vol7d_in%ana),1,SIZE(vol7d_in%level)))
2462  ENDIF
2463 
2464  ENDIF
2465 
2466  DO inetwork = 1, SIZE(vol7d_in%network)
2467  DO ivar = 1, SIZE(vol7d_in%dativar%r)
2468  DO itimerange = 1, SIZE(vol7d_in%timerange)
2469  DO itime = 1, SIZE(vol7d_in%time)
2470 
2471 ! dirty trick to make voldatir look like a 2d-array of shape (nana,1)
2472  IF (c_e(lvar_coord_vol)) THEN
2473  IF (c_e(spos)) THEN ! compute difference wrt surface coordinate
2474  IF (spos == 0) THEN ! error condition, set all to missing and goodnight
2475  coord_3d_in(:,:,levshift+1:levshift+levused) = rmiss
2476  ELSE
2477  DO ilevel = levshift+1, levshift+levused
2478  WHERE(c_e(vol7d_in%voldatir(:,itime:itime,ilevel,itimerange,lvar_coord_vol,inetwork)) .AND. &
2479  c_e(vol7d_in%voldatir(:,itime:itime,spos,itimerange,lvar_coord_vol,inetwork)))
2480  coord_3d_in(:,:,ilevel) = vol7d_in%voldatir(:,itime:itime,ilevel,itimerange,lvar_coord_vol,inetwork) - &
2481  vol7d_in%voldatir(:,itime:itime,spos,itimerange,lvar_coord_vol,inetwork)
2482  ELSEWHERE
2483  coord_3d_in(:,:,ilevel) = rmiss
2484  END WHERE
2485  ENDDO
2486  ENDIF
2487  CALL compute(this, &
2488  vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), &
2489  vol7d_out%voldatir(:,itime:itime,:,itimerange,ivar,inetwork), &
2490  var=vol7d_in%dativar%r(ivar), &
2491  coord_3d_in=coord_3d_in)
2492  ELSE
2493  CALL compute(this, &
2494  vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), &
2495  vol7d_out%voldatir(:,itime:itime,:,itimerange,ivar,inetwork), &
2496  var=vol7d_in%dativar%r(ivar), &
2497  coord_3d_in=vol7d_in%voldatir(:,itime:itime,:,itimerange, &
2498  lvar_coord_vol,inetwork))
2499  ENDIF
2500  ELSE
2501  CALL compute(this, &
2502  vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), &
2503  vol7d_out%voldatir(:,itime:itime,:,itimerange,ivar,inetwork), &
2504  var=vol7d_in%dativar%r(ivar))
2505  ENDIF
2506  ENDDO
2507  ENDDO
2508  ENDDO
2509  ENDDO
2510 
2511 ENDIF
2512 
2513 END SUBROUTINE v7d_v7d_transform_compute
2514 
2515 
2523 SUBROUTINE v7d_v7d_transform(this, vol7d_in, vol7d_out, v7d, lev_out, vol7d_coord_in, &
2524  categoryappend)
2525 TYPE(transform_def),INTENT(in) :: this
2526 type(vol7d),INTENT(inout) :: vol7d_in
2527 type(vol7d),INTENT(out) :: vol7d_out
2528 type(vol7d),INTENT(in),OPTIONAL :: v7d
2529 type(vol7d_level),INTENT(in),OPTIONAL,TARGET :: lev_out(:)
2530 type(vol7d),INTENT(in),OPTIONAL :: vol7d_coord_in
2531 CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend
2532 
2533 INTEGER :: nvar, inetwork
2534 TYPE(grid_transform) :: grid_trans
2535 TYPE(vol7d_level),POINTER :: llev_out(:)
2536 TYPE(vol7d_level) :: input_levtype, output_levtype
2537 TYPE(vol7d_var) :: vcoord_var
2538 REAL,ALLOCATABLE :: coord_3d_in(:,:,:)
2539 INTEGER :: var_coord_in, var_coord_vol, i, k, ulstart, ulend, spos
2540 INTEGER,ALLOCATABLE :: point_index(:)
2541 TYPE(vol7d) :: v7d_locana, vol7d_tmpana
2542 CHARACTER(len=80) :: trans_type
2543 LOGICAL,ALLOCATABLE :: mask_in(:), point_mask(:)
2544 
2545 CALL vol7d_alloc_vol(vol7d_in) ! be safe
2546 nvar=0
2547 IF (ASSOCIATED(vol7d_in%dativar%r)) nvar=SIZE(vol7d_in%dativar%r)
2548 
2549 CALL init(v7d_locana)
2550 IF (present(v7d)) v7d_locana = v7d
2551 CALL init(vol7d_out, time_definition=vol7d_in%time_definition)
2552 
2553 CALL get_val(this, trans_type=trans_type)
2554 
2555 var_coord_vol = imiss
2556 IF (trans_type == 'vertint') THEN
2557 
2558  IF (present(lev_out)) THEN
2559 
2560 ! if vol7d_coord_in provided and allocated, check that it fits
2561  var_coord_in = -1
2562  IF (present(vol7d_coord_in)) THEN
2563  IF (ASSOCIATED(vol7d_coord_in%voldatir) .AND. &
2564  ASSOCIATED(vol7d_coord_in%dativar%r)) THEN
2565 
2566 ! strictly 1 time, 1 timerange and 1 network
2567  IF (SIZE(vol7d_coord_in%voldatir,2) /= 1 .OR. &
2568  SIZE(vol7d_coord_in%voldatir,4) /= 1 .OR. &
2569  SIZE(vol7d_coord_in%voldatir,6) /= 1) THEN
2570  CALL l4f_log(l4f_error, &
2571  'volume providing constant input vertical coordinate must have &
2572  &only 1 time, 1 timerange and 1 network')
2573  CALL raise_error()
2574  RETURN
2575  ENDIF
2576 
2577 ! search for variable providing vertical coordinate
2578  CALL get_val(this, output_levtype=output_levtype)
2579  vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
2580  IF (.NOT.c_e(vcoord_var)) THEN
2581  CALL l4f_log(l4f_error, &
2582  'requested output level type '//t2c(output_levtype%level1)// &
2583  ' does not correspond to any known physical variable for &
2584  &providing vertical coordinate')
2585  CALL raise_error()
2586  RETURN
2587  ENDIF
2588 
2589  var_coord_in = index(vol7d_coord_in%dativar%r, vcoord_var)
2590 
2591  IF (var_coord_in <= 0) THEN
2592  CALL l4f_log(l4f_error, &
2593  'volume providing constant input vertical coordinate contains no &
2594  &real variables matching output level type '//t2c(output_levtype%level1))
2595  CALL raise_error()
2596  RETURN
2597  ENDIF
2598  CALL l4f_log(l4f_info, &
2599  'Coordinate for vertint found in coord volume at position '// &
2600  t2c(var_coord_in))
2601 
2602 ! check vertical coordinate system
2603  CALL get_val(this, input_levtype=input_levtype)
2604  mask_in = & ! implicit allocation
2605  (vol7d_coord_in%level(:)%level1 == input_levtype%level1) .AND. &
2606  (vol7d_coord_in%level(:)%level2 == input_levtype%level2)
2607  ulstart = firsttrue(mask_in)
2608  ulend = lasttrue(mask_in)
2609  IF (ulstart == 0 .OR. ulend == 0) THEN
2610  CALL l4f_log(l4f_error, &
2611  'coordinate file does not contain levels of type '// &
2612  t2c(input_levtype%level1)//'/'//t2c(input_levtype%level2)// &
2613  ' specified for input data')
2614  CALL raise_error()
2615  RETURN
2616  ENDIF
2617 
2618  coord_3d_in = vol7d_coord_in%voldatir(:,1:1,ulstart:ulend,1,var_coord_in,1) ! dirty 1:1, implicit allocation
2619 ! special case
2620  IF (output_levtype%level1 == 103 &
2621  .OR. output_levtype%level1 == 108) THEN ! surface coordinate needed
2622  spos = firsttrue(vol7d_coord_in%level(:) == vol7d_level_new(1))
2623  IF (spos == 0) THEN
2624  CALL l4f_log(l4f_error, &
2625  'output level '//t2c(output_levtype%level1)// &
2626  ' requested, but height/press of surface not provided in coordinate file')
2627  CALL raise_error()
2628  RETURN
2629  ENDIF
2630  DO k = 1, SIZE(coord_3d_in,3)
2631  WHERE(c_e(coord_3d_in(:,:,k)) .AND. &
2632  c_e(vol7d_coord_in%voldatir(:,1:1,spos,1,var_coord_in,1)))
2633  coord_3d_in(:,:,k) = coord_3d_in(:,:,k) - &
2634  vol7d_coord_in%voldatir(:,1:1,spos,1,var_coord_in,1)
2635  ELSEWHERE
2636  coord_3d_in(:,:,k) = rmiss
2637  END WHERE
2638  ENDDO
2639  ENDIF
2640 
2641  ENDIF
2642  ENDIF
2643 
2644  IF (var_coord_in <= 0) THEN ! search for coordinate within volume
2645 ! search for variable providing vertical coordinate
2646  CALL get_val(this, output_levtype=output_levtype)
2647  vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
2648  IF (c_e(vcoord_var)) THEN
2649  DO i = 1, SIZE(vol7d_in%dativar%r)
2650  IF (vol7d_in%dativar%r(i) == vcoord_var) THEN
2651  var_coord_vol = i
2652  EXIT
2653  ENDIF
2654  ENDDO
2655 
2656  IF (c_e(var_coord_vol)) THEN
2657  CALL l4f_log(l4f_info, &
2658  'Coordinate for vertint found in input volume at position '// &
2659  t2c(var_coord_vol))
2660  ENDIF
2661 
2662  ENDIF
2663  ENDIF
2664 
2665  IF (var_coord_in > 0) THEN
2666  CALL init(grid_trans, this, lev_in=vol7d_in%level, lev_out=lev_out, &
2667  coord_3d_in=coord_3d_in, categoryappend=categoryappend)
2668  ELSE
2669  CALL init(grid_trans, this, lev_in=vol7d_in%level, lev_out=lev_out, &
2670  categoryappend=categoryappend)
2671  ENDIF
2672 
2673  CALL get_val(grid_trans, output_level_auto=llev_out) ! get levels if auto-generated
2674  IF (.NOT.associated(llev_out)) llev_out => lev_out
2675 
2676  IF (c_e(grid_trans)) then! .AND. nvar > 0) THEN
2677 
2678  CALL vol7d_alloc(vol7d_out, nana=SIZE(vol7d_in%ana), &
2679  ntime=SIZE(vol7d_in%time), ntimerange=SIZE(vol7d_in%timerange), &
2680  nlevel=SIZE(llev_out), nnetwork=SIZE(vol7d_in%network), ndativarr=nvar)
2681  vol7d_out%ana(:) = vol7d_in%ana(:)
2682 
2683  CALL vol7d_alloc_vol(vol7d_out)
2684 
2685 ! no need to check c_e(var_coord_vol) here since the presence of
2686 ! this%coord_3d_in (external) has precedence over coord_3d_in internal
2687 ! in grid_transform_compute
2688  CALL compute(grid_trans, vol7d_in, vol7d_out, llev_out, &
2689  var_coord_vol=var_coord_vol)
2690  ELSE
2691  CALL l4f_log(l4f_error, 'v7d_v7d_transform: transformation not valid')
2692  CALL raise_error()
2693  ENDIF
2694  ELSE
2695  CALL l4f_log(l4f_error, &
2696  'v7d_v7d_transform: vertint requested but lev_out not provided')
2697  CALL raise_error()
2698  ENDIF
2699 
2700 ELSE
2701 
2702  CALL init(grid_trans, this, vol7d_in, v7d_locana, &
2703  categoryappend=categoryappend)
2704 ! if this init is successful, I am sure that v7d_locana%ana is associated
2705 
2706  IF (c_e(grid_trans)) then! .AND. nvar > 0) THEN
2707 
2708  CALL vol7d_alloc(vol7d_out, nana=SIZE(v7d_locana%ana), &
2709  ntime=SIZE(vol7d_in%time), ntimerange=SIZE(vol7d_in%timerange), &
2710  nlevel=SIZE(vol7d_in%level), nnetwork=SIZE(vol7d_in%network), ndativarr=nvar)
2711  vol7d_out%ana = v7d_locana%ana
2712 
2713  CALL get_val(grid_trans, point_mask=point_mask, output_point_index=point_index)
2714 
2715  IF (ALLOCATED(point_index)) THEN
2716  CALL vol7d_alloc(vol7d_out, nanavari=1)
2717  CALL init(vol7d_out%anavar%i(1), 'B01192')
2718  ENDIF
2719 
2720  CALL vol7d_alloc_vol(vol7d_out)
2721 
2722  IF (ALLOCATED(point_index)) THEN
2723  DO inetwork = 1, SIZE(vol7d_in%network)
2724  vol7d_out%volanai(:,1,inetwork) = point_index(:)
2725  ENDDO
2726  ENDIF
2727  CALL compute(grid_trans, vol7d_in, vol7d_out)
2728 
2729  IF (ALLOCATED(point_mask)) THEN ! keep full ana
2730  IF (SIZE(point_mask) /= SIZE(vol7d_in%ana)) THEN
2731  CALL l4f_log(l4f_warn, &
2732  'v7d_v7d_transform: inconsistency in point size: '//t2c(SIZE(point_mask)) &
2733  //':'//t2c(SIZE(vol7d_in%ana)))
2734  ELSE
2735 #ifdef DEBUG
2736  CALL l4f_log(l4f_debug, 'v7d_v7d_transform: merging ana from in to out')
2737 #endif
2738  CALL vol7d_copy(vol7d_in, vol7d_tmpana, &
2739  lana=point_mask, lnetwork=(/.true./), &
2740  ltime=(/.false./), ltimerange=(/.false./), llevel=(/.false./))
2741  CALL vol7d_append(vol7d_out, vol7d_tmpana)
2742  ENDIF
2743  ENDIF
2744 
2745  ELSE
2746  CALL l4f_log(l4f_error, 'v7d_v7d_transform: transformation not valid')
2747  CALL raise_error()
2748  ENDIF
2749 
2750 ENDIF
2751 
2752 CALL delete(grid_trans)
2753 IF (.NOT. present(v7d)) CALL delete(v7d_locana)
2754 
2755 END SUBROUTINE v7d_v7d_transform
2756 
2757 
2765 subroutine vg6d_wind_unrot(this)
2766 type(volgrid6d) :: this
2767 
2768 integer :: component_flag
2769 
2770 call get_val(this%griddim,component_flag=component_flag)
2771 
2772 if (component_flag == 1) then
2773  call l4f_category_log(this%category,l4f_info, &
2774  "unrotating vector components")
2775  call vg6d_wind__un_rot(this,.false.)
2776  call set_val(this%griddim,component_flag=0)
2777 else
2778  call l4f_category_log(this%category,l4f_info, &
2779  "no need to unrotate vector components")
2780 end if
2781 
2782 end subroutine vg6d_wind_unrot
2783 
2784 
2790 subroutine vg6d_wind_rot(this)
2791 type(volgrid6d) :: this
2792 
2793 integer :: component_flag
2794 
2795 call get_val(this%griddim,component_flag=component_flag)
2796 
2797 if (component_flag == 0) then
2798  call l4f_category_log(this%category,l4f_info, &
2799  "rotating vector components")
2800  call vg6d_wind__un_rot(this,.true.)
2801  call set_val(this%griddim,component_flag=1)
2802 else
2803  call l4f_category_log(this%category,l4f_info, &
2804  "no need to rotate vector components")
2805 end if
2806 
2807 end subroutine vg6d_wind_rot
2808 
2809 
2810 ! Generic UnRotate the wind components.
2811 SUBROUTINE vg6d_wind__un_rot(this,rot)
2812 TYPE(volgrid6d) :: this ! object containing wind to be (un)rotated
2813 LOGICAL :: rot ! if .true. rotate else unrotate
2814 
2815 INTEGER :: i, j, k, l, a11, a12, a21, a22, stallo
2816 double precision,pointer :: rot_mat(:,:,:)
2817 real,allocatable :: tmp_arr(:,:)
2818 REAL,POINTER :: voldatiu(:,:), voldativ(:,:)
2819 INTEGER,POINTER :: iu(:), iv(:)
2820 
2821 IF (.NOT.ASSOCIATED(this%var)) THEN
2822  CALL l4f_category_log(this%category, l4f_error, &
2823  "trying to unrotate an incomplete volgrid6d object")
2824  CALL raise_fatal_error()
2825 ! RETURN
2826 ENDIF
2827 
2828 CALL volgrid6d_var_hor_comp_index(this%var, iu, iv)
2829 IF (.NOT.ASSOCIATED(iu)) THEN
2830  CALL l4f_category_log(this%category,l4f_error, &
2831  "unrotation impossible")
2832  CALL raise_fatal_error()
2833 ! RETURN
2834 ENDIF
2835 
2836 ! Temporary workspace
2837 ALLOCATE(tmp_arr(this%griddim%dim%nx, this%griddim%dim%ny),stat=stallo)
2838 IF (stallo /= 0) THEN
2839  CALL l4f_category_log(this%category, l4f_fatal, "allocating memory")
2840  CALL raise_fatal_error()
2841 ENDIF
2842 ! allocate once for speed
2843 IF (.NOT.ASSOCIATED(this%voldati)) THEN
2844  ALLOCATE(voldatiu(this%griddim%dim%nx, this%griddim%dim%ny), &
2845  voldativ(this%griddim%dim%nx, this%griddim%dim%ny))
2846 ENDIF
2847 
2848 CALL griddim_unproj(this%griddim)
2849 CALL wind_unrot(this%griddim, rot_mat)
2850 
2851 a11=1
2852 if (rot)then
2853  a12=2
2854  a21=3
2855 else
2856  a12=3
2857  a21=2
2858 end if
2859 a22=4
2860 
2861 DO l = 1, SIZE(iu)
2862  DO k = 1, SIZE(this%timerange)
2863  DO j = 1, SIZE(this%time)
2864  DO i = 1, SIZE(this%level)
2865 ! get data
2866  CALL volgrid_get_vol_2d(this, i, j, k, iu(l), voldatiu)
2867  CALL volgrid_get_vol_2d(this, i, j, k, iv(l), voldativ)
2868 ! convert units forward
2869 ! CALL compute(conv_fwd(iu(l)), voldatiu)
2870 ! CALL compute(conv_fwd(iv(l)), voldativ)
2871 
2872 ! multiply wind components by rotation matrix
2873  WHERE(voldatiu /= rmiss .AND. voldativ /= rmiss)
2874  tmp_arr(:,:) = real(voldatiu(:,:)*rot_mat(:,:,a11) + &
2875  voldativ(:,:)*rot_mat(:,:,a12))
2876  voldativ(:,:) = real(voldatiu(:,:)*rot_mat(:,:,a21) + &
2877  voldativ(:,:)*rot_mat(:,:,a22))
2878  voldatiu(:,:) = tmp_arr(:,:)
2879  END WHERE
2880 ! convert units backward
2881 ! CALL uncompute(conv_fwd(iu(l)), voldatiu)
2882 ! CALL uncompute(conv_fwd(iv(l)), voldativ)
2883 ! put data
2884  CALL volgrid_set_vol_2d(this, i, j, k, iu(l), voldatiu)
2885  CALL volgrid_set_vol_2d(this, i, j, k, iv(l), voldativ)
2886  ENDDO
2887  ENDDO
2888  ENDDO
2889 ENDDO
2890 
2891 IF (.NOT.ASSOCIATED(this%voldati)) THEN
2892  DEALLOCATE(voldatiu, voldativ)
2893 ENDIF
2894 DEALLOCATE(rot_mat, tmp_arr, iu, iv)
2895 
2896 END SUBROUTINE vg6d_wind__un_rot
2897 
2898 
2899 !!$ try to understand the problem:
2900 !!$
2901 !!$ case:
2902 !!$
2903 !!$ 1) we have only one volume: we have to provide the direction of shift
2904 !!$ compute H and traslate on it
2905 !!$ 2) we have two volumes:
2906 !!$ 1) volume U and volume V: compute H and traslate on it
2907 !!$ 2) volume U/V and volume H : translate U/V on H
2908 !!$ 3) we have tree volumes: translate U and V on H
2909 !!$
2910 !!$ strange cases:
2911 !!$ 1) do not have U in volume U
2912 !!$ 2) do not have V in volume V
2913 !!$ 3) we have others variables more than U and V in volumes U e V
2914 !!$
2915 !!$
2916 !!$ so the steps are:
2917 !!$ 1) find the volumes
2918 !!$ 2) define or compute H grid
2919 !!$ 3) trasform the volumes in H
2920 
2921 !!$ N.B.
2922 !!$ case 1) for only one vg6d (U or V) is not managed, but
2923 !!$ the not pubblic subroutines will work but you have to know what you want to do
2924 
2925 
2942 subroutine vg6d_c2a (this)
2943 
2944 TYPE(volgrid6d),INTENT(inout) :: this(:)
2945 
2946 integer :: ngrid,igrid,jgrid,ugrid,vgrid,tgrid
2947 doubleprecision :: xmin, xmax, ymin, ymax
2948 doubleprecision :: xmin_t, xmax_t, ymin_t, ymax_t
2949 doubleprecision :: step_lon_t,step_lat_t
2950 character(len=80) :: type_t,type
2951 TYPE(griddim_def):: griddim_t
2952 
2953 ngrid=size(this)
2954 
2955 do igrid=1,ngrid
2956 
2957  call init(griddim_t)
2958 
2959  call get_val(this(igrid)%griddim,xmin=xmin_t, xmax=xmax_t, ymin=ymin_t, ymax=ymax_t,proj_type=type_t)
2960  step_lon_t=(xmax_t-xmin_t)/dble(this(igrid)%griddim%dim%nx-1)
2961  step_lat_t=(ymax_t-ymin_t)/dble(this(igrid)%griddim%dim%ny-1)
2962 
2963  do jgrid=1,ngrid
2964 
2965  ugrid=imiss
2966  vgrid=imiss
2967  tgrid=imiss
2968 
2969 #ifdef DEBUG
2970  call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: search U/V/T points:"//to_char(igrid)//to_char(jgrid))
2971  call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: test delta: "//to_char(step_lon_t)//to_char(step_lat_t))
2972 #endif
2973 
2974  if (this(igrid)%griddim == this(jgrid)%griddim ) cycle
2975 
2976  if (this(igrid)%griddim%dim%nx == this(jgrid)%griddim%dim%nx .and. &
2977  this(igrid)%griddim%dim%ny == this(jgrid)%griddim%dim%ny ) then
2978 
2979  call get_val(this(jgrid)%griddim,xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,proj_type=type)
2980 
2981  if (type_t /= type )cycle
2982 
2983 #ifdef DEBUG
2984  call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: test U "//&
2985  to_char(xmin)//to_char(xmax)//to_char(ymin)//to_char(ymax))
2986 
2987  call l4f_category_log(this(igrid)%category,l4f_debug,"diff coordinate lon"//&
2988  to_char(abs(xmin - (xmin_t+step_lon_t/2.d0)))//&
2989  to_char(abs(xmax - (xmax_t+step_lon_t/2.d0))))
2990  call l4f_category_log(this(igrid)%category,l4f_debug,"diff coordinate lat"//&
2991  to_char(abs(ymin - (ymin_t+step_lat_t/2.d0)))//&
2992  to_char(abs(ymax - (ymax_t+step_lat_t/2.d0))))
2993 #endif
2994 
2995  if ( abs(xmin - (xmin_t+step_lon_t/2.d0)) < 1.d-3 .and. abs(xmax - (xmax_t+step_lon_t/2.d0)) < 1.d-3 ) then
2996  if ( abs(ymin - ymin_t) < 1.d-3 .and. abs(ymax - ymax_t) < 1.d-3 ) then
2997 
2998 #ifdef DEBUG
2999  call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: found U")
3000 #endif
3001  ugrid=jgrid
3002  tgrid=igrid
3003 
3004  end if
3005  end if
3006 
3007 #ifdef DEBUG
3008  call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: test V "//&
3009  to_char(xmin)//to_char(xmax)//to_char(ymin)//to_char(ymax))
3010 #endif
3011 
3012  if ( abs(ymin - (ymin_t+step_lat_t/2.d0)) < 1.d-3 .and. abs(ymax - (ymax_t+step_lat_t/2.d0)) < 1.d-3 ) then
3013  if ( abs(xmin - xmin_t) < 1.d-3 .and. abs(xmax - xmax_t) < 1.d-3 ) then
3014 
3015 #ifdef DEBUG
3016  call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: found V")
3017 #endif
3018  vgrid=jgrid
3019  tgrid=igrid
3020 
3021  end if
3022  end if
3023 
3024 
3025  ! test if we have U and V
3026 
3027 #ifdef DEBUG
3028  call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: test U and V"//&
3029  to_char(xmin_t)//to_char(xmax_t)//to_char(ymin_t)//to_char(ymax_t)//&
3030  to_char(xmin)//to_char(xmax)//to_char(ymin)//to_char(ymax))
3031 
3032  call l4f_category_log(this(igrid)%category,l4f_debug,"UV diff coordinate lon"//&
3033  to_char(abs(xmin_t - xmin)-step_lon_t/2.d0)//&
3034  to_char(abs(xmax_t - xmax)-step_lon_t/2.d0))
3035  call l4f_category_log(this(igrid)%category,l4f_debug,"UV diff coordinate lat"//&
3036  to_char(abs(ymin_t - ymin) -step_lat_t/2.d0)//&
3037  to_char(abs(ymax_t - ymax)-step_lat_t/2.d0))
3038 #endif
3039 
3040  if ( abs(ymin - (ymin_t+step_lat_t/2.d0)) < 2.d-3 .and. abs(ymax - (ymax_t+step_lat_t/2.d0)) < 2.d-3 ) then
3041  if ( abs(xmin_t - (xmin+step_lon_t/2.d0)) < 2.d-3 .and. abs(xmax_t - (xmax+step_lon_t/2.d0)) < 2.d-3 ) then
3042 
3043 #ifdef DEBUG
3044  call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: found U and V case up and right")
3045 #endif
3046 
3047  vgrid=jgrid
3048  ugrid=igrid
3049 
3050  call init(griddim_t,xmin=xmin, xmax=xmax, ymin=ymin_t, ymax=ymax_t)
3051 
3052  end if
3053  end if
3054  end if
3055 
3056  ! abbiamo almeno due volumi: riportiamo U e/o V su T
3057  if (c_e(ugrid)) then
3058 #ifdef DEBUG
3059  call l4f_category_log(this(igrid)%category,l4f_debug,"C grid U points: interpolate U on T "//to_char(tgrid)//to_char(ugrid))
3060 #endif
3061  if (c_e(tgrid))then
3062  call vg6d_c2a_grid(this(ugrid),this(tgrid)%griddim,cgrid=1)
3063  else
3064  call vg6d_c2a_grid(this(ugrid),griddim_t,cgrid=1)
3065  end if
3066  call vg6d_c2a_mat(this(ugrid),cgrid=1)
3067  end if
3068 
3069  if (c_e(vgrid)) then
3070 #ifdef DEBUG
3071  call l4f_category_log(this(igrid)%category,l4f_debug,"C grid V points: interpolate V on T "//to_char(tgrid)//to_char(vgrid))
3072 #endif
3073  if (c_e(tgrid))then
3074  call vg6d_c2a_grid(this(vgrid),this(tgrid)%griddim,cgrid=2)
3075  else
3076  call vg6d_c2a_grid(this(vgrid),griddim_t,cgrid=2)
3077  end if
3078  call vg6d_c2a_mat(this(vgrid),cgrid=2)
3079  end if
3080 
3081  end do
3082 
3083  call delete(griddim_t)
3084 
3085 end do
3086 
3087 
3088 end subroutine vg6d_c2a
3089 
3090 
3091 ! Convert C grid to A grid
3092 subroutine vg6d_c2a_grid(this,griddim_t,cgrid)
3093 
3094 type(volgrid6d),intent(inout) :: this ! object containing fields to be translated (U or V points)
3095 type(griddim_def),intent(in),optional :: griddim_t ! object containing grid of T points
3096 integer,intent(in) :: cgrid ! in C grid (Arakawa) we have 0=T,1=U,2=V points
3097 
3098 doubleprecision :: xmin, xmax, ymin, ymax
3099 doubleprecision :: step_lon,step_lat
3100 
3101 
3102 if (present(griddim_t)) then
3103 
3104  call get_val(griddim_t,xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax)
3105  call set_val(this%griddim,xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax)
3106 ! improve grid definition precision
3107  CALL griddim_setsteps(this%griddim)
3108 
3109 else
3110 
3111  select case (cgrid)
3112 
3113  case(0)
3114 
3115 #ifdef DEBUG
3116  call l4f_category_log(this%category,l4f_debug,"C grid: T points, nothing to do")
3117 #endif
3118  return
3119 
3120  case (1)
3121 
3122 #ifdef DEBUG
3123  call l4f_category_log(this%category,l4f_debug,"C grid: U points, we need interpolation")
3124 #endif
3125 
3126  call get_val(this%griddim, xmin=xmin, xmax=xmax)
3127  step_lon=(xmax-xmin)/dble(this%griddim%dim%nx-1)
3128  xmin=xmin-step_lon/2.d0
3129  xmax=xmax-step_lon/2.d0
3130  call set_val(this%griddim, xmin=xmin, xmax=xmax)
3131 ! improve grid definition precision
3132  CALL griddim_setsteps(this%griddim)
3133 
3134  case (2)
3135 
3136 #ifdef DEBUG
3137  call l4f_category_log(this%category,l4f_debug,"C grid: V points, we need interpolation")
3138 #endif
3139 
3140  call get_val(this%griddim, ymin=ymin, ymax=ymax)
3141  step_lat=(ymax-ymin)/dble(this%griddim%dim%ny-1)
3142  ymin=ymin-step_lat/2.d0
3143  ymax=ymax-step_lat/2.d0
3144  call set_val(this%griddim, ymin=ymin, ymax=ymax)
3145 ! improve grid definition precision
3146  CALL griddim_setsteps(this%griddim)
3147 
3148  case default
3149 
3150  call l4f_category_log(this%category,l4f_fatal,"C grid type not known")
3151  call raise_fatal_error()
3152 
3153  end select
3154 
3155 end if
3156 
3157 
3158 call griddim_unproj(this%griddim)
3159 
3160 
3161 end subroutine vg6d_c2a_grid
3162 
3163 ! Convert C grid to A grid
3164 subroutine vg6d_c2a_mat(this,cgrid)
3165 
3166 type(volgrid6d),intent(inout) :: this ! object containing fields to be translated
3167 integer,intent(in) :: cgrid ! in C grid (Arakawa) we have 0=T,1=U,2=V points
3168 
3169 INTEGER :: i, j, k, iv, stallo
3170 REAL,ALLOCATABLE :: tmp_arr(:,:)
3171 REAL,POINTER :: voldatiuv(:,:)
3172 
3173 
3174 IF (cgrid == 0) RETURN ! nothing to do
3175 IF (cgrid /= 1 .AND. cgrid /= 2) THEN ! wrong cgrid
3176  CALL l4f_category_log(this%category,l4f_fatal,"C grid type "// &
3177  trim(to_char(cgrid))//" not known")
3178  CALL raise_error()
3179  RETURN
3180 ENDIF
3181 
3182 ! Temporary workspace
3183 ALLOCATE(tmp_arr(this%griddim%dim%nx, this%griddim%dim%ny),stat=stallo)
3184 if (stallo /=0)then
3185  call l4f_log(l4f_fatal,"allocating memory")
3186  call raise_fatal_error()
3187 end if
3188 
3189 ! allocate once for speed
3190 IF (.NOT.ASSOCIATED(this%voldati)) THEN
3191  ALLOCATE(voldatiuv(this%griddim%dim%nx, this%griddim%dim%ny), stat=stallo)
3192  IF (stallo /= 0) THEN
3193  CALL l4f_log(l4f_fatal,"allocating memory")
3194  CALL raise_fatal_error()
3195  ENDIF
3196 ENDIF
3197 
3198 IF (cgrid == 1) THEN ! U points to H points
3199  DO iv = 1, SIZE(this%var)
3200  DO k = 1, SIZE(this%timerange)
3201  DO j = 1, SIZE(this%time)
3202  DO i = 1, SIZE(this%level)
3203  tmp_arr(:,:) = rmiss
3204  CALL volgrid_get_vol_2d(this, i, j, k, iv, voldatiuv)
3205 
3206 ! West boundary extrapolation
3207  WHERE(voldatiuv(1,:) /= rmiss .AND. voldatiuv(2,:) /= rmiss)
3208  tmp_arr(1,:) = voldatiuv(1,:) - (voldatiuv(2,:) - voldatiuv(1,:)) / 2.
3209  END WHERE
3210 
3211 ! Rest of the field
3212  WHERE(voldatiuv(1:this%griddim%dim%nx-1,:) /= rmiss .AND. &
3213  voldatiuv(2:this%griddim%dim%nx,:) /= rmiss)
3214  tmp_arr(2:this%griddim%dim%nx,:) = &
3215  (voldatiuv(1:this%griddim%dim%nx-1,:) + &
3216  voldatiuv(2:this%griddim%dim%nx,:)) / 2.
3217  END WHERE
3218 
3219  voldatiuv(:,:) = tmp_arr
3220  CALL volgrid_set_vol_2d(this, i, j, k, iv, voldatiuv)
3221  ENDDO
3222  ENDDO
3223  ENDDO
3224  ENDDO
3225 
3226 ELSE IF (cgrid == 2) THEN ! V points to H points
3227  DO iv = 1, SIZE(this%var)
3228  DO k = 1, SIZE(this%timerange)
3229  DO j = 1, SIZE(this%time)
3230  DO i = 1, SIZE(this%level)
3231  tmp_arr(:,:) = rmiss
3232  CALL volgrid_get_vol_2d(this, i, j, k, iv, voldatiuv)
3233 
3234 ! South boundary extrapolation
3235  WHERE(voldatiuv(:,1) /= rmiss .AND. voldatiuv(:,2) /= rmiss)
3236  tmp_arr(:,1) = voldatiuv(:,1) - (voldatiuv(:,2) - voldatiuv(:,1)) / 2.
3237  END WHERE
3238 
3239 ! Rest of the field
3240  WHERE(voldatiuv(:,1:this%griddim%dim%ny-1) /= rmiss .AND. &
3241  voldatiuv(:,2:this%griddim%dim%ny) /= rmiss)
3242  tmp_arr(:,2:this%griddim%dim%ny) = &
3243  (voldatiuv(:,1:this%griddim%dim%ny-1) + &
3244  voldatiuv(:,2:this%griddim%dim%ny)) / 2.
3245  END WHERE
3246 
3247  voldatiuv(:,:) = tmp_arr
3248  CALL volgrid_set_vol_2d(this, i, j, k, iv, voldatiuv)
3249  ENDDO
3250  ENDDO
3251  ENDDO
3252  ENDDO
3253 ENDIF ! tertium non datur
3254 
3255 IF (.NOT.ASSOCIATED(this%voldati)) THEN
3256  DEALLOCATE(voldatiuv)
3257 ENDIF
3258 DEALLOCATE (tmp_arr)
3259 
3260 end subroutine vg6d_c2a_mat
3261 
3262 
3266 subroutine display_volgrid6d (this)
3267 
3268 TYPE(volgrid6d),intent(in) :: this
3269 integer :: i
3270 
3271 #ifdef DEBUG
3272 call l4f_category_log(this%category,l4f_debug,"ora mostro gridinfo " )
3273 #endif
3274 
3275 print*,"----------------------- volgrid6d display ---------------------"
3276 call display(this%griddim)
3277 
3278 IF (ASSOCIATED(this%time))then
3279  print*,"---- time vector ----"
3280  print*,"elements=",size(this%time)
3281  do i=1, size(this%time)
3282  call display(this%time(i))
3283  end do
3284 end IF
3285 
3286 IF (ASSOCIATED(this%timerange))then
3287  print*,"---- timerange vector ----"
3288  print*,"elements=",size(this%timerange)
3289  do i=1, size(this%timerange)
3290  call display(this%timerange(i))
3291  end do
3292 end IF
3293 
3294 IF (ASSOCIATED(this%level))then
3295  print*,"---- level vector ----"
3296  print*,"elements=",size(this%level)
3297  do i=1, size(this%level)
3298  call display(this%level(i))
3299  end do
3300 end IF
3301 
3302 IF (ASSOCIATED(this%var))then
3303  print*,"---- var vector ----"
3304  print*,"elements=",size(this%var)
3305  do i=1, size(this%var)
3306  call display(this%var(i))
3307  end do
3308 end IF
3309 
3310 IF (ASSOCIATED(this%gaid))then
3311  print*,"---- gaid vector (present mask only) ----"
3312  print*,"elements=",shape(this%gaid)
3313  print*,c_e(reshape(this%gaid,(/SIZE(this%gaid)/)))
3314 end IF
3315 
3316 print*,"--------------------------------------------------------------"
3317 
3318 
3319 end subroutine display_volgrid6d
3320 
3321 
3325 subroutine display_volgrid6dv (this)
3326 
3327 TYPE(volgrid6d),intent(in) :: this(:)
3328 integer :: i
3329 
3330 print*,"----------------------- volgrid6d vector ---------------------"
3331 
3332 print*,"elements=",size(this)
3333 
3334 do i=1, size(this)
3335 
3336 #ifdef DEBUG
3337  call l4f_category_log(this(i)%category,l4f_debug,"ora mostro il vettore volgrid6d" )
3338 #endif
3339 
3340  call display(this(i))
3341 
3342 end do
3343 print*,"--------------------------------------------------------------"
3344 
3345 end subroutine display_volgrid6dv
3346 
3347 
3350 subroutine vg6dv_rounding(vg6din,vg6dout,level,timerange,nostatproc,merge)
3351 type(volgrid6d),intent(in) :: vg6din(:)
3352 type(volgrid6d),intent(out),pointer :: vg6dout(:) !> output volume
3353 type(vol7d_level),intent(in),optional :: level(:)
3354 type(vol7d_timerange),intent(in),optional :: timerange(:)
3355 logical,intent(in),optional :: merge
3356 logical,intent(in),optional :: nostatproc
3357 
3358 integer :: i
3359 
3360 allocate(vg6dout(size(vg6din)))
3361 
3362 do i = 1, size(vg6din)
3363  call vg6d_rounding(vg6din(i),vg6dout(i),level,timerange,nostatproc,merge)
3364 end do
3365 
3366 end subroutine vg6dv_rounding
3367 
3379 subroutine vg6d_rounding(vg6din,vg6dout,level,timerange,nostatproc,merge)
3380 type(volgrid6d),intent(in) :: vg6din
3381 type(volgrid6d),intent(out) :: vg6dout !> output volume
3382 type(vol7d_level),intent(in),optional :: level(:)
3383 type(vol7d_timerange),intent(in),optional :: timerange(:)
3384 logical,intent(in),optional :: merge
3386 logical,intent(in),optional :: nostatproc
3387 
3388 integer :: ilevel,itimerange
3389 type(vol7d_level) :: roundlevel(size(vg6din%level))
3390 type(vol7d_timerange) :: roundtimerange(size(vg6din%timerange))
3391 
3392 roundlevel=vg6din%level
3393 
3394 if (present(level))then
3395  do ilevel = 1, size(vg6din%level)
3396  if ((any(vg6din%level(ilevel) .almosteq. level))) then
3397  roundlevel(ilevel)=level(1)
3398  end if
3399  end do
3400 end if
3401 
3402 roundtimerange=vg6din%timerange
3403 
3404 if (present(timerange))then
3405  do itimerange = 1, size(vg6din%timerange)
3406  if ((any(vg6din%timerange(itimerange) .almosteq. timerange))) then
3407  roundtimerange(itimerange)=timerange(1)
3408  end if
3409  end do
3410 end if
3411 
3412 !set istantaneous values everywere
3413 !preserve p1 for forecast time
3414 if (optio_log(nostatproc)) then
3415  roundtimerange(:)%timerange=254
3416  roundtimerange(:)%p2=imiss
3417 end if
3418 
3419 
3420 call vg6d_reduce(vg6din,vg6dout,roundlevel,roundtimerange,merge)
3421 
3422 end subroutine vg6d_rounding
3423 
3432 subroutine vg6d_reduce(vg6din,vg6dout,roundlevel,roundtimerange,merge)
3433 type(volgrid6d),intent(in) :: vg6din
3434 type(volgrid6d),intent(out) :: vg6dout
3435 type(vol7d_level),intent(in) :: roundlevel(:)
3436 type(vol7d_timerange),intent(in) :: roundtimerange(:)
3437 logical,intent(in),optional :: merge
3438 
3439 integer :: nlevel,ntime,ntimerange,nvar,ilevel,itimerange,ivar,indl,indt,itime,nx,ny
3440 real,allocatable :: vol2d(:,:)
3441 
3442 nx=vg6din%griddim%dim%nx
3443 ny=vg6din%griddim%dim%ny
3444 nlevel=count_distinct(roundlevel,back=.true.)
3445 ntime=size(vg6din%time)
3446 ntimerange=count_distinct(roundtimerange,back=.true.)
3447 nvar=size(vg6din%var)
3448 
3449 call init(vg6dout, vg6din%griddim, vg6din%time_definition, categoryappend="generated by vg6d_reduce")
3450 call volgrid6d_alloc(vg6dout, vg6din%griddim%dim, ntime, nlevel, ntimerange, nvar)
3451 
3452 if ( ASSOCIATED(vg6din%voldati) .or. optio_log(merge)) then
3453  call volgrid6d_alloc_vol(vg6dout,inivol=.true.,decode=.true.)
3454  allocate(vol2d(nx,ny))
3455 else
3456  call volgrid6d_alloc_vol(vg6dout,inivol=.true.,decode=.false.)
3457 end if
3458 
3459 vg6dout%time=vg6din%time
3460 vg6dout%var=vg6din%var
3461 vg6dout%timerange=pack_distinct(roundtimerange,ntimerange,back=.true.)
3462 vg6dout%level=pack_distinct(roundlevel,nlevel,back=.true.)
3463 ! sort modified dimensions
3464 CALL sort(vg6dout%timerange)
3465 CALL sort(vg6dout%level)
3466 
3467 do ilevel=1,size(vg6din%level)
3468  indl=index(vg6dout%level,roundlevel(ilevel))
3469  do itimerange=1,size(vg6din%timerange)
3470  indt=index(vg6dout%timerange,roundtimerange(itimerange))
3471  do ivar=1, nvar
3472  do itime=1,ntime
3473 
3474  if ( ASSOCIATED(vg6din%voldati)) then
3475  vol2d=vg6din%voldati(:,:,ilevel,itime,itimerange,ivar)
3476  end if
3477 
3478  if (optio_log(merge)) then
3479 
3480  if ( .not. ASSOCIATED(vg6din%voldati)) then
3481  CALL grid_id_decode_data(vg6din%gaid(ilevel,itime,itimerange,ivar), vol2d)
3482  end if
3483 
3484  !! merge present data point by point
3485  where (.not. c_e(vg6dout%voldati(:,:,indl,itime,indt,ivar)))
3486 
3487  vg6dout%voldati(:,:,indl,itime,indt,ivar)=vol2d
3488 
3489  end where
3490  else if ( ASSOCIATED(vg6din%voldati)) then
3491  if (.not. any(c_e(vg6dout%voldati(:,:,indl,itime,indt,ivar))))then
3492  vg6dout%voldati(:,:,indl,itime,indt,ivar)=vol2d
3493  end if
3494  end if
3495 
3496  if (c_e(vg6din%gaid(ilevel,itime,itimerange,ivar)).and. .not. c_e(vg6dout%gaid(indl,itime,indt,ivar)))then
3497  call copy(vg6din%gaid(ilevel,itime,itimerange,ivar), vg6dout%gaid(indl,itime,indt,ivar))
3498  end if
3499  end do
3500  end do
3501  end do
3502 end do
3503 
3504 if ( ASSOCIATED(vg6din%voldati) .or. optio_log(merge)) then
3505  deallocate(vol2d)
3506 end if
3507 
3508 end subroutine vg6d_reduce
3509 
3510 
3511 end module volgrid6d_class
3512 
3513 
3514 
3519 
3526 
3529 
3532 
3535 
Classi per la gestione delle coordinate temporali.
Derived type defining a dynamically extensible array of TYPE(gridinfo_def) elements.
Functions that return a trimmed CHARACTER representation of the input variable.
Compute forward coordinate transformation from geographical system to projected system.
Legge un oggetto datetime/timedelta o un vettore di oggetti datetime/timedelta da un file FORMATTED o...
Constructor, it creates a new instance of the object.
Derived type associated to a block/message/record/band of gridded data coming from a file-like object...
Definition of a physical variable in grib coding style.
Operatore di valore assoluto di un intervallo.
This module defines an abstract interface to different drivers for access to files containing gridded...
Destructor, it releases every information and memory buffer associated with the object.
Object describing a rectangular, homogeneous gridded dataset.
Module for describing geographically referenced regular grids.
Definition: grid_class.F90:243
Transform between any combination of volgrid6d and vol7d objects by means of a transform_def object d...
Apply the conversion function this to values.
Decode and return the data array from a grid_id object associated to a gridinfo object.
Module for defining the extension and coordinates of a rectangular georeferenced grid.
Encode a data array into a grid_id object associated to a gridinfo object.
Module for defining transformations between rectangular georeferenced grids and between grids and spa...
This object completely describes a grid on a geographic projection.
Definition: grid_class.F90:276
Classe per la gestione dei livelli verticali in osservazioni meteo e affini.
Method for inserting elements of the array at a desired position.
Scrive un oggetto datetime/timedelta o un vettore di oggetti datetime/timedelta su un file FORMATTED ...
Import an object dirctly from a native file, from a gridinfo object or from a supported file format t...
Classe per la gestione di un volume completo di dati osservati.
This module defines objects and methods for managing data volumes on rectangular georeferenced grids...
This object fully defines a transformation between a couple of particular griddim_def or vol7d object...
Export an object dirctly to a native file, to a gridinfo object or to a supported file format through...
Index method.
Restituiscono il valore dell&#39;oggetto in forma di stringa stampabile.
Classe per la gestione degli intervalli temporali di osservazioni meteo e affini. ...
Definisce un oggetto contenente i volumi anagrafica e dati e tutti i descrittori delle loro dimension...
Copy an object, creating a fully new instance.
classe per import ed export di volumi da e in DB-All.e
Display on standard output a description of the volgrid6d object provided.
Method for returning the contents of the object.
Object describing a single gridded message/band.
Class for managing physical variables in a grib 1/2 fashion.
Clone the object, creating a new independent instance of the object exactly equal to the starting one...
This object defines in an abstract way the type of transformation to be applied.
classe per la gestione del logging
Definisce una variabile meteorologica osservata o un suo attributo.
Module for quickly interpreting the OPTIONAL parameters passed to a subprogram.
Method for setting the contents of the object.
Class defining a real conversion function between units.
Class for managing information about a single gridded georeferenced field, typically imported from an...
Emit log message for a category with specific priority.
Derived type describing the extension of a grid and the geographical coordinates of each point...
This module defines usefull general purpose function and subroutine.
Reduce some dimensions (level and timerage) for semplification (rounding).

Generated with Doxygen.