libsim Versione 7.2.0
vol7d_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
85MODULE vol7d_class
86USE kinds
92USE io_units
99IMPLICIT NONE
100
101
102INTEGER, PARAMETER :: vol7d_maxdim_a = 3, vol7d_maxdim_aa = 4, &
103 vol7d_maxdim_d = 6, vol7d_maxdim_ad = 7
104
105INTEGER, PARAMETER :: vol7d_ana_a=1
106INTEGER, PARAMETER :: vol7d_var_a=2
107INTEGER, PARAMETER :: vol7d_network_a=3
108INTEGER, PARAMETER :: vol7d_attr_a=4
109INTEGER, PARAMETER :: vol7d_ana_d=1
110INTEGER, PARAMETER :: vol7d_time_d=2
111INTEGER, PARAMETER :: vol7d_level_d=3
112INTEGER, PARAMETER :: vol7d_timerange_d=4
113INTEGER, PARAMETER :: vol7d_var_d=5
114INTEGER, PARAMETER :: vol7d_network_d=6
115INTEGER, PARAMETER :: vol7d_attr_d=7
116INTEGER, PARAMETER :: vol7d_cdatalen=32
117
118TYPE vol7d_varmap
119 INTEGER :: r, d, i, b, c
120END TYPE vol7d_varmap
121
124TYPE vol7d
126 TYPE(vol7d_ana),POINTER :: ana(:)
128 TYPE(datetime),POINTER :: time(:)
130 TYPE(vol7d_level),POINTER :: level(:)
132 TYPE(vol7d_timerange),POINTER :: timerange(:)
134 TYPE(vol7d_network),POINTER :: network(:)
136 TYPE(vol7d_varvect) :: anavar
138 TYPE(vol7d_varvect) :: anaattr
140 TYPE(vol7d_varvect) :: anavarattr
142 TYPE(vol7d_varvect) :: dativar
144 TYPE(vol7d_varvect) :: datiattr
146 TYPE(vol7d_varvect) :: dativarattr
147
149 REAL,POINTER :: volanar(:,:,:)
151 DOUBLE PRECISION,POINTER :: volanad(:,:,:)
153 INTEGER,POINTER :: volanai(:,:,:)
155 INTEGER(kind=int_b),POINTER :: volanab(:,:,:)
157 CHARACTER(len=vol7d_cdatalen),POINTER :: volanac(:,:,:)
158
160 REAL,POINTER :: volanaattrr(:,:,:,:)
162 DOUBLE PRECISION,POINTER :: volanaattrd(:,:,:,:)
164 INTEGER,POINTER :: volanaattri(:,:,:,:)
166 INTEGER(kind=int_b),POINTER :: volanaattrb(:,:,:,:)
168 CHARACTER(len=vol7d_cdatalen),POINTER :: volanaattrc(:,:,:,:)
169
171 REAL,POINTER :: voldatir(:,:,:,:,:,:) ! sono i dati
173 DOUBLE PRECISION,POINTER :: voldatid(:,:,:,:,:,:)
175 INTEGER,POINTER :: voldatii(:,:,:,:,:,:)
177 INTEGER(kind=int_b),POINTER :: voldatib(:,:,:,:,:,:)
179 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatic(:,:,:,:,:,:)
180
182 REAL,POINTER :: voldatiattrr(:,:,:,:,:,:,:)
184 DOUBLE PRECISION,POINTER :: voldatiattrd(:,:,:,:,:,:,:)
186 INTEGER,POINTER :: voldatiattri(:,:,:,:,:,:,:)
188 INTEGER(kind=int_b),POINTER :: voldatiattrb(:,:,:,:,:,:,:)
190 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatiattrc(:,:,:,:,:,:,:)
191
193 integer :: time_definition
194
195END TYPE vol7d
196
200INTERFACE init
201 MODULE PROCEDURE vol7d_init
202END INTERFACE
203
205INTERFACE delete
206 MODULE PROCEDURE vol7d_delete
207END INTERFACE
208
210INTERFACE export
211 MODULE PROCEDURE vol7d_write_on_file
212END INTERFACE
213
215INTERFACE import
216 MODULE PROCEDURE vol7d_read_from_file
217END INTERFACE
218
220INTERFACE display
221 MODULE PROCEDURE vol7d_display, dat_display, dat_vect_display
222END INTERFACE
223
225INTERFACE to_char
226 MODULE PROCEDURE to_char_dat
227END INTERFACE
228
230INTERFACE doubledat
231 MODULE PROCEDURE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
232END INTERFACE
233
235INTERFACE realdat
236 MODULE PROCEDURE realdatd,realdatr,realdati,realdatb,realdatc
237END INTERFACE
238
240INTERFACE integerdat
241 MODULE PROCEDURE integerdatd,integerdatr,integerdati,integerdatb,integerdatc
242END INTERFACE
243
245INTERFACE copy
246 MODULE PROCEDURE vol7d_copy
247END INTERFACE
248
250INTERFACE c_e
251 MODULE PROCEDURE vol7d_c_e
252END INTERFACE
253
257INTERFACE check
258 MODULE PROCEDURE vol7d_check
259END INTERFACE
260
274INTERFACE rounding
275 MODULE PROCEDURE v7d_rounding
276END INTERFACE
277
278!!$INTERFACE get_volana
279!!$ MODULE PROCEDURE vol7d_get_volanar, vol7d_get_volanad, vol7d_get_volanai, &
280!!$ vol7d_get_volanab, vol7d_get_volanac
281!!$END INTERFACE
282!!$
283!!$INTERFACE get_voldati
284!!$ MODULE PROCEDURE vol7d_get_voldatir, vol7d_get_voldatid, vol7d_get_voldatii, &
285!!$ vol7d_get_voldatib, vol7d_get_voldatic
286!!$END INTERFACE
287!!$
288!!$INTERFACE get_volanaattr
289!!$ MODULE PROCEDURE vol7d_get_volanaattrr, vol7d_get_volanaattrd, &
290!!$ vol7d_get_volanaattri, vol7d_get_volanaattrb, vol7d_get_volanaattrc
291!!$END INTERFACE
292!!$
293!!$INTERFACE get_voldatiattr
294!!$ MODULE PROCEDURE vol7d_get_voldatiattrr, vol7d_get_voldatiattrd, &
295!!$ vol7d_get_voldatiattri, vol7d_get_voldatiattrb, vol7d_get_voldatiattrc
296!!$END INTERFACE
298PRIVATE vol7d_get_volr, vol7d_get_vold, vol7d_get_voli, vol7d_get_volb, &
299 vol7d_get_volc, &
300 volptr1dr, volptr2dr, volptr3dr, volptr4dr, volptr5dr, volptr6dr, volptr7dr, &
301 volptr1dd, volptr2dd, volptr3dd, volptr4dd, volptr5dd, volptr6dd, volptr7dd, &
302 volptr1di, volptr2di, volptr3di, volptr4di, volptr5di, volptr6di, volptr7di, &
303 volptr1db, volptr2db, volptr3db, volptr4db, volptr5db, volptr6db, volptr7db, &
304 volptr1dc, volptr2dc, volptr3dc, volptr4dc, volptr5dc, volptr6dc, volptr7dc, &
305 vol7d_nullifyr, vol7d_nullifyd, vol7d_nullifyi, vol7d_nullifyb, vol7d_nullifyc, &
306 vol7d_init, vol7d_delete, vol7d_write_on_file, vol7d_read_from_file, &
307 vol7d_check_alloc_ana, vol7d_force_alloc_ana, &
308 vol7d_check_alloc_dati, vol7d_force_alloc_dati, vol7d_force_alloc, &
309 vol7d_display, dat_display, dat_vect_display, &
310 to_char_dat, vol7d_check
311
312PRIVATE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
313
314PRIVATE vol7d_c_e
315
316CONTAINS
317
323SUBROUTINE vol7d_init(this,time_definition)
324TYPE(vol7d),intent(out) :: this
325integer,INTENT(IN),OPTIONAL :: time_definition
327CALL init(this%anavar)
328CALL init(this%anaattr)
329CALL init(this%anavarattr)
330CALL init(this%dativar)
331CALL init(this%datiattr)
332CALL init(this%dativarattr)
333CALL vol7d_var_features_init() ! initialise var features table once
335NULLIFY(this%ana, this%time, this%level, this%timerange, this%network)
336
337NULLIFY(this%volanar, this%volanaattrr, this%voldatir, this%voldatiattrr)
338NULLIFY(this%volanad, this%volanaattrd, this%voldatid, this%voldatiattrd)
339NULLIFY(this%volanai, this%volanaattri, this%voldatii, this%voldatiattri)
340NULLIFY(this%volanab, this%volanaattrb, this%voldatib, this%voldatiattrb)
341NULLIFY(this%volanac, this%volanaattrc, this%voldatic, this%voldatiattrc)
342
343if(present(time_definition)) then
344 this%time_definition=time_definition
345else
346 this%time_definition=1 !default to validity time
347end if
349END SUBROUTINE vol7d_init
351
355ELEMENTAL SUBROUTINE vol7d_delete(this, dataonly)
356TYPE(vol7d),intent(inout) :: this
357LOGICAL, INTENT(in), OPTIONAL :: dataonly
358
360IF (.NOT. optio_log(dataonly)) THEN
361 IF (ASSOCIATED(this%volanar)) DEALLOCATE(this%volanar)
362 IF (ASSOCIATED(this%volanad)) DEALLOCATE(this%volanad)
363 IF (ASSOCIATED(this%volanai)) DEALLOCATE(this%volanai)
364 IF (ASSOCIATED(this%volanab)) DEALLOCATE(this%volanab)
365 IF (ASSOCIATED(this%volanac)) DEALLOCATE(this%volanac)
366 IF (ASSOCIATED(this%volanaattrr)) DEALLOCATE(this%volanaattrr)
367 IF (ASSOCIATED(this%volanaattrd)) DEALLOCATE(this%volanaattrd)
368 IF (ASSOCIATED(this%volanaattri)) DEALLOCATE(this%volanaattri)
369 IF (ASSOCIATED(this%volanaattrb)) DEALLOCATE(this%volanaattrb)
370 IF (ASSOCIATED(this%volanaattrc)) DEALLOCATE(this%volanaattrc)
371ENDIF
372IF (ASSOCIATED(this%voldatir)) DEALLOCATE(this%voldatir)
373IF (ASSOCIATED(this%voldatid)) DEALLOCATE(this%voldatid)
374IF (ASSOCIATED(this%voldatii)) DEALLOCATE(this%voldatii)
375IF (ASSOCIATED(this%voldatib)) DEALLOCATE(this%voldatib)
376IF (ASSOCIATED(this%voldatic)) DEALLOCATE(this%voldatic)
377IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
378IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
379IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
380IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
381IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
382
383IF (.NOT. optio_log(dataonly)) THEN
384 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
385 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
386ENDIF
387IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
388IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
389IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
390
391IF (.NOT. optio_log(dataonly)) THEN
392 CALL delete(this%anavar)
393 CALL delete(this%anaattr)
394 CALL delete(this%anavarattr)
395ENDIF
396CALL delete(this%dativar)
397CALL delete(this%datiattr)
398CALL delete(this%dativarattr)
399
400END SUBROUTINE vol7d_delete
401
402
404integer function vol7d_check(this)
405TYPE(vol7d),intent(in) :: this
406integer :: i,j,k,l,m,n
407
408vol7d_check=0
409
410if (associated(this%voldatii)) then
411do i = 1,size(this%voldatii,1)
412 do j = 1,size(this%voldatii,2)
413 do k = 1,size(this%voldatii,3)
414 do l = 1,size(this%voldatii,4)
415 do m = 1,size(this%voldatii,5)
416 do n = 1,size(this%voldatii,6)
417 if (this%voldatii(i,j,k,l,m,n) /= this%voldatii(i,j,k,l,m,n) ) then
418 CALL l4f_log(l4f_warn,"check: abnormal value at voldatii("&
419 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
420 vol7d_check=1
421 end if
422 end do
423 end do
424 end do
425 end do
426 end do
427end do
428end if
429
430
431if (associated(this%voldatir)) then
432do i = 1,size(this%voldatir,1)
433 do j = 1,size(this%voldatir,2)
434 do k = 1,size(this%voldatir,3)
435 do l = 1,size(this%voldatir,4)
436 do m = 1,size(this%voldatir,5)
437 do n = 1,size(this%voldatir,6)
438 if (this%voldatir(i,j,k,l,m,n) /= this%voldatir(i,j,k,l,m,n) ) then
439 CALL l4f_log(l4f_warn,"check: abnormal value at voldatir("&
440 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
441 vol7d_check=2
442 end if
443 end do
444 end do
445 end do
446 end do
447 end do
448end do
449end if
450
451if (associated(this%voldatid)) then
452do i = 1,size(this%voldatid,1)
453 do j = 1,size(this%voldatid,2)
454 do k = 1,size(this%voldatid,3)
455 do l = 1,size(this%voldatid,4)
456 do m = 1,size(this%voldatid,5)
457 do n = 1,size(this%voldatid,6)
458 if (this%voldatid(i,j,k,l,m,n) /= this%voldatid(i,j,k,l,m,n) ) then
459 CALL l4f_log(l4f_warn,"check: abnormal value at voldatid("&
460 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
461 vol7d_check=3
462 end if
463 end do
464 end do
465 end do
466 end do
467 end do
468end do
469end if
470
471if (associated(this%voldatib)) then
472do i = 1,size(this%voldatib,1)
473 do j = 1,size(this%voldatib,2)
474 do k = 1,size(this%voldatib,3)
475 do l = 1,size(this%voldatib,4)
476 do m = 1,size(this%voldatib,5)
477 do n = 1,size(this%voldatib,6)
478 if (this%voldatib(i,j,k,l,m,n) /= this%voldatib(i,j,k,l,m,n) ) then
479 CALL l4f_log(l4f_warn,"check: abnormal value at voldatib("&
480 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
481 vol7d_check=4
482 end if
483 end do
484 end do
485 end do
486 end do
487 end do
488end do
489end if
490
491end function vol7d_check
492
493
494
495!TODO da completare ! aborta se i volumi sono allocati a dimensione 0
497SUBROUTINE vol7d_display(this)
498TYPE(vol7d),intent(in) :: this
499integer :: i
500
501REAL :: rdat
502DOUBLE PRECISION :: ddat
503INTEGER :: idat
504INTEGER(kind=int_b) :: bdat
505CHARACTER(len=vol7d_cdatalen) :: cdat
506
507
508print*,"<<<<<<<<<<<<<<<<<<< vol7d object >>>>>>>>>>>>>>>>>>>>"
509if (this%time_definition == 0) then
510 print*,"TIME DEFINITION: time is reference time"
511else if (this%time_definition == 1) then
512 print*,"TIME DEFINITION: time is validity time"
513else
514 print*,"Time definition have a wrong walue:", this%time_definition
515end if
516
517IF (ASSOCIATED(this%network))then
518 print*,"---- network vector ----"
519 print*,"elements=",size(this%network)
520 do i=1, size(this%network)
521 call display(this%network(i))
522 end do
523end IF
524
525IF (ASSOCIATED(this%ana))then
526 print*,"---- ana vector ----"
527 print*,"elements=",size(this%ana)
528 do i=1, size(this%ana)
529 call display(this%ana(i))
530 end do
531end IF
532
533IF (ASSOCIATED(this%time))then
534 print*,"---- time vector ----"
535 print*,"elements=",size(this%time)
536 do i=1, size(this%time)
537 call display(this%time(i))
538 end do
539end if
540
541IF (ASSOCIATED(this%level)) then
542 print*,"---- level vector ----"
543 print*,"elements=",size(this%level)
544 do i =1,size(this%level)
545 call display(this%level(i))
546 end do
547end if
548
549IF (ASSOCIATED(this%timerange))then
550 print*,"---- timerange vector ----"
551 print*,"elements=",size(this%timerange)
552 do i =1,size(this%timerange)
553 call display(this%timerange(i))
554 end do
555end if
556
557
558print*,"---- ana vector ----"
559print*,""
560print*,"->>>>>>>>> anavar -"
561call display(this%anavar)
562print*,""
563print*,"->>>>>>>>> anaattr -"
564call display(this%anaattr)
565print*,""
566print*,"->>>>>>>>> anavarattr -"
567call display(this%anavarattr)
568
569print*,"-- ana data section (first point) --"
570
571idat=imiss
572rdat=rmiss
573ddat=dmiss
574bdat=ibmiss
575cdat=cmiss
576
577!ntime = MIN(SIZE(this%time),nprint)
578!ntimerange = MIN(SIZE(this%timerange),nprint)
579!nlevel = MIN(SIZE(this%level),nprint)
580!nnetwork = MIN(SIZE(this%network),nprint)
581!nana = MIN(SIZE(this%ana),nprint)
582
583IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0) THEN
584if (associated(this%volanai)) then
585 do i=1,size(this%anavar%i)
586 idat=this%volanai(1,i,1)
587 if (associated(this%anavar%i)) call display(this%anavar%i(i),idat,rdat,ddat,bdat,cdat)
588 end do
589end if
590idat=imiss
591
592if (associated(this%volanar)) then
593 do i=1,size(this%anavar%r)
594 rdat=this%volanar(1,i,1)
595 if (associated(this%anavar%r)) call display(this%anavar%r(i),idat,rdat,ddat,bdat,cdat)
596 end do
597end if
598rdat=rmiss
599
600if (associated(this%volanad)) then
601 do i=1,size(this%anavar%d)
602 ddat=this%volanad(1,i,1)
603 if (associated(this%anavar%d)) call display(this%anavar%d(i),idat,rdat,ddat,bdat,cdat)
604 end do
605end if
606ddat=dmiss
607
608if (associated(this%volanab)) then
609 do i=1,size(this%anavar%b)
610 bdat=this%volanab(1,i,1)
611 if (associated(this%anavar%b)) call display(this%anavar%b(i),idat,rdat,ddat,bdat,cdat)
612 end do
613end if
614bdat=ibmiss
615
616if (associated(this%volanac)) then
617 do i=1,size(this%anavar%c)
618 cdat=this%volanac(1,i,1)
619 if (associated(this%anavar%c)) call display(this%anavar%c(i),idat,rdat,ddat,bdat,cdat)
620 end do
621end if
622cdat=cmiss
623ENDIF
624
625print*,"---- data vector ----"
626print*,""
627print*,"->>>>>>>>> dativar -"
628call display(this%dativar)
629print*,""
630print*,"->>>>>>>>> datiattr -"
631call display(this%datiattr)
632print*,""
633print*,"->>>>>>>>> dativarattr -"
634call display(this%dativarattr)
635
636print*,"-- data data section (first point) --"
637
638idat=imiss
639rdat=rmiss
640ddat=dmiss
641bdat=ibmiss
642cdat=cmiss
643
644IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0 .AND. size(this%time) > 0 &
645 .AND. size(this%level) > 0 .AND. size(this%timerange) > 0) THEN
646if (associated(this%voldatii)) then
647 do i=1,size(this%dativar%i)
648 idat=this%voldatii(1,1,1,1,i,1)
649 if (associated(this%dativar%i)) call display(this%dativar%i(i),idat,rdat,ddat,bdat,cdat)
650 end do
651end if
652idat=imiss
653
654if (associated(this%voldatir)) then
655 do i=1,size(this%dativar%r)
656 rdat=this%voldatir(1,1,1,1,i,1)
657 if (associated(this%dativar%r)) call display(this%dativar%r(i),idat,rdat,ddat,bdat,cdat)
658 end do
659end if
660rdat=rmiss
661
662if (associated(this%voldatid)) then
663 do i=1,size(this%dativar%d)
664 ddat=this%voldatid(1,1,1,1,i,1)
665 if (associated(this%dativar%d)) call display(this%dativar%d(i),idat,rdat,ddat,bdat,cdat)
666 end do
667end if
668ddat=dmiss
669
670if (associated(this%voldatib)) then
671 do i=1,size(this%dativar%b)
672 bdat=this%voldatib(1,1,1,1,i,1)
673 if (associated(this%dativar%b)) call display(this%dativar%b(i),idat,rdat,ddat,bdat,cdat)
674 end do
675end if
676bdat=ibmiss
677
678if (associated(this%voldatic)) then
679 do i=1,size(this%dativar%c)
680 cdat=this%voldatic(1,1,1,1,i,1)
681 if (associated(this%dativar%c)) call display(this%dativar%c(i),idat,rdat,ddat,bdat,cdat)
682 end do
683end if
684cdat=cmiss
685ENDIF
686
687print*,"<<<<<<<<<<<<<<<<<<< END vol7d object >>>>>>>>>>>>>>>>>>>>"
688
689END SUBROUTINE vol7d_display
690
691
693SUBROUTINE dat_display(this,idat,rdat,ddat,bdat,cdat)
694TYPE(vol7d_var),intent(in) :: this
696REAL :: rdat
698DOUBLE PRECISION :: ddat
700INTEGER :: idat
702INTEGER(kind=int_b) :: bdat
704CHARACTER(len=*) :: cdat
705
706print *, to_char_dat(this,idat,rdat,ddat,bdat,cdat)
707
708end SUBROUTINE dat_display
709
711SUBROUTINE dat_vect_display(this,idat,rdat,ddat,bdat,cdat)
712
713TYPE(vol7d_var),intent(in) :: this(:)
715REAL :: rdat(:)
717DOUBLE PRECISION :: ddat(:)
719INTEGER :: idat(:)
721INTEGER(kind=int_b) :: bdat(:)
723CHARACTER(len=*):: cdat(:)
724
725integer :: i
726
727do i =1,size(this)
728 call display(this(i),idat(i),rdat(i),ddat(i),bdat(i),cdat(i))
729end do
730
731end SUBROUTINE dat_vect_display
732
733
734FUNCTION to_char_dat(this,idat,rdat,ddat,bdat,cdat)
735#ifdef HAVE_DBALLE
736USE dballef
737#endif
738TYPE(vol7d_var),INTENT(in) :: this
740REAL :: rdat
742DOUBLE PRECISION :: ddat
744INTEGER :: idat
746INTEGER(kind=int_b) :: bdat
748CHARACTER(len=*) :: cdat
749CHARACTER(len=80) :: to_char_dat
750
751CHARACTER(len=LEN(to_char_dat)) :: to_char_tmp
752
753
754#ifdef HAVE_DBALLE
755INTEGER :: handle, ier
756
757handle = 0
758to_char_dat="VALUE: "
759
760if (c_e(idat)) to_char_dat=trim(to_char_dat)//" ;int> "//trim(to_char(idat))
761if (c_e(rdat)) to_char_dat=trim(to_char_dat)//" ;real> "//trim(to_char(rdat))
762if (c_e(ddat)) to_char_dat=trim(to_char_dat)//" ;double> "//trim(to_char(ddat))
763if (c_e(bdat)) to_char_dat=trim(to_char_dat)//" ;byte> "//trim(to_char(bdat))
764
765if ( c_e(cdat))then
766 ier = idba_messaggi(handle,"/dev/null", "w", "BUFR")
767 ier = idba_spiegab(handle,this%btable,cdat,to_char_tmp)
768 ier = idba_fatto(handle)
769 to_char_dat=trim(to_char_dat)//" ;char> "//trim(to_char_tmp)
770endif
771
772#else
773
774to_char_dat="VALUE: "
775if (c_e(idat)) to_char_dat=trim(to_char_dat)//" ;int> "//trim(to_char(idat))
776if (c_e(rdat)) to_char_dat=trim(to_char_dat)//" ;real> "//trim(to_char(rdat))
777if (c_e(ddat)) to_char_dat=trim(to_char_dat)//" ;double> "//trim(to_char(ddat))
778if (c_e(bdat)) to_char_dat=trim(to_char_dat)//" ;byte> "//trim(to_char(bdat))
779if (c_e(cdat)) to_char_dat=trim(to_char_dat)//" ;char> "//trim(cdat)
780
781#endif
782
783END FUNCTION to_char_dat
784
785
788FUNCTION vol7d_c_e(this) RESULT(c_e)
789TYPE(vol7d), INTENT(in) :: this
790
791LOGICAL :: c_e
792
793c_e = ASSOCIATED(this%ana) .OR. ASSOCIATED(this%time) .OR. &
794 ASSOCIATED(this%level) .OR. ASSOCIATED(this%timerange) .OR. &
795 ASSOCIATED(this%network) .OR. &
796 ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
797 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
798 ASSOCIATED(this%anavar%c) .OR. &
799 ASSOCIATED(this%anaattr%r) .OR. ASSOCIATED(this%anaattr%d) .OR. &
800 ASSOCIATED(this%anaattr%i) .OR. ASSOCIATED(this%anaattr%b) .OR. &
801 ASSOCIATED(this%anaattr%c) .OR. &
802 ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
803 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
804 ASSOCIATED(this%dativar%c) .OR. &
805 ASSOCIATED(this%datiattr%r) .OR. ASSOCIATED(this%datiattr%d) .OR. &
806 ASSOCIATED(this%datiattr%i) .OR. ASSOCIATED(this%datiattr%b) .OR. &
807 ASSOCIATED(this%datiattr%c)
808
809END FUNCTION vol7d_c_e
810
811
850SUBROUTINE vol7d_alloc(this, nana, ntime, nlevel, ntimerange, nnetwork, &
851 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
852 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
853 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
854 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
855 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
856 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc, &
857 ini)
858TYPE(vol7d),INTENT(inout) :: this
859INTEGER,INTENT(in),OPTIONAL :: nana
860INTEGER,INTENT(in),OPTIONAL :: ntime
861INTEGER,INTENT(in),OPTIONAL :: nlevel
862INTEGER,INTENT(in),OPTIONAL :: ntimerange
863INTEGER,INTENT(in),OPTIONAL :: nnetwork
865INTEGER,INTENT(in),OPTIONAL :: &
866 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
867 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
868 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
869 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
870 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
871 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc
872LOGICAL,INTENT(in),OPTIONAL :: ini
873
874INTEGER :: i
875LOGICAL :: linit
876
877IF (PRESENT(ini)) THEN
878 linit = ini
879ELSE
880 linit = .false.
881ENDIF
882
883! Dimensioni principali
884IF (PRESENT(nana)) THEN
885 IF (nana >= 0) THEN
886 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
887 ALLOCATE(this%ana(nana))
888 IF (linit) THEN
889 DO i = 1, nana
890 CALL init(this%ana(i))
891 ENDDO
892 ENDIF
893 ENDIF
894ENDIF
895IF (PRESENT(ntime)) THEN
896 IF (ntime >= 0) THEN
897 IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
898 ALLOCATE(this%time(ntime))
899 IF (linit) THEN
900 DO i = 1, ntime
901 CALL init(this%time(i))
902 ENDDO
903 ENDIF
904 ENDIF
905ENDIF
906IF (PRESENT(nlevel)) THEN
907 IF (nlevel >= 0) THEN
908 IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
909 ALLOCATE(this%level(nlevel))
910 IF (linit) THEN
911 DO i = 1, nlevel
912 CALL init(this%level(i))
913 ENDDO
914 ENDIF
915 ENDIF
916ENDIF
917IF (PRESENT(ntimerange)) THEN
918 IF (ntimerange >= 0) THEN
919 IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
920 ALLOCATE(this%timerange(ntimerange))
921 IF (linit) THEN
922 DO i = 1, ntimerange
923 CALL init(this%timerange(i))
924 ENDDO
925 ENDIF
926 ENDIF
927ENDIF
928IF (PRESENT(nnetwork)) THEN
929 IF (nnetwork >= 0) THEN
930 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
931 ALLOCATE(this%network(nnetwork))
932 IF (linit) THEN
933 DO i = 1, nnetwork
934 CALL init(this%network(i))
935 ENDDO
936 ENDIF
937 ENDIF
938ENDIF
939! Dimensioni dei tipi delle variabili
940CALL vol7d_varvect_alloc(this%anavar, nanavarr, nanavard, &
941 nanavari, nanavarb, nanavarc, ini)
942CALL vol7d_varvect_alloc(this%anaattr, nanaattrr, nanaattrd, &
943 nanaattri, nanaattrb, nanaattrc, ini)
944CALL vol7d_varvect_alloc(this%anavarattr, nanavarattrr, nanavarattrd, &
945 nanavarattri, nanavarattrb, nanavarattrc, ini)
946CALL vol7d_varvect_alloc(this%dativar, ndativarr, ndativard, &
947 ndativari, ndativarb, ndativarc, ini)
948CALL vol7d_varvect_alloc(this%datiattr, ndatiattrr, ndatiattrd, &
949 ndatiattri, ndatiattrb, ndatiattrc, ini)
950CALL vol7d_varvect_alloc(this%dativarattr, ndativarattrr, ndativarattrd, &
951 ndativarattri, ndativarattrb, ndativarattrc, ini)
952
953END SUBROUTINE vol7d_alloc
954
955
956FUNCTION vol7d_check_alloc_ana(this)
957TYPE(vol7d),INTENT(in) :: this
958LOGICAL :: vol7d_check_alloc_ana
959
960vol7d_check_alloc_ana = ASSOCIATED(this%ana) .AND. ASSOCIATED(this%network)
961
962END FUNCTION vol7d_check_alloc_ana
963
964SUBROUTINE vol7d_force_alloc_ana(this, ini)
965TYPE(vol7d),INTENT(inout) :: this
966LOGICAL,INTENT(in),OPTIONAL :: ini
968! Alloco i descrittori minimi per avere un volume di anagrafica
969IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=1, ini=ini)
970IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=1, ini=ini)
971
972END SUBROUTINE vol7d_force_alloc_ana
973
974
975FUNCTION vol7d_check_alloc_dati(this)
976TYPE(vol7d),INTENT(in) :: this
977LOGICAL :: vol7d_check_alloc_dati
978
979vol7d_check_alloc_dati = vol7d_check_alloc_ana(this) .AND. &
980 ASSOCIATED(this%time) .AND. ASSOCIATED(this%level) .AND. &
981 ASSOCIATED(this%timerange)
982
983END FUNCTION vol7d_check_alloc_dati
984
985SUBROUTINE vol7d_force_alloc_dati(this, ini)
986TYPE(vol7d),INTENT(inout) :: this
987LOGICAL,INTENT(in),OPTIONAL :: ini
988
989! Alloco i descrittori minimi per avere un volume di dati
990CALL vol7d_force_alloc_ana(this, ini)
991IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=1, ini=ini)
992IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=1, ini=ini)
993IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=1, ini=ini)
994
995END SUBROUTINE vol7d_force_alloc_dati
996
997
998SUBROUTINE vol7d_force_alloc(this)
999TYPE(vol7d),INTENT(inout) :: this
1000
1001! If anything really not allocated yet, allocate with size 0
1002IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=0)
1003IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=0)
1004IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=0)
1005IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=0)
1006IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=0)
1007
1008END SUBROUTINE vol7d_force_alloc
1009
1010
1011FUNCTION vol7d_check_vol(this)
1012TYPE(vol7d),INTENT(in) :: this
1013LOGICAL :: vol7d_check_vol
1014
1015vol7d_check_vol = c_e(this)
1016
1017! Anagrafica
1018IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
1019 vol7d_check_vol = .false.
1020ENDIF
1021
1022IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
1023 vol7d_check_vol = .false.
1024ENDIF
1025
1026IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
1027 vol7d_check_vol = .false.
1028ENDIF
1030IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
1031 vol7d_check_vol = .false.
1032ENDIF
1033
1034IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
1035 vol7d_check_vol = .false.
1036ENDIF
1037IF (ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
1038 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
1039 ASSOCIATED(this%anavar%c)) THEN
1040 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_ana(this)
1041ENDIF
1042
1043! Attributi dell'anagrafica
1044IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
1045 .NOT.ASSOCIATED(this%volanaattrr)) THEN
1046 vol7d_check_vol = .false.
1047ENDIF
1048
1049IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
1050 .NOT.ASSOCIATED(this%volanaattrd)) THEN
1051 vol7d_check_vol = .false.
1052ENDIF
1053
1054IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
1055 .NOT.ASSOCIATED(this%volanaattri)) THEN
1056 vol7d_check_vol = .false.
1057ENDIF
1058
1059IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
1060 .NOT.ASSOCIATED(this%volanaattrb)) THEN
1061 vol7d_check_vol = .false.
1062ENDIF
1063
1064IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
1065 .NOT.ASSOCIATED(this%volanaattrc)) THEN
1066 vol7d_check_vol = .false.
1067ENDIF
1068
1069! Dati
1070IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
1071 vol7d_check_vol = .false.
1072ENDIF
1073
1074IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
1075 vol7d_check_vol = .false.
1076ENDIF
1077
1078IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
1079 vol7d_check_vol = .false.
1080ENDIF
1081
1082IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
1083 vol7d_check_vol = .false.
1084ENDIF
1085
1086IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
1087 vol7d_check_vol = .false.
1088ENDIF
1089
1090! Attributi dei dati
1091IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
1092 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
1093 vol7d_check_vol = .false.
1094ENDIF
1095
1096IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
1097 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
1098 vol7d_check_vol = .false.
1099ENDIF
1100
1101IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
1102 .NOT.ASSOCIATED(this%voldatiattri)) THEN
1103 vol7d_check_vol = .false.
1104ENDIF
1105
1106IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
1107 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
1108 vol7d_check_vol = .false.
1109ENDIF
1110
1111IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
1112 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
1113 vol7d_check_vol = .false.
1114ENDIF
1115IF (ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
1116 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
1117 ASSOCIATED(this%dativar%c)) THEN
1118 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_dati(this)
1119ENDIF
1120
1121END FUNCTION vol7d_check_vol
1122
1123
1138SUBROUTINE vol7d_alloc_vol(this, ini, inivol)
1139TYPE(vol7d),INTENT(inout) :: this
1140LOGICAL,INTENT(in),OPTIONAL :: ini
1141LOGICAL,INTENT(in),OPTIONAL :: inivol
1142
1143LOGICAL :: linivol
1144
1145IF (PRESENT(inivol)) THEN
1146 linivol = inivol
1147ELSE
1148 linivol = .true.
1149ENDIF
1150
1151! Anagrafica
1152IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
1153 CALL vol7d_force_alloc_ana(this, ini)
1154 ALLOCATE(this%volanar(SIZE(this%ana), SIZE(this%anavar%r), SIZE(this%network)))
1155 IF (linivol) this%volanar(:,:,:) = rmiss
1156ENDIF
1157
1158IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
1159 CALL vol7d_force_alloc_ana(this, ini)
1160 ALLOCATE(this%volanad(SIZE(this%ana), SIZE(this%anavar%d), SIZE(this%network)))
1161 IF (linivol) this%volanad(:,:,:) = rdmiss
1162ENDIF
1163
1164IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
1165 CALL vol7d_force_alloc_ana(this, ini)
1166 ALLOCATE(this%volanai(SIZE(this%ana), SIZE(this%anavar%i), SIZE(this%network)))
1167 IF (linivol) this%volanai(:,:,:) = imiss
1168ENDIF
1169
1170IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
1171 CALL vol7d_force_alloc_ana(this, ini)
1172 ALLOCATE(this%volanab(SIZE(this%ana), SIZE(this%anavar%b), SIZE(this%network)))
1173 IF (linivol) this%volanab(:,:,:) = ibmiss
1174ENDIF
1175
1176IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
1177 CALL vol7d_force_alloc_ana(this, ini)
1178 ALLOCATE(this%volanac(SIZE(this%ana), SIZE(this%anavar%c), SIZE(this%network)))
1179 IF (linivol) this%volanac(:,:,:) = cmiss
1180ENDIF
1181
1182! Attributi dell'anagrafica
1183IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
1184 .NOT.ASSOCIATED(this%volanaattrr)) THEN
1185 CALL vol7d_force_alloc_ana(this, ini)
1186 ALLOCATE(this%volanaattrr(SIZE(this%ana), SIZE(this%anavarattr%r), &
1187 SIZE(this%network), SIZE(this%anaattr%r)))
1188 IF (linivol) this%volanaattrr(:,:,:,:) = rmiss
1189ENDIF
1191IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
1192 .NOT.ASSOCIATED(this%volanaattrd)) THEN
1193 CALL vol7d_force_alloc_ana(this, ini)
1194 ALLOCATE(this%volanaattrd(SIZE(this%ana), SIZE(this%anavarattr%d), &
1195 SIZE(this%network), SIZE(this%anaattr%d)))
1196 IF (linivol) this%volanaattrd(:,:,:,:) = rdmiss
1197ENDIF
1198
1199IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
1200 .NOT.ASSOCIATED(this%volanaattri)) THEN
1201 CALL vol7d_force_alloc_ana(this, ini)
1202 ALLOCATE(this%volanaattri(SIZE(this%ana), SIZE(this%anavarattr%i), &
1203 SIZE(this%network), SIZE(this%anaattr%i)))
1204 IF (linivol) this%volanaattri(:,:,:,:) = imiss
1205ENDIF
1206
1207IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
1208 .NOT.ASSOCIATED(this%volanaattrb)) THEN
1209 CALL vol7d_force_alloc_ana(this, ini)
1210 ALLOCATE(this%volanaattrb(SIZE(this%ana), SIZE(this%anavarattr%b), &
1211 SIZE(this%network), SIZE(this%anaattr%b)))
1212 IF (linivol) this%volanaattrb(:,:,:,:) = ibmiss
1213ENDIF
1214
1215IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
1216 .NOT.ASSOCIATED(this%volanaattrc)) THEN
1217 CALL vol7d_force_alloc_ana(this, ini)
1218 ALLOCATE(this%volanaattrc(SIZE(this%ana), SIZE(this%anavarattr%c), &
1219 SIZE(this%network), SIZE(this%anaattr%c)))
1220 IF (linivol) this%volanaattrc(:,:,:,:) = cmiss
1221ENDIF
1222
1223! Dati
1224IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
1225 CALL vol7d_force_alloc_dati(this, ini)
1226 ALLOCATE(this%voldatir(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
1227 SIZE(this%timerange), SIZE(this%dativar%r), SIZE(this%network)))
1228 IF (linivol) this%voldatir(:,:,:,:,:,:) = rmiss
1229ENDIF
1230
1231IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
1232 CALL vol7d_force_alloc_dati(this, ini)
1233 ALLOCATE(this%voldatid(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
1234 SIZE(this%timerange), SIZE(this%dativar%d), SIZE(this%network)))
1235 IF (linivol) this%voldatid(:,:,:,:,:,:) = rdmiss
1236ENDIF
1237
1238IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
1239 CALL vol7d_force_alloc_dati(this, ini)
1240 ALLOCATE(this%voldatii(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
1241 SIZE(this%timerange), SIZE(this%dativar%i), SIZE(this%network)))
1242 IF (linivol) this%voldatii(:,:,:,:,:,:) = imiss
1243ENDIF
1244
1245IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
1246 CALL vol7d_force_alloc_dati(this, ini)
1247 ALLOCATE(this%voldatib(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
1248 SIZE(this%timerange), SIZE(this%dativar%b), SIZE(this%network)))
1249 IF (linivol) this%voldatib(:,:,:,:,:,:) = ibmiss
1250ENDIF
1251
1252IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
1253 CALL vol7d_force_alloc_dati(this, ini)
1254 ALLOCATE(this%voldatic(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
1255 SIZE(this%timerange), SIZE(this%dativar%c), SIZE(this%network)))
1256 IF (linivol) this%voldatic(:,:,:,:,:,:) = cmiss
1257ENDIF
1258
1259! Attributi dei dati
1260IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
1261 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
1262 CALL vol7d_force_alloc_dati(this, ini)
1263 ALLOCATE(this%voldatiattrr(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
1264 SIZE(this%timerange), SIZE(this%dativarattr%r), SIZE(this%network), &
1265 SIZE(this%datiattr%r)))
1266 IF (linivol) this%voldatiattrr(:,:,:,:,:,:,:) = rmiss
1267ENDIF
1268
1269IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
1270 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
1271 CALL vol7d_force_alloc_dati(this, ini)
1272 ALLOCATE(this%voldatiattrd(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
1273 SIZE(this%timerange), SIZE(this%dativarattr%d), SIZE(this%network), &
1274 SIZE(this%datiattr%d)))
1275 IF (linivol) this%voldatiattrd(:,:,:,:,:,:,:) = rdmiss
1276ENDIF
1277
1278IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
1279 .NOT.ASSOCIATED(this%voldatiattri)) THEN
1280 CALL vol7d_force_alloc_dati(this, ini)
1281 ALLOCATE(this%voldatiattri(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
1282 SIZE(this%timerange), SIZE(this%dativarattr%i), SIZE(this%network), &
1283 SIZE(this%datiattr%i)))
1284 IF (linivol) this%voldatiattri(:,:,:,:,:,:,:) = imiss
1285ENDIF
1286
1287IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
1288 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
1289 CALL vol7d_force_alloc_dati(this, ini)
1290 ALLOCATE(this%voldatiattrb(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
1291 SIZE(this%timerange), SIZE(this%dativarattr%b), SIZE(this%network), &
1292 SIZE(this%datiattr%b)))
1293 IF (linivol) this%voldatiattrb(:,:,:,:,:,:,:) = ibmiss
1294ENDIF
1295
1296IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
1297 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
1298 CALL vol7d_force_alloc_dati(this, ini)
1299 ALLOCATE(this%voldatiattrc(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
1300 SIZE(this%timerange), SIZE(this%dativarattr%c), SIZE(this%network), &
1301 SIZE(this%datiattr%c)))
1302 IF (linivol) this%voldatiattrc(:,:,:,:,:,:,:) = cmiss
1303ENDIF
1304
1305! Catch-all method
1306CALL vol7d_force_alloc(this)
1307
1308! Creo gli indici var-attr
1309
1310#ifdef DEBUG
1311CALL l4f_log(l4f_debug,"calling: vol7d_set_attr_ind")
1312#endif
1313
1314CALL vol7d_set_attr_ind(this)
1315
1316
1318END SUBROUTINE vol7d_alloc_vol
1319
1320
1327SUBROUTINE vol7d_set_attr_ind(this)
1328TYPE(vol7d),INTENT(inout) :: this
1329
1330INTEGER :: i
1331
1332! real
1333IF (ASSOCIATED(this%dativar%r)) THEN
1334 IF (ASSOCIATED(this%dativarattr%r)) THEN
1335 DO i = 1, SIZE(this%dativar%r)
1336 this%dativar%r(i)%r = &
1337 firsttrue(this%dativar%r(i)%btable == this%dativarattr%r(:)%btable)
1338 ENDDO
1339 ENDIF
1340
1341 IF (ASSOCIATED(this%dativarattr%d)) THEN
1342 DO i = 1, SIZE(this%dativar%r)
1343 this%dativar%r(i)%d = &
1344 firsttrue(this%dativar%r(i)%btable == this%dativarattr%d(:)%btable)
1345 ENDDO
1346 ENDIF
1347
1348 IF (ASSOCIATED(this%dativarattr%i)) THEN
1349 DO i = 1, SIZE(this%dativar%r)
1350 this%dativar%r(i)%i = &
1351 firsttrue(this%dativar%r(i)%btable == this%dativarattr%i(:)%btable)
1352 ENDDO
1353 ENDIF
1354
1355 IF (ASSOCIATED(this%dativarattr%b)) THEN
1356 DO i = 1, SIZE(this%dativar%r)
1357 this%dativar%r(i)%b = &
1358 firsttrue(this%dativar%r(i)%btable == this%dativarattr%b(:)%btable)
1359 ENDDO
1360 ENDIF
1361
1362 IF (ASSOCIATED(this%dativarattr%c)) THEN
1363 DO i = 1, SIZE(this%dativar%r)
1364 this%dativar%r(i)%c = &
1365 firsttrue(this%dativar%r(i)%btable == this%dativarattr%c(:)%btable)
1366 ENDDO
1367 ENDIF
1368ENDIF
1369! double
1370IF (ASSOCIATED(this%dativar%d)) THEN
1371 IF (ASSOCIATED(this%dativarattr%r)) THEN
1372 DO i = 1, SIZE(this%dativar%d)
1373 this%dativar%d(i)%r = &
1374 firsttrue(this%dativar%d(i)%btable == this%dativarattr%r(:)%btable)
1375 ENDDO
1376 ENDIF
1377
1378 IF (ASSOCIATED(this%dativarattr%d)) THEN
1379 DO i = 1, SIZE(this%dativar%d)
1380 this%dativar%d(i)%d = &
1381 firsttrue(this%dativar%d(i)%btable == this%dativarattr%d(:)%btable)
1382 ENDDO
1383 ENDIF
1384
1385 IF (ASSOCIATED(this%dativarattr%i)) THEN
1386 DO i = 1, SIZE(this%dativar%d)
1387 this%dativar%d(i)%i = &
1388 firsttrue(this%dativar%d(i)%btable == this%dativarattr%i(:)%btable)
1389 ENDDO
1390 ENDIF
1391
1392 IF (ASSOCIATED(this%dativarattr%b)) THEN
1393 DO i = 1, SIZE(this%dativar%d)
1394 this%dativar%d(i)%b = &
1395 firsttrue(this%dativar%d(i)%btable == this%dativarattr%b(:)%btable)
1396 ENDDO
1397 ENDIF
1398
1399 IF (ASSOCIATED(this%dativarattr%c)) THEN
1400 DO i = 1, SIZE(this%dativar%d)
1401 this%dativar%d(i)%c = &
1402 firsttrue(this%dativar%d(i)%btable == this%dativarattr%c(:)%btable)
1403 ENDDO
1404 ENDIF
1405ENDIF
1406! integer
1407IF (ASSOCIATED(this%dativar%i)) THEN
1408 IF (ASSOCIATED(this%dativarattr%r)) THEN
1409 DO i = 1, SIZE(this%dativar%i)
1410 this%dativar%i(i)%r = &
1411 firsttrue(this%dativar%i(i)%btable == this%dativarattr%r(:)%btable)
1412 ENDDO
1413 ENDIF
1414
1415 IF (ASSOCIATED(this%dativarattr%d)) THEN
1416 DO i = 1, SIZE(this%dativar%i)
1417 this%dativar%i(i)%d = &
1418 firsttrue(this%dativar%i(i)%btable == this%dativarattr%d(:)%btable)
1419 ENDDO
1420 ENDIF
1421
1422 IF (ASSOCIATED(this%dativarattr%i)) THEN
1423 DO i = 1, SIZE(this%dativar%i)
1424 this%dativar%i(i)%i = &
1425 firsttrue(this%dativar%i(i)%btable == this%dativarattr%i(:)%btable)
1426 ENDDO
1427 ENDIF
1428
1429 IF (ASSOCIATED(this%dativarattr%b)) THEN
1430 DO i = 1, SIZE(this%dativar%i)
1431 this%dativar%i(i)%b = &
1432 firsttrue(this%dativar%i(i)%btable == this%dativarattr%b(:)%btable)
1433 ENDDO
1434 ENDIF
1435
1436 IF (ASSOCIATED(this%dativarattr%c)) THEN
1437 DO i = 1, SIZE(this%dativar%i)
1438 this%dativar%i(i)%c = &
1439 firsttrue(this%dativar%i(i)%btable == this%dativarattr%c(:)%btable)
1440 ENDDO
1441 ENDIF
1442ENDIF
1443! byte
1444IF (ASSOCIATED(this%dativar%b)) THEN
1445 IF (ASSOCIATED(this%dativarattr%r)) THEN
1446 DO i = 1, SIZE(this%dativar%b)
1447 this%dativar%b(i)%r = &
1448 firsttrue(this%dativar%b(i)%btable == this%dativarattr%r(:)%btable)
1449 ENDDO
1450 ENDIF
1451
1452 IF (ASSOCIATED(this%dativarattr%d)) THEN
1453 DO i = 1, SIZE(this%dativar%b)
1454 this%dativar%b(i)%d = &
1455 firsttrue(this%dativar%b(i)%btable == this%dativarattr%d(:)%btable)
1456 ENDDO
1457 ENDIF
1458
1459 IF (ASSOCIATED(this%dativarattr%i)) THEN
1460 DO i = 1, SIZE(this%dativar%b)
1461 this%dativar%b(i)%i = &
1462 firsttrue(this%dativar%b(i)%btable == this%dativarattr%i(:)%btable)
1463 ENDDO
1464 ENDIF
1465
1466 IF (ASSOCIATED(this%dativarattr%b)) THEN
1467 DO i = 1, SIZE(this%dativar%b)
1468 this%dativar%b(i)%b = &
1469 firsttrue(this%dativar%b(i)%btable == this%dativarattr%b(:)%btable)
1470 ENDDO
1471 ENDIF
1472
1473 IF (ASSOCIATED(this%dativarattr%c)) THEN
1474 DO i = 1, SIZE(this%dativar%b)
1475 this%dativar%b(i)%c = &
1476 firsttrue(this%dativar%b(i)%btable == this%dativarattr%c(:)%btable)
1477 ENDDO
1478 ENDIF
1479ENDIF
1480! character
1481IF (ASSOCIATED(this%dativar%c)) THEN
1482 IF (ASSOCIATED(this%dativarattr%r)) THEN
1483 DO i = 1, SIZE(this%dativar%c)
1484 this%dativar%c(i)%r = &
1485 firsttrue(this%dativar%c(i)%btable == this%dativarattr%r(:)%btable)
1486 ENDDO
1487 ENDIF
1488
1489 IF (ASSOCIATED(this%dativarattr%d)) THEN
1490 DO i = 1, SIZE(this%dativar%c)
1491 this%dativar%c(i)%d = &
1492 firsttrue(this%dativar%c(i)%btable == this%dativarattr%d(:)%btable)
1493 ENDDO
1494 ENDIF
1495
1496 IF (ASSOCIATED(this%dativarattr%i)) THEN
1497 DO i = 1, SIZE(this%dativar%c)
1498 this%dativar%c(i)%i = &
1499 firsttrue(this%dativar%c(i)%btable == this%dativarattr%i(:)%btable)
1500 ENDDO
1501 ENDIF
1502
1503 IF (ASSOCIATED(this%dativarattr%b)) THEN
1504 DO i = 1, SIZE(this%dativar%c)
1505 this%dativar%c(i)%b = &
1506 firsttrue(this%dativar%c(i)%btable == this%dativarattr%b(:)%btable)
1507 ENDDO
1508 ENDIF
1509
1510 IF (ASSOCIATED(this%dativarattr%c)) THEN
1511 DO i = 1, SIZE(this%dativar%c)
1512 this%dativar%c(i)%c = &
1513 firsttrue(this%dativar%c(i)%btable == this%dativarattr%c(:)%btable)
1514 ENDDO
1515 ENDIF
1516ENDIF
1517
1518END SUBROUTINE vol7d_set_attr_ind
1519
1520
1525SUBROUTINE vol7d_merge(this, that, sort, bestdata, &
1526 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
1527TYPE(vol7d),INTENT(INOUT) :: this
1528TYPE(vol7d),INTENT(INOUT) :: that
1529LOGICAL,INTENT(IN),OPTIONAL :: sort
1530LOGICAL,INTENT(in),OPTIONAL :: bestdata
1531LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple ! experimental, please do not use outside the library now
1532
1533TYPE(vol7d) :: v7d_clean
1534
1535
1536IF (.NOT.c_e(this)) THEN ! speedup
1537 this = that
1538 CALL init(v7d_clean)
1539 that = v7d_clean ! destroy that without deallocating
1540ELSE ! Append that to this and destroy that
1541 CALL vol7d_append(this, that, sort, bestdata, &
1542 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
1543 CALL delete(that)
1544ENDIF
1545
1546END SUBROUTINE vol7d_merge
1547
1548
1577SUBROUTINE vol7d_append(this, that, sort, bestdata, &
1578 ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple)
1579TYPE(vol7d),INTENT(INOUT) :: this
1580TYPE(vol7d),INTENT(IN) :: that
1581LOGICAL,INTENT(IN),OPTIONAL :: sort
1582! experimental, please do not use outside the library now, they force the use
1583! of a simplified mapping algorithm which is valid only whene the dimension
1584! content is the same in both volumes , or when one of them is empty
1585LOGICAL,INTENT(in),OPTIONAL :: bestdata
1586LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple
1587
1588
1589TYPE(vol7d) :: v7dtmp
1590LOGICAL :: lsort, lbestdata
1591INTEGER,POINTER :: remapt1(:), remapt2(:), remaptr1(:), remaptr2(:), &
1592 remapl1(:), remapl2(:), remapa1(:), remapa2(:), remapn1(:), remapn2(:)
1593
1594IF (.NOT.c_e(that)) RETURN ! speedup, nothing to do
1595IF (.NOT.vol7d_check_vol(that)) RETURN ! be safe
1596IF (.NOT.c_e(this)) THEN ! this case is like a vol7d_copy, more efficient to copy?
1597 CALL vol7d_copy(that, this, sort=sort)
1598 RETURN
1599ENDIF
1600
1601IF (this%time_definition /= that%time_definition) THEN
1602 CALL l4f_log(l4f_fatal, &
1603 'in vol7d_append, cannot append volumes with different &
1604 &time definition')
1605 CALL raise_fatal_error()
1606ENDIF
1607
1608! Completo l'allocazione per avere volumi a norma
1609CALL vol7d_alloc_vol(this)
1610
1611CALL init(v7dtmp, time_definition=this%time_definition)
1612CALL optio(sort, lsort)
1613CALL optio(bestdata, lbestdata)
1614
1615! Calcolo le mappature tra volumi vecchi e volume nuovo
1616! I puntatori remap* vengono tutti o allocati o nullificati
1617IF (optio_log(ltimesimple)) THEN
1618 CALL vol7d_remap2simple_datetime(this%time, that%time, v7dtmp%time, &
1619 lsort, remapt1, remapt2)
1620ELSE
1621 CALL vol7d_remap2_datetime(this%time, that%time, v7dtmp%time, &
1622 lsort, remapt1, remapt2)
1623ENDIF
1624IF (optio_log(ltimerangesimple)) THEN
1625 CALL vol7d_remap2simple_vol7d_timerange(this%timerange, that%timerange, &
1626 v7dtmp%timerange, lsort, remaptr1, remaptr2)
1627ELSE
1628 CALL vol7d_remap2_vol7d_timerange(this%timerange, that%timerange, &
1629 v7dtmp%timerange, lsort, remaptr1, remaptr2)
1630ENDIF
1631IF (optio_log(llevelsimple)) THEN
1632 CALL vol7d_remap2simple_vol7d_level(this%level, that%level, v7dtmp%level, &
1633 lsort, remapl1, remapl2)
1634ELSE
1635 CALL vol7d_remap2_vol7d_level(this%level, that%level, v7dtmp%level, &
1636 lsort, remapl1, remapl2)
1637ENDIF
1638IF (optio_log(lanasimple)) THEN
1639 CALL vol7d_remap2simple_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
1640 .false., remapa1, remapa2)
1641ELSE
1642 CALL vol7d_remap2_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
1643 .false., remapa1, remapa2)
1644ENDIF
1645IF (optio_log(lnetworksimple)) THEN
1646 CALL vol7d_remap2simple_vol7d_network(this%network, that%network, v7dtmp%network, &
1647 .false., remapn1, remapn2)
1648ELSE
1649 CALL vol7d_remap2_vol7d_network(this%network, that%network, v7dtmp%network, &
1650 .false., remapn1, remapn2)
1651ENDIF
1652
1653! Faccio la fusione fisica dei volumi
1654CALL vol7d_merge_finalr(this, that, v7dtmp, &
1655 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
1656 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
1657CALL vol7d_merge_finald(this, that, v7dtmp, &
1658 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
1659 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
1660CALL vol7d_merge_finali(this, that, v7dtmp, &
1661 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
1662 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
1663CALL vol7d_merge_finalb(this, that, v7dtmp, &
1664 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
1665 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
1666CALL vol7d_merge_finalc(this, that, v7dtmp, &
1667 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
1668 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
1669
1670! Dealloco i vettori di rimappatura
1671IF (ASSOCIATED(remapt1)) DEALLOCATE(remapt1)
1672IF (ASSOCIATED(remapt2)) DEALLOCATE(remapt2)
1673IF (ASSOCIATED(remaptr1)) DEALLOCATE(remaptr1)
1674IF (ASSOCIATED(remaptr2)) DEALLOCATE(remaptr2)
1675IF (ASSOCIATED(remapl1)) DEALLOCATE(remapl1)
1676IF (ASSOCIATED(remapl2)) DEALLOCATE(remapl2)
1677IF (ASSOCIATED(remapa1)) DEALLOCATE(remapa1)
1678IF (ASSOCIATED(remapa2)) DEALLOCATE(remapa2)
1679IF (ASSOCIATED(remapn1)) DEALLOCATE(remapn1)
1680IF (ASSOCIATED(remapn2)) DEALLOCATE(remapn2)
1681
1682! Distruggo il vecchio volume e assegno il nuovo a this
1683CALL delete(this)
1684this = v7dtmp
1685! Ricreo gli indici var-attr
1686CALL vol7d_set_attr_ind(this)
1687
1688END SUBROUTINE vol7d_append
1689
1690
1723SUBROUTINE vol7d_copy(this, that, sort, unique, miss, &
1724 lsort_time, lsort_timerange, lsort_level, &
1725 ltime, ltimerange, llevel, lana, lnetwork, &
1726 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
1727 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
1728 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
1729 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
1730 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
1731 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
1732TYPE(vol7d),INTENT(IN) :: this
1733TYPE(vol7d),INTENT(INOUT) :: that
1734LOGICAL,INTENT(IN),OPTIONAL :: sort
1735LOGICAL,INTENT(IN),OPTIONAL :: unique
1736LOGICAL,INTENT(IN),OPTIONAL :: miss
1737LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
1738LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
1739LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
1747LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
1749LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
1751LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
1753LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
1755LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
1757LOGICAL,INTENT(in),OPTIONAL :: &
1758 lanavarr(:), lanavard(:), lanavari(:), lanavarb(:), lanavarc(:), &
1759 lanaattrr(:), lanaattrd(:), lanaattri(:), lanaattrb(:), lanaattrc(:), &
1760 lanavarattrr(:), lanavarattrd(:), lanavarattri(:), lanavarattrb(:), lanavarattrc(:), &
1761 ldativarr(:), ldativard(:), ldativari(:), ldativarb(:), ldativarc(:), &
1762 ldatiattrr(:), ldatiattrd(:), ldatiattri(:), ldatiattrb(:), ldatiattrc(:), &
1763 ldativarattrr(:), ldativarattrd(:), ldativarattri(:), ldativarattrb(:), ldativarattrc(:)
1764
1765LOGICAL :: lsort, lunique, lmiss
1766INTEGER,POINTER :: remapt(:), remaptr(:), remapl(:), remapa(:), remapn(:)
1767
1768CALL init(that)
1769IF (.NOT.c_e(this)) RETURN ! speedup, nothing to do
1770IF (.NOT.vol7d_check_vol(this)) RETURN ! be safe
1771
1772CALL optio(sort, lsort)
1773CALL optio(unique, lunique)
1774CALL optio(miss, lmiss)
1775
1776! Calcolo le mappature tra volume vecchio e volume nuovo
1777! I puntatori remap* vengono tutti o allocati o nullificati
1778CALL vol7d_remap1_datetime(this%time, that%time, &
1779 lsort.OR.optio_log(lsort_time), lunique, lmiss, remapt, ltime)
1780CALL vol7d_remap1_vol7d_timerange(this%timerange, that%timerange, &
1781 lsort.OR.optio_log(lsort_timerange), lunique, lmiss, remaptr, ltimerange)
1782CALL vol7d_remap1_vol7d_level(this%level, that%level, &
1783 lsort.OR.optio_log(lsort_level), lunique, lmiss, remapl, llevel)
1784CALL vol7d_remap1_vol7d_ana(this%ana, that%ana, &
1785 lsort, lunique, lmiss, remapa, lana)
1786CALL vol7d_remap1_vol7d_network(this%network, that%network, &
1787 lsort, lunique, lmiss, remapn, lnetwork)
1788
1789! lanavari, lanavarb, lanavarc, &
1790! lanaattri, lanaattrb, lanaattrc, &
1791! lanavarattri, lanavarattrb, lanavarattrc, &
1792! ldativari, ldativarb, ldativarc, &
1793! ldatiattri, ldatiattrb, ldatiattrc, &
1794! ldativarattri, ldativarattrb, ldativarattrc
1795! Faccio la riforma fisica dei volumi
1796CALL vol7d_reform_finalr(this, that, &
1797 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
1798 lanavarr, lanaattrr, lanavarattrr, ldativarr, ldatiattrr, ldativarattrr)
1799CALL vol7d_reform_finald(this, that, &
1800 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
1801 lanavard, lanaattrd, lanavarattrd, ldativard, ldatiattrd, ldativarattrd)
1802CALL vol7d_reform_finali(this, that, &
1803 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
1804 lanavari, lanaattri, lanavarattri, ldativari, ldatiattri, ldativarattri)
1805CALL vol7d_reform_finalb(this, that, &
1806 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
1807 lanavarb, lanaattrb, lanavarattrb, ldativarb, ldatiattrb, ldativarattrb)
1808CALL vol7d_reform_finalc(this, that, &
1809 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
1810 lanavarc, lanaattrc, lanavarattrc, ldativarc, ldatiattrc, ldativarattrc)
1811
1812! Dealloco i vettori di rimappatura
1813IF (ASSOCIATED(remapt)) DEALLOCATE(remapt)
1814IF (ASSOCIATED(remaptr)) DEALLOCATE(remaptr)
1815IF (ASSOCIATED(remapl)) DEALLOCATE(remapl)
1816IF (ASSOCIATED(remapa)) DEALLOCATE(remapa)
1817IF (ASSOCIATED(remapn)) DEALLOCATE(remapn)
1818
1819! Ricreo gli indici var-attr
1820CALL vol7d_set_attr_ind(that)
1821that%time_definition = this%time_definition
1822
1823END SUBROUTINE vol7d_copy
1824
1825
1836SUBROUTINE vol7d_reform(this, sort, unique, miss, &
1837 lsort_time, lsort_timerange, lsort_level, &
1838 ltime, ltimerange, llevel, lana, lnetwork, &
1839 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
1840 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
1841 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
1842 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
1843 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
1844 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc&
1845 ,purgeana)
1846TYPE(vol7d),INTENT(INOUT) :: this
1847LOGICAL,INTENT(IN),OPTIONAL :: sort
1848LOGICAL,INTENT(IN),OPTIONAL :: unique
1849LOGICAL,INTENT(IN),OPTIONAL :: miss
1850LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
1851LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
1852LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
1860LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
1861LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
1862LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
1863LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
1864LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
1866LOGICAL,INTENT(in),OPTIONAL :: &
1867 lanavarr(:), lanavard(:), lanavari(:), lanavarb(:), lanavarc(:), &
1868 lanaattrr(:), lanaattrd(:), lanaattri(:), lanaattrb(:), lanaattrc(:), &
1869 lanavarattrr(:), lanavarattrd(:), lanavarattri(:), lanavarattrb(:), lanavarattrc(:), &
1870 ldativarr(:), ldativard(:), ldativari(:), ldativarb(:), ldativarc(:), &
1871 ldatiattrr(:), ldatiattrd(:), ldatiattri(:), ldatiattrb(:), ldatiattrc(:), &
1872 ldativarattrr(:), ldativarattrd(:), ldativarattri(:), ldativarattrb(:), ldativarattrc(:)
1873LOGICAL,INTENT(IN),OPTIONAL :: purgeana
1874
1875TYPE(vol7d) :: v7dtmp
1876logical,allocatable :: llana(:)
1877integer :: i
1878
1879CALL vol7d_copy(this, v7dtmp, sort, unique, miss, &
1880 lsort_time, lsort_timerange, lsort_level, &
1881 ltime, ltimerange, llevel, lana, lnetwork, &
1882 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
1883 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
1884 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
1885 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
1886 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
1887 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
1888
1889! destroy old volume
1890CALL delete(this)
1891
1892if (optio_log(purgeana)) then
1893 allocate(llana(size(v7dtmp%ana)))
1894 llana =.false.
1895 do i =1,size(v7dtmp%ana)
1896 if (associated(v7dtmp%voldatii)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatii(i,:,:,:,:,:)))
1897 if (associated(v7dtmp%voldatir)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatir(i,:,:,:,:,:)))
1898 if (associated(v7dtmp%voldatid)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatid(i,:,:,:,:,:)))
1899 if (associated(v7dtmp%voldatib)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatib(i,:,:,:,:,:)))
1900 if (associated(v7dtmp%voldatic)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatic(i,:,:,:,:,:)))
1901 end do
1902 CALL vol7d_copy(v7dtmp, this,lana=llana)
1903 CALL delete(v7dtmp)
1904 deallocate(llana)
1905else
1906 this=v7dtmp
1907end if
1908
1909END SUBROUTINE vol7d_reform
1910
1911
1919SUBROUTINE vol7d_smart_sort(this, lsort_time, lsort_timerange, lsort_level)
1920TYPE(vol7d),INTENT(INOUT) :: this
1921LOGICAL,OPTIONAL,INTENT(in) :: lsort_time
1922LOGICAL,OPTIONAL,INTENT(in) :: lsort_timerange
1923LOGICAL,OPTIONAL,INTENT(in) :: lsort_level
1924
1925INTEGER :: i
1926LOGICAL :: to_be_sorted
1927
1928to_be_sorted = .false.
1929CALL vol7d_alloc_vol(this) ! usual safety check
1930
1931IF (optio_log(lsort_time)) THEN
1932 DO i = 2, SIZE(this%time)
1933 IF (this%time(i) < this%time(i-1)) THEN
1934 to_be_sorted = .true.
1935 EXIT
1936 ENDIF
1937 ENDDO
1938ENDIF
1939IF (optio_log(lsort_timerange)) THEN
1940 DO i = 2, SIZE(this%timerange)
1941 IF (this%timerange(i) < this%timerange(i-1)) THEN
1942 to_be_sorted = .true.
1943 EXIT
1944 ENDIF
1945 ENDDO
1946ENDIF
1947IF (optio_log(lsort_level)) THEN
1948 DO i = 2, SIZE(this%level)
1949 IF (this%level(i) < this%level(i-1)) THEN
1950 to_be_sorted = .true.
1951 EXIT
1952 ENDIF
1953 ENDDO
1954ENDIF
1955
1956IF (to_be_sorted) CALL vol7d_reform(this, &
1957 lsort_time=lsort_time, lsort_timerange=lsort_timerange, lsort_level=lsort_level )
1958
1959END SUBROUTINE vol7d_smart_sort
1960
1968SUBROUTINE vol7d_filter(this, avl, vl, nl, s_d, e_d)
1969TYPE(vol7d),INTENT(inout) :: this
1970CHARACTER(len=*),INTENT(in),OPTIONAL :: avl(:)
1971CHARACTER(len=*),INTENT(in),OPTIONAL :: vl(:)
1972TYPE(vol7d_network),OPTIONAL :: nl(:)
1973TYPE(datetime),INTENT(in),OPTIONAL :: s_d
1974TYPE(datetime),INTENT(in),OPTIONAL :: e_d
1975
1976INTEGER :: i
1977
1978IF (PRESENT(avl)) THEN
1979 IF (SIZE(avl) > 0) THEN
1980
1981 IF (ASSOCIATED(this%anavar%r)) THEN
1982 DO i = 1, SIZE(this%anavar%r)
1983 IF (all(this%anavar%r(i)%btable /= avl)) this%anavar%r(i) = vol7d_var_miss
1984 ENDDO
1985 ENDIF
1986
1987 IF (ASSOCIATED(this%anavar%i)) THEN
1988 DO i = 1, SIZE(this%anavar%i)
1989 IF (all(this%anavar%i(i)%btable /= avl)) this%anavar%i(i) = vol7d_var_miss
1990 ENDDO
1991 ENDIF
1992
1993 IF (ASSOCIATED(this%anavar%b)) THEN
1994 DO i = 1, SIZE(this%anavar%b)
1995 IF (all(this%anavar%b(i)%btable /= avl)) this%anavar%b(i) = vol7d_var_miss
1996 ENDDO
1997 ENDIF
1998
1999 IF (ASSOCIATED(this%anavar%d)) THEN
2000 DO i = 1, SIZE(this%anavar%d)
2001 IF (all(this%anavar%d(i)%btable /= avl)) this%anavar%d(i) = vol7d_var_miss
2002 ENDDO
2003 ENDIF
2004
2005 IF (ASSOCIATED(this%anavar%c)) THEN
2006 DO i = 1, SIZE(this%anavar%c)
2007 IF (all(this%anavar%c(i)%btable /= avl)) this%anavar%c(i) = vol7d_var_miss
2008 ENDDO
2009 ENDIF
2010
2011 ENDIF
2012ENDIF
2013
2014
2015IF (PRESENT(vl)) THEN
2016 IF (size(vl) > 0) THEN
2017 IF (ASSOCIATED(this%dativar%r)) THEN
2018 DO i = 1, SIZE(this%dativar%r)
2019 IF (all(this%dativar%r(i)%btable /= vl)) this%dativar%r(i) = vol7d_var_miss
2020 ENDDO
2021 ENDIF
2022
2023 IF (ASSOCIATED(this%dativar%i)) THEN
2024 DO i = 1, SIZE(this%dativar%i)
2025 IF (all(this%dativar%i(i)%btable /= vl)) this%dativar%i(i) = vol7d_var_miss
2026 ENDDO
2027 ENDIF
2028
2029 IF (ASSOCIATED(this%dativar%b)) THEN
2030 DO i = 1, SIZE(this%dativar%b)
2031 IF (all(this%dativar%b(i)%btable /= vl)) this%dativar%b(i) = vol7d_var_miss
2032 ENDDO
2033 ENDIF
2034
2035 IF (ASSOCIATED(this%dativar%d)) THEN
2036 DO i = 1, SIZE(this%dativar%d)
2037 IF (all(this%dativar%d(i)%btable /= vl)) this%dativar%d(i) = vol7d_var_miss
2038 ENDDO
2039 ENDIF
2040
2041 IF (ASSOCIATED(this%dativar%c)) THEN
2042 DO i = 1, SIZE(this%dativar%c)
2043 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
2044 ENDDO
2045 ENDIF
2046
2047 IF (ASSOCIATED(this%dativar%c)) THEN
2048 DO i = 1, SIZE(this%dativar%c)
2049 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
2050 ENDDO
2051 ENDIF
2052
2053 ENDIF
2054ENDIF
2055
2056IF (PRESENT(nl)) THEN
2057 IF (SIZE(nl) > 0) THEN
2058 DO i = 1, SIZE(this%network)
2059 IF (all(this%network(i) /= nl)) this%network(i) = vol7d_network_miss
2060 ENDDO
2061 ENDIF
2062ENDIF
2063
2064IF (PRESENT(s_d)) THEN
2065 IF (c_e(s_d)) THEN
2066 WHERE (this%time < s_d)
2067 this%time = datetime_miss
2068 END WHERE
2069 ENDIF
2070ENDIF
2071
2072IF (PRESENT(e_d)) THEN
2073 IF (c_e(e_d)) THEN
2074 WHERE (this%time > e_d)
2075 this%time = datetime_miss
2076 END WHERE
2077 ENDIF
2078ENDIF
2079
2080CALL vol7d_reform(this, miss=.true.)
2081
2082END SUBROUTINE vol7d_filter
2083
2084
2091SUBROUTINE vol7d_convr(this, that, anaconv)
2092TYPE(vol7d),INTENT(IN) :: this
2093TYPE(vol7d),INTENT(INOUT) :: that
2094LOGICAL,OPTIONAL,INTENT(in) :: anaconv
2095INTEGER :: i
2096LOGICAL :: fv(1)=(/.false./), tv(1)=(/.true./), acp(1), acn(1)
2097TYPE(vol7d) :: v7d_tmp
2099IF (optio_log(anaconv)) THEN
2100 acp=fv
2101 acn=tv
2102ELSE
2103 acp=tv
2104 acn=fv
2105ENDIF
2106
2107! Volume con solo i dati reali e tutti gli attributi
2108! l'anagrafica e` copiata interamente se necessario
2109CALL vol7d_copy(this, that, &
2110 lanavarr=tv, lanavard=acp, lanavari=acp, lanavarb=acp, lanavarc=acp, &
2111 ldativarr=tv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=fv)
2112
2113! Volume solo di dati double
2114CALL vol7d_copy(this, v7d_tmp, &
2115 lanavarr=fv, lanavard=acn, lanavari=fv, lanavarb=fv, lanavarc=fv, &
2116 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
2117 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
2118 ldativarr=fv, ldativard=tv, ldativari=fv, ldativarb=fv, ldativarc=fv, &
2119 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
2120 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
2121
2122! converto a dati reali
2123IF (ASSOCIATED(v7d_tmp%anavar%d) .OR. ASSOCIATED(v7d_tmp%dativar%d)) THEN
2124
2125 IF (ASSOCIATED(v7d_tmp%anavar%d)) THEN
2126! alloco i dati reali e vi trasferisco i double
2127 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanad, 1), SIZE(v7d_tmp%volanad, 2), &
2128 SIZE(v7d_tmp%volanad, 3)))
2129 DO i = 1, SIZE(v7d_tmp%anavar%d)
2130 v7d_tmp%volanar(:,i,:) = &
2131 realdat(v7d_tmp%volanad(:,i,:), v7d_tmp%anavar%d(i))
2132 ENDDO
2133 DEALLOCATE(v7d_tmp%volanad)
2134! trasferisco le variabili
2135 v7d_tmp%anavar%r => v7d_tmp%anavar%d
2136 NULLIFY(v7d_tmp%anavar%d)
2137 ENDIF
2138
2139 IF (ASSOCIATED(v7d_tmp%dativar%d)) THEN
2140! alloco i dati reali e vi trasferisco i double
2141 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatid, 1), SIZE(v7d_tmp%voldatid, 2), &
2142 SIZE(v7d_tmp%voldatid, 3), SIZE(v7d_tmp%voldatid, 4), SIZE(v7d_tmp%voldatid, 5), &
2143 SIZE(v7d_tmp%voldatid, 6)))
2144 DO i = 1, SIZE(v7d_tmp%dativar%d)
2145 v7d_tmp%voldatir(:,:,:,:,i,:) = &
2146 realdat(v7d_tmp%voldatid(:,:,:,:,i,:), v7d_tmp%dativar%d(i))
2147 ENDDO
2148 DEALLOCATE(v7d_tmp%voldatid)
2149! trasferisco le variabili
2150 v7d_tmp%dativar%r => v7d_tmp%dativar%d
2151 NULLIFY(v7d_tmp%dativar%d)
2152 ENDIF
2153
2154! fondo con il volume definitivo
2155 CALL vol7d_merge(that, v7d_tmp)
2156ELSE
2157 CALL delete(v7d_tmp)
2158ENDIF
2159
2160
2161! Volume solo di dati interi
2162CALL vol7d_copy(this, v7d_tmp, &
2163 lanavarr=fv, lanavard=fv, lanavari=acn, lanavarb=fv, lanavarc=fv, &
2164 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
2165 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
2166 ldativarr=fv, ldativard=fv, ldativari=tv, ldativarb=fv, ldativarc=fv, &
2167 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
2168 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
2169
2170! converto a dati reali
2171IF (ASSOCIATED(v7d_tmp%anavar%i) .OR. ASSOCIATED(v7d_tmp%dativar%i)) THEN
2172
2173 IF (ASSOCIATED(v7d_tmp%anavar%i)) THEN
2174! alloco i dati reali e vi trasferisco gli interi
2175 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanai, 1), SIZE(v7d_tmp%volanai, 2), &
2176 SIZE(v7d_tmp%volanai, 3)))
2177 DO i = 1, SIZE(v7d_tmp%anavar%i)
2178 v7d_tmp%volanar(:,i,:) = &
2179 realdat(v7d_tmp%volanai(:,i,:), v7d_tmp%anavar%i(i))
2180 ENDDO
2181 DEALLOCATE(v7d_tmp%volanai)
2182! trasferisco le variabili
2183 v7d_tmp%anavar%r => v7d_tmp%anavar%i
2184 NULLIFY(v7d_tmp%anavar%i)
2185 ENDIF
2186
2187 IF (ASSOCIATED(v7d_tmp%dativar%i)) THEN
2188! alloco i dati reali e vi trasferisco gli interi
2189 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatii, 1), SIZE(v7d_tmp%voldatii, 2), &
2190 SIZE(v7d_tmp%voldatii, 3), SIZE(v7d_tmp%voldatii, 4), SIZE(v7d_tmp%voldatii, 5), &
2191 SIZE(v7d_tmp%voldatii, 6)))
2192 DO i = 1, SIZE(v7d_tmp%dativar%i)
2193 v7d_tmp%voldatir(:,:,:,:,i,:) = &
2194 realdat(v7d_tmp%voldatii(:,:,:,:,i,:), v7d_tmp%dativar%i(i))
2195 ENDDO
2196 DEALLOCATE(v7d_tmp%voldatii)
2197! trasferisco le variabili
2198 v7d_tmp%dativar%r => v7d_tmp%dativar%i
2199 NULLIFY(v7d_tmp%dativar%i)
2200 ENDIF
2201
2202! fondo con il volume definitivo
2203 CALL vol7d_merge(that, v7d_tmp)
2204ELSE
2205 CALL delete(v7d_tmp)
2206ENDIF
2207
2208
2209! Volume solo di dati byte
2210CALL vol7d_copy(this, v7d_tmp, &
2211 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=acn, lanavarc=fv, &
2212 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
2213 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
2214 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=tv, ldativarc=fv, &
2215 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
2216 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
2217
2218! converto a dati reali
2219IF (ASSOCIATED(v7d_tmp%anavar%b) .OR. ASSOCIATED(v7d_tmp%dativar%b)) THEN
2220
2221 IF (ASSOCIATED(v7d_tmp%anavar%b)) THEN
2222! alloco i dati reali e vi trasferisco i byte
2223 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanab, 1), SIZE(v7d_tmp%volanab, 2), &
2224 SIZE(v7d_tmp%volanab, 3)))
2225 DO i = 1, SIZE(v7d_tmp%anavar%b)
2226 v7d_tmp%volanar(:,i,:) = &
2227 realdat(v7d_tmp%volanab(:,i,:), v7d_tmp%anavar%b(i))
2228 ENDDO
2229 DEALLOCATE(v7d_tmp%volanab)
2230! trasferisco le variabili
2231 v7d_tmp%anavar%r => v7d_tmp%anavar%b
2232 NULLIFY(v7d_tmp%anavar%b)
2233 ENDIF
2234
2235 IF (ASSOCIATED(v7d_tmp%dativar%b)) THEN
2236! alloco i dati reali e vi trasferisco i byte
2237 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatib, 1), SIZE(v7d_tmp%voldatib, 2), &
2238 SIZE(v7d_tmp%voldatib, 3), SIZE(v7d_tmp%voldatib, 4), SIZE(v7d_tmp%voldatib, 5), &
2239 SIZE(v7d_tmp%voldatib, 6)))
2240 DO i = 1, SIZE(v7d_tmp%dativar%b)
2241 v7d_tmp%voldatir(:,:,:,:,i,:) = &
2242 realdat(v7d_tmp%voldatib(:,:,:,:,i,:), v7d_tmp%dativar%b(i))
2243 ENDDO
2244 DEALLOCATE(v7d_tmp%voldatib)
2245! trasferisco le variabili
2246 v7d_tmp%dativar%r => v7d_tmp%dativar%b
2247 NULLIFY(v7d_tmp%dativar%b)
2248 ENDIF
2249
2250! fondo con il volume definitivo
2251 CALL vol7d_merge(that, v7d_tmp)
2252ELSE
2253 CALL delete(v7d_tmp)
2254ENDIF
2255
2256
2257! Volume solo di dati character
2258CALL vol7d_copy(this, v7d_tmp, &
2259 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=fv, lanavarc=acn, &
2260 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
2261 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
2262 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=tv, &
2263 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
2264 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
2265
2266! converto a dati reali
2267IF (ASSOCIATED(v7d_tmp%anavar%c) .OR. ASSOCIATED(v7d_tmp%dativar%c)) THEN
2268
2269 IF (ASSOCIATED(v7d_tmp%anavar%c)) THEN
2270! alloco i dati reali e vi trasferisco i character
2271 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanac, 1), SIZE(v7d_tmp%volanac, 2), &
2272 SIZE(v7d_tmp%volanac, 3)))
2273 DO i = 1, SIZE(v7d_tmp%anavar%c)
2274 v7d_tmp%volanar(:,i,:) = &
2275 realdat(v7d_tmp%volanac(:,i,:), v7d_tmp%anavar%c(i))
2276 ENDDO
2277 DEALLOCATE(v7d_tmp%volanac)
2278! trasferisco le variabili
2279 v7d_tmp%anavar%r => v7d_tmp%anavar%c
2280 NULLIFY(v7d_tmp%anavar%c)
2281 ENDIF
2282
2283 IF (ASSOCIATED(v7d_tmp%dativar%c)) THEN
2284! alloco i dati reali e vi trasferisco i character
2285 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatic, 1), SIZE(v7d_tmp%voldatic, 2), &
2286 SIZE(v7d_tmp%voldatic, 3), SIZE(v7d_tmp%voldatic, 4), SIZE(v7d_tmp%voldatic, 5), &
2287 SIZE(v7d_tmp%voldatic, 6)))
2288 DO i = 1, SIZE(v7d_tmp%dativar%c)
2289 v7d_tmp%voldatir(:,:,:,:,i,:) = &
2290 realdat(v7d_tmp%voldatic(:,:,:,:,i,:), v7d_tmp%dativar%c(i))
2291 ENDDO
2292 DEALLOCATE(v7d_tmp%voldatic)
2293! trasferisco le variabili
2294 v7d_tmp%dativar%r => v7d_tmp%dativar%c
2295 NULLIFY(v7d_tmp%dativar%c)
2296 ENDIF
2297
2298! fondo con il volume definitivo
2299 CALL vol7d_merge(that, v7d_tmp)
2300ELSE
2301 CALL delete(v7d_tmp)
2302ENDIF
2303
2304END SUBROUTINE vol7d_convr
2305
2306
2310SUBROUTINE vol7d_diff_only (this, that, data_only,ana)
2311TYPE(vol7d),INTENT(IN) :: this
2312TYPE(vol7d),INTENT(OUT) :: that
2313logical , optional, intent(in) :: data_only
2314logical , optional, intent(in) :: ana
2315logical :: ldata_only,lana
2316
2317IF (PRESENT(data_only)) THEN
2318 ldata_only = data_only
2319ELSE
2320 ldata_only = .false.
2321ENDIF
2322
2323IF (PRESENT(ana)) THEN
2324 lana = ana
2325ELSE
2326 lana = .false.
2327ENDIF
2328
2329
2330#undef VOL7D_POLY_ARRAY
2331#define VOL7D_POLY_ARRAY voldati
2332#include "vol7d_class_diff.F90"
2333#undef VOL7D_POLY_ARRAY
2334#define VOL7D_POLY_ARRAY voldatiattr
2335#include "vol7d_class_diff.F90"
2336#undef VOL7D_POLY_ARRAY
2337
2338if ( .not. ldata_only) then
2339
2340#define VOL7D_POLY_ARRAY volana
2341#include "vol7d_class_diff.F90"
2342#undef VOL7D_POLY_ARRAY
2343#define VOL7D_POLY_ARRAY volanaattr
2344#include "vol7d_class_diff.F90"
2345#undef VOL7D_POLY_ARRAY
2346
2347 if(lana)then
2348 where ( this%ana == that%ana )
2349 that%ana = vol7d_ana_miss
2350 end where
2351 end if
2352
2353end if
2354
2355
2356
2357END SUBROUTINE vol7d_diff_only
2358
2359
2360
2361! Creo le routine da ripetere per i vari tipi di dati di v7d
2362! tramite un template e il preprocessore
2363#undef VOL7D_POLY_TYPE
2364#undef VOL7D_POLY_TYPES
2365#define VOL7D_POLY_TYPE REAL
2366#define VOL7D_POLY_TYPES r
2367#include "vol7d_class_type_templ.F90"
2368#undef VOL7D_POLY_TYPE
2369#undef VOL7D_POLY_TYPES
2370#define VOL7D_POLY_TYPE DOUBLE PRECISION
2371#define VOL7D_POLY_TYPES d
2372#include "vol7d_class_type_templ.F90"
2373#undef VOL7D_POLY_TYPE
2374#undef VOL7D_POLY_TYPES
2375#define VOL7D_POLY_TYPE INTEGER
2376#define VOL7D_POLY_TYPES i
2377#include "vol7d_class_type_templ.F90"
2378#undef VOL7D_POLY_TYPE
2379#undef VOL7D_POLY_TYPES
2380#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
2381#define VOL7D_POLY_TYPES b
2382#include "vol7d_class_type_templ.F90"
2383#undef VOL7D_POLY_TYPE
2384#undef VOL7D_POLY_TYPES
2385#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
2386#define VOL7D_POLY_TYPES c
2387#include "vol7d_class_type_templ.F90"
2388
2389! Creo le routine da ripetere per i vari descrittori di dimensioni di v7d
2390! tramite un template e il preprocessore
2391#define VOL7D_SORT
2392#undef VOL7D_NO_ZERO_ALLOC
2393#undef VOL7D_POLY_TYPE
2394#define VOL7D_POLY_TYPE datetime
2395#include "vol7d_class_desc_templ.F90"
2396#undef VOL7D_POLY_TYPE
2397#define VOL7D_POLY_TYPE vol7d_timerange
2398#include "vol7d_class_desc_templ.F90"
2399#undef VOL7D_POLY_TYPE
2400#define VOL7D_POLY_TYPE vol7d_level
2401#include "vol7d_class_desc_templ.F90"
2402#undef VOL7D_SORT
2403#undef VOL7D_POLY_TYPE
2404#define VOL7D_POLY_TYPE vol7d_network
2405#include "vol7d_class_desc_templ.F90"
2406#undef VOL7D_POLY_TYPE
2407#define VOL7D_POLY_TYPE vol7d_ana
2408#include "vol7d_class_desc_templ.F90"
2409#define VOL7D_NO_ZERO_ALLOC
2410#undef VOL7D_POLY_TYPE
2411#define VOL7D_POLY_TYPE vol7d_var
2412#include "vol7d_class_desc_templ.F90"
2413
2423subroutine vol7d_write_on_file (this,unit,description,filename,filename_auto)
2424
2425TYPE(vol7d),INTENT(IN) :: this
2426integer,optional,intent(inout) :: unit
2427character(len=*),intent(in),optional :: filename
2428character(len=*),intent(out),optional :: filename_auto
2429character(len=*),INTENT(IN),optional :: description
2430
2431integer :: lunit
2432character(len=254) :: ldescription,arg,lfilename
2433integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
2434 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
2435 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
2436 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
2437 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
2438 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
2439 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
2440!integer :: im,id,iy
2441integer :: tarray(8)
2442logical :: opened,exist
2443
2444 nana=0
2445 ntime=0
2446 ntimerange=0
2447 nlevel=0
2448 nnetwork=0
2449 ndativarr=0
2450 ndativari=0
2451 ndativarb=0
2452 ndativard=0
2453 ndativarc=0
2454 ndatiattrr=0
2455 ndatiattri=0
2456 ndatiattrb=0
2457 ndatiattrd=0
2458 ndatiattrc=0
2459 ndativarattrr=0
2460 ndativarattri=0
2461 ndativarattrb=0
2462 ndativarattrd=0
2463 ndativarattrc=0
2464 nanavarr=0
2465 nanavari=0
2466 nanavarb=0
2467 nanavard=0
2468 nanavarc=0
2469 nanaattrr=0
2470 nanaattri=0
2471 nanaattrb=0
2472 nanaattrd=0
2473 nanaattrc=0
2474 nanavarattrr=0
2475 nanavarattri=0
2476 nanavarattrb=0
2477 nanavarattrd=0
2478 nanavarattrc=0
2479
2480
2481!call idate(im,id,iy)
2482call date_and_time(values=tarray)
2483call getarg(0,arg)
2484
2485if (present(description))then
2486 ldescription=description
2487else
2488 ldescription="Vol7d generated by: "//trim(arg)
2489end if
2490
2491if (.not. present(unit))then
2492 lunit=getunit()
2493else
2494 if (unit==0)then
2495 lunit=getunit()
2496 unit=lunit
2497 else
2498 lunit=unit
2499 end if
2500end if
2501
2502lfilename=trim(arg)//".v7d"
2503if (index(arg,'/',back=.true.) > 0) lfilename=lfilename(index(arg,'/',back=.true.)+1 : )
2504
2505if (present(filename))then
2506 if (filename /= "")then
2507 lfilename=filename
2508 end if
2509end if
2510
2511if (present(filename_auto))filename_auto=lfilename
2512
2513
2514inquire(unit=lunit,opened=opened)
2515if (.not. opened) then
2516! inquire(file=lfilename, EXIST=exist)
2517! IF (exist) THEN
2518! CALL l4f_log(L4F_FATAL, &
2519! 'in vol7d_write_on_file, file exists, cannot open file '//TRIM(lfilename))
2520! CALL raise_fatal_error()
2521! ENDIF
2522 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM')
2523 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
2524end if
2525
2526if (associated(this%ana)) nana=size(this%ana)
2527if (associated(this%time)) ntime=size(this%time)
2528if (associated(this%timerange)) ntimerange=size(this%timerange)
2529if (associated(this%level)) nlevel=size(this%level)
2530if (associated(this%network)) nnetwork=size(this%network)
2531
2532if (associated(this%dativar%r)) ndativarr=size(this%dativar%r)
2533if (associated(this%dativar%i)) ndativari=size(this%dativar%i)
2534if (associated(this%dativar%b)) ndativarb=size(this%dativar%b)
2535if (associated(this%dativar%d)) ndativard=size(this%dativar%d)
2536if (associated(this%dativar%c)) ndativarc=size(this%dativar%c)
2537
2538if (associated(this%datiattr%r)) ndatiattrr=size(this%datiattr%r)
2539if (associated(this%datiattr%i)) ndatiattri=size(this%datiattr%i)
2540if (associated(this%datiattr%b)) ndatiattrb=size(this%datiattr%b)
2541if (associated(this%datiattr%d)) ndatiattrd=size(this%datiattr%d)
2542if (associated(this%datiattr%c)) ndatiattrc=size(this%datiattr%c)
2543
2544if (associated(this%dativarattr%r)) ndativarattrr=size(this%dativarattr%r)
2545if (associated(this%dativarattr%i)) ndativarattri=size(this%dativarattr%i)
2546if (associated(this%dativarattr%b)) ndativarattrb=size(this%dativarattr%b)
2547if (associated(this%dativarattr%d)) ndativarattrd=size(this%dativarattr%d)
2548if (associated(this%dativarattr%c)) ndativarattrc=size(this%dativarattr%c)
2549
2550if (associated(this%anavar%r)) nanavarr=size(this%anavar%r)
2551if (associated(this%anavar%i)) nanavari=size(this%anavar%i)
2552if (associated(this%anavar%b)) nanavarb=size(this%anavar%b)
2553if (associated(this%anavar%d)) nanavard=size(this%anavar%d)
2554if (associated(this%anavar%c)) nanavarc=size(this%anavar%c)
2555
2556if (associated(this%anaattr%r)) nanaattrr=size(this%anaattr%r)
2557if (associated(this%anaattr%i)) nanaattri=size(this%anaattr%i)
2558if (associated(this%anaattr%b)) nanaattrb=size(this%anaattr%b)
2559if (associated(this%anaattr%d)) nanaattrd=size(this%anaattr%d)
2560if (associated(this%anaattr%c)) nanaattrc=size(this%anaattr%c)
2561
2562if (associated(this%anavarattr%r)) nanavarattrr=size(this%anavarattr%r)
2563if (associated(this%anavarattr%i)) nanavarattri=size(this%anavarattr%i)
2564if (associated(this%anavarattr%b)) nanavarattrb=size(this%anavarattr%b)
2565if (associated(this%anavarattr%d)) nanavarattrd=size(this%anavarattr%d)
2566if (associated(this%anavarattr%c)) nanavarattrc=size(this%anavarattr%c)
2567
2568write(unit=lunit)ldescription
2569write(unit=lunit)tarray
2570
2571write(unit=lunit)&
2572 nana, ntime, ntimerange, nlevel, nnetwork, &
2573 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
2574 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
2575 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
2576 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
2577 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
2578 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
2579 this%time_definition
2580
2581
2582!write(unit=lunit)this
2583
2584
2585!! prime 5 dimensioni
2586if (associated(this%ana)) call write_unit(this%ana, lunit)
2587if (associated(this%time)) call write_unit(this%time, lunit)
2588if (associated(this%level)) write(unit=lunit)this%level
2589if (associated(this%timerange)) write(unit=lunit)this%timerange
2590if (associated(this%network)) write(unit=lunit)this%network
2591
2592 !! 6a dimensione: variabile dell'anagrafica e dei dati
2593 !! con relativi attributi e in 5 tipi diversi
2594
2595if (associated(this%anavar%r)) write(unit=lunit)this%anavar%r
2596if (associated(this%anavar%i)) write(unit=lunit)this%anavar%i
2597if (associated(this%anavar%b)) write(unit=lunit)this%anavar%b
2598if (associated(this%anavar%d)) write(unit=lunit)this%anavar%d
2599if (associated(this%anavar%c)) write(unit=lunit)this%anavar%c
2600
2601if (associated(this%anaattr%r)) write(unit=lunit)this%anaattr%r
2602if (associated(this%anaattr%i)) write(unit=lunit)this%anaattr%i
2603if (associated(this%anaattr%b)) write(unit=lunit)this%anaattr%b
2604if (associated(this%anaattr%d)) write(unit=lunit)this%anaattr%d
2605if (associated(this%anaattr%c)) write(unit=lunit)this%anaattr%c
2606
2607if (associated(this%anavarattr%r)) write(unit=lunit)this%anavarattr%r
2608if (associated(this%anavarattr%i)) write(unit=lunit)this%anavarattr%i
2609if (associated(this%anavarattr%b)) write(unit=lunit)this%anavarattr%b
2610if (associated(this%anavarattr%d)) write(unit=lunit)this%anavarattr%d
2611if (associated(this%anavarattr%c)) write(unit=lunit)this%anavarattr%c
2612
2613if (associated(this%dativar%r)) write(unit=lunit)this%dativar%r
2614if (associated(this%dativar%i)) write(unit=lunit)this%dativar%i
2615if (associated(this%dativar%b)) write(unit=lunit)this%dativar%b
2616if (associated(this%dativar%d)) write(unit=lunit)this%dativar%d
2617if (associated(this%dativar%c)) write(unit=lunit)this%dativar%c
2618
2619if (associated(this%datiattr%r)) write(unit=lunit)this%datiattr%r
2620if (associated(this%datiattr%i)) write(unit=lunit)this%datiattr%i
2621if (associated(this%datiattr%b)) write(unit=lunit)this%datiattr%b
2622if (associated(this%datiattr%d)) write(unit=lunit)this%datiattr%d
2623if (associated(this%datiattr%c)) write(unit=lunit)this%datiattr%c
2624
2625if (associated(this%dativarattr%r)) write(unit=lunit)this%dativarattr%r
2626if (associated(this%dativarattr%i)) write(unit=lunit)this%dativarattr%i
2627if (associated(this%dativarattr%b)) write(unit=lunit)this%dativarattr%b
2628if (associated(this%dativarattr%d)) write(unit=lunit)this%dativarattr%d
2629if (associated(this%dativarattr%c)) write(unit=lunit)this%dativarattr%c
2630
2631!! Volumi di valori e attributi per anagrafica e dati
2632
2633if (associated(this%volanar)) write(unit=lunit)this%volanar
2634if (associated(this%volanaattrr)) write(unit=lunit)this%volanaattrr
2635if (associated(this%voldatir)) write(unit=lunit)this%voldatir
2636if (associated(this%voldatiattrr)) write(unit=lunit)this%voldatiattrr
2637
2638if (associated(this%volanai)) write(unit=lunit)this%volanai
2639if (associated(this%volanaattri)) write(unit=lunit)this%volanaattri
2640if (associated(this%voldatii)) write(unit=lunit)this%voldatii
2641if (associated(this%voldatiattri)) write(unit=lunit)this%voldatiattri
2642
2643if (associated(this%volanab)) write(unit=lunit)this%volanab
2644if (associated(this%volanaattrb)) write(unit=lunit)this%volanaattrb
2645if (associated(this%voldatib)) write(unit=lunit)this%voldatib
2646if (associated(this%voldatiattrb)) write(unit=lunit)this%voldatiattrb
2647
2648if (associated(this%volanad)) write(unit=lunit)this%volanad
2649if (associated(this%volanaattrd)) write(unit=lunit)this%volanaattrd
2650if (associated(this%voldatid)) write(unit=lunit)this%voldatid
2651if (associated(this%voldatiattrd)) write(unit=lunit)this%voldatiattrd
2652
2653if (associated(this%volanac)) write(unit=lunit)this%volanac
2654if (associated(this%volanaattrc)) write(unit=lunit)this%volanaattrc
2655if (associated(this%voldatic)) write(unit=lunit)this%voldatic
2656if (associated(this%voldatiattrc)) write(unit=lunit)this%voldatiattrc
2657
2658if (.not. present(unit)) close(unit=lunit)
2659
2660end subroutine vol7d_write_on_file
2661
2662
2669
2670
2671subroutine vol7d_read_from_file (this,unit,filename,description,tarray,filename_auto)
2672
2673TYPE(vol7d),INTENT(OUT) :: this
2674integer,intent(inout),optional :: unit
2675character(len=*),INTENT(in),optional :: filename
2676character(len=*),intent(out),optional :: filename_auto
2677character(len=*),INTENT(out),optional :: description
2678integer,intent(out),optional :: tarray(8)
2679
2680
2681integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
2682 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
2683 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
2684 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
2685 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
2686 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
2687 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
2688
2689character(len=254) :: ldescription,lfilename,arg
2690integer :: ltarray(8),lunit,ios
2691logical :: opened,exist
2692
2693
2694call getarg(0,arg)
2695
2696if (.not. present(unit))then
2697 lunit=getunit()
2698else
2699 if (unit==0)then
2700 lunit=getunit()
2701 unit=lunit
2702 else
2703 lunit=unit
2704 end if
2705end if
2706
2707lfilename=trim(arg)//".v7d"
2708if (index(arg,'/',back=.true.) > 0) lfilename=lfilename(index(arg,'/',back=.true.)+1 : )
2709
2710if (present(filename))then
2711 if (filename /= "")then
2712 lfilename=filename
2713 end if
2714end if
2715
2716if (present(filename_auto))filename_auto=lfilename
2717
2718
2719inquire(unit=lunit,opened=opened)
2720IF (.NOT. opened) THEN
2721 inquire(file=lfilename,exist=exist)
2722 IF (.NOT.exist) THEN
2723 CALL l4f_log(l4f_fatal, &
2724 'in vol7d_read_from_file, file does not exists, cannot open')
2725 CALL raise_fatal_error()
2726 ENDIF
2727 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM', &
2728 status='OLD', action='READ')
2729 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
2730end if
2731
2732
2733call init(this)
2734read(unit=lunit,iostat=ios)ldescription
2735
2736if (ios < 0) then ! A negative value indicates that the End of File or End of Record
2737 call vol7d_alloc (this)
2738 call vol7d_alloc_vol (this)
2739 if (present(description))description=ldescription
2740 if (present(tarray))tarray=ltarray
2741 if (.not. present(unit)) close(unit=lunit)
2742end if
2743
2744read(unit=lunit)ltarray
2745
2746CALL l4f_log(l4f_info, 'Reading vol7d from file')
2747CALL l4f_log(l4f_info, 'description: '//trim(ldescription))
2748CALL l4f_log(l4f_info, 'written on '//trim(to_char(ltarray(1)))//' '// &
2749 trim(to_char(ltarray(2)))//' '//trim(to_char(ltarray(3))))
2750
2751if (present(description))description=ldescription
2752if (present(tarray))tarray=ltarray
2753
2754read(unit=lunit)&
2755 nana, ntime, ntimerange, nlevel, nnetwork, &
2756 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
2757 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
2758 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
2759 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
2760 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
2761 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
2762 this%time_definition
2763
2764call vol7d_alloc (this, &
2765 nana=nana, ntime=ntime, ntimerange=ntimerange, nlevel=nlevel, nnetwork=nnetwork,&
2766 ndativarr=ndativarr, ndativari=ndativari, ndativarb=ndativarb,&
2767 ndativard=ndativard, ndativarc=ndativarc,&
2768 ndatiattrr=ndatiattrr, ndatiattri=ndatiattri, ndatiattrb=ndatiattrb,&
2769 ndatiattrd=ndatiattrd, ndatiattrc=ndatiattrc,&
2770 ndativarattrr=ndativarattrr, ndativarattri=ndativarattri, ndativarattrb=ndativarattrb, &
2771 ndativarattrd=ndativarattrd, ndativarattrc=ndativarattrc,&
2772 nanavarr=nanavarr, nanavari=nanavari, nanavarb=nanavarb, &
2773 nanavard=nanavard, nanavarc=nanavarc,&
2774 nanaattrr=nanaattrr, nanaattri=nanaattri, nanaattrb=nanaattrb,&
2775 nanaattrd=nanaattrd, nanaattrc=nanaattrc,&
2776 nanavarattrr=nanavarattrr, nanavarattri=nanavarattri, nanavarattrb=nanavarattrb, &
2777 nanavarattrd=nanavarattrd, nanavarattrc=nanavarattrc)
2778
2779
2780if (associated(this%ana)) call read_unit(this%ana, lunit)
2781if (associated(this%time)) call read_unit(this%time, lunit)
2782if (associated(this%level)) read(unit=lunit)this%level
2783if (associated(this%timerange)) read(unit=lunit)this%timerange
2784if (associated(this%network)) read(unit=lunit)this%network
2785
2786if (associated(this%anavar%r)) read(unit=lunit)this%anavar%r
2787if (associated(this%anavar%i)) read(unit=lunit)this%anavar%i
2788if (associated(this%anavar%b)) read(unit=lunit)this%anavar%b
2789if (associated(this%anavar%d)) read(unit=lunit)this%anavar%d
2790if (associated(this%anavar%c)) read(unit=lunit)this%anavar%c
2791
2792if (associated(this%anaattr%r)) read(unit=lunit)this%anaattr%r
2793if (associated(this%anaattr%i)) read(unit=lunit)this%anaattr%i
2794if (associated(this%anaattr%b)) read(unit=lunit)this%anaattr%b
2795if (associated(this%anaattr%d)) read(unit=lunit)this%anaattr%d
2796if (associated(this%anaattr%c)) read(unit=lunit)this%anaattr%c
2797
2798if (associated(this%anavarattr%r)) read(unit=lunit)this%anavarattr%r
2799if (associated(this%anavarattr%i)) read(unit=lunit)this%anavarattr%i
2800if (associated(this%anavarattr%b)) read(unit=lunit)this%anavarattr%b
2801if (associated(this%anavarattr%d)) read(unit=lunit)this%anavarattr%d
2802if (associated(this%anavarattr%c)) read(unit=lunit)this%anavarattr%c
2803
2804if (associated(this%dativar%r)) read(unit=lunit)this%dativar%r
2805if (associated(this%dativar%i)) read(unit=lunit)this%dativar%i
2806if (associated(this%dativar%b)) read(unit=lunit)this%dativar%b
2807if (associated(this%dativar%d)) read(unit=lunit)this%dativar%d
2808if (associated(this%dativar%c)) read(unit=lunit)this%dativar%c
2809
2810if (associated(this%datiattr%r)) read(unit=lunit)this%datiattr%r
2811if (associated(this%datiattr%i)) read(unit=lunit)this%datiattr%i
2812if (associated(this%datiattr%b)) read(unit=lunit)this%datiattr%b
2813if (associated(this%datiattr%d)) read(unit=lunit)this%datiattr%d
2814if (associated(this%datiattr%c)) read(unit=lunit)this%datiattr%c
2815
2816if (associated(this%dativarattr%r)) read(unit=lunit)this%dativarattr%r
2817if (associated(this%dativarattr%i)) read(unit=lunit)this%dativarattr%i
2818if (associated(this%dativarattr%b)) read(unit=lunit)this%dativarattr%b
2819if (associated(this%dativarattr%d)) read(unit=lunit)this%dativarattr%d
2820if (associated(this%dativarattr%c)) read(unit=lunit)this%dativarattr%c
2821
2822call vol7d_alloc_vol (this)
2823
2824!! Volumi di valori e attributi per anagrafica e dati
2825
2826if (associated(this%volanar)) read(unit=lunit)this%volanar
2827if (associated(this%volanaattrr)) read(unit=lunit)this%volanaattrr
2828if (associated(this%voldatir)) read(unit=lunit)this%voldatir
2829if (associated(this%voldatiattrr)) read(unit=lunit)this%voldatiattrr
2830
2831if (associated(this%volanai)) read(unit=lunit)this%volanai
2832if (associated(this%volanaattri)) read(unit=lunit)this%volanaattri
2833if (associated(this%voldatii)) read(unit=lunit)this%voldatii
2834if (associated(this%voldatiattri)) read(unit=lunit)this%voldatiattri
2835
2836if (associated(this%volanab)) read(unit=lunit)this%volanab
2837if (associated(this%volanaattrb)) read(unit=lunit)this%volanaattrb
2838if (associated(this%voldatib)) read(unit=lunit)this%voldatib
2839if (associated(this%voldatiattrb)) read(unit=lunit)this%voldatiattrb
2840
2841if (associated(this%volanad)) read(unit=lunit)this%volanad
2842if (associated(this%volanaattrd)) read(unit=lunit)this%volanaattrd
2843if (associated(this%voldatid)) read(unit=lunit)this%voldatid
2844if (associated(this%voldatiattrd)) read(unit=lunit)this%voldatiattrd
2845
2846if (associated(this%volanac)) read(unit=lunit)this%volanac
2847if (associated(this%volanaattrc)) read(unit=lunit)this%volanaattrc
2848if (associated(this%voldatic)) read(unit=lunit)this%voldatic
2849if (associated(this%voldatiattrc)) read(unit=lunit)this%voldatiattrc
2850
2851if (.not. present(unit)) close(unit=lunit)
2852
2853end subroutine vol7d_read_from_file
2854
2855
2856! to double precision
2857elemental doubleprecision function doubledatd(voldat,var)
2858doubleprecision,intent(in) :: voldat
2859type(vol7d_var),intent(in) :: var
2860
2861doubledatd=voldat
2862
2863end function doubledatd
2864
2865
2866elemental doubleprecision function doubledatr(voldat,var)
2867real,intent(in) :: voldat
2868type(vol7d_var),intent(in) :: var
2869
2870if (c_e(voldat))then
2871 doubledatr=dble(voldat)
2872else
2873 doubledatr=dmiss
2874end if
2875
2876end function doubledatr
2877
2878
2879elemental doubleprecision function doubledati(voldat,var)
2880integer,intent(in) :: voldat
2881type(vol7d_var),intent(in) :: var
2882
2883if (c_e(voldat)) then
2884 if (c_e(var%scalefactor))then
2885 doubledati=dble(voldat)/10.d0**var%scalefactor
2886 else
2887 doubledati=dble(voldat)
2888 endif
2889else
2890 doubledati=dmiss
2891end if
2892
2893end function doubledati
2894
2895
2896elemental doubleprecision function doubledatb(voldat,var)
2897integer(kind=int_b),intent(in) :: voldat
2898type(vol7d_var),intent(in) :: var
2899
2900if (c_e(voldat)) then
2901 if (c_e(var%scalefactor))then
2902 doubledatb=dble(voldat)/10.d0**var%scalefactor
2903 else
2904 doubledatb=dble(voldat)
2905 endif
2906else
2907 doubledatb=dmiss
2908end if
2909
2910end function doubledatb
2911
2912
2913elemental doubleprecision function doubledatc(voldat,var)
2914CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
2915type(vol7d_var),intent(in) :: var
2916
2917doubledatc = c2d(voldat)
2918if (c_e(doubledatc) .and. c_e(var%scalefactor))then
2919 doubledatc=doubledatc/10.d0**var%scalefactor
2920end if
2921
2922end function doubledatc
2923
2924
2925! to integer
2926elemental integer function integerdatd(voldat,var)
2927doubleprecision,intent(in) :: voldat
2928type(vol7d_var),intent(in) :: var
2929
2930if (c_e(voldat))then
2931 if (c_e(var%scalefactor)) then
2932 integerdatd=nint(voldat*10d0**var%scalefactor)
2933 else
2934 integerdatd=nint(voldat)
2935 endif
2936else
2937 integerdatd=imiss
2938end if
2939
2940end function integerdatd
2941
2942
2943elemental integer function integerdatr(voldat,var)
2944real,intent(in) :: voldat
2945type(vol7d_var),intent(in) :: var
2946
2947if (c_e(voldat))then
2948 if (c_e(var%scalefactor)) then
2949 integerdatr=nint(voldat*10d0**var%scalefactor)
2950 else
2951 integerdatr=nint(voldat)
2952 endif
2953else
2954 integerdatr=imiss
2955end if
2956
2957end function integerdatr
2958
2959
2960elemental integer function integerdati(voldat,var)
2961integer,intent(in) :: voldat
2962type(vol7d_var),intent(in) :: var
2963
2964integerdati=voldat
2965
2966end function integerdati
2967
2968
2969elemental integer function integerdatb(voldat,var)
2970integer(kind=int_b),intent(in) :: voldat
2971type(vol7d_var),intent(in) :: var
2972
2973if (c_e(voldat))then
2974 integerdatb=voldat
2975else
2976 integerdatb=imiss
2977end if
2978
2979end function integerdatb
2980
2981
2982elemental integer function integerdatc(voldat,var)
2983CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
2984type(vol7d_var),intent(in) :: var
2985
2986integerdatc=c2i(voldat)
2987
2988end function integerdatc
2989
2990
2991! to real
2992elemental real function realdatd(voldat,var)
2993doubleprecision,intent(in) :: voldat
2994type(vol7d_var),intent(in) :: var
2995
2996if (c_e(voldat))then
2997 realdatd=real(voldat)
2998else
2999 realdatd=rmiss
3000end if
3001
3002end function realdatd
3003
3004
3005elemental real function realdatr(voldat,var)
3006real,intent(in) :: voldat
3007type(vol7d_var),intent(in) :: var
3008
3009realdatr=voldat
3010
3011end function realdatr
3012
3013
3014elemental real function realdati(voldat,var)
3015integer,intent(in) :: voldat
3016type(vol7d_var),intent(in) :: var
3017
3018if (c_e(voldat)) then
3019 if (c_e(var%scalefactor))then
3020 realdati=float(voldat)/10.**var%scalefactor
3021 else
3022 realdati=float(voldat)
3023 endif
3024else
3025 realdati=rmiss
3026end if
3027
3028end function realdati
3029
3030
3031elemental real function realdatb(voldat,var)
3032integer(kind=int_b),intent(in) :: voldat
3033type(vol7d_var),intent(in) :: var
3034
3035if (c_e(voldat)) then
3036 if (c_e(var%scalefactor))then
3037 realdatb=float(voldat)/10**var%scalefactor
3038 else
3039 realdatb=float(voldat)
3040 endif
3041else
3042 realdatb=rmiss
3043end if
3044
3045end function realdatb
3046
3047
3048elemental real function realdatc(voldat,var)
3049CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
3050type(vol7d_var),intent(in) :: var
3051
3052realdatc=c2r(voldat)
3053if (c_e(realdatc) .and. c_e(var%scalefactor))then
3054 realdatc=realdatc/10.**var%scalefactor
3055end if
3056
3057end function realdatc
3058
3059
3065FUNCTION realanavol(this, var) RESULT(vol)
3066TYPE(vol7d),INTENT(in) :: this
3067TYPE(vol7d_var),INTENT(in) :: var
3068REAL :: vol(SIZE(this%ana),size(this%network))
3069
3070CHARACTER(len=1) :: dtype
3071INTEGER :: indvar
3072
3073dtype = cmiss
3074indvar = index(this%anavar, var, type=dtype)
3075
3076IF (indvar > 0) THEN
3077 SELECT CASE (dtype)
3078 CASE("d")
3079 vol = realdat(this%volanad(:,indvar,:), var)
3080 CASE("r")
3081 vol = this%volanar(:,indvar,:)
3082 CASE("i")
3083 vol = realdat(this%volanai(:,indvar,:), var)
3084 CASE("b")
3085 vol = realdat(this%volanab(:,indvar,:), var)
3086 CASE("c")
3087 vol = realdat(this%volanac(:,indvar,:), var)
3088 CASE default
3089 vol = rmiss
3090 END SELECT
3091ELSE
3092 vol = rmiss
3093ENDIF
3094
3095END FUNCTION realanavol
3096
3097
3103FUNCTION integeranavol(this, var) RESULT(vol)
3104TYPE(vol7d),INTENT(in) :: this
3105TYPE(vol7d_var),INTENT(in) :: var
3106INTEGER :: vol(SIZE(this%ana),size(this%network))
3107
3108CHARACTER(len=1) :: dtype
3109INTEGER :: indvar
3110
3111dtype = cmiss
3112indvar = index(this%anavar, var, type=dtype)
3113
3114IF (indvar > 0) THEN
3115 SELECT CASE (dtype)
3116 CASE("d")
3117 vol = integerdat(this%volanad(:,indvar,:), var)
3118 CASE("r")
3119 vol = integerdat(this%volanar(:,indvar,:), var)
3120 CASE("i")
3121 vol = this%volanai(:,indvar,:)
3122 CASE("b")
3123 vol = integerdat(this%volanab(:,indvar,:), var)
3124 CASE("c")
3125 vol = integerdat(this%volanac(:,indvar,:), var)
3126 CASE default
3127 vol = imiss
3128 END SELECT
3129ELSE
3130 vol = imiss
3131ENDIF
3132
3133END FUNCTION integeranavol
3134
3135
3141subroutine move_datac (v7d,&
3142 indana,indtime,indlevel,indtimerange,indnetwork,&
3143 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
3144
3145TYPE(vol7d),intent(inout) :: v7d
3146
3147integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
3148integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
3149integer :: inddativar,inddativarattr
3150
3151
3152do inddativar=1,size(v7d%dativar%c)
3153
3154 if (c_e(v7d%voldatic(indana,indtime,indlevel,indtimerange,inddativar,indnetwork)) .and. &
3155 .not. c_e(v7d%voldatic(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
3156 ) then
3157
3158 ! dati
3159 v7d%voldatic &
3160 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
3161 v7d%voldatic &
3162 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
3163
3164
3165 ! attributi
3166 if (associated (v7d%dativarattr%i)) then
3167 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%c(inddativar))
3168 if (inddativarattr > 0 ) then
3169 v7d%voldatiattri &
3170 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
3171 v7d%voldatiattri &
3172 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
3173 end if
3174 end if
3175
3176 if (associated (v7d%dativarattr%r)) then
3177 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%c(inddativar))
3178 if (inddativarattr > 0 ) then
3179 v7d%voldatiattrr &
3180 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
3181 v7d%voldatiattrr &
3182 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
3183 end if
3184 end if
3185
3186 if (associated (v7d%dativarattr%d)) then
3187 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%c(inddativar))
3188 if (inddativarattr > 0 ) then
3189 v7d%voldatiattrd &
3190 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
3191 v7d%voldatiattrd &
3192 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
3193 end if
3194 end if
3195
3196 if (associated (v7d%dativarattr%b)) then
3197 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%c(inddativar))
3198 if (inddativarattr > 0 ) then
3199 v7d%voldatiattrb &
3200 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
3201 v7d%voldatiattrb &
3202 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
3203 end if
3204 end if
3205
3206 if (associated (v7d%dativarattr%c)) then
3207 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%c(inddativar))
3208 if (inddativarattr > 0 ) then
3209 v7d%voldatiattrc &
3210 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
3211 v7d%voldatiattrc &
3212 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
3213 end if
3214 end if
3215
3216 end if
3217
3218end do
3219
3220end subroutine move_datac
3221
3227subroutine move_datar (v7d,&
3228 indana,indtime,indlevel,indtimerange,indnetwork,&
3229 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
3230
3231TYPE(vol7d),intent(inout) :: v7d
3232
3233integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
3234integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
3235integer :: inddativar,inddativarattr
3236
3237
3238do inddativar=1,size(v7d%dativar%r)
3239
3240 if (c_e(v7d%voldatir(indana,indtime,indlevel,indtimerange,inddativar,indnetwork)) .and. &
3241 .not. c_e(v7d%voldatir(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
3242 ) then
3243
3244 ! dati
3245 v7d%voldatir &
3246 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
3247 v7d%voldatir &
3248 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
3249
3250
3251 ! attributi
3252 if (associated (v7d%dativarattr%i)) then
3253 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%r(inddativar))
3254 if (inddativarattr > 0 ) then
3255 v7d%voldatiattri &
3256 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
3257 v7d%voldatiattri &
3258 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
3259 end if
3260 end if
3261
3262 if (associated (v7d%dativarattr%r)) then
3263 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%r(inddativar))
3264 if (inddativarattr > 0 ) then
3265 v7d%voldatiattrr &
3266 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
3267 v7d%voldatiattrr &
3268 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
3269 end if
3270 end if
3271
3272 if (associated (v7d%dativarattr%d)) then
3273 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%r(inddativar))
3274 if (inddativarattr > 0 ) then
3275 v7d%voldatiattrd &
3276 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
3277 v7d%voldatiattrd &
3278 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
3279 end if
3280 end if
3281
3282 if (associated (v7d%dativarattr%b)) then
3283 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%r(inddativar))
3284 if (inddativarattr > 0 ) then
3285 v7d%voldatiattrb &
3286 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
3287 v7d%voldatiattrb &
3288 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
3289 end if
3290 end if
3291
3292 if (associated (v7d%dativarattr%c)) then
3293 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%r(inddativar))
3294 if (inddativarattr > 0 ) then
3295 v7d%voldatiattrc &
3296 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
3297 v7d%voldatiattrc &
3298 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
3299 end if
3300 end if
3301
3302 end if
3303
3304end do
3305
3306end subroutine move_datar
3307
3308
3322subroutine v7d_rounding(v7din,v7dout,level,timerange,nostatproc)
3323type(vol7d),intent(inout) :: v7din
3324type(vol7d),intent(out) :: v7dout
3325type(vol7d_level),intent(in),optional :: level(:)
3326type(vol7d_timerange),intent(in),optional :: timerange(:)
3327!logical,intent(in),optional :: merge !< if there are data on more then one almost equal levels or timeranges
3328!! will be merged POINT BY POINT with priority for the fird data found ordered by icreasing var index
3329logical,intent(in),optional :: nostatproc
3330
3331integer :: nana,nlevel,ntime,ntimerange,nnetwork,nbin
3332integer :: iana,ilevel,itimerange,indl,indt,itime,inetwork
3333type(vol7d_level) :: roundlevel(size(v7din%level))
3334type(vol7d_timerange) :: roundtimerange(size(v7din%timerange))
3335type(vol7d) :: v7d_tmp
3336
3337
3338nbin=0
3339
3340if (associated(v7din%dativar%r)) nbin = nbin + size(v7din%dativar%r)
3341if (associated(v7din%dativar%i)) nbin = nbin + size(v7din%dativar%i)
3342if (associated(v7din%dativar%d)) nbin = nbin + size(v7din%dativar%d)
3343if (associated(v7din%dativar%b)) nbin = nbin + size(v7din%dativar%b)
3344
3345call init(v7d_tmp)
3346
3347roundlevel=v7din%level
3348
3349if (present(level))then
3350 do ilevel = 1, size(v7din%level)
3351 if ((any(v7din%level(ilevel) .almosteq. level))) then
3352 roundlevel(ilevel)=level(1)
3353 end if
3354 end do
3355end if
3356
3357roundtimerange=v7din%timerange
3358
3359if (present(timerange))then
3360 do itimerange = 1, size(v7din%timerange)
3361 if ((any(v7din%timerange(itimerange) .almosteq. timerange))) then
3362 roundtimerange(itimerange)=timerange(1)
3363 end if
3364 end do
3365end if
3366
3367!set istantaneous values everywere
3368!preserve p1 for forecast time
3369if (optio_log(nostatproc)) then
3370 roundtimerange(:)%timerange=254
3371 roundtimerange(:)%p2=0
3372end if
3373
3374
3375nana=size(v7din%ana)
3376nlevel=count_distinct(roundlevel,back=.true.)
3377ntime=size(v7din%time)
3378ntimerange=count_distinct(roundtimerange,back=.true.)
3379nnetwork=size(v7din%network)
3380
3381call init(v7d_tmp)
3382
3383if (nbin == 0) then
3384 call copy(v7din,v7d_tmp)
3385else
3386 call vol7d_convr(v7din,v7d_tmp)
3387end if
3388
3389v7d_tmp%level=roundlevel
3390v7d_tmp%timerange=roundtimerange
3391
3392do ilevel=1, size(v7d_tmp%level)
3393 indl=index(v7d_tmp%level,roundlevel(ilevel))
3394 do itimerange=1,size(v7d_tmp%timerange)
3395 indt=index(v7d_tmp%timerange,roundtimerange(itimerange))
3396
3397 if (indl /= ilevel .or. indt /= itimerange) then
3398
3399 do iana=1, nana
3400 do itime=1,ntime
3401 do inetwork=1,nnetwork
3402
3403 if (nbin > 0) then
3404 call move_datar (v7d_tmp,&
3405 iana,itime,ilevel,itimerange,inetwork,&
3406 iana,itime,indl,indt,inetwork)
3407 else
3408 call move_datac (v7d_tmp,&
3409 iana,itime,ilevel,itimerange,inetwork,&
3410 iana,itime,indl,indt,inetwork)
3411 end if
3412
3413 end do
3414 end do
3415 end do
3416
3417 end if
3418
3419 end do
3420end do
3421
3422! set to missing level and time > nlevel
3423do ilevel=nlevel+1,size(v7d_tmp%level)
3424 call init (v7d_tmp%level(ilevel))
3425end do
3426
3427do itimerange=ntimerange+1,size(v7d_tmp%timerange)
3428 call init (v7d_tmp%timerange(itimerange))
3429end do
3430
3431!copy with remove
3432CALL copy(v7d_tmp,v7dout,miss=.true.,lsort_timerange=.true.,lsort_level=.true.)
3433CALL delete(v7d_tmp)
3434
3435!call display(v7dout)
3436
3437end subroutine v7d_rounding
3438
3439
3440END MODULE vol7d_class
3441
3447
3448
Set of functions that return a trimmed CHARACTER representation of the input variable.
Legge un oggetto datetime/timedelta o un vettore di oggetti datetime/timedelta da un file FORMATTED o...
Scrive un oggetto datetime/timedelta o un vettore di oggetti datetime/timedelta su un file FORMATTED ...
Index method.
Generic subroutine for checking OPTIONAL parameters.
Test for a missing volume.
Check for problems return 0 if all check passed print diagnostics with log4f.
Distruttore per la classe vol7d.
doubleprecision data conversion
Scrittura su file.
Costruttore per la classe vol7d.
integer data conversion
real data conversion
Reduce some dimensions (level and timerage) for semplification (rounding).
Represent data in a pretty string.
This module defines usefull general purpose function and subroutine.
Utilities for CHARACTER variables.
Classi per la gestione delle coordinate temporali.
Gestione degli errori.
Definition of constants related to I/O units.
Definition: io_units.F90:225
Definition of constants to be used for declaring variables of a desired type.
Definition: kinds.F90:245
classe per la gestione del logging
Module for quickly interpreting the OPTIONAL parameters passed to a subprogram.
Classe per la gestione dell'anagrafica di stazioni meteo e affini.
Classe per la gestione di un volume completo di dati osservati.
Classe per la gestione dei livelli verticali in osservazioni meteo e affini.
Classe per la gestione delle reti di stazioni per osservazioni meteo e affini.
Classe per la gestione degli intervalli temporali di osservazioni meteo e affini.
Classe per gestire un vettore di oggetti di tipo vol7d_var_class::vol7d_var.
Class for expressing an absolute time value.
Definisce l'anagrafica di una stazione.
Definisce un oggetto contenente i volumi anagrafica e dati e tutti i descrittori delle loro dimension...
Definisce il livello verticale di un'osservazione.
Definisce la rete a cui appartiene una stazione.
Definisce l'intervallo temporale di un'osservazione meteo.
Definisce un vettore di vol7d_var_class::vol7d_var per ogni tipo di dato supportato.

Generated with Doxygen.