libsim Versione 7.1.11
|
◆ move_datar()
Move data for all variables from one coordinate in the real volume to other. Only not missing data will be copyed and all attributes will be moved together. Usefull to colapse data spread in more indices (level or time or ....). After the move is possible to set to missing some descriptor and make a copy with miss=.true. to obtain a new vol7d with less data shape.
Definizione alla linea 9278 del file vol7d_class.F90. 9281! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
9282! authors:
9283! Davide Cesari <dcesari@arpa.emr.it>
9284! Paolo Patruno <ppatruno@arpa.emr.it>
9285
9286! This program is free software; you can redistribute it and/or
9287! modify it under the terms of the GNU General Public License as
9288! published by the Free Software Foundation; either version 2 of
9289! the License, or (at your option) any later version.
9290
9291! This program is distributed in the hope that it will be useful,
9292! but WITHOUT ANY WARRANTY; without even the implied warranty of
9293! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
9294! GNU General Public License for more details.
9295
9296! You should have received a copy of the GNU General Public License
9297! along with this program. If not, see <http://www.gnu.org/licenses/>.
9298#include "config.h"
9299
9311
9379IMPLICIT NONE
9380
9381
9382INTEGER, PARAMETER :: vol7d_maxdim_a = 3, vol7d_maxdim_aa = 4, &
9383 vol7d_maxdim_d = 6, vol7d_maxdim_ad = 7
9384
9385INTEGER, PARAMETER :: vol7d_ana_a=1
9386INTEGER, PARAMETER :: vol7d_var_a=2
9387INTEGER, PARAMETER :: vol7d_network_a=3
9388INTEGER, PARAMETER :: vol7d_attr_a=4
9389INTEGER, PARAMETER :: vol7d_ana_d=1
9390INTEGER, PARAMETER :: vol7d_time_d=2
9391INTEGER, PARAMETER :: vol7d_level_d=3
9392INTEGER, PARAMETER :: vol7d_timerange_d=4
9393INTEGER, PARAMETER :: vol7d_var_d=5
9394INTEGER, PARAMETER :: vol7d_network_d=6
9395INTEGER, PARAMETER :: vol7d_attr_d=7
9396INTEGER, PARAMETER :: vol7d_cdatalen=32
9397
9398TYPE vol7d_varmap
9399 INTEGER :: r, d, i, b, c
9400END TYPE vol7d_varmap
9401
9406 TYPE(vol7d_ana),POINTER :: ana(:)
9408 TYPE(datetime),POINTER :: time(:)
9410 TYPE(vol7d_level),POINTER :: level(:)
9412 TYPE(vol7d_timerange),POINTER :: timerange(:)
9414 TYPE(vol7d_network),POINTER :: network(:)
9416 TYPE(vol7d_varvect) :: anavar
9418 TYPE(vol7d_varvect) :: anaattr
9420 TYPE(vol7d_varvect) :: anavarattr
9422 TYPE(vol7d_varvect) :: dativar
9424 TYPE(vol7d_varvect) :: datiattr
9426 TYPE(vol7d_varvect) :: dativarattr
9427
9429 REAL,POINTER :: volanar(:,:,:)
9431 DOUBLE PRECISION,POINTER :: volanad(:,:,:)
9433 INTEGER,POINTER :: volanai(:,:,:)
9435 INTEGER(kind=int_b),POINTER :: volanab(:,:,:)
9437 CHARACTER(len=vol7d_cdatalen),POINTER :: volanac(:,:,:)
9438
9440 REAL,POINTER :: volanaattrr(:,:,:,:)
9442 DOUBLE PRECISION,POINTER :: volanaattrd(:,:,:,:)
9444 INTEGER,POINTER :: volanaattri(:,:,:,:)
9446 INTEGER(kind=int_b),POINTER :: volanaattrb(:,:,:,:)
9448 CHARACTER(len=vol7d_cdatalen),POINTER :: volanaattrc(:,:,:,:)
9449
9451 REAL,POINTER :: voldatir(:,:,:,:,:,:) ! sono i dati
9453 DOUBLE PRECISION,POINTER :: voldatid(:,:,:,:,:,:)
9455 INTEGER,POINTER :: voldatii(:,:,:,:,:,:)
9457 INTEGER(kind=int_b),POINTER :: voldatib(:,:,:,:,:,:)
9459 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatic(:,:,:,:,:,:)
9460
9462 REAL,POINTER :: voldatiattrr(:,:,:,:,:,:,:)
9464 DOUBLE PRECISION,POINTER :: voldatiattrd(:,:,:,:,:,:,:)
9466 INTEGER,POINTER :: voldatiattri(:,:,:,:,:,:,:)
9468 INTEGER(kind=int_b),POINTER :: voldatiattrb(:,:,:,:,:,:,:)
9470 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatiattrc(:,:,:,:,:,:,:)
9471
9473 integer :: time_definition
9474
9476
9481 MODULE PROCEDURE vol7d_init
9482END INTERFACE
9483
9486 MODULE PROCEDURE vol7d_delete
9487END INTERFACE
9488
9491 MODULE PROCEDURE vol7d_write_on_file
9492END INTERFACE
9493
9495INTERFACE import
9496 MODULE PROCEDURE vol7d_read_from_file
9497END INTERFACE
9498
9501 MODULE PROCEDURE vol7d_display, dat_display, dat_vect_display
9502END INTERFACE
9503
9506 MODULE PROCEDURE to_char_dat
9507END INTERFACE
9508
9511 MODULE PROCEDURE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
9512END INTERFACE
9513
9516 MODULE PROCEDURE realdatd,realdatr,realdati,realdatb,realdatc
9517END INTERFACE
9518
9521 MODULE PROCEDURE integerdatd,integerdatr,integerdati,integerdatb,integerdatc
9522END INTERFACE
9523
9526 MODULE PROCEDURE vol7d_copy
9527END INTERFACE
9528
9531 MODULE PROCEDURE vol7d_c_e
9532END INTERFACE
9533
9538 MODULE PROCEDURE vol7d_check
9539END INTERFACE
9540
9555 MODULE PROCEDURE v7d_rounding
9556END INTERFACE
9557
9558!!$INTERFACE get_volana
9559!!$ MODULE PROCEDURE vol7d_get_volanar, vol7d_get_volanad, vol7d_get_volanai, &
9560!!$ vol7d_get_volanab, vol7d_get_volanac
9561!!$END INTERFACE
9562!!$
9563!!$INTERFACE get_voldati
9564!!$ MODULE PROCEDURE vol7d_get_voldatir, vol7d_get_voldatid, vol7d_get_voldatii, &
9565!!$ vol7d_get_voldatib, vol7d_get_voldatic
9566!!$END INTERFACE
9567!!$
9568!!$INTERFACE get_volanaattr
9569!!$ MODULE PROCEDURE vol7d_get_volanaattrr, vol7d_get_volanaattrd, &
9570!!$ vol7d_get_volanaattri, vol7d_get_volanaattrb, vol7d_get_volanaattrc
9571!!$END INTERFACE
9572!!$
9573!!$INTERFACE get_voldatiattr
9574!!$ MODULE PROCEDURE vol7d_get_voldatiattrr, vol7d_get_voldatiattrd, &
9575!!$ vol7d_get_voldatiattri, vol7d_get_voldatiattrb, vol7d_get_voldatiattrc
9576!!$END INTERFACE
9577
9578PRIVATE vol7d_get_volr, vol7d_get_vold, vol7d_get_voli, vol7d_get_volb, &
9579 vol7d_get_volc, &
9580 volptr1dr, volptr2dr, volptr3dr, volptr4dr, volptr5dr, volptr6dr, volptr7dr, &
9581 volptr1dd, volptr2dd, volptr3dd, volptr4dd, volptr5dd, volptr6dd, volptr7dd, &
9582 volptr1di, volptr2di, volptr3di, volptr4di, volptr5di, volptr6di, volptr7di, &
9583 volptr1db, volptr2db, volptr3db, volptr4db, volptr5db, volptr6db, volptr7db, &
9584 volptr1dc, volptr2dc, volptr3dc, volptr4dc, volptr5dc, volptr6dc, volptr7dc, &
9585 vol7d_nullifyr, vol7d_nullifyd, vol7d_nullifyi, vol7d_nullifyb, vol7d_nullifyc, &
9586 vol7d_init, vol7d_delete, vol7d_write_on_file, vol7d_read_from_file, &
9587 vol7d_check_alloc_ana, vol7d_force_alloc_ana, &
9588 vol7d_check_alloc_dati, vol7d_force_alloc_dati, vol7d_force_alloc, &
9589 vol7d_display, dat_display, dat_vect_display, &
9590 to_char_dat, vol7d_check
9591
9592PRIVATE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
9593
9594PRIVATE vol7d_c_e
9595
9596CONTAINS
9597
9598
9603SUBROUTINE vol7d_init(this,time_definition)
9604TYPE(vol7d),intent(out) :: this
9605integer,INTENT(IN),OPTIONAL :: time_definition
9606
9613CALL vol7d_var_features_init() ! initialise var features table once
9614
9615NULLIFY(this%ana, this%time, this%level, this%timerange, this%network)
9616
9617NULLIFY(this%volanar, this%volanaattrr, this%voldatir, this%voldatiattrr)
9618NULLIFY(this%volanad, this%volanaattrd, this%voldatid, this%voldatiattrd)
9619NULLIFY(this%volanai, this%volanaattri, this%voldatii, this%voldatiattri)
9620NULLIFY(this%volanab, this%volanaattrb, this%voldatib, this%voldatiattrb)
9621NULLIFY(this%volanac, this%volanaattrc, this%voldatic, this%voldatiattrc)
9622
9623if(present(time_definition)) then
9624 this%time_definition=time_definition
9625else
9626 this%time_definition=1 !default to validity time
9627end if
9628
9629END SUBROUTINE vol7d_init
9630
9631
9635ELEMENTAL SUBROUTINE vol7d_delete(this, dataonly)
9636TYPE(vol7d),intent(inout) :: this
9637LOGICAL, INTENT(in), OPTIONAL :: dataonly
9638
9639
9640IF (.NOT. optio_log(dataonly)) THEN
9641 IF (ASSOCIATED(this%volanar)) DEALLOCATE(this%volanar)
9642 IF (ASSOCIATED(this%volanad)) DEALLOCATE(this%volanad)
9643 IF (ASSOCIATED(this%volanai)) DEALLOCATE(this%volanai)
9644 IF (ASSOCIATED(this%volanab)) DEALLOCATE(this%volanab)
9645 IF (ASSOCIATED(this%volanac)) DEALLOCATE(this%volanac)
9646 IF (ASSOCIATED(this%volanaattrr)) DEALLOCATE(this%volanaattrr)
9647 IF (ASSOCIATED(this%volanaattrd)) DEALLOCATE(this%volanaattrd)
9648 IF (ASSOCIATED(this%volanaattri)) DEALLOCATE(this%volanaattri)
9649 IF (ASSOCIATED(this%volanaattrb)) DEALLOCATE(this%volanaattrb)
9650 IF (ASSOCIATED(this%volanaattrc)) DEALLOCATE(this%volanaattrc)
9651ENDIF
9652IF (ASSOCIATED(this%voldatir)) DEALLOCATE(this%voldatir)
9653IF (ASSOCIATED(this%voldatid)) DEALLOCATE(this%voldatid)
9654IF (ASSOCIATED(this%voldatii)) DEALLOCATE(this%voldatii)
9655IF (ASSOCIATED(this%voldatib)) DEALLOCATE(this%voldatib)
9656IF (ASSOCIATED(this%voldatic)) DEALLOCATE(this%voldatic)
9657IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
9658IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
9659IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
9660IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
9661IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
9662
9663IF (.NOT. optio_log(dataonly)) THEN
9664 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
9665 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
9666ENDIF
9667IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
9668IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
9669IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
9670
9671IF (.NOT. optio_log(dataonly)) THEN
9675ENDIF
9679
9680END SUBROUTINE vol7d_delete
9681
9682
9683
9684integer function vol7d_check(this)
9685TYPE(vol7d),intent(in) :: this
9686integer :: i,j,k,l,m,n
9687
9688vol7d_check=0
9689
9690if (associated(this%voldatii)) then
9691do i = 1,size(this%voldatii,1)
9692 do j = 1,size(this%voldatii,2)
9693 do k = 1,size(this%voldatii,3)
9694 do l = 1,size(this%voldatii,4)
9695 do m = 1,size(this%voldatii,5)
9696 do n = 1,size(this%voldatii,6)
9697 if (this%voldatii(i,j,k,l,m,n) /= this%voldatii(i,j,k,l,m,n) ) then
9698 CALL l4f_log(l4f_warn,"check: abnormal value at voldatii("&
9700 vol7d_check=1
9701 end if
9702 end do
9703 end do
9704 end do
9705 end do
9706 end do
9707end do
9708end if
9709
9710
9711if (associated(this%voldatir)) then
9712do i = 1,size(this%voldatir,1)
9713 do j = 1,size(this%voldatir,2)
9714 do k = 1,size(this%voldatir,3)
9715 do l = 1,size(this%voldatir,4)
9716 do m = 1,size(this%voldatir,5)
9717 do n = 1,size(this%voldatir,6)
9718 if (this%voldatir(i,j,k,l,m,n) /= this%voldatir(i,j,k,l,m,n) ) then
9719 CALL l4f_log(l4f_warn,"check: abnormal value at voldatir("&
9721 vol7d_check=2
9722 end if
9723 end do
9724 end do
9725 end do
9726 end do
9727 end do
9728end do
9729end if
9730
9731if (associated(this%voldatid)) then
9732do i = 1,size(this%voldatid,1)
9733 do j = 1,size(this%voldatid,2)
9734 do k = 1,size(this%voldatid,3)
9735 do l = 1,size(this%voldatid,4)
9736 do m = 1,size(this%voldatid,5)
9737 do n = 1,size(this%voldatid,6)
9738 if (this%voldatid(i,j,k,l,m,n) /= this%voldatid(i,j,k,l,m,n) ) then
9739 CALL l4f_log(l4f_warn,"check: abnormal value at voldatid("&
9741 vol7d_check=3
9742 end if
9743 end do
9744 end do
9745 end do
9746 end do
9747 end do
9748end do
9749end if
9750
9751if (associated(this%voldatib)) then
9752do i = 1,size(this%voldatib,1)
9753 do j = 1,size(this%voldatib,2)
9754 do k = 1,size(this%voldatib,3)
9755 do l = 1,size(this%voldatib,4)
9756 do m = 1,size(this%voldatib,5)
9757 do n = 1,size(this%voldatib,6)
9758 if (this%voldatib(i,j,k,l,m,n) /= this%voldatib(i,j,k,l,m,n) ) then
9759 CALL l4f_log(l4f_warn,"check: abnormal value at voldatib("&
9761 vol7d_check=4
9762 end if
9763 end do
9764 end do
9765 end do
9766 end do
9767 end do
9768end do
9769end if
9770
9771end function vol7d_check
9772
9773
9774
9775!TODO da completare ! aborta se i volumi sono allocati a dimensione 0
9777SUBROUTINE vol7d_display(this)
9778TYPE(vol7d),intent(in) :: this
9779integer :: i
9780
9781REAL :: rdat
9782DOUBLE PRECISION :: ddat
9783INTEGER :: idat
9784INTEGER(kind=int_b) :: bdat
9785CHARACTER(len=vol7d_cdatalen) :: cdat
9786
9787
9788print*,"<<<<<<<<<<<<<<<<<<< vol7d object >>>>>>>>>>>>>>>>>>>>"
9789if (this%time_definition == 0) then
9790 print*,"TIME DEFINITION: time is reference time"
9791else if (this%time_definition == 1) then
9792 print*,"TIME DEFINITION: time is validity time"
9793else
9794 print*,"Time definition have a wrong walue:", this%time_definition
9795end if
9796
9797IF (ASSOCIATED(this%network))then
9798 print*,"---- network vector ----"
9799 print*,"elements=",size(this%network)
9800 do i=1, size(this%network)
9802 end do
9803end IF
9804
9805IF (ASSOCIATED(this%ana))then
9806 print*,"---- ana vector ----"
9807 print*,"elements=",size(this%ana)
9808 do i=1, size(this%ana)
9810 end do
9811end IF
9812
9813IF (ASSOCIATED(this%time))then
9814 print*,"---- time vector ----"
9815 print*,"elements=",size(this%time)
9816 do i=1, size(this%time)
9818 end do
9819end if
9820
9821IF (ASSOCIATED(this%level)) then
9822 print*,"---- level vector ----"
9823 print*,"elements=",size(this%level)
9824 do i =1,size(this%level)
9826 end do
9827end if
9828
9829IF (ASSOCIATED(this%timerange))then
9830 print*,"---- timerange vector ----"
9831 print*,"elements=",size(this%timerange)
9832 do i =1,size(this%timerange)
9834 end do
9835end if
9836
9837
9838print*,"---- ana vector ----"
9839print*,""
9840print*,"->>>>>>>>> anavar -"
9842print*,""
9843print*,"->>>>>>>>> anaattr -"
9845print*,""
9846print*,"->>>>>>>>> anavarattr -"
9848
9849print*,"-- ana data section (first point) --"
9850
9851idat=imiss
9852rdat=rmiss
9853ddat=dmiss
9854bdat=ibmiss
9855cdat=cmiss
9856
9857!ntime = MIN(SIZE(this%time),nprint)
9858!ntimerange = MIN(SIZE(this%timerange),nprint)
9859!nlevel = MIN(SIZE(this%level),nprint)
9860!nnetwork = MIN(SIZE(this%network),nprint)
9861!nana = MIN(SIZE(this%ana),nprint)
9862
9863IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0) THEN
9864if (associated(this%volanai)) then
9865 do i=1,size(this%anavar%i)
9866 idat=this%volanai(1,i,1)
9868 end do
9869end if
9870idat=imiss
9871
9872if (associated(this%volanar)) then
9873 do i=1,size(this%anavar%r)
9874 rdat=this%volanar(1,i,1)
9876 end do
9877end if
9878rdat=rmiss
9879
9880if (associated(this%volanad)) then
9881 do i=1,size(this%anavar%d)
9882 ddat=this%volanad(1,i,1)
9884 end do
9885end if
9886ddat=dmiss
9887
9888if (associated(this%volanab)) then
9889 do i=1,size(this%anavar%b)
9890 bdat=this%volanab(1,i,1)
9892 end do
9893end if
9894bdat=ibmiss
9895
9896if (associated(this%volanac)) then
9897 do i=1,size(this%anavar%c)
9898 cdat=this%volanac(1,i,1)
9900 end do
9901end if
9902cdat=cmiss
9903ENDIF
9904
9905print*,"---- data vector ----"
9906print*,""
9907print*,"->>>>>>>>> dativar -"
9909print*,""
9910print*,"->>>>>>>>> datiattr -"
9912print*,""
9913print*,"->>>>>>>>> dativarattr -"
9915
9916print*,"-- data data section (first point) --"
9917
9918idat=imiss
9919rdat=rmiss
9920ddat=dmiss
9921bdat=ibmiss
9922cdat=cmiss
9923
9924IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0 .AND. size(this%time) > 0 &
9925 .AND. size(this%level) > 0 .AND. size(this%timerange) > 0) THEN
9926if (associated(this%voldatii)) then
9927 do i=1,size(this%dativar%i)
9928 idat=this%voldatii(1,1,1,1,i,1)
9930 end do
9931end if
9932idat=imiss
9933
9934if (associated(this%voldatir)) then
9935 do i=1,size(this%dativar%r)
9936 rdat=this%voldatir(1,1,1,1,i,1)
9938 end do
9939end if
9940rdat=rmiss
9941
9942if (associated(this%voldatid)) then
9943 do i=1,size(this%dativar%d)
9944 ddat=this%voldatid(1,1,1,1,i,1)
9946 end do
9947end if
9948ddat=dmiss
9949
9950if (associated(this%voldatib)) then
9951 do i=1,size(this%dativar%b)
9952 bdat=this%voldatib(1,1,1,1,i,1)
9954 end do
9955end if
9956bdat=ibmiss
9957
9958if (associated(this%voldatic)) then
9959 do i=1,size(this%dativar%c)
9960 cdat=this%voldatic(1,1,1,1,i,1)
9962 end do
9963end if
9964cdat=cmiss
9965ENDIF
9966
9967print*,"<<<<<<<<<<<<<<<<<<< END vol7d object >>>>>>>>>>>>>>>>>>>>"
9968
9969END SUBROUTINE vol7d_display
9970
9971
9973SUBROUTINE dat_display(this,idat,rdat,ddat,bdat,cdat)
9974TYPE(vol7d_var),intent(in) :: this
9976REAL :: rdat
9978DOUBLE PRECISION :: ddat
9980INTEGER :: idat
9982INTEGER(kind=int_b) :: bdat
9984CHARACTER(len=*) :: cdat
9985
9986print *, to_char_dat(this,idat,rdat,ddat,bdat,cdat)
9987
9988end SUBROUTINE dat_display
9989
9991SUBROUTINE dat_vect_display(this,idat,rdat,ddat,bdat,cdat)
9992
9993TYPE(vol7d_var),intent(in) :: this(:)
9995REAL :: rdat(:)
9997DOUBLE PRECISION :: ddat(:)
9999INTEGER :: idat(:)
10001INTEGER(kind=int_b) :: bdat(:)
10003CHARACTER(len=*):: cdat(:)
10004
10005integer :: i
10006
10007do i =1,size(this)
10009end do
10010
10011end SUBROUTINE dat_vect_display
10012
10013
10014FUNCTION to_char_dat(this,idat,rdat,ddat,bdat,cdat)
10015#ifdef HAVE_DBALLE
10016USE dballef
10017#endif
10018TYPE(vol7d_var),INTENT(in) :: this
10020REAL :: rdat
10022DOUBLE PRECISION :: ddat
10024INTEGER :: idat
10026INTEGER(kind=int_b) :: bdat
10028CHARACTER(len=*) :: cdat
10029CHARACTER(len=80) :: to_char_dat
10030
10031CHARACTER(len=LEN(to_char_dat)) :: to_char_tmp
10032
10033
10034#ifdef HAVE_DBALLE
10035INTEGER :: handle, ier
10036
10037handle = 0
10038to_char_dat="VALUE: "
10039
10044
10046 ier = idba_messaggi(handle,"/dev/null", "w", "BUFR")
10047 ier = idba_spiegab(handle,this%btable,cdat,to_char_tmp)
10048 ier = idba_fatto(handle)
10049 to_char_dat=trim(to_char_dat)//" ;char> "//trim(to_char_tmp)
10050endif
10051
10052#else
10053
10054to_char_dat="VALUE: "
10060
10061#endif
10062
10063END FUNCTION to_char_dat
10064
10065
10068FUNCTION vol7d_c_e(this) RESULT(c_e)
10069TYPE(vol7d), INTENT(in) :: this
10070
10071LOGICAL :: c_e
10072
10074 ASSOCIATED(this%level) .OR. ASSOCIATED(this%timerange) .OR. &
10075 ASSOCIATED(this%network) .OR. &
10076 ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
10077 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
10078 ASSOCIATED(this%anavar%c) .OR. &
10079 ASSOCIATED(this%anaattr%r) .OR. ASSOCIATED(this%anaattr%d) .OR. &
10080 ASSOCIATED(this%anaattr%i) .OR. ASSOCIATED(this%anaattr%b) .OR. &
10081 ASSOCIATED(this%anaattr%c) .OR. &
10082 ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
10083 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
10084 ASSOCIATED(this%dativar%c) .OR. &
10085 ASSOCIATED(this%datiattr%r) .OR. ASSOCIATED(this%datiattr%d) .OR. &
10086 ASSOCIATED(this%datiattr%i) .OR. ASSOCIATED(this%datiattr%b) .OR. &
10087 ASSOCIATED(this%datiattr%c)
10088
10089END FUNCTION vol7d_c_e
10090
10091
10130SUBROUTINE vol7d_alloc(this, nana, ntime, nlevel, ntimerange, nnetwork, &
10131 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
10132 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
10133 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
10134 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
10135 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
10136 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc, &
10137 ini)
10138TYPE(vol7d),INTENT(inout) :: this
10139INTEGER,INTENT(in),OPTIONAL :: nana
10140INTEGER,INTENT(in),OPTIONAL :: ntime
10141INTEGER,INTENT(in),OPTIONAL :: nlevel
10142INTEGER,INTENT(in),OPTIONAL :: ntimerange
10143INTEGER,INTENT(in),OPTIONAL :: nnetwork
10145INTEGER,INTENT(in),OPTIONAL :: &
10146 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
10147 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
10148 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
10149 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
10150 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
10151 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc
10152LOGICAL,INTENT(in),OPTIONAL :: ini
10153
10154INTEGER :: i
10155LOGICAL :: linit
10156
10157IF (PRESENT(ini)) THEN
10158 linit = ini
10159ELSE
10160 linit = .false.
10161ENDIF
10162
10163! Dimensioni principali
10164IF (PRESENT(nana)) THEN
10165 IF (nana >= 0) THEN
10166 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
10167 ALLOCATE(this%ana(nana))
10168 IF (linit) THEN
10169 DO i = 1, nana
10171 ENDDO
10172 ENDIF
10173 ENDIF
10174ENDIF
10175IF (PRESENT(ntime)) THEN
10176 IF (ntime >= 0) THEN
10177 IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
10178 ALLOCATE(this%time(ntime))
10179 IF (linit) THEN
10180 DO i = 1, ntime
10182 ENDDO
10183 ENDIF
10184 ENDIF
10185ENDIF
10186IF (PRESENT(nlevel)) THEN
10187 IF (nlevel >= 0) THEN
10188 IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
10189 ALLOCATE(this%level(nlevel))
10190 IF (linit) THEN
10191 DO i = 1, nlevel
10193 ENDDO
10194 ENDIF
10195 ENDIF
10196ENDIF
10197IF (PRESENT(ntimerange)) THEN
10198 IF (ntimerange >= 0) THEN
10199 IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
10200 ALLOCATE(this%timerange(ntimerange))
10201 IF (linit) THEN
10202 DO i = 1, ntimerange
10204 ENDDO
10205 ENDIF
10206 ENDIF
10207ENDIF
10208IF (PRESENT(nnetwork)) THEN
10209 IF (nnetwork >= 0) THEN
10210 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
10211 ALLOCATE(this%network(nnetwork))
10212 IF (linit) THEN
10213 DO i = 1, nnetwork
10215 ENDDO
10216 ENDIF
10217 ENDIF
10218ENDIF
10219! Dimensioni dei tipi delle variabili
10220CALL vol7d_varvect_alloc(this%anavar, nanavarr, nanavard, &
10221 nanavari, nanavarb, nanavarc, ini)
10222CALL vol7d_varvect_alloc(this%anaattr, nanaattrr, nanaattrd, &
10223 nanaattri, nanaattrb, nanaattrc, ini)
10224CALL vol7d_varvect_alloc(this%anavarattr, nanavarattrr, nanavarattrd, &
10225 nanavarattri, nanavarattrb, nanavarattrc, ini)
10226CALL vol7d_varvect_alloc(this%dativar, ndativarr, ndativard, &
10227 ndativari, ndativarb, ndativarc, ini)
10228CALL vol7d_varvect_alloc(this%datiattr, ndatiattrr, ndatiattrd, &
10229 ndatiattri, ndatiattrb, ndatiattrc, ini)
10230CALL vol7d_varvect_alloc(this%dativarattr, ndativarattrr, ndativarattrd, &
10231 ndativarattri, ndativarattrb, ndativarattrc, ini)
10232
10233END SUBROUTINE vol7d_alloc
10234
10235
10236FUNCTION vol7d_check_alloc_ana(this)
10237TYPE(vol7d),INTENT(in) :: this
10238LOGICAL :: vol7d_check_alloc_ana
10239
10240vol7d_check_alloc_ana = ASSOCIATED(this%ana) .AND. ASSOCIATED(this%network)
10241
10242END FUNCTION vol7d_check_alloc_ana
10243
10244SUBROUTINE vol7d_force_alloc_ana(this, ini)
10245TYPE(vol7d),INTENT(inout) :: this
10246LOGICAL,INTENT(in),OPTIONAL :: ini
10247
10248! Alloco i descrittori minimi per avere un volume di anagrafica
10249IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=1, ini=ini)
10250IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=1, ini=ini)
10251
10252END SUBROUTINE vol7d_force_alloc_ana
10253
10254
10255FUNCTION vol7d_check_alloc_dati(this)
10256TYPE(vol7d),INTENT(in) :: this
10257LOGICAL :: vol7d_check_alloc_dati
10258
10259vol7d_check_alloc_dati = vol7d_check_alloc_ana(this) .AND. &
10260 ASSOCIATED(this%time) .AND. ASSOCIATED(this%level) .AND. &
10261 ASSOCIATED(this%timerange)
10262
10263END FUNCTION vol7d_check_alloc_dati
10264
10265SUBROUTINE vol7d_force_alloc_dati(this, ini)
10266TYPE(vol7d),INTENT(inout) :: this
10267LOGICAL,INTENT(in),OPTIONAL :: ini
10268
10269! Alloco i descrittori minimi per avere un volume di dati
10270CALL vol7d_force_alloc_ana(this, ini)
10271IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=1, ini=ini)
10272IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=1, ini=ini)
10273IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=1, ini=ini)
10274
10275END SUBROUTINE vol7d_force_alloc_dati
10276
10277
10278SUBROUTINE vol7d_force_alloc(this)
10279TYPE(vol7d),INTENT(inout) :: this
10280
10281! If anything really not allocated yet, allocate with size 0
10282IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=0)
10283IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=0)
10284IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=0)
10285IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=0)
10286IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=0)
10287
10288END SUBROUTINE vol7d_force_alloc
10289
10290
10291FUNCTION vol7d_check_vol(this)
10292TYPE(vol7d),INTENT(in) :: this
10293LOGICAL :: vol7d_check_vol
10294
10295vol7d_check_vol = c_e(this)
10296
10297! Anagrafica
10298IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
10299 vol7d_check_vol = .false.
10300ENDIF
10301
10302IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
10303 vol7d_check_vol = .false.
10304ENDIF
10305
10306IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
10307 vol7d_check_vol = .false.
10308ENDIF
10309
10310IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
10311 vol7d_check_vol = .false.
10312ENDIF
10313
10314IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
10315 vol7d_check_vol = .false.
10316ENDIF
10317IF (ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
10318 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
10319 ASSOCIATED(this%anavar%c)) THEN
10320 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_ana(this)
10321ENDIF
10322
10323! Attributi dell'anagrafica
10324IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
10325 .NOT.ASSOCIATED(this%volanaattrr)) THEN
10326 vol7d_check_vol = .false.
10327ENDIF
10328
10329IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
10330 .NOT.ASSOCIATED(this%volanaattrd)) THEN
10331 vol7d_check_vol = .false.
10332ENDIF
10333
10334IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
10335 .NOT.ASSOCIATED(this%volanaattri)) THEN
10336 vol7d_check_vol = .false.
10337ENDIF
10338
10339IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
10340 .NOT.ASSOCIATED(this%volanaattrb)) THEN
10341 vol7d_check_vol = .false.
10342ENDIF
10343
10344IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
10345 .NOT.ASSOCIATED(this%volanaattrc)) THEN
10346 vol7d_check_vol = .false.
10347ENDIF
10348
10349! Dati
10350IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
10351 vol7d_check_vol = .false.
10352ENDIF
10353
10354IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
10355 vol7d_check_vol = .false.
10356ENDIF
10357
10358IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
10359 vol7d_check_vol = .false.
10360ENDIF
10361
10362IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
10363 vol7d_check_vol = .false.
10364ENDIF
10365
10366IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
10367 vol7d_check_vol = .false.
10368ENDIF
10369
10370! Attributi dei dati
10371IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
10372 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
10373 vol7d_check_vol = .false.
10374ENDIF
10375
10376IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
10377 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
10378 vol7d_check_vol = .false.
10379ENDIF
10380
10381IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
10382 .NOT.ASSOCIATED(this%voldatiattri)) THEN
10383 vol7d_check_vol = .false.
10384ENDIF
10385
10386IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
10387 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
10388 vol7d_check_vol = .false.
10389ENDIF
10390
10391IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
10392 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
10393 vol7d_check_vol = .false.
10394ENDIF
10395IF (ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
10396 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
10397 ASSOCIATED(this%dativar%c)) THEN
10398 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_dati(this)
10399ENDIF
10400
10401END FUNCTION vol7d_check_vol
10402
10403
10418SUBROUTINE vol7d_alloc_vol(this, ini, inivol)
10419TYPE(vol7d),INTENT(inout) :: this
10420LOGICAL,INTENT(in),OPTIONAL :: ini
10421LOGICAL,INTENT(in),OPTIONAL :: inivol
10422
10423LOGICAL :: linivol
10424
10425IF (PRESENT(inivol)) THEN
10426 linivol = inivol
10427ELSE
10428 linivol = .true.
10429ENDIF
10430
10431! Anagrafica
10432IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
10433 CALL vol7d_force_alloc_ana(this, ini)
10434 ALLOCATE(this%volanar(SIZE(this%ana), SIZE(this%anavar%r), SIZE(this%network)))
10435 IF (linivol) this%volanar(:,:,:) = rmiss
10436ENDIF
10437
10438IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
10439 CALL vol7d_force_alloc_ana(this, ini)
10440 ALLOCATE(this%volanad(SIZE(this%ana), SIZE(this%anavar%d), SIZE(this%network)))
10441 IF (linivol) this%volanad(:,:,:) = rdmiss
10442ENDIF
10443
10444IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
10445 CALL vol7d_force_alloc_ana(this, ini)
10446 ALLOCATE(this%volanai(SIZE(this%ana), SIZE(this%anavar%i), SIZE(this%network)))
10447 IF (linivol) this%volanai(:,:,:) = imiss
10448ENDIF
10449
10450IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
10451 CALL vol7d_force_alloc_ana(this, ini)
10452 ALLOCATE(this%volanab(SIZE(this%ana), SIZE(this%anavar%b), SIZE(this%network)))
10453 IF (linivol) this%volanab(:,:,:) = ibmiss
10454ENDIF
10455
10456IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
10457 CALL vol7d_force_alloc_ana(this, ini)
10458 ALLOCATE(this%volanac(SIZE(this%ana), SIZE(this%anavar%c), SIZE(this%network)))
10459 IF (linivol) this%volanac(:,:,:) = cmiss
10460ENDIF
10461
10462! Attributi dell'anagrafica
10463IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
10464 .NOT.ASSOCIATED(this%volanaattrr)) THEN
10465 CALL vol7d_force_alloc_ana(this, ini)
10466 ALLOCATE(this%volanaattrr(SIZE(this%ana), SIZE(this%anavarattr%r), &
10467 SIZE(this%network), SIZE(this%anaattr%r)))
10468 IF (linivol) this%volanaattrr(:,:,:,:) = rmiss
10469ENDIF
10470
10471IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
10472 .NOT.ASSOCIATED(this%volanaattrd)) THEN
10473 CALL vol7d_force_alloc_ana(this, ini)
10474 ALLOCATE(this%volanaattrd(SIZE(this%ana), SIZE(this%anavarattr%d), &
10475 SIZE(this%network), SIZE(this%anaattr%d)))
10476 IF (linivol) this%volanaattrd(:,:,:,:) = rdmiss
10477ENDIF
10478
10479IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
10480 .NOT.ASSOCIATED(this%volanaattri)) THEN
10481 CALL vol7d_force_alloc_ana(this, ini)
10482 ALLOCATE(this%volanaattri(SIZE(this%ana), SIZE(this%anavarattr%i), &
10483 SIZE(this%network), SIZE(this%anaattr%i)))
10484 IF (linivol) this%volanaattri(:,:,:,:) = imiss
10485ENDIF
10486
10487IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
10488 .NOT.ASSOCIATED(this%volanaattrb)) THEN
10489 CALL vol7d_force_alloc_ana(this, ini)
10490 ALLOCATE(this%volanaattrb(SIZE(this%ana), SIZE(this%anavarattr%b), &
10491 SIZE(this%network), SIZE(this%anaattr%b)))
10492 IF (linivol) this%volanaattrb(:,:,:,:) = ibmiss
10493ENDIF
10494
10495IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
10496 .NOT.ASSOCIATED(this%volanaattrc)) THEN
10497 CALL vol7d_force_alloc_ana(this, ini)
10498 ALLOCATE(this%volanaattrc(SIZE(this%ana), SIZE(this%anavarattr%c), &
10499 SIZE(this%network), SIZE(this%anaattr%c)))
10500 IF (linivol) this%volanaattrc(:,:,:,:) = cmiss
10501ENDIF
10502
10503! Dati
10504IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
10505 CALL vol7d_force_alloc_dati(this, ini)
10506 ALLOCATE(this%voldatir(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10507 SIZE(this%timerange), SIZE(this%dativar%r), SIZE(this%network)))
10508 IF (linivol) this%voldatir(:,:,:,:,:,:) = rmiss
10509ENDIF
10510
10511IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
10512 CALL vol7d_force_alloc_dati(this, ini)
10513 ALLOCATE(this%voldatid(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10514 SIZE(this%timerange), SIZE(this%dativar%d), SIZE(this%network)))
10515 IF (linivol) this%voldatid(:,:,:,:,:,:) = rdmiss
10516ENDIF
10517
10518IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
10519 CALL vol7d_force_alloc_dati(this, ini)
10520 ALLOCATE(this%voldatii(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10521 SIZE(this%timerange), SIZE(this%dativar%i), SIZE(this%network)))
10522 IF (linivol) this%voldatii(:,:,:,:,:,:) = imiss
10523ENDIF
10524
10525IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
10526 CALL vol7d_force_alloc_dati(this, ini)
10527 ALLOCATE(this%voldatib(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10528 SIZE(this%timerange), SIZE(this%dativar%b), SIZE(this%network)))
10529 IF (linivol) this%voldatib(:,:,:,:,:,:) = ibmiss
10530ENDIF
10531
10532IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
10533 CALL vol7d_force_alloc_dati(this, ini)
10534 ALLOCATE(this%voldatic(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10535 SIZE(this%timerange), SIZE(this%dativar%c), SIZE(this%network)))
10536 IF (linivol) this%voldatic(:,:,:,:,:,:) = cmiss
10537ENDIF
10538
10539! Attributi dei dati
10540IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
10541 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
10542 CALL vol7d_force_alloc_dati(this, ini)
10543 ALLOCATE(this%voldatiattrr(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10544 SIZE(this%timerange), SIZE(this%dativarattr%r), SIZE(this%network), &
10545 SIZE(this%datiattr%r)))
10546 IF (linivol) this%voldatiattrr(:,:,:,:,:,:,:) = rmiss
10547ENDIF
10548
10549IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
10550 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
10551 CALL vol7d_force_alloc_dati(this, ini)
10552 ALLOCATE(this%voldatiattrd(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10553 SIZE(this%timerange), SIZE(this%dativarattr%d), SIZE(this%network), &
10554 SIZE(this%datiattr%d)))
10555 IF (linivol) this%voldatiattrd(:,:,:,:,:,:,:) = rdmiss
10556ENDIF
10557
10558IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
10559 .NOT.ASSOCIATED(this%voldatiattri)) THEN
10560 CALL vol7d_force_alloc_dati(this, ini)
10561 ALLOCATE(this%voldatiattri(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10562 SIZE(this%timerange), SIZE(this%dativarattr%i), SIZE(this%network), &
10563 SIZE(this%datiattr%i)))
10564 IF (linivol) this%voldatiattri(:,:,:,:,:,:,:) = imiss
10565ENDIF
10566
10567IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
10568 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
10569 CALL vol7d_force_alloc_dati(this, ini)
10570 ALLOCATE(this%voldatiattrb(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10571 SIZE(this%timerange), SIZE(this%dativarattr%b), SIZE(this%network), &
10572 SIZE(this%datiattr%b)))
10573 IF (linivol) this%voldatiattrb(:,:,:,:,:,:,:) = ibmiss
10574ENDIF
10575
10576IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
10577 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
10578 CALL vol7d_force_alloc_dati(this, ini)
10579 ALLOCATE(this%voldatiattrc(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10580 SIZE(this%timerange), SIZE(this%dativarattr%c), SIZE(this%network), &
10581 SIZE(this%datiattr%c)))
10582 IF (linivol) this%voldatiattrc(:,:,:,:,:,:,:) = cmiss
10583ENDIF
10584
10585! Catch-all method
10586CALL vol7d_force_alloc(this)
10587
10588! Creo gli indici var-attr
10589
10590#ifdef DEBUG
10591CALL l4f_log(l4f_debug,"calling: vol7d_set_attr_ind")
10592#endif
10593
10594CALL vol7d_set_attr_ind(this)
10595
10596
10597
10598END SUBROUTINE vol7d_alloc_vol
10599
10600
10607SUBROUTINE vol7d_set_attr_ind(this)
10608TYPE(vol7d),INTENT(inout) :: this
10609
10610INTEGER :: i
10611
10612! real
10613IF (ASSOCIATED(this%dativar%r)) THEN
10614 IF (ASSOCIATED(this%dativarattr%r)) THEN
10615 DO i = 1, SIZE(this%dativar%r)
10616 this%dativar%r(i)%r = &
10617 firsttrue(this%dativar%r(i)%btable == this%dativarattr%r(:)%btable)
10618 ENDDO
10619 ENDIF
10620
10621 IF (ASSOCIATED(this%dativarattr%d)) THEN
10622 DO i = 1, SIZE(this%dativar%r)
10623 this%dativar%r(i)%d = &
10624 firsttrue(this%dativar%r(i)%btable == this%dativarattr%d(:)%btable)
10625 ENDDO
10626 ENDIF
10627
10628 IF (ASSOCIATED(this%dativarattr%i)) THEN
10629 DO i = 1, SIZE(this%dativar%r)
10630 this%dativar%r(i)%i = &
10631 firsttrue(this%dativar%r(i)%btable == this%dativarattr%i(:)%btable)
10632 ENDDO
10633 ENDIF
10634
10635 IF (ASSOCIATED(this%dativarattr%b)) THEN
10636 DO i = 1, SIZE(this%dativar%r)
10637 this%dativar%r(i)%b = &
10638 firsttrue(this%dativar%r(i)%btable == this%dativarattr%b(:)%btable)
10639 ENDDO
10640 ENDIF
10641
10642 IF (ASSOCIATED(this%dativarattr%c)) THEN
10643 DO i = 1, SIZE(this%dativar%r)
10644 this%dativar%r(i)%c = &
10645 firsttrue(this%dativar%r(i)%btable == this%dativarattr%c(:)%btable)
10646 ENDDO
10647 ENDIF
10648ENDIF
10649! double
10650IF (ASSOCIATED(this%dativar%d)) THEN
10651 IF (ASSOCIATED(this%dativarattr%r)) THEN
10652 DO i = 1, SIZE(this%dativar%d)
10653 this%dativar%d(i)%r = &
10654 firsttrue(this%dativar%d(i)%btable == this%dativarattr%r(:)%btable)
10655 ENDDO
10656 ENDIF
10657
10658 IF (ASSOCIATED(this%dativarattr%d)) THEN
10659 DO i = 1, SIZE(this%dativar%d)
10660 this%dativar%d(i)%d = &
10661 firsttrue(this%dativar%d(i)%btable == this%dativarattr%d(:)%btable)
10662 ENDDO
10663 ENDIF
10664
10665 IF (ASSOCIATED(this%dativarattr%i)) THEN
10666 DO i = 1, SIZE(this%dativar%d)
10667 this%dativar%d(i)%i = &
10668 firsttrue(this%dativar%d(i)%btable == this%dativarattr%i(:)%btable)
10669 ENDDO
10670 ENDIF
10671
10672 IF (ASSOCIATED(this%dativarattr%b)) THEN
10673 DO i = 1, SIZE(this%dativar%d)
10674 this%dativar%d(i)%b = &
10675 firsttrue(this%dativar%d(i)%btable == this%dativarattr%b(:)%btable)
10676 ENDDO
10677 ENDIF
10678
10679 IF (ASSOCIATED(this%dativarattr%c)) THEN
10680 DO i = 1, SIZE(this%dativar%d)
10681 this%dativar%d(i)%c = &
10682 firsttrue(this%dativar%d(i)%btable == this%dativarattr%c(:)%btable)
10683 ENDDO
10684 ENDIF
10685ENDIF
10686! integer
10687IF (ASSOCIATED(this%dativar%i)) THEN
10688 IF (ASSOCIATED(this%dativarattr%r)) THEN
10689 DO i = 1, SIZE(this%dativar%i)
10690 this%dativar%i(i)%r = &
10691 firsttrue(this%dativar%i(i)%btable == this%dativarattr%r(:)%btable)
10692 ENDDO
10693 ENDIF
10694
10695 IF (ASSOCIATED(this%dativarattr%d)) THEN
10696 DO i = 1, SIZE(this%dativar%i)
10697 this%dativar%i(i)%d = &
10698 firsttrue(this%dativar%i(i)%btable == this%dativarattr%d(:)%btable)
10699 ENDDO
10700 ENDIF
10701
10702 IF (ASSOCIATED(this%dativarattr%i)) THEN
10703 DO i = 1, SIZE(this%dativar%i)
10704 this%dativar%i(i)%i = &
10705 firsttrue(this%dativar%i(i)%btable == this%dativarattr%i(:)%btable)
10706 ENDDO
10707 ENDIF
10708
10709 IF (ASSOCIATED(this%dativarattr%b)) THEN
10710 DO i = 1, SIZE(this%dativar%i)
10711 this%dativar%i(i)%b = &
10712 firsttrue(this%dativar%i(i)%btable == this%dativarattr%b(:)%btable)
10713 ENDDO
10714 ENDIF
10715
10716 IF (ASSOCIATED(this%dativarattr%c)) THEN
10717 DO i = 1, SIZE(this%dativar%i)
10718 this%dativar%i(i)%c = &
10719 firsttrue(this%dativar%i(i)%btable == this%dativarattr%c(:)%btable)
10720 ENDDO
10721 ENDIF
10722ENDIF
10723! byte
10724IF (ASSOCIATED(this%dativar%b)) THEN
10725 IF (ASSOCIATED(this%dativarattr%r)) THEN
10726 DO i = 1, SIZE(this%dativar%b)
10727 this%dativar%b(i)%r = &
10728 firsttrue(this%dativar%b(i)%btable == this%dativarattr%r(:)%btable)
10729 ENDDO
10730 ENDIF
10731
10732 IF (ASSOCIATED(this%dativarattr%d)) THEN
10733 DO i = 1, SIZE(this%dativar%b)
10734 this%dativar%b(i)%d = &
10735 firsttrue(this%dativar%b(i)%btable == this%dativarattr%d(:)%btable)
10736 ENDDO
10737 ENDIF
10738
10739 IF (ASSOCIATED(this%dativarattr%i)) THEN
10740 DO i = 1, SIZE(this%dativar%b)
10741 this%dativar%b(i)%i = &
10742 firsttrue(this%dativar%b(i)%btable == this%dativarattr%i(:)%btable)
10743 ENDDO
10744 ENDIF
10745
10746 IF (ASSOCIATED(this%dativarattr%b)) THEN
10747 DO i = 1, SIZE(this%dativar%b)
10748 this%dativar%b(i)%b = &
10749 firsttrue(this%dativar%b(i)%btable == this%dativarattr%b(:)%btable)
10750 ENDDO
10751 ENDIF
10752
10753 IF (ASSOCIATED(this%dativarattr%c)) THEN
10754 DO i = 1, SIZE(this%dativar%b)
10755 this%dativar%b(i)%c = &
10756 firsttrue(this%dativar%b(i)%btable == this%dativarattr%c(:)%btable)
10757 ENDDO
10758 ENDIF
10759ENDIF
10760! character
10761IF (ASSOCIATED(this%dativar%c)) THEN
10762 IF (ASSOCIATED(this%dativarattr%r)) THEN
10763 DO i = 1, SIZE(this%dativar%c)
10764 this%dativar%c(i)%r = &
10765 firsttrue(this%dativar%c(i)%btable == this%dativarattr%r(:)%btable)
10766 ENDDO
10767 ENDIF
10768
10769 IF (ASSOCIATED(this%dativarattr%d)) THEN
10770 DO i = 1, SIZE(this%dativar%c)
10771 this%dativar%c(i)%d = &
10772 firsttrue(this%dativar%c(i)%btable == this%dativarattr%d(:)%btable)
10773 ENDDO
10774 ENDIF
10775
10776 IF (ASSOCIATED(this%dativarattr%i)) THEN
10777 DO i = 1, SIZE(this%dativar%c)
10778 this%dativar%c(i)%i = &
10779 firsttrue(this%dativar%c(i)%btable == this%dativarattr%i(:)%btable)
10780 ENDDO
10781 ENDIF
10782
10783 IF (ASSOCIATED(this%dativarattr%b)) THEN
10784 DO i = 1, SIZE(this%dativar%c)
10785 this%dativar%c(i)%b = &
10786 firsttrue(this%dativar%c(i)%btable == this%dativarattr%b(:)%btable)
10787 ENDDO
10788 ENDIF
10789
10790 IF (ASSOCIATED(this%dativarattr%c)) THEN
10791 DO i = 1, SIZE(this%dativar%c)
10792 this%dativar%c(i)%c = &
10793 firsttrue(this%dativar%c(i)%btable == this%dativarattr%c(:)%btable)
10794 ENDDO
10795 ENDIF
10796ENDIF
10797
10798END SUBROUTINE vol7d_set_attr_ind
10799
10800
10805SUBROUTINE vol7d_merge(this, that, sort, bestdata, &
10806 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
10807TYPE(vol7d),INTENT(INOUT) :: this
10808TYPE(vol7d),INTENT(INOUT) :: that
10809LOGICAL,INTENT(IN),OPTIONAL :: sort
10810LOGICAL,INTENT(in),OPTIONAL :: bestdata
10811LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple ! experimental, please do not use outside the library now
10812
10813TYPE(vol7d) :: v7d_clean
10814
10815
10817 this = that
10819 that = v7d_clean ! destroy that without deallocating
10820ELSE ! Append that to this and destroy that
10822 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
10824ENDIF
10825
10826END SUBROUTINE vol7d_merge
10827
10828
10857SUBROUTINE vol7d_append(this, that, sort, bestdata, &
10858 ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple)
10859TYPE(vol7d),INTENT(INOUT) :: this
10860TYPE(vol7d),INTENT(IN) :: that
10861LOGICAL,INTENT(IN),OPTIONAL :: sort
10862! experimental, please do not use outside the library now, they force the use
10863! of a simplified mapping algorithm which is valid only whene the dimension
10864! content is the same in both volumes , or when one of them is empty
10865LOGICAL,INTENT(in),OPTIONAL :: bestdata
10866LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple
10867
10868
10869TYPE(vol7d) :: v7dtmp
10870LOGICAL :: lsort, lbestdata
10871INTEGER,POINTER :: remapt1(:), remapt2(:), remaptr1(:), remaptr2(:), &
10872 remapl1(:), remapl2(:), remapa1(:), remapa2(:), remapn1(:), remapn2(:)
10873
10875IF (.NOT.vol7d_check_vol(that)) RETURN ! be safe
10878 RETURN
10879ENDIF
10880
10881IF (this%time_definition /= that%time_definition) THEN
10882 CALL l4f_log(l4f_fatal, &
10883 'in vol7d_append, cannot append volumes with different &
10884 &time definition')
10885 CALL raise_fatal_error()
10886ENDIF
10887
10888! Completo l'allocazione per avere volumi a norma
10889CALL vol7d_alloc_vol(this)
10890
10894
10895! Calcolo le mappature tra volumi vecchi e volume nuovo
10896! I puntatori remap* vengono tutti o allocati o nullificati
10897IF (optio_log(ltimesimple)) THEN
10898 CALL vol7d_remap2simple_datetime(this%time, that%time, v7dtmp%time, &
10899 lsort, remapt1, remapt2)
10900ELSE
10901 CALL vol7d_remap2_datetime(this%time, that%time, v7dtmp%time, &
10902 lsort, remapt1, remapt2)
10903ENDIF
10904IF (optio_log(ltimerangesimple)) THEN
10905 CALL vol7d_remap2simple_vol7d_timerange(this%timerange, that%timerange, &
10906 v7dtmp%timerange, lsort, remaptr1, remaptr2)
10907ELSE
10908 CALL vol7d_remap2_vol7d_timerange(this%timerange, that%timerange, &
10909 v7dtmp%timerange, lsort, remaptr1, remaptr2)
10910ENDIF
10911IF (optio_log(llevelsimple)) THEN
10912 CALL vol7d_remap2simple_vol7d_level(this%level, that%level, v7dtmp%level, &
10913 lsort, remapl1, remapl2)
10914ELSE
10915 CALL vol7d_remap2_vol7d_level(this%level, that%level, v7dtmp%level, &
10916 lsort, remapl1, remapl2)
10917ENDIF
10918IF (optio_log(lanasimple)) THEN
10919 CALL vol7d_remap2simple_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
10920 .false., remapa1, remapa2)
10921ELSE
10922 CALL vol7d_remap2_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
10923 .false., remapa1, remapa2)
10924ENDIF
10925IF (optio_log(lnetworksimple)) THEN
10926 CALL vol7d_remap2simple_vol7d_network(this%network, that%network, v7dtmp%network, &
10927 .false., remapn1, remapn2)
10928ELSE
10929 CALL vol7d_remap2_vol7d_network(this%network, that%network, v7dtmp%network, &
10930 .false., remapn1, remapn2)
10931ENDIF
10932
10933! Faccio la fusione fisica dei volumi
10934CALL vol7d_merge_finalr(this, that, v7dtmp, &
10935 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10936 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10937CALL vol7d_merge_finald(this, that, v7dtmp, &
10938 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10939 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10940CALL vol7d_merge_finali(this, that, v7dtmp, &
10941 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10942 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10943CALL vol7d_merge_finalb(this, that, v7dtmp, &
10944 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10945 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10946CALL vol7d_merge_finalc(this, that, v7dtmp, &
10947 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10948 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10949
10950! Dealloco i vettori di rimappatura
10951IF (ASSOCIATED(remapt1)) DEALLOCATE(remapt1)
10952IF (ASSOCIATED(remapt2)) DEALLOCATE(remapt2)
10953IF (ASSOCIATED(remaptr1)) DEALLOCATE(remaptr1)
10954IF (ASSOCIATED(remaptr2)) DEALLOCATE(remaptr2)
10955IF (ASSOCIATED(remapl1)) DEALLOCATE(remapl1)
10956IF (ASSOCIATED(remapl2)) DEALLOCATE(remapl2)
10957IF (ASSOCIATED(remapa1)) DEALLOCATE(remapa1)
10958IF (ASSOCIATED(remapa2)) DEALLOCATE(remapa2)
10959IF (ASSOCIATED(remapn1)) DEALLOCATE(remapn1)
10960IF (ASSOCIATED(remapn2)) DEALLOCATE(remapn2)
10961
10962! Distruggo il vecchio volume e assegno il nuovo a this
10964this = v7dtmp
10965! Ricreo gli indici var-attr
10966CALL vol7d_set_attr_ind(this)
10967
10968END SUBROUTINE vol7d_append
10969
10970
11003SUBROUTINE vol7d_copy(this, that, sort, unique, miss, &
11004 lsort_time, lsort_timerange, lsort_level, &
11005 ltime, ltimerange, llevel, lana, lnetwork, &
11006 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
11007 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
11008 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
11009 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
11010 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
11011 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
11012TYPE(vol7d),INTENT(IN) :: this
11013TYPE(vol7d),INTENT(INOUT) :: that
11014LOGICAL,INTENT(IN),OPTIONAL :: sort
11015LOGICAL,INTENT(IN),OPTIONAL :: unique
11016LOGICAL,INTENT(IN),OPTIONAL :: miss
11017LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
11018LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
11019LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
11027LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
11029LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
11031LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
11033LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
11035LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
11037LOGICAL,INTENT(in),OPTIONAL :: &
11038 lanavarr(:), lanavard(:), lanavari(:), lanavarb(:), lanavarc(:), &
11039 lanaattrr(:), lanaattrd(:), lanaattri(:), lanaattrb(:), lanaattrc(:), &
11040 lanavarattrr(:), lanavarattrd(:), lanavarattri(:), lanavarattrb(:), lanavarattrc(:), &
11041 ldativarr(:), ldativard(:), ldativari(:), ldativarb(:), ldativarc(:), &
11042 ldatiattrr(:), ldatiattrd(:), ldatiattri(:), ldatiattrb(:), ldatiattrc(:), &
11043 ldativarattrr(:), ldativarattrd(:), ldativarattri(:), ldativarattrb(:), ldativarattrc(:)
11044
11045LOGICAL :: lsort, lunique, lmiss
11046INTEGER,POINTER :: remapt(:), remaptr(:), remapl(:), remapa(:), remapn(:)
11047
11050IF (.NOT.vol7d_check_vol(this)) RETURN ! be safe
11051
11055
11056! Calcolo le mappature tra volume vecchio e volume nuovo
11057! I puntatori remap* vengono tutti o allocati o nullificati
11058CALL vol7d_remap1_datetime(this%time, that%time, &
11059 lsort.OR.optio_log(lsort_time), lunique, lmiss, remapt, ltime)
11060CALL vol7d_remap1_vol7d_timerange(this%timerange, that%timerange, &
11061 lsort.OR.optio_log(lsort_timerange), lunique, lmiss, remaptr, ltimerange)
11062CALL vol7d_remap1_vol7d_level(this%level, that%level, &
11063 lsort.OR.optio_log(lsort_level), lunique, lmiss, remapl, llevel)
11064CALL vol7d_remap1_vol7d_ana(this%ana, that%ana, &
11065 lsort, lunique, lmiss, remapa, lana)
11066CALL vol7d_remap1_vol7d_network(this%network, that%network, &
11067 lsort, lunique, lmiss, remapn, lnetwork)
11068
11069! lanavari, lanavarb, lanavarc, &
11070! lanaattri, lanaattrb, lanaattrc, &
11071! lanavarattri, lanavarattrb, lanavarattrc, &
11072! ldativari, ldativarb, ldativarc, &
11073! ldatiattri, ldatiattrb, ldatiattrc, &
11074! ldativarattri, ldativarattrb, ldativarattrc
11075! Faccio la riforma fisica dei volumi
11076CALL vol7d_reform_finalr(this, that, &
11077 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
11078 lanavarr, lanaattrr, lanavarattrr, ldativarr, ldatiattrr, ldativarattrr)
11079CALL vol7d_reform_finald(this, that, &
11080 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
11081 lanavard, lanaattrd, lanavarattrd, ldativard, ldatiattrd, ldativarattrd)
11082CALL vol7d_reform_finali(this, that, &
11083 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
11084 lanavari, lanaattri, lanavarattri, ldativari, ldatiattri, ldativarattri)
11085CALL vol7d_reform_finalb(this, that, &
11086 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
11087 lanavarb, lanaattrb, lanavarattrb, ldativarb, ldatiattrb, ldativarattrb)
11088CALL vol7d_reform_finalc(this, that, &
11089 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
11090 lanavarc, lanaattrc, lanavarattrc, ldativarc, ldatiattrc, ldativarattrc)
11091
11092! Dealloco i vettori di rimappatura
11093IF (ASSOCIATED(remapt)) DEALLOCATE(remapt)
11094IF (ASSOCIATED(remaptr)) DEALLOCATE(remaptr)
11095IF (ASSOCIATED(remapl)) DEALLOCATE(remapl)
11096IF (ASSOCIATED(remapa)) DEALLOCATE(remapa)
11097IF (ASSOCIATED(remapn)) DEALLOCATE(remapn)
11098
11099! Ricreo gli indici var-attr
11100CALL vol7d_set_attr_ind(that)
11101that%time_definition = this%time_definition
11102
11103END SUBROUTINE vol7d_copy
11104
11105
11116SUBROUTINE vol7d_reform(this, sort, unique, miss, &
11117 lsort_time, lsort_timerange, lsort_level, &
11118 ltime, ltimerange, llevel, lana, lnetwork, &
11119 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
11120 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
11121 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
11122 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
11123 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
11124 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc&
11125 ,purgeana)
11126TYPE(vol7d),INTENT(INOUT) :: this
11127LOGICAL,INTENT(IN),OPTIONAL :: sort
11128LOGICAL,INTENT(IN),OPTIONAL :: unique
11129LOGICAL,INTENT(IN),OPTIONAL :: miss
11130LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
11131LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
11132LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
11140LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
11141LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
11142LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
11143LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
11144LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
11146LOGICAL,INTENT(in),OPTIONAL :: &
11147 lanavarr(:), lanavard(:), lanavari(:), lanavarb(:), lanavarc(:), &
11148 lanaattrr(:), lanaattrd(:), lanaattri(:), lanaattrb(:), lanaattrc(:), &
11149 lanavarattrr(:), lanavarattrd(:), lanavarattri(:), lanavarattrb(:), lanavarattrc(:), &
11150 ldativarr(:), ldativard(:), ldativari(:), ldativarb(:), ldativarc(:), &
11151 ldatiattrr(:), ldatiattrd(:), ldatiattri(:), ldatiattrb(:), ldatiattrc(:), &
11152 ldativarattrr(:), ldativarattrd(:), ldativarattri(:), ldativarattrb(:), ldativarattrc(:)
11153LOGICAL,INTENT(IN),OPTIONAL :: purgeana
11154
11155TYPE(vol7d) :: v7dtmp
11156logical,allocatable :: llana(:)
11157integer :: i
11158
11160 lsort_time, lsort_timerange, lsort_level, &
11161 ltime, ltimerange, llevel, lana, lnetwork, &
11162 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
11163 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
11164 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
11165 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
11166 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
11167 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
11168
11169! destroy old volume
11171
11172if (optio_log(purgeana)) then
11173 allocate(llana(size(v7dtmp%ana)))
11174 llana =.false.
11175 do i =1,size(v7dtmp%ana)
11176 if (associated(v7dtmp%voldatii)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatii(i,:,:,:,:,:)))
11177 if (associated(v7dtmp%voldatir)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatir(i,:,:,:,:,:)))
11178 if (associated(v7dtmp%voldatid)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatid(i,:,:,:,:,:)))
11179 if (associated(v7dtmp%voldatib)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatib(i,:,:,:,:,:)))
11180 if (associated(v7dtmp%voldatic)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatic(i,:,:,:,:,:)))
11181 end do
11182 CALL vol7d_copy(v7dtmp, this,lana=llana)
11184 deallocate(llana)
11185else
11186 this=v7dtmp
11187end if
11188
11189END SUBROUTINE vol7d_reform
11190
11191
11199SUBROUTINE vol7d_smart_sort(this, lsort_time, lsort_timerange, lsort_level)
11200TYPE(vol7d),INTENT(INOUT) :: this
11201LOGICAL,OPTIONAL,INTENT(in) :: lsort_time
11202LOGICAL,OPTIONAL,INTENT(in) :: lsort_timerange
11203LOGICAL,OPTIONAL,INTENT(in) :: lsort_level
11204
11205INTEGER :: i
11206LOGICAL :: to_be_sorted
11207
11208to_be_sorted = .false.
11209CALL vol7d_alloc_vol(this) ! usual safety check
11210
11211IF (optio_log(lsort_time)) THEN
11212 DO i = 2, SIZE(this%time)
11213 IF (this%time(i) < this%time(i-1)) THEN
11214 to_be_sorted = .true.
11215 EXIT
11216 ENDIF
11217 ENDDO
11218ENDIF
11219IF (optio_log(lsort_timerange)) THEN
11220 DO i = 2, SIZE(this%timerange)
11221 IF (this%timerange(i) < this%timerange(i-1)) THEN
11222 to_be_sorted = .true.
11223 EXIT
11224 ENDIF
11225 ENDDO
11226ENDIF
11227IF (optio_log(lsort_level)) THEN
11228 DO i = 2, SIZE(this%level)
11229 IF (this%level(i) < this%level(i-1)) THEN
11230 to_be_sorted = .true.
11231 EXIT
11232 ENDIF
11233 ENDDO
11234ENDIF
11235
11236IF (to_be_sorted) CALL vol7d_reform(this, &
11237 lsort_time=lsort_time, lsort_timerange=lsort_timerange, lsort_level=lsort_level )
11238
11239END SUBROUTINE vol7d_smart_sort
11240
11248SUBROUTINE vol7d_filter(this, avl, vl, nl, s_d, e_d)
11249TYPE(vol7d),INTENT(inout) :: this
11250CHARACTER(len=*),INTENT(in),OPTIONAL :: avl(:)
11251CHARACTER(len=*),INTENT(in),OPTIONAL :: vl(:)
11252TYPE(vol7d_network),OPTIONAL :: nl(:)
11253TYPE(datetime),INTENT(in),OPTIONAL :: s_d
11254TYPE(datetime),INTENT(in),OPTIONAL :: e_d
11255
11256INTEGER :: i
11257
11258IF (PRESENT(avl)) THEN
11259 IF (SIZE(avl) > 0) THEN
11260
11261 IF (ASSOCIATED(this%anavar%r)) THEN
11262 DO i = 1, SIZE(this%anavar%r)
11263 IF (all(this%anavar%r(i)%btable /= avl)) this%anavar%r(i) = vol7d_var_miss
11264 ENDDO
11265 ENDIF
11266
11267 IF (ASSOCIATED(this%anavar%i)) THEN
11268 DO i = 1, SIZE(this%anavar%i)
11269 IF (all(this%anavar%i(i)%btable /= avl)) this%anavar%i(i) = vol7d_var_miss
11270 ENDDO
11271 ENDIF
11272
11273 IF (ASSOCIATED(this%anavar%b)) THEN
11274 DO i = 1, SIZE(this%anavar%b)
11275 IF (all(this%anavar%b(i)%btable /= avl)) this%anavar%b(i) = vol7d_var_miss
11276 ENDDO
11277 ENDIF
11278
11279 IF (ASSOCIATED(this%anavar%d)) THEN
11280 DO i = 1, SIZE(this%anavar%d)
11281 IF (all(this%anavar%d(i)%btable /= avl)) this%anavar%d(i) = vol7d_var_miss
11282 ENDDO
11283 ENDIF
11284
11285 IF (ASSOCIATED(this%anavar%c)) THEN
11286 DO i = 1, SIZE(this%anavar%c)
11287 IF (all(this%anavar%c(i)%btable /= avl)) this%anavar%c(i) = vol7d_var_miss
11288 ENDDO
11289 ENDIF
11290
11291 ENDIF
11292ENDIF
11293
11294
11295IF (PRESENT(vl)) THEN
11296 IF (size(vl) > 0) THEN
11297 IF (ASSOCIATED(this%dativar%r)) THEN
11298 DO i = 1, SIZE(this%dativar%r)
11299 IF (all(this%dativar%r(i)%btable /= vl)) this%dativar%r(i) = vol7d_var_miss
11300 ENDDO
11301 ENDIF
11302
11303 IF (ASSOCIATED(this%dativar%i)) THEN
11304 DO i = 1, SIZE(this%dativar%i)
11305 IF (all(this%dativar%i(i)%btable /= vl)) this%dativar%i(i) = vol7d_var_miss
11306 ENDDO
11307 ENDIF
11308
11309 IF (ASSOCIATED(this%dativar%b)) THEN
11310 DO i = 1, SIZE(this%dativar%b)
11311 IF (all(this%dativar%b(i)%btable /= vl)) this%dativar%b(i) = vol7d_var_miss
11312 ENDDO
11313 ENDIF
11314
11315 IF (ASSOCIATED(this%dativar%d)) THEN
11316 DO i = 1, SIZE(this%dativar%d)
11317 IF (all(this%dativar%d(i)%btable /= vl)) this%dativar%d(i) = vol7d_var_miss
11318 ENDDO
11319 ENDIF
11320
11321 IF (ASSOCIATED(this%dativar%c)) THEN
11322 DO i = 1, SIZE(this%dativar%c)
11323 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
11324 ENDDO
11325 ENDIF
11326
11327 IF (ASSOCIATED(this%dativar%c)) THEN
11328 DO i = 1, SIZE(this%dativar%c)
11329 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
11330 ENDDO
11331 ENDIF
11332
11333 ENDIF
11334ENDIF
11335
11336IF (PRESENT(nl)) THEN
11337 IF (SIZE(nl) > 0) THEN
11338 DO i = 1, SIZE(this%network)
11339 IF (all(this%network(i) /= nl)) this%network(i) = vol7d_network_miss
11340 ENDDO
11341 ENDIF
11342ENDIF
11343
11344IF (PRESENT(s_d)) THEN
11346 WHERE (this%time < s_d)
11347 this%time = datetime_miss
11348 END WHERE
11349 ENDIF
11350ENDIF
11351
11352IF (PRESENT(e_d)) THEN
11354 WHERE (this%time > e_d)
11355 this%time = datetime_miss
11356 END WHERE
11357 ENDIF
11358ENDIF
11359
11360CALL vol7d_reform(this, miss=.true.)
11361
11362END SUBROUTINE vol7d_filter
11363
11364
11371SUBROUTINE vol7d_convr(this, that, anaconv)
11372TYPE(vol7d),INTENT(IN) :: this
11373TYPE(vol7d),INTENT(INOUT) :: that
11374LOGICAL,OPTIONAL,INTENT(in) :: anaconv
11375INTEGER :: i
11376LOGICAL :: fv(1)=(/.false./), tv(1)=(/.true./), acp(1), acn(1)
11377TYPE(vol7d) :: v7d_tmp
11378
11379IF (optio_log(anaconv)) THEN
11380 acp=fv
11381 acn=tv
11382ELSE
11383 acp=tv
11384 acn=fv
11385ENDIF
11386
11387! Volume con solo i dati reali e tutti gli attributi
11388! l'anagrafica e` copiata interamente se necessario
11389CALL vol7d_copy(this, that, &
11390 lanavarr=tv, lanavard=acp, lanavari=acp, lanavarb=acp, lanavarc=acp, &
11391 ldativarr=tv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=fv)
11392
11393! Volume solo di dati double
11394CALL vol7d_copy(this, v7d_tmp, &
11395 lanavarr=fv, lanavard=acn, lanavari=fv, lanavarb=fv, lanavarc=fv, &
11396 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11397 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11398 ldativarr=fv, ldativard=tv, ldativari=fv, ldativarb=fv, ldativarc=fv, &
11399 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11400 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11401
11402! converto a dati reali
11403IF (ASSOCIATED(v7d_tmp%anavar%d) .OR. ASSOCIATED(v7d_tmp%dativar%d)) THEN
11404
11405 IF (ASSOCIATED(v7d_tmp%anavar%d)) THEN
11406! alloco i dati reali e vi trasferisco i double
11407 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanad, 1), SIZE(v7d_tmp%volanad, 2), &
11408 SIZE(v7d_tmp%volanad, 3)))
11409 DO i = 1, SIZE(v7d_tmp%anavar%d)
11410 v7d_tmp%volanar(:,i,:) = &
11411 realdat(v7d_tmp%volanad(:,i,:), v7d_tmp%anavar%d(i))
11412 ENDDO
11413 DEALLOCATE(v7d_tmp%volanad)
11414! trasferisco le variabili
11415 v7d_tmp%anavar%r => v7d_tmp%anavar%d
11416 NULLIFY(v7d_tmp%anavar%d)
11417 ENDIF
11418
11419 IF (ASSOCIATED(v7d_tmp%dativar%d)) THEN
11420! alloco i dati reali e vi trasferisco i double
11421 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatid, 1), SIZE(v7d_tmp%voldatid, 2), &
11422 SIZE(v7d_tmp%voldatid, 3), SIZE(v7d_tmp%voldatid, 4), SIZE(v7d_tmp%voldatid, 5), &
11423 SIZE(v7d_tmp%voldatid, 6)))
11424 DO i = 1, SIZE(v7d_tmp%dativar%d)
11425 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11426 realdat(v7d_tmp%voldatid(:,:,:,:,i,:), v7d_tmp%dativar%d(i))
11427 ENDDO
11428 DEALLOCATE(v7d_tmp%voldatid)
11429! trasferisco le variabili
11430 v7d_tmp%dativar%r => v7d_tmp%dativar%d
11431 NULLIFY(v7d_tmp%dativar%d)
11432 ENDIF
11433
11434! fondo con il volume definitivo
11435 CALL vol7d_merge(that, v7d_tmp)
11436ELSE
11438ENDIF
11439
11440
11441! Volume solo di dati interi
11442CALL vol7d_copy(this, v7d_tmp, &
11443 lanavarr=fv, lanavard=fv, lanavari=acn, lanavarb=fv, lanavarc=fv, &
11444 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11445 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11446 ldativarr=fv, ldativard=fv, ldativari=tv, ldativarb=fv, ldativarc=fv, &
11447 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11448 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11449
11450! converto a dati reali
11451IF (ASSOCIATED(v7d_tmp%anavar%i) .OR. ASSOCIATED(v7d_tmp%dativar%i)) THEN
11452
11453 IF (ASSOCIATED(v7d_tmp%anavar%i)) THEN
11454! alloco i dati reali e vi trasferisco gli interi
11455 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanai, 1), SIZE(v7d_tmp%volanai, 2), &
11456 SIZE(v7d_tmp%volanai, 3)))
11457 DO i = 1, SIZE(v7d_tmp%anavar%i)
11458 v7d_tmp%volanar(:,i,:) = &
11459 realdat(v7d_tmp%volanai(:,i,:), v7d_tmp%anavar%i(i))
11460 ENDDO
11461 DEALLOCATE(v7d_tmp%volanai)
11462! trasferisco le variabili
11463 v7d_tmp%anavar%r => v7d_tmp%anavar%i
11464 NULLIFY(v7d_tmp%anavar%i)
11465 ENDIF
11466
11467 IF (ASSOCIATED(v7d_tmp%dativar%i)) THEN
11468! alloco i dati reali e vi trasferisco gli interi
11469 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatii, 1), SIZE(v7d_tmp%voldatii, 2), &
11470 SIZE(v7d_tmp%voldatii, 3), SIZE(v7d_tmp%voldatii, 4), SIZE(v7d_tmp%voldatii, 5), &
11471 SIZE(v7d_tmp%voldatii, 6)))
11472 DO i = 1, SIZE(v7d_tmp%dativar%i)
11473 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11474 realdat(v7d_tmp%voldatii(:,:,:,:,i,:), v7d_tmp%dativar%i(i))
11475 ENDDO
11476 DEALLOCATE(v7d_tmp%voldatii)
11477! trasferisco le variabili
11478 v7d_tmp%dativar%r => v7d_tmp%dativar%i
11479 NULLIFY(v7d_tmp%dativar%i)
11480 ENDIF
11481
11482! fondo con il volume definitivo
11483 CALL vol7d_merge(that, v7d_tmp)
11484ELSE
11486ENDIF
11487
11488
11489! Volume solo di dati byte
11490CALL vol7d_copy(this, v7d_tmp, &
11491 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=acn, lanavarc=fv, &
11492 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11493 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11494 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=tv, ldativarc=fv, &
11495 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11496 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11497
11498! converto a dati reali
11499IF (ASSOCIATED(v7d_tmp%anavar%b) .OR. ASSOCIATED(v7d_tmp%dativar%b)) THEN
11500
11501 IF (ASSOCIATED(v7d_tmp%anavar%b)) THEN
11502! alloco i dati reali e vi trasferisco i byte
11503 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanab, 1), SIZE(v7d_tmp%volanab, 2), &
11504 SIZE(v7d_tmp%volanab, 3)))
11505 DO i = 1, SIZE(v7d_tmp%anavar%b)
11506 v7d_tmp%volanar(:,i,:) = &
11507 realdat(v7d_tmp%volanab(:,i,:), v7d_tmp%anavar%b(i))
11508 ENDDO
11509 DEALLOCATE(v7d_tmp%volanab)
11510! trasferisco le variabili
11511 v7d_tmp%anavar%r => v7d_tmp%anavar%b
11512 NULLIFY(v7d_tmp%anavar%b)
11513 ENDIF
11514
11515 IF (ASSOCIATED(v7d_tmp%dativar%b)) THEN
11516! alloco i dati reali e vi trasferisco i byte
11517 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatib, 1), SIZE(v7d_tmp%voldatib, 2), &
11518 SIZE(v7d_tmp%voldatib, 3), SIZE(v7d_tmp%voldatib, 4), SIZE(v7d_tmp%voldatib, 5), &
11519 SIZE(v7d_tmp%voldatib, 6)))
11520 DO i = 1, SIZE(v7d_tmp%dativar%b)
11521 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11522 realdat(v7d_tmp%voldatib(:,:,:,:,i,:), v7d_tmp%dativar%b(i))
11523 ENDDO
11524 DEALLOCATE(v7d_tmp%voldatib)
11525! trasferisco le variabili
11526 v7d_tmp%dativar%r => v7d_tmp%dativar%b
11527 NULLIFY(v7d_tmp%dativar%b)
11528 ENDIF
11529
11530! fondo con il volume definitivo
11531 CALL vol7d_merge(that, v7d_tmp)
11532ELSE
11534ENDIF
11535
11536
11537! Volume solo di dati character
11538CALL vol7d_copy(this, v7d_tmp, &
11539 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=fv, lanavarc=acn, &
11540 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11541 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11542 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=tv, &
11543 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11544 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11545
11546! converto a dati reali
11547IF (ASSOCIATED(v7d_tmp%anavar%c) .OR. ASSOCIATED(v7d_tmp%dativar%c)) THEN
11548
11549 IF (ASSOCIATED(v7d_tmp%anavar%c)) THEN
11550! alloco i dati reali e vi trasferisco i character
11551 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanac, 1), SIZE(v7d_tmp%volanac, 2), &
11552 SIZE(v7d_tmp%volanac, 3)))
11553 DO i = 1, SIZE(v7d_tmp%anavar%c)
11554 v7d_tmp%volanar(:,i,:) = &
11555 realdat(v7d_tmp%volanac(:,i,:), v7d_tmp%anavar%c(i))
11556 ENDDO
11557 DEALLOCATE(v7d_tmp%volanac)
11558! trasferisco le variabili
11559 v7d_tmp%anavar%r => v7d_tmp%anavar%c
11560 NULLIFY(v7d_tmp%anavar%c)
11561 ENDIF
11562
11563 IF (ASSOCIATED(v7d_tmp%dativar%c)) THEN
11564! alloco i dati reali e vi trasferisco i character
11565 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatic, 1), SIZE(v7d_tmp%voldatic, 2), &
11566 SIZE(v7d_tmp%voldatic, 3), SIZE(v7d_tmp%voldatic, 4), SIZE(v7d_tmp%voldatic, 5), &
11567 SIZE(v7d_tmp%voldatic, 6)))
11568 DO i = 1, SIZE(v7d_tmp%dativar%c)
11569 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11570 realdat(v7d_tmp%voldatic(:,:,:,:,i,:), v7d_tmp%dativar%c(i))
11571 ENDDO
11572 DEALLOCATE(v7d_tmp%voldatic)
11573! trasferisco le variabili
11574 v7d_tmp%dativar%r => v7d_tmp%dativar%c
11575 NULLIFY(v7d_tmp%dativar%c)
11576 ENDIF
11577
11578! fondo con il volume definitivo
11579 CALL vol7d_merge(that, v7d_tmp)
11580ELSE
11582ENDIF
11583
11584END SUBROUTINE vol7d_convr
11585
11586
11590SUBROUTINE vol7d_diff_only (this, that, data_only,ana)
11591TYPE(vol7d),INTENT(IN) :: this
11592TYPE(vol7d),INTENT(OUT) :: that
11593logical , optional, intent(in) :: data_only
11594logical , optional, intent(in) :: ana
11595logical :: ldata_only,lana
11596
11597IF (PRESENT(data_only)) THEN
11598 ldata_only = data_only
11599ELSE
11600 ldata_only = .false.
11601ENDIF
11602
11603IF (PRESENT(ana)) THEN
11604 lana = ana
11605ELSE
11606 lana = .false.
11607ENDIF
11608
11609
11610#undef VOL7D_POLY_ARRAY
11611#define VOL7D_POLY_ARRAY voldati
11612#include "vol7d_class_diff.F90"
11613#undef VOL7D_POLY_ARRAY
11614#define VOL7D_POLY_ARRAY voldatiattr
11615#include "vol7d_class_diff.F90"
11616#undef VOL7D_POLY_ARRAY
11617
11618if ( .not. ldata_only) then
11619
11620#define VOL7D_POLY_ARRAY volana
11621#include "vol7d_class_diff.F90"
11622#undef VOL7D_POLY_ARRAY
11623#define VOL7D_POLY_ARRAY volanaattr
11624#include "vol7d_class_diff.F90"
11625#undef VOL7D_POLY_ARRAY
11626
11627 if(lana)then
11628 where ( this%ana == that%ana )
11629 that%ana = vol7d_ana_miss
11630 end where
11631 end if
11632
11633end if
11634
11635
11636
11637END SUBROUTINE vol7d_diff_only
11638
11639
11640
11641! Creo le routine da ripetere per i vari tipi di dati di v7d
11642! tramite un template e il preprocessore
11643#undef VOL7D_POLY_TYPE
11644#undef VOL7D_POLY_TYPES
11645#define VOL7D_POLY_TYPE REAL
11646#define VOL7D_POLY_TYPES r
11647#include "vol7d_class_type_templ.F90"
11648#undef VOL7D_POLY_TYPE
11649#undef VOL7D_POLY_TYPES
11650#define VOL7D_POLY_TYPE DOUBLE PRECISION
11651#define VOL7D_POLY_TYPES d
11652#include "vol7d_class_type_templ.F90"
11653#undef VOL7D_POLY_TYPE
11654#undef VOL7D_POLY_TYPES
11655#define VOL7D_POLY_TYPE INTEGER
11656#define VOL7D_POLY_TYPES i
11657#include "vol7d_class_type_templ.F90"
11658#undef VOL7D_POLY_TYPE
11659#undef VOL7D_POLY_TYPES
11660#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
11661#define VOL7D_POLY_TYPES b
11662#include "vol7d_class_type_templ.F90"
11663#undef VOL7D_POLY_TYPE
11664#undef VOL7D_POLY_TYPES
11665#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
11666#define VOL7D_POLY_TYPES c
11667#include "vol7d_class_type_templ.F90"
11668
11669! Creo le routine da ripetere per i vari descrittori di dimensioni di v7d
11670! tramite un template e il preprocessore
11671#define VOL7D_SORT
11672#undef VOL7D_NO_ZERO_ALLOC
11673#undef VOL7D_POLY_TYPE
11674#define VOL7D_POLY_TYPE datetime
11675#include "vol7d_class_desc_templ.F90"
11676#undef VOL7D_POLY_TYPE
11677#define VOL7D_POLY_TYPE vol7d_timerange
11678#include "vol7d_class_desc_templ.F90"
11679#undef VOL7D_POLY_TYPE
11680#define VOL7D_POLY_TYPE vol7d_level
11681#include "vol7d_class_desc_templ.F90"
11682#undef VOL7D_SORT
11683#undef VOL7D_POLY_TYPE
11684#define VOL7D_POLY_TYPE vol7d_network
11685#include "vol7d_class_desc_templ.F90"
11686#undef VOL7D_POLY_TYPE
11687#define VOL7D_POLY_TYPE vol7d_ana
11688#include "vol7d_class_desc_templ.F90"
11689#define VOL7D_NO_ZERO_ALLOC
11690#undef VOL7D_POLY_TYPE
11691#define VOL7D_POLY_TYPE vol7d_var
11692#include "vol7d_class_desc_templ.F90"
11693
11703subroutine vol7d_write_on_file (this,unit,description,filename,filename_auto)
11704
11705TYPE(vol7d),INTENT(IN) :: this
11706integer,optional,intent(inout) :: unit
11707character(len=*),intent(in),optional :: filename
11708character(len=*),intent(out),optional :: filename_auto
11709character(len=*),INTENT(IN),optional :: description
11710
11711integer :: lunit
11712character(len=254) :: ldescription,arg,lfilename
11713integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
11714 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11715 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11716 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11717 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11718 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11719 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
11720!integer :: im,id,iy
11721integer :: tarray(8)
11722logical :: opened,exist
11723
11724 nana=0
11725 ntime=0
11726 ntimerange=0
11727 nlevel=0
11728 nnetwork=0
11729 ndativarr=0
11730 ndativari=0
11731 ndativarb=0
11732 ndativard=0
11733 ndativarc=0
11734 ndatiattrr=0
11735 ndatiattri=0
11736 ndatiattrb=0
11737 ndatiattrd=0
11738 ndatiattrc=0
11739 ndativarattrr=0
11740 ndativarattri=0
11741 ndativarattrb=0
11742 ndativarattrd=0
11743 ndativarattrc=0
11744 nanavarr=0
11745 nanavari=0
11746 nanavarb=0
11747 nanavard=0
11748 nanavarc=0
11749 nanaattrr=0
11750 nanaattri=0
11751 nanaattrb=0
11752 nanaattrd=0
11753 nanaattrc=0
11754 nanavarattrr=0
11755 nanavarattri=0
11756 nanavarattrb=0
11757 nanavarattrd=0
11758 nanavarattrc=0
11759
11760
11761!call idate(im,id,iy)
11762call date_and_time(values=tarray)
11763call getarg(0,arg)
11764
11765if (present(description))then
11766 ldescription=description
11767else
11768 ldescription="Vol7d generated by: "//trim(arg)
11769end if
11770
11771if (.not. present(unit))then
11772 lunit=getunit()
11773else
11774 if (unit==0)then
11775 lunit=getunit()
11776 unit=lunit
11777 else
11778 lunit=unit
11779 end if
11780end if
11781
11782lfilename=trim(arg)//".v7d"
11784
11785if (present(filename))then
11786 if (filename /= "")then
11787 lfilename=filename
11788 end if
11789end if
11790
11791if (present(filename_auto))filename_auto=lfilename
11792
11793
11794inquire(unit=lunit,opened=opened)
11795if (.not. opened) then
11796! inquire(file=lfilename, EXIST=exist)
11797! IF (exist) THEN
11798! CALL l4f_log(L4F_FATAL, &
11799! 'in vol7d_write_on_file, file exists, cannot open file '//TRIM(lfilename))
11800! CALL raise_fatal_error()
11801! ENDIF
11802 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM')
11803 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
11804end if
11805
11806if (associated(this%ana)) nana=size(this%ana)
11807if (associated(this%time)) ntime=size(this%time)
11808if (associated(this%timerange)) ntimerange=size(this%timerange)
11809if (associated(this%level)) nlevel=size(this%level)
11810if (associated(this%network)) nnetwork=size(this%network)
11811
11812if (associated(this%dativar%r)) ndativarr=size(this%dativar%r)
11813if (associated(this%dativar%i)) ndativari=size(this%dativar%i)
11814if (associated(this%dativar%b)) ndativarb=size(this%dativar%b)
11815if (associated(this%dativar%d)) ndativard=size(this%dativar%d)
11816if (associated(this%dativar%c)) ndativarc=size(this%dativar%c)
11817
11818if (associated(this%datiattr%r)) ndatiattrr=size(this%datiattr%r)
11819if (associated(this%datiattr%i)) ndatiattri=size(this%datiattr%i)
11820if (associated(this%datiattr%b)) ndatiattrb=size(this%datiattr%b)
11821if (associated(this%datiattr%d)) ndatiattrd=size(this%datiattr%d)
11822if (associated(this%datiattr%c)) ndatiattrc=size(this%datiattr%c)
11823
11824if (associated(this%dativarattr%r)) ndativarattrr=size(this%dativarattr%r)
11825if (associated(this%dativarattr%i)) ndativarattri=size(this%dativarattr%i)
11826if (associated(this%dativarattr%b)) ndativarattrb=size(this%dativarattr%b)
11827if (associated(this%dativarattr%d)) ndativarattrd=size(this%dativarattr%d)
11828if (associated(this%dativarattr%c)) ndativarattrc=size(this%dativarattr%c)
11829
11830if (associated(this%anavar%r)) nanavarr=size(this%anavar%r)
11831if (associated(this%anavar%i)) nanavari=size(this%anavar%i)
11832if (associated(this%anavar%b)) nanavarb=size(this%anavar%b)
11833if (associated(this%anavar%d)) nanavard=size(this%anavar%d)
11834if (associated(this%anavar%c)) nanavarc=size(this%anavar%c)
11835
11836if (associated(this%anaattr%r)) nanaattrr=size(this%anaattr%r)
11837if (associated(this%anaattr%i)) nanaattri=size(this%anaattr%i)
11838if (associated(this%anaattr%b)) nanaattrb=size(this%anaattr%b)
11839if (associated(this%anaattr%d)) nanaattrd=size(this%anaattr%d)
11840if (associated(this%anaattr%c)) nanaattrc=size(this%anaattr%c)
11841
11842if (associated(this%anavarattr%r)) nanavarattrr=size(this%anavarattr%r)
11843if (associated(this%anavarattr%i)) nanavarattri=size(this%anavarattr%i)
11844if (associated(this%anavarattr%b)) nanavarattrb=size(this%anavarattr%b)
11845if (associated(this%anavarattr%d)) nanavarattrd=size(this%anavarattr%d)
11846if (associated(this%anavarattr%c)) nanavarattrc=size(this%anavarattr%c)
11847
11848write(unit=lunit)ldescription
11849write(unit=lunit)tarray
11850
11851write(unit=lunit)&
11852 nana, ntime, ntimerange, nlevel, nnetwork, &
11853 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11854 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11855 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11856 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11857 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11858 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
11859 this%time_definition
11860
11861
11862!write(unit=lunit)this
11863
11864
11865!! prime 5 dimensioni
11868if (associated(this%level)) write(unit=lunit)this%level
11869if (associated(this%timerange)) write(unit=lunit)this%timerange
11870if (associated(this%network)) write(unit=lunit)this%network
11871
11872 !! 6a dimensione: variabile dell'anagrafica e dei dati
11873 !! con relativi attributi e in 5 tipi diversi
11874
11875if (associated(this%anavar%r)) write(unit=lunit)this%anavar%r
11876if (associated(this%anavar%i)) write(unit=lunit)this%anavar%i
11877if (associated(this%anavar%b)) write(unit=lunit)this%anavar%b
11878if (associated(this%anavar%d)) write(unit=lunit)this%anavar%d
11879if (associated(this%anavar%c)) write(unit=lunit)this%anavar%c
11880
11881if (associated(this%anaattr%r)) write(unit=lunit)this%anaattr%r
11882if (associated(this%anaattr%i)) write(unit=lunit)this%anaattr%i
11883if (associated(this%anaattr%b)) write(unit=lunit)this%anaattr%b
11884if (associated(this%anaattr%d)) write(unit=lunit)this%anaattr%d
11885if (associated(this%anaattr%c)) write(unit=lunit)this%anaattr%c
11886
11887if (associated(this%anavarattr%r)) write(unit=lunit)this%anavarattr%r
11888if (associated(this%anavarattr%i)) write(unit=lunit)this%anavarattr%i
11889if (associated(this%anavarattr%b)) write(unit=lunit)this%anavarattr%b
11890if (associated(this%anavarattr%d)) write(unit=lunit)this%anavarattr%d
11891if (associated(this%anavarattr%c)) write(unit=lunit)this%anavarattr%c
11892
11893if (associated(this%dativar%r)) write(unit=lunit)this%dativar%r
11894if (associated(this%dativar%i)) write(unit=lunit)this%dativar%i
11895if (associated(this%dativar%b)) write(unit=lunit)this%dativar%b
11896if (associated(this%dativar%d)) write(unit=lunit)this%dativar%d
11897if (associated(this%dativar%c)) write(unit=lunit)this%dativar%c
11898
11899if (associated(this%datiattr%r)) write(unit=lunit)this%datiattr%r
11900if (associated(this%datiattr%i)) write(unit=lunit)this%datiattr%i
11901if (associated(this%datiattr%b)) write(unit=lunit)this%datiattr%b
11902if (associated(this%datiattr%d)) write(unit=lunit)this%datiattr%d
11903if (associated(this%datiattr%c)) write(unit=lunit)this%datiattr%c
11904
11905if (associated(this%dativarattr%r)) write(unit=lunit)this%dativarattr%r
11906if (associated(this%dativarattr%i)) write(unit=lunit)this%dativarattr%i
11907if (associated(this%dativarattr%b)) write(unit=lunit)this%dativarattr%b
11908if (associated(this%dativarattr%d)) write(unit=lunit)this%dativarattr%d
11909if (associated(this%dativarattr%c)) write(unit=lunit)this%dativarattr%c
11910
11911!! Volumi di valori e attributi per anagrafica e dati
11912
11913if (associated(this%volanar)) write(unit=lunit)this%volanar
11914if (associated(this%volanaattrr)) write(unit=lunit)this%volanaattrr
11915if (associated(this%voldatir)) write(unit=lunit)this%voldatir
11916if (associated(this%voldatiattrr)) write(unit=lunit)this%voldatiattrr
11917
11918if (associated(this%volanai)) write(unit=lunit)this%volanai
11919if (associated(this%volanaattri)) write(unit=lunit)this%volanaattri
11920if (associated(this%voldatii)) write(unit=lunit)this%voldatii
11921if (associated(this%voldatiattri)) write(unit=lunit)this%voldatiattri
11922
11923if (associated(this%volanab)) write(unit=lunit)this%volanab
11924if (associated(this%volanaattrb)) write(unit=lunit)this%volanaattrb
11925if (associated(this%voldatib)) write(unit=lunit)this%voldatib
11926if (associated(this%voldatiattrb)) write(unit=lunit)this%voldatiattrb
11927
11928if (associated(this%volanad)) write(unit=lunit)this%volanad
11929if (associated(this%volanaattrd)) write(unit=lunit)this%volanaattrd
11930if (associated(this%voldatid)) write(unit=lunit)this%voldatid
11931if (associated(this%voldatiattrd)) write(unit=lunit)this%voldatiattrd
11932
11933if (associated(this%volanac)) write(unit=lunit)this%volanac
11934if (associated(this%volanaattrc)) write(unit=lunit)this%volanaattrc
11935if (associated(this%voldatic)) write(unit=lunit)this%voldatic
11936if (associated(this%voldatiattrc)) write(unit=lunit)this%voldatiattrc
11937
11938if (.not. present(unit)) close(unit=lunit)
11939
11940end subroutine vol7d_write_on_file
11941
11942
11949
11950
11951subroutine vol7d_read_from_file (this,unit,filename,description,tarray,filename_auto)
11952
11953TYPE(vol7d),INTENT(OUT) :: this
11954integer,intent(inout),optional :: unit
11955character(len=*),INTENT(in),optional :: filename
11956character(len=*),intent(out),optional :: filename_auto
11957character(len=*),INTENT(out),optional :: description
11958integer,intent(out),optional :: tarray(8)
11959
11960
11961integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
11962 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11963 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11964 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11965 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11966 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11967 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
11968
11969character(len=254) :: ldescription,lfilename,arg
11970integer :: ltarray(8),lunit,ios
11971logical :: opened,exist
11972
11973
11974call getarg(0,arg)
11975
11976if (.not. present(unit))then
11977 lunit=getunit()
11978else
11979 if (unit==0)then
11980 lunit=getunit()
11981 unit=lunit
11982 else
11983 lunit=unit
11984 end if
11985end if
11986
11987lfilename=trim(arg)//".v7d"
11989
11990if (present(filename))then
11991 if (filename /= "")then
11992 lfilename=filename
11993 end if
11994end if
11995
11996if (present(filename_auto))filename_auto=lfilename
11997
11998
11999inquire(unit=lunit,opened=opened)
12000IF (.NOT. opened) THEN
12001 inquire(file=lfilename,exist=exist)
12002 IF (.NOT.exist) THEN
12003 CALL l4f_log(l4f_fatal, &
12004 'in vol7d_read_from_file, file does not exists, cannot open')
12005 CALL raise_fatal_error()
12006 ENDIF
12007 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM', &
12008 status='OLD', action='READ')
12009 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
12010end if
12011
12012
12014read(unit=lunit,iostat=ios)ldescription
12015
12016if (ios < 0) then ! A negative value indicates that the End of File or End of Record
12017 call vol7d_alloc (this)
12018 call vol7d_alloc_vol (this)
12019 if (present(description))description=ldescription
12020 if (present(tarray))tarray=ltarray
12021 if (.not. present(unit)) close(unit=lunit)
12022end if
12023
12024read(unit=lunit)ltarray
12025
12026CALL l4f_log(l4f_info, 'Reading vol7d from file')
12027CALL l4f_log(l4f_info, 'description: '//trim(ldescription))
12030
12031if (present(description))description=ldescription
12032if (present(tarray))tarray=ltarray
12033
12034read(unit=lunit)&
12035 nana, ntime, ntimerange, nlevel, nnetwork, &
12036 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
12037 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
12038 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
12039 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
12040 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
12041 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
12042 this%time_definition
12043
12044call vol7d_alloc (this, &
12045 nana=nana, ntime=ntime, ntimerange=ntimerange, nlevel=nlevel, nnetwork=nnetwork,&
12046 ndativarr=ndativarr, ndativari=ndativari, ndativarb=ndativarb,&
12047 ndativard=ndativard, ndativarc=ndativarc,&
12048 ndatiattrr=ndatiattrr, ndatiattri=ndatiattri, ndatiattrb=ndatiattrb,&
12049 ndatiattrd=ndatiattrd, ndatiattrc=ndatiattrc,&
12050 ndativarattrr=ndativarattrr, ndativarattri=ndativarattri, ndativarattrb=ndativarattrb, &
12051 ndativarattrd=ndativarattrd, ndativarattrc=ndativarattrc,&
12052 nanavarr=nanavarr, nanavari=nanavari, nanavarb=nanavarb, &
12053 nanavard=nanavard, nanavarc=nanavarc,&
12054 nanaattrr=nanaattrr, nanaattri=nanaattri, nanaattrb=nanaattrb,&
12055 nanaattrd=nanaattrd, nanaattrc=nanaattrc,&
12056 nanavarattrr=nanavarattrr, nanavarattri=nanavarattri, nanavarattrb=nanavarattrb, &
12057 nanavarattrd=nanavarattrd, nanavarattrc=nanavarattrc)
12058
12059
12062if (associated(this%level)) read(unit=lunit)this%level
12063if (associated(this%timerange)) read(unit=lunit)this%timerange
12064if (associated(this%network)) read(unit=lunit)this%network
12065
12066if (associated(this%anavar%r)) read(unit=lunit)this%anavar%r
12067if (associated(this%anavar%i)) read(unit=lunit)this%anavar%i
12068if (associated(this%anavar%b)) read(unit=lunit)this%anavar%b
12069if (associated(this%anavar%d)) read(unit=lunit)this%anavar%d
12070if (associated(this%anavar%c)) read(unit=lunit)this%anavar%c
12071
12072if (associated(this%anaattr%r)) read(unit=lunit)this%anaattr%r
12073if (associated(this%anaattr%i)) read(unit=lunit)this%anaattr%i
12074if (associated(this%anaattr%b)) read(unit=lunit)this%anaattr%b
12075if (associated(this%anaattr%d)) read(unit=lunit)this%anaattr%d
12076if (associated(this%anaattr%c)) read(unit=lunit)this%anaattr%c
12077
12078if (associated(this%anavarattr%r)) read(unit=lunit)this%anavarattr%r
12079if (associated(this%anavarattr%i)) read(unit=lunit)this%anavarattr%i
12080if (associated(this%anavarattr%b)) read(unit=lunit)this%anavarattr%b
12081if (associated(this%anavarattr%d)) read(unit=lunit)this%anavarattr%d
12082if (associated(this%anavarattr%c)) read(unit=lunit)this%anavarattr%c
12083
12084if (associated(this%dativar%r)) read(unit=lunit)this%dativar%r
12085if (associated(this%dativar%i)) read(unit=lunit)this%dativar%i
12086if (associated(this%dativar%b)) read(unit=lunit)this%dativar%b
12087if (associated(this%dativar%d)) read(unit=lunit)this%dativar%d
12088if (associated(this%dativar%c)) read(unit=lunit)this%dativar%c
12089
12090if (associated(this%datiattr%r)) read(unit=lunit)this%datiattr%r
12091if (associated(this%datiattr%i)) read(unit=lunit)this%datiattr%i
12092if (associated(this%datiattr%b)) read(unit=lunit)this%datiattr%b
12093if (associated(this%datiattr%d)) read(unit=lunit)this%datiattr%d
12094if (associated(this%datiattr%c)) read(unit=lunit)this%datiattr%c
12095
12096if (associated(this%dativarattr%r)) read(unit=lunit)this%dativarattr%r
12097if (associated(this%dativarattr%i)) read(unit=lunit)this%dativarattr%i
12098if (associated(this%dativarattr%b)) read(unit=lunit)this%dativarattr%b
12099if (associated(this%dativarattr%d)) read(unit=lunit)this%dativarattr%d
12100if (associated(this%dativarattr%c)) read(unit=lunit)this%dativarattr%c
12101
12102call vol7d_alloc_vol (this)
12103
12104!! Volumi di valori e attributi per anagrafica e dati
12105
12106if (associated(this%volanar)) read(unit=lunit)this%volanar
12107if (associated(this%volanaattrr)) read(unit=lunit)this%volanaattrr
12108if (associated(this%voldatir)) read(unit=lunit)this%voldatir
12109if (associated(this%voldatiattrr)) read(unit=lunit)this%voldatiattrr
12110
12111if (associated(this%volanai)) read(unit=lunit)this%volanai
12112if (associated(this%volanaattri)) read(unit=lunit)this%volanaattri
12113if (associated(this%voldatii)) read(unit=lunit)this%voldatii
12114if (associated(this%voldatiattri)) read(unit=lunit)this%voldatiattri
12115
12116if (associated(this%volanab)) read(unit=lunit)this%volanab
12117if (associated(this%volanaattrb)) read(unit=lunit)this%volanaattrb
12118if (associated(this%voldatib)) read(unit=lunit)this%voldatib
12119if (associated(this%voldatiattrb)) read(unit=lunit)this%voldatiattrb
12120
12121if (associated(this%volanad)) read(unit=lunit)this%volanad
12122if (associated(this%volanaattrd)) read(unit=lunit)this%volanaattrd
12123if (associated(this%voldatid)) read(unit=lunit)this%voldatid
12124if (associated(this%voldatiattrd)) read(unit=lunit)this%voldatiattrd
12125
12126if (associated(this%volanac)) read(unit=lunit)this%volanac
12127if (associated(this%volanaattrc)) read(unit=lunit)this%volanaattrc
12128if (associated(this%voldatic)) read(unit=lunit)this%voldatic
12129if (associated(this%voldatiattrc)) read(unit=lunit)this%voldatiattrc
12130
12131if (.not. present(unit)) close(unit=lunit)
12132
12133end subroutine vol7d_read_from_file
12134
12135
12136! to double precision
12137elemental doubleprecision function doubledatd(voldat,var)
12138doubleprecision,intent(in) :: voldat
12139type(vol7d_var),intent(in) :: var
12140
12141doubledatd=voldat
12142
12143end function doubledatd
12144
12145
12146elemental doubleprecision function doubledatr(voldat,var)
12147real,intent(in) :: voldat
12148type(vol7d_var),intent(in) :: var
12149
12151 doubledatr=dble(voldat)
12152else
12153 doubledatr=dmiss
12154end if
12155
12156end function doubledatr
12157
12158
12159elemental doubleprecision function doubledati(voldat,var)
12160integer,intent(in) :: voldat
12161type(vol7d_var),intent(in) :: var
12162
12165 doubledati=dble(voldat)/10.d0**var%scalefactor
12166 else
12167 doubledati=dble(voldat)
12168 endif
12169else
12170 doubledati=dmiss
12171end if
12172
12173end function doubledati
12174
12175
12176elemental doubleprecision function doubledatb(voldat,var)
12177integer(kind=int_b),intent(in) :: voldat
12178type(vol7d_var),intent(in) :: var
12179
12182 doubledatb=dble(voldat)/10.d0**var%scalefactor
12183 else
12184 doubledatb=dble(voldat)
12185 endif
12186else
12187 doubledatb=dmiss
12188end if
12189
12190end function doubledatb
12191
12192
12193elemental doubleprecision function doubledatc(voldat,var)
12194CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
12195type(vol7d_var),intent(in) :: var
12196
12197doubledatc = c2d(voldat)
12199 doubledatc=doubledatc/10.d0**var%scalefactor
12200end if
12201
12202end function doubledatc
12203
12204
12205! to integer
12206elemental integer function integerdatd(voldat,var)
12207doubleprecision,intent(in) :: voldat
12208type(vol7d_var),intent(in) :: var
12209
12212 integerdatd=nint(voldat*10d0**var%scalefactor)
12213 else
12214 integerdatd=nint(voldat)
12215 endif
12216else
12217 integerdatd=imiss
12218end if
12219
12220end function integerdatd
12221
12222
12223elemental integer function integerdatr(voldat,var)
12224real,intent(in) :: voldat
12225type(vol7d_var),intent(in) :: var
12226
12229 integerdatr=nint(voldat*10d0**var%scalefactor)
12230 else
12231 integerdatr=nint(voldat)
12232 endif
12233else
12234 integerdatr=imiss
12235end if
12236
12237end function integerdatr
12238
12239
12240elemental integer function integerdati(voldat,var)
12241integer,intent(in) :: voldat
12242type(vol7d_var),intent(in) :: var
12243
12244integerdati=voldat
12245
12246end function integerdati
12247
12248
12249elemental integer function integerdatb(voldat,var)
12250integer(kind=int_b),intent(in) :: voldat
12251type(vol7d_var),intent(in) :: var
12252
12254 integerdatb=voldat
12255else
12256 integerdatb=imiss
12257end if
12258
12259end function integerdatb
12260
12261
12262elemental integer function integerdatc(voldat,var)
12263CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
12264type(vol7d_var),intent(in) :: var
12265
12266integerdatc=c2i(voldat)
12267
12268end function integerdatc
12269
12270
12271! to real
12272elemental real function realdatd(voldat,var)
12273doubleprecision,intent(in) :: voldat
12274type(vol7d_var),intent(in) :: var
12275
12277 realdatd=real(voldat)
12278else
12279 realdatd=rmiss
12280end if
12281
12282end function realdatd
12283
12284
12285elemental real function realdatr(voldat,var)
12286real,intent(in) :: voldat
12287type(vol7d_var),intent(in) :: var
12288
12289realdatr=voldat
12290
12291end function realdatr
12292
12293
12294elemental real function realdati(voldat,var)
12295integer,intent(in) :: voldat
12296type(vol7d_var),intent(in) :: var
12297
12300 realdati=float(voldat)/10.**var%scalefactor
12301 else
12302 realdati=float(voldat)
12303 endif
12304else
12305 realdati=rmiss
12306end if
12307
12308end function realdati
12309
12310
12311elemental real function realdatb(voldat,var)
12312integer(kind=int_b),intent(in) :: voldat
12313type(vol7d_var),intent(in) :: var
12314
12317 realdatb=float(voldat)/10**var%scalefactor
12318 else
12319 realdatb=float(voldat)
12320 endif
12321else
12322 realdatb=rmiss
12323end if
12324
12325end function realdatb
12326
12327
12328elemental real function realdatc(voldat,var)
12329CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
12330type(vol7d_var),intent(in) :: var
12331
12332realdatc=c2r(voldat)
12334 realdatc=realdatc/10.**var%scalefactor
12335end if
12336
12337end function realdatc
12338
12339
12345FUNCTION realanavol(this, var) RESULT(vol)
12346TYPE(vol7d),INTENT(in) :: this
12347TYPE(vol7d_var),INTENT(in) :: var
12348REAL :: vol(SIZE(this%ana),size(this%network))
12349
12350CHARACTER(len=1) :: dtype
12351INTEGER :: indvar
12352
12353dtype = cmiss
12354indvar = index(this%anavar, var, type=dtype)
12355
12356IF (indvar > 0) THEN
12357 SELECT CASE (dtype)
12358 CASE("d")
12359 vol = realdat(this%volanad(:,indvar,:), var)
12360 CASE("r")
12361 vol = this%volanar(:,indvar,:)
12362 CASE("i")
12363 vol = realdat(this%volanai(:,indvar,:), var)
12364 CASE("b")
12365 vol = realdat(this%volanab(:,indvar,:), var)
12366 CASE("c")
12367 vol = realdat(this%volanac(:,indvar,:), var)
12368 CASE default
12369 vol = rmiss
12370 END SELECT
12371ELSE
12372 vol = rmiss
12373ENDIF
12374
12375END FUNCTION realanavol
12376
12377
12383FUNCTION integeranavol(this, var) RESULT(vol)
12384TYPE(vol7d),INTENT(in) :: this
12385TYPE(vol7d_var),INTENT(in) :: var
12386INTEGER :: vol(SIZE(this%ana),size(this%network))
12387
12388CHARACTER(len=1) :: dtype
12389INTEGER :: indvar
12390
12391dtype = cmiss
12392indvar = index(this%anavar, var, type=dtype)
12393
12394IF (indvar > 0) THEN
12395 SELECT CASE (dtype)
12396 CASE("d")
12397 vol = integerdat(this%volanad(:,indvar,:), var)
12398 CASE("r")
12399 vol = integerdat(this%volanar(:,indvar,:), var)
12400 CASE("i")
12401 vol = this%volanai(:,indvar,:)
12402 CASE("b")
12403 vol = integerdat(this%volanab(:,indvar,:), var)
12404 CASE("c")
12405 vol = integerdat(this%volanac(:,indvar,:), var)
12406 CASE default
12407 vol = imiss
12408 END SELECT
12409ELSE
12410 vol = imiss
12411ENDIF
12412
12413END FUNCTION integeranavol
12414
12415
12421subroutine move_datac (v7d,&
12422 indana,indtime,indlevel,indtimerange,indnetwork,&
12423 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
12424
12425TYPE(vol7d),intent(inout) :: v7d
12426
12427integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
12428integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
12429integer :: inddativar,inddativarattr
12430
12431
12432do inddativar=1,size(v7d%dativar%c)
12433
12435 .not. c_e(v7d%voldatic(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
12436 ) then
12437
12438 ! dati
12439 v7d%voldatic &
12440 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
12441 v7d%voldatic &
12442 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
12443
12444
12445 ! attributi
12446 if (associated (v7d%dativarattr%i)) then
12447 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%c(inddativar))
12448 if (inddativarattr > 0 ) then
12449 v7d%voldatiattri &
12450 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12451 v7d%voldatiattri &
12452 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12453 end if
12454 end if
12455
12456 if (associated (v7d%dativarattr%r)) then
12457 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%c(inddativar))
12458 if (inddativarattr > 0 ) then
12459 v7d%voldatiattrr &
12460 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12461 v7d%voldatiattrr &
12462 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12463 end if
12464 end if
12465
12466 if (associated (v7d%dativarattr%d)) then
12467 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%c(inddativar))
12468 if (inddativarattr > 0 ) then
12469 v7d%voldatiattrd &
12470 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12471 v7d%voldatiattrd &
12472 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12473 end if
12474 end if
12475
12476 if (associated (v7d%dativarattr%b)) then
12477 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%c(inddativar))
12478 if (inddativarattr > 0 ) then
12479 v7d%voldatiattrb &
12480 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12481 v7d%voldatiattrb &
12482 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12483 end if
12484 end if
12485
12486 if (associated (v7d%dativarattr%c)) then
12487 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%c(inddativar))
12488 if (inddativarattr > 0 ) then
12489 v7d%voldatiattrc &
12490 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12491 v7d%voldatiattrc &
12492 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12493 end if
12494 end if
12495
12496 end if
12497
12498end do
12499
12500end subroutine move_datac
12501
12507subroutine move_datar (v7d,&
12508 indana,indtime,indlevel,indtimerange,indnetwork,&
12509 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
12510
12511TYPE(vol7d),intent(inout) :: v7d
12512
12513integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
12514integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
12515integer :: inddativar,inddativarattr
12516
12517
12518do inddativar=1,size(v7d%dativar%r)
12519
12521 .not. c_e(v7d%voldatir(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
12522 ) then
12523
12524 ! dati
12525 v7d%voldatir &
12526 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
12527 v7d%voldatir &
12528 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
12529
12530
12531 ! attributi
12532 if (associated (v7d%dativarattr%i)) then
12533 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%r(inddativar))
12534 if (inddativarattr > 0 ) then
12535 v7d%voldatiattri &
12536 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12537 v7d%voldatiattri &
12538 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12539 end if
12540 end if
12541
12542 if (associated (v7d%dativarattr%r)) then
12543 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%r(inddativar))
12544 if (inddativarattr > 0 ) then
12545 v7d%voldatiattrr &
12546 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12547 v7d%voldatiattrr &
12548 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12549 end if
12550 end if
12551
12552 if (associated (v7d%dativarattr%d)) then
12553 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%r(inddativar))
12554 if (inddativarattr > 0 ) then
12555 v7d%voldatiattrd &
12556 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12557 v7d%voldatiattrd &
12558 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12559 end if
12560 end if
12561
12562 if (associated (v7d%dativarattr%b)) then
12563 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%r(inddativar))
12564 if (inddativarattr > 0 ) then
12565 v7d%voldatiattrb &
12566 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12567 v7d%voldatiattrb &
12568 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12569 end if
12570 end if
12571
12572 if (associated (v7d%dativarattr%c)) then
12573 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%r(inddativar))
12574 if (inddativarattr > 0 ) then
12575 v7d%voldatiattrc &
12576 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12577 v7d%voldatiattrc &
12578 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12579 end if
12580 end if
12581
12582 end if
12583
12584end do
12585
12586end subroutine move_datar
12587
12588
12602subroutine v7d_rounding(v7din,v7dout,level,timerange,nostatproc)
12603type(vol7d),intent(inout) :: v7din
12604type(vol7d),intent(out) :: v7dout
12605type(vol7d_level),intent(in),optional :: level(:)
12606type(vol7d_timerange),intent(in),optional :: timerange(:)
12607!logical,intent(in),optional :: merge !< if there are data on more then one almost equal levels or timeranges
12608!! will be merged POINT BY POINT with priority for the fird data found ordered by icreasing var index
12609logical,intent(in),optional :: nostatproc
12610
12611integer :: nana,nlevel,ntime,ntimerange,nnetwork,nbin
12612integer :: iana,ilevel,itimerange,indl,indt,itime,inetwork
12613type(vol7d_level) :: roundlevel(size(v7din%level))
12614type(vol7d_timerange) :: roundtimerange(size(v7din%timerange))
12615type(vol7d) :: v7d_tmp
12616
12617
12618nbin=0
12619
12620if (associated(v7din%dativar%r)) nbin = nbin + size(v7din%dativar%r)
12621if (associated(v7din%dativar%i)) nbin = nbin + size(v7din%dativar%i)
12622if (associated(v7din%dativar%d)) nbin = nbin + size(v7din%dativar%d)
12623if (associated(v7din%dativar%b)) nbin = nbin + size(v7din%dativar%b)
12624
12626
12627roundlevel=v7din%level
12628
12629if (present(level))then
12630 do ilevel = 1, size(v7din%level)
12631 if ((any(v7din%level(ilevel) .almosteq. level))) then
12632 roundlevel(ilevel)=level(1)
12633 end if
12634 end do
12635end if
12636
12637roundtimerange=v7din%timerange
12638
12639if (present(timerange))then
12640 do itimerange = 1, size(v7din%timerange)
12641 if ((any(v7din%timerange(itimerange) .almosteq. timerange))) then
12642 roundtimerange(itimerange)=timerange(1)
12643 end if
12644 end do
12645end if
12646
12647!set istantaneous values everywere
12648!preserve p1 for forecast time
12649if (optio_log(nostatproc)) then
12650 roundtimerange(:)%timerange=254
12651 roundtimerange(:)%p2=0
12652end if
12653
12654
12655nana=size(v7din%ana)
12656nlevel=count_distinct(roundlevel,back=.true.)
12657ntime=size(v7din%time)
12658ntimerange=count_distinct(roundtimerange,back=.true.)
12659nnetwork=size(v7din%network)
12660
12662
12663if (nbin == 0) then
12665else
12666 call vol7d_convr(v7din,v7d_tmp)
12667end if
12668
12669v7d_tmp%level=roundlevel
12670v7d_tmp%timerange=roundtimerange
12671
12672do ilevel=1, size(v7d_tmp%level)
12673 indl=index(v7d_tmp%level,roundlevel(ilevel))
12674 do itimerange=1,size(v7d_tmp%timerange)
12675 indt=index(v7d_tmp%timerange,roundtimerange(itimerange))
12676
12677 if (indl /= ilevel .or. indt /= itimerange) then
12678
12679 do iana=1, nana
12680 do itime=1,ntime
12681 do inetwork=1,nnetwork
12682
12683 if (nbin > 0) then
12684 call move_datar (v7d_tmp,&
12685 iana,itime,ilevel,itimerange,inetwork,&
12686 iana,itime,indl,indt,inetwork)
12687 else
12688 call move_datac (v7d_tmp,&
12689 iana,itime,ilevel,itimerange,inetwork,&
12690 iana,itime,indl,indt,inetwork)
12691 end if
12692
12693 end do
12694 end do
12695 end do
12696
12697 end if
12698
12699 end do
12700end do
12701
12702! set to missing level and time > nlevel
12703do ilevel=nlevel+1,size(v7d_tmp%level)
12705end do
12706
12707do itimerange=ntimerange+1,size(v7d_tmp%timerange)
12709end do
12710
12711!copy with remove
12714
12715!call display(v7dout)
12716
12717end subroutine v7d_rounding
12718
12719
12721
12727
12728
Set of functions that return a trimmed CHARACTER representation of the input variable. Definition: char_utilities.F90:284 Legge un oggetto datetime/timedelta o un vettore di oggetti datetime/timedelta da un file FORMATTED o... Definition: datetime_class.F90:484 Scrive un oggetto datetime/timedelta o un vettore di oggetti datetime/timedelta su un file FORMATTED ... Definition: datetime_class.F90:491 Generic subroutine for checking OPTIONAL parameters. Definition: optional_values.f90:36 Check for problems return 0 if all check passed print diagnostics with log4f. Definition: vol7d_class.F90:451 Reduce some dimensions (level and timerage) for semplification (rounding). Definition: vol7d_class.F90:468 This module defines usefull general purpose function and subroutine. Definition: array_utilities.F90:218 Definition of constants to be used for declaring variables of a desired type. Definition: kinds.F90:251 Module for quickly interpreting the OPTIONAL parameters passed to a subprogram. Definition: optional_values.f90:28 Classe per la gestione dell'anagrafica di stazioni meteo e affini. Definition: vol7d_ana_class.F90:218 Classe per la gestione di un volume completo di dati osservati. Definition: vol7d_class.F90:279 Classe per la gestione dei livelli verticali in osservazioni meteo e affini. Definition: vol7d_level_class.F90:219 Classe per la gestione delle reti di stazioni per osservazioni meteo e affini. Definition: vol7d_network_class.F90:220 Classe per la gestione degli intervalli temporali di osservazioni meteo e affini. Definition: vol7d_timerange_class.F90:221 Classe per gestire un vettore di oggetti di tipo vol7d_var_class::vol7d_var. Definition: vol7d_varvect_class.f90:22 Definisce un oggetto contenente i volumi anagrafica e dati e tutti i descrittori delle loro dimension... Definition: vol7d_class.F90:318 |