libsim Versione 7.1.11
|
◆ integeranavol()
Return an ana volume of a requested variable as integer data. It returns a 2-d array of the proper shape (ana x network) for the ana variable requested, converted to integer type. If the conversion fails or if the variable is not contained in the ana volume, missing data are returned.
Definizione alla linea 9154 del file vol7d_class.F90. 9155! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
9156! authors:
9157! Davide Cesari <dcesari@arpa.emr.it>
9158! Paolo Patruno <ppatruno@arpa.emr.it>
9159
9160! This program is free software; you can redistribute it and/or
9161! modify it under the terms of the GNU General Public License as
9162! published by the Free Software Foundation; either version 2 of
9163! the License, or (at your option) any later version.
9164
9165! This program is distributed in the hope that it will be useful,
9166! but WITHOUT ANY WARRANTY; without even the implied warranty of
9167! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
9168! GNU General Public License for more details.
9169
9170! You should have received a copy of the GNU General Public License
9171! along with this program. If not, see <http://www.gnu.org/licenses/>.
9172#include "config.h"
9173
9185
9253IMPLICIT NONE
9254
9255
9256INTEGER, PARAMETER :: vol7d_maxdim_a = 3, vol7d_maxdim_aa = 4, &
9257 vol7d_maxdim_d = 6, vol7d_maxdim_ad = 7
9258
9259INTEGER, PARAMETER :: vol7d_ana_a=1
9260INTEGER, PARAMETER :: vol7d_var_a=2
9261INTEGER, PARAMETER :: vol7d_network_a=3
9262INTEGER, PARAMETER :: vol7d_attr_a=4
9263INTEGER, PARAMETER :: vol7d_ana_d=1
9264INTEGER, PARAMETER :: vol7d_time_d=2
9265INTEGER, PARAMETER :: vol7d_level_d=3
9266INTEGER, PARAMETER :: vol7d_timerange_d=4
9267INTEGER, PARAMETER :: vol7d_var_d=5
9268INTEGER, PARAMETER :: vol7d_network_d=6
9269INTEGER, PARAMETER :: vol7d_attr_d=7
9270INTEGER, PARAMETER :: vol7d_cdatalen=32
9271
9272TYPE vol7d_varmap
9273 INTEGER :: r, d, i, b, c
9274END TYPE vol7d_varmap
9275
9280 TYPE(vol7d_ana),POINTER :: ana(:)
9282 TYPE(datetime),POINTER :: time(:)
9284 TYPE(vol7d_level),POINTER :: level(:)
9286 TYPE(vol7d_timerange),POINTER :: timerange(:)
9288 TYPE(vol7d_network),POINTER :: network(:)
9290 TYPE(vol7d_varvect) :: anavar
9292 TYPE(vol7d_varvect) :: anaattr
9294 TYPE(vol7d_varvect) :: anavarattr
9296 TYPE(vol7d_varvect) :: dativar
9298 TYPE(vol7d_varvect) :: datiattr
9300 TYPE(vol7d_varvect) :: dativarattr
9301
9303 REAL,POINTER :: volanar(:,:,:)
9305 DOUBLE PRECISION,POINTER :: volanad(:,:,:)
9307 INTEGER,POINTER :: volanai(:,:,:)
9309 INTEGER(kind=int_b),POINTER :: volanab(:,:,:)
9311 CHARACTER(len=vol7d_cdatalen),POINTER :: volanac(:,:,:)
9312
9314 REAL,POINTER :: volanaattrr(:,:,:,:)
9316 DOUBLE PRECISION,POINTER :: volanaattrd(:,:,:,:)
9318 INTEGER,POINTER :: volanaattri(:,:,:,:)
9320 INTEGER(kind=int_b),POINTER :: volanaattrb(:,:,:,:)
9322 CHARACTER(len=vol7d_cdatalen),POINTER :: volanaattrc(:,:,:,:)
9323
9325 REAL,POINTER :: voldatir(:,:,:,:,:,:) ! sono i dati
9327 DOUBLE PRECISION,POINTER :: voldatid(:,:,:,:,:,:)
9329 INTEGER,POINTER :: voldatii(:,:,:,:,:,:)
9331 INTEGER(kind=int_b),POINTER :: voldatib(:,:,:,:,:,:)
9333 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatic(:,:,:,:,:,:)
9334
9336 REAL,POINTER :: voldatiattrr(:,:,:,:,:,:,:)
9338 DOUBLE PRECISION,POINTER :: voldatiattrd(:,:,:,:,:,:,:)
9340 INTEGER,POINTER :: voldatiattri(:,:,:,:,:,:,:)
9342 INTEGER(kind=int_b),POINTER :: voldatiattrb(:,:,:,:,:,:,:)
9344 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatiattrc(:,:,:,:,:,:,:)
9345
9347 integer :: time_definition
9348
9350
9355 MODULE PROCEDURE vol7d_init
9356END INTERFACE
9357
9360 MODULE PROCEDURE vol7d_delete
9361END INTERFACE
9362
9365 MODULE PROCEDURE vol7d_write_on_file
9366END INTERFACE
9367
9369INTERFACE import
9370 MODULE PROCEDURE vol7d_read_from_file
9371END INTERFACE
9372
9375 MODULE PROCEDURE vol7d_display, dat_display, dat_vect_display
9376END INTERFACE
9377
9380 MODULE PROCEDURE to_char_dat
9381END INTERFACE
9382
9385 MODULE PROCEDURE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
9386END INTERFACE
9387
9390 MODULE PROCEDURE realdatd,realdatr,realdati,realdatb,realdatc
9391END INTERFACE
9392
9395 MODULE PROCEDURE integerdatd,integerdatr,integerdati,integerdatb,integerdatc
9396END INTERFACE
9397
9400 MODULE PROCEDURE vol7d_copy
9401END INTERFACE
9402
9405 MODULE PROCEDURE vol7d_c_e
9406END INTERFACE
9407
9412 MODULE PROCEDURE vol7d_check
9413END INTERFACE
9414
9429 MODULE PROCEDURE v7d_rounding
9430END INTERFACE
9431
9432!!$INTERFACE get_volana
9433!!$ MODULE PROCEDURE vol7d_get_volanar, vol7d_get_volanad, vol7d_get_volanai, &
9434!!$ vol7d_get_volanab, vol7d_get_volanac
9435!!$END INTERFACE
9436!!$
9437!!$INTERFACE get_voldati
9438!!$ MODULE PROCEDURE vol7d_get_voldatir, vol7d_get_voldatid, vol7d_get_voldatii, &
9439!!$ vol7d_get_voldatib, vol7d_get_voldatic
9440!!$END INTERFACE
9441!!$
9442!!$INTERFACE get_volanaattr
9443!!$ MODULE PROCEDURE vol7d_get_volanaattrr, vol7d_get_volanaattrd, &
9444!!$ vol7d_get_volanaattri, vol7d_get_volanaattrb, vol7d_get_volanaattrc
9445!!$END INTERFACE
9446!!$
9447!!$INTERFACE get_voldatiattr
9448!!$ MODULE PROCEDURE vol7d_get_voldatiattrr, vol7d_get_voldatiattrd, &
9449!!$ vol7d_get_voldatiattri, vol7d_get_voldatiattrb, vol7d_get_voldatiattrc
9450!!$END INTERFACE
9451
9452PRIVATE vol7d_get_volr, vol7d_get_vold, vol7d_get_voli, vol7d_get_volb, &
9453 vol7d_get_volc, &
9454 volptr1dr, volptr2dr, volptr3dr, volptr4dr, volptr5dr, volptr6dr, volptr7dr, &
9455 volptr1dd, volptr2dd, volptr3dd, volptr4dd, volptr5dd, volptr6dd, volptr7dd, &
9456 volptr1di, volptr2di, volptr3di, volptr4di, volptr5di, volptr6di, volptr7di, &
9457 volptr1db, volptr2db, volptr3db, volptr4db, volptr5db, volptr6db, volptr7db, &
9458 volptr1dc, volptr2dc, volptr3dc, volptr4dc, volptr5dc, volptr6dc, volptr7dc, &
9459 vol7d_nullifyr, vol7d_nullifyd, vol7d_nullifyi, vol7d_nullifyb, vol7d_nullifyc, &
9460 vol7d_init, vol7d_delete, vol7d_write_on_file, vol7d_read_from_file, &
9461 vol7d_check_alloc_ana, vol7d_force_alloc_ana, &
9462 vol7d_check_alloc_dati, vol7d_force_alloc_dati, vol7d_force_alloc, &
9463 vol7d_display, dat_display, dat_vect_display, &
9464 to_char_dat, vol7d_check
9465
9466PRIVATE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
9467
9468PRIVATE vol7d_c_e
9469
9470CONTAINS
9471
9472
9477SUBROUTINE vol7d_init(this,time_definition)
9478TYPE(vol7d),intent(out) :: this
9479integer,INTENT(IN),OPTIONAL :: time_definition
9480
9487CALL vol7d_var_features_init() ! initialise var features table once
9488
9489NULLIFY(this%ana, this%time, this%level, this%timerange, this%network)
9490
9491NULLIFY(this%volanar, this%volanaattrr, this%voldatir, this%voldatiattrr)
9492NULLIFY(this%volanad, this%volanaattrd, this%voldatid, this%voldatiattrd)
9493NULLIFY(this%volanai, this%volanaattri, this%voldatii, this%voldatiattri)
9494NULLIFY(this%volanab, this%volanaattrb, this%voldatib, this%voldatiattrb)
9495NULLIFY(this%volanac, this%volanaattrc, this%voldatic, this%voldatiattrc)
9496
9497if(present(time_definition)) then
9498 this%time_definition=time_definition
9499else
9500 this%time_definition=1 !default to validity time
9501end if
9502
9503END SUBROUTINE vol7d_init
9504
9505
9509ELEMENTAL SUBROUTINE vol7d_delete(this, dataonly)
9510TYPE(vol7d),intent(inout) :: this
9511LOGICAL, INTENT(in), OPTIONAL :: dataonly
9512
9513
9514IF (.NOT. optio_log(dataonly)) THEN
9515 IF (ASSOCIATED(this%volanar)) DEALLOCATE(this%volanar)
9516 IF (ASSOCIATED(this%volanad)) DEALLOCATE(this%volanad)
9517 IF (ASSOCIATED(this%volanai)) DEALLOCATE(this%volanai)
9518 IF (ASSOCIATED(this%volanab)) DEALLOCATE(this%volanab)
9519 IF (ASSOCIATED(this%volanac)) DEALLOCATE(this%volanac)
9520 IF (ASSOCIATED(this%volanaattrr)) DEALLOCATE(this%volanaattrr)
9521 IF (ASSOCIATED(this%volanaattrd)) DEALLOCATE(this%volanaattrd)
9522 IF (ASSOCIATED(this%volanaattri)) DEALLOCATE(this%volanaattri)
9523 IF (ASSOCIATED(this%volanaattrb)) DEALLOCATE(this%volanaattrb)
9524 IF (ASSOCIATED(this%volanaattrc)) DEALLOCATE(this%volanaattrc)
9525ENDIF
9526IF (ASSOCIATED(this%voldatir)) DEALLOCATE(this%voldatir)
9527IF (ASSOCIATED(this%voldatid)) DEALLOCATE(this%voldatid)
9528IF (ASSOCIATED(this%voldatii)) DEALLOCATE(this%voldatii)
9529IF (ASSOCIATED(this%voldatib)) DEALLOCATE(this%voldatib)
9530IF (ASSOCIATED(this%voldatic)) DEALLOCATE(this%voldatic)
9531IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
9532IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
9533IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
9534IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
9535IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
9536
9537IF (.NOT. optio_log(dataonly)) THEN
9538 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
9539 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
9540ENDIF
9541IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
9542IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
9543IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
9544
9545IF (.NOT. optio_log(dataonly)) THEN
9549ENDIF
9553
9554END SUBROUTINE vol7d_delete
9555
9556
9557
9558integer function vol7d_check(this)
9559TYPE(vol7d),intent(in) :: this
9560integer :: i,j,k,l,m,n
9561
9562vol7d_check=0
9563
9564if (associated(this%voldatii)) then
9565do i = 1,size(this%voldatii,1)
9566 do j = 1,size(this%voldatii,2)
9567 do k = 1,size(this%voldatii,3)
9568 do l = 1,size(this%voldatii,4)
9569 do m = 1,size(this%voldatii,5)
9570 do n = 1,size(this%voldatii,6)
9571 if (this%voldatii(i,j,k,l,m,n) /= this%voldatii(i,j,k,l,m,n) ) then
9572 CALL l4f_log(l4f_warn,"check: abnormal value at voldatii("&
9574 vol7d_check=1
9575 end if
9576 end do
9577 end do
9578 end do
9579 end do
9580 end do
9581end do
9582end if
9583
9584
9585if (associated(this%voldatir)) then
9586do i = 1,size(this%voldatir,1)
9587 do j = 1,size(this%voldatir,2)
9588 do k = 1,size(this%voldatir,3)
9589 do l = 1,size(this%voldatir,4)
9590 do m = 1,size(this%voldatir,5)
9591 do n = 1,size(this%voldatir,6)
9592 if (this%voldatir(i,j,k,l,m,n) /= this%voldatir(i,j,k,l,m,n) ) then
9593 CALL l4f_log(l4f_warn,"check: abnormal value at voldatir("&
9595 vol7d_check=2
9596 end if
9597 end do
9598 end do
9599 end do
9600 end do
9601 end do
9602end do
9603end if
9604
9605if (associated(this%voldatid)) then
9606do i = 1,size(this%voldatid,1)
9607 do j = 1,size(this%voldatid,2)
9608 do k = 1,size(this%voldatid,3)
9609 do l = 1,size(this%voldatid,4)
9610 do m = 1,size(this%voldatid,5)
9611 do n = 1,size(this%voldatid,6)
9612 if (this%voldatid(i,j,k,l,m,n) /= this%voldatid(i,j,k,l,m,n) ) then
9613 CALL l4f_log(l4f_warn,"check: abnormal value at voldatid("&
9615 vol7d_check=3
9616 end if
9617 end do
9618 end do
9619 end do
9620 end do
9621 end do
9622end do
9623end if
9624
9625if (associated(this%voldatib)) then
9626do i = 1,size(this%voldatib,1)
9627 do j = 1,size(this%voldatib,2)
9628 do k = 1,size(this%voldatib,3)
9629 do l = 1,size(this%voldatib,4)
9630 do m = 1,size(this%voldatib,5)
9631 do n = 1,size(this%voldatib,6)
9632 if (this%voldatib(i,j,k,l,m,n) /= this%voldatib(i,j,k,l,m,n) ) then
9633 CALL l4f_log(l4f_warn,"check: abnormal value at voldatib("&
9635 vol7d_check=4
9636 end if
9637 end do
9638 end do
9639 end do
9640 end do
9641 end do
9642end do
9643end if
9644
9645end function vol7d_check
9646
9647
9648
9649!TODO da completare ! aborta se i volumi sono allocati a dimensione 0
9651SUBROUTINE vol7d_display(this)
9652TYPE(vol7d),intent(in) :: this
9653integer :: i
9654
9655REAL :: rdat
9656DOUBLE PRECISION :: ddat
9657INTEGER :: idat
9658INTEGER(kind=int_b) :: bdat
9659CHARACTER(len=vol7d_cdatalen) :: cdat
9660
9661
9662print*,"<<<<<<<<<<<<<<<<<<< vol7d object >>>>>>>>>>>>>>>>>>>>"
9663if (this%time_definition == 0) then
9664 print*,"TIME DEFINITION: time is reference time"
9665else if (this%time_definition == 1) then
9666 print*,"TIME DEFINITION: time is validity time"
9667else
9668 print*,"Time definition have a wrong walue:", this%time_definition
9669end if
9670
9671IF (ASSOCIATED(this%network))then
9672 print*,"---- network vector ----"
9673 print*,"elements=",size(this%network)
9674 do i=1, size(this%network)
9676 end do
9677end IF
9678
9679IF (ASSOCIATED(this%ana))then
9680 print*,"---- ana vector ----"
9681 print*,"elements=",size(this%ana)
9682 do i=1, size(this%ana)
9684 end do
9685end IF
9686
9687IF (ASSOCIATED(this%time))then
9688 print*,"---- time vector ----"
9689 print*,"elements=",size(this%time)
9690 do i=1, size(this%time)
9692 end do
9693end if
9694
9695IF (ASSOCIATED(this%level)) then
9696 print*,"---- level vector ----"
9697 print*,"elements=",size(this%level)
9698 do i =1,size(this%level)
9700 end do
9701end if
9702
9703IF (ASSOCIATED(this%timerange))then
9704 print*,"---- timerange vector ----"
9705 print*,"elements=",size(this%timerange)
9706 do i =1,size(this%timerange)
9708 end do
9709end if
9710
9711
9712print*,"---- ana vector ----"
9713print*,""
9714print*,"->>>>>>>>> anavar -"
9716print*,""
9717print*,"->>>>>>>>> anaattr -"
9719print*,""
9720print*,"->>>>>>>>> anavarattr -"
9722
9723print*,"-- ana data section (first point) --"
9724
9725idat=imiss
9726rdat=rmiss
9727ddat=dmiss
9728bdat=ibmiss
9729cdat=cmiss
9730
9731!ntime = MIN(SIZE(this%time),nprint)
9732!ntimerange = MIN(SIZE(this%timerange),nprint)
9733!nlevel = MIN(SIZE(this%level),nprint)
9734!nnetwork = MIN(SIZE(this%network),nprint)
9735!nana = MIN(SIZE(this%ana),nprint)
9736
9737IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0) THEN
9738if (associated(this%volanai)) then
9739 do i=1,size(this%anavar%i)
9740 idat=this%volanai(1,i,1)
9742 end do
9743end if
9744idat=imiss
9745
9746if (associated(this%volanar)) then
9747 do i=1,size(this%anavar%r)
9748 rdat=this%volanar(1,i,1)
9750 end do
9751end if
9752rdat=rmiss
9753
9754if (associated(this%volanad)) then
9755 do i=1,size(this%anavar%d)
9756 ddat=this%volanad(1,i,1)
9758 end do
9759end if
9760ddat=dmiss
9761
9762if (associated(this%volanab)) then
9763 do i=1,size(this%anavar%b)
9764 bdat=this%volanab(1,i,1)
9766 end do
9767end if
9768bdat=ibmiss
9769
9770if (associated(this%volanac)) then
9771 do i=1,size(this%anavar%c)
9772 cdat=this%volanac(1,i,1)
9774 end do
9775end if
9776cdat=cmiss
9777ENDIF
9778
9779print*,"---- data vector ----"
9780print*,""
9781print*,"->>>>>>>>> dativar -"
9783print*,""
9784print*,"->>>>>>>>> datiattr -"
9786print*,""
9787print*,"->>>>>>>>> dativarattr -"
9789
9790print*,"-- data data section (first point) --"
9791
9792idat=imiss
9793rdat=rmiss
9794ddat=dmiss
9795bdat=ibmiss
9796cdat=cmiss
9797
9798IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0 .AND. size(this%time) > 0 &
9799 .AND. size(this%level) > 0 .AND. size(this%timerange) > 0) THEN
9800if (associated(this%voldatii)) then
9801 do i=1,size(this%dativar%i)
9802 idat=this%voldatii(1,1,1,1,i,1)
9804 end do
9805end if
9806idat=imiss
9807
9808if (associated(this%voldatir)) then
9809 do i=1,size(this%dativar%r)
9810 rdat=this%voldatir(1,1,1,1,i,1)
9812 end do
9813end if
9814rdat=rmiss
9815
9816if (associated(this%voldatid)) then
9817 do i=1,size(this%dativar%d)
9818 ddat=this%voldatid(1,1,1,1,i,1)
9820 end do
9821end if
9822ddat=dmiss
9823
9824if (associated(this%voldatib)) then
9825 do i=1,size(this%dativar%b)
9826 bdat=this%voldatib(1,1,1,1,i,1)
9828 end do
9829end if
9830bdat=ibmiss
9831
9832if (associated(this%voldatic)) then
9833 do i=1,size(this%dativar%c)
9834 cdat=this%voldatic(1,1,1,1,i,1)
9836 end do
9837end if
9838cdat=cmiss
9839ENDIF
9840
9841print*,"<<<<<<<<<<<<<<<<<<< END vol7d object >>>>>>>>>>>>>>>>>>>>"
9842
9843END SUBROUTINE vol7d_display
9844
9845
9847SUBROUTINE dat_display(this,idat,rdat,ddat,bdat,cdat)
9848TYPE(vol7d_var),intent(in) :: this
9850REAL :: rdat
9852DOUBLE PRECISION :: ddat
9854INTEGER :: idat
9856INTEGER(kind=int_b) :: bdat
9858CHARACTER(len=*) :: cdat
9859
9860print *, to_char_dat(this,idat,rdat,ddat,bdat,cdat)
9861
9862end SUBROUTINE dat_display
9863
9865SUBROUTINE dat_vect_display(this,idat,rdat,ddat,bdat,cdat)
9866
9867TYPE(vol7d_var),intent(in) :: this(:)
9869REAL :: rdat(:)
9871DOUBLE PRECISION :: ddat(:)
9873INTEGER :: idat(:)
9875INTEGER(kind=int_b) :: bdat(:)
9877CHARACTER(len=*):: cdat(:)
9878
9879integer :: i
9880
9881do i =1,size(this)
9883end do
9884
9885end SUBROUTINE dat_vect_display
9886
9887
9888FUNCTION to_char_dat(this,idat,rdat,ddat,bdat,cdat)
9889#ifdef HAVE_DBALLE
9890USE dballef
9891#endif
9892TYPE(vol7d_var),INTENT(in) :: this
9894REAL :: rdat
9896DOUBLE PRECISION :: ddat
9898INTEGER :: idat
9900INTEGER(kind=int_b) :: bdat
9902CHARACTER(len=*) :: cdat
9903CHARACTER(len=80) :: to_char_dat
9904
9905CHARACTER(len=LEN(to_char_dat)) :: to_char_tmp
9906
9907
9908#ifdef HAVE_DBALLE
9909INTEGER :: handle, ier
9910
9911handle = 0
9912to_char_dat="VALUE: "
9913
9918
9920 ier = idba_messaggi(handle,"/dev/null", "w", "BUFR")
9921 ier = idba_spiegab(handle,this%btable,cdat,to_char_tmp)
9922 ier = idba_fatto(handle)
9923 to_char_dat=trim(to_char_dat)//" ;char> "//trim(to_char_tmp)
9924endif
9925
9926#else
9927
9928to_char_dat="VALUE: "
9934
9935#endif
9936
9937END FUNCTION to_char_dat
9938
9939
9942FUNCTION vol7d_c_e(this) RESULT(c_e)
9943TYPE(vol7d), INTENT(in) :: this
9944
9945LOGICAL :: c_e
9946
9948 ASSOCIATED(this%level) .OR. ASSOCIATED(this%timerange) .OR. &
9949 ASSOCIATED(this%network) .OR. &
9950 ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
9951 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
9952 ASSOCIATED(this%anavar%c) .OR. &
9953 ASSOCIATED(this%anaattr%r) .OR. ASSOCIATED(this%anaattr%d) .OR. &
9954 ASSOCIATED(this%anaattr%i) .OR. ASSOCIATED(this%anaattr%b) .OR. &
9955 ASSOCIATED(this%anaattr%c) .OR. &
9956 ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
9957 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
9958 ASSOCIATED(this%dativar%c) .OR. &
9959 ASSOCIATED(this%datiattr%r) .OR. ASSOCIATED(this%datiattr%d) .OR. &
9960 ASSOCIATED(this%datiattr%i) .OR. ASSOCIATED(this%datiattr%b) .OR. &
9961 ASSOCIATED(this%datiattr%c)
9962
9963END FUNCTION vol7d_c_e
9964
9965
10004SUBROUTINE vol7d_alloc(this, nana, ntime, nlevel, ntimerange, nnetwork, &
10005 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
10006 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
10007 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
10008 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
10009 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
10010 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc, &
10011 ini)
10012TYPE(vol7d),INTENT(inout) :: this
10013INTEGER,INTENT(in),OPTIONAL :: nana
10014INTEGER,INTENT(in),OPTIONAL :: ntime
10015INTEGER,INTENT(in),OPTIONAL :: nlevel
10016INTEGER,INTENT(in),OPTIONAL :: ntimerange
10017INTEGER,INTENT(in),OPTIONAL :: nnetwork
10019INTEGER,INTENT(in),OPTIONAL :: &
10020 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
10021 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
10022 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
10023 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
10024 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
10025 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc
10026LOGICAL,INTENT(in),OPTIONAL :: ini
10027
10028INTEGER :: i
10029LOGICAL :: linit
10030
10031IF (PRESENT(ini)) THEN
10032 linit = ini
10033ELSE
10034 linit = .false.
10035ENDIF
10036
10037! Dimensioni principali
10038IF (PRESENT(nana)) THEN
10039 IF (nana >= 0) THEN
10040 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
10041 ALLOCATE(this%ana(nana))
10042 IF (linit) THEN
10043 DO i = 1, nana
10045 ENDDO
10046 ENDIF
10047 ENDIF
10048ENDIF
10049IF (PRESENT(ntime)) THEN
10050 IF (ntime >= 0) THEN
10051 IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
10052 ALLOCATE(this%time(ntime))
10053 IF (linit) THEN
10054 DO i = 1, ntime
10056 ENDDO
10057 ENDIF
10058 ENDIF
10059ENDIF
10060IF (PRESENT(nlevel)) THEN
10061 IF (nlevel >= 0) THEN
10062 IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
10063 ALLOCATE(this%level(nlevel))
10064 IF (linit) THEN
10065 DO i = 1, nlevel
10067 ENDDO
10068 ENDIF
10069 ENDIF
10070ENDIF
10071IF (PRESENT(ntimerange)) THEN
10072 IF (ntimerange >= 0) THEN
10073 IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
10074 ALLOCATE(this%timerange(ntimerange))
10075 IF (linit) THEN
10076 DO i = 1, ntimerange
10078 ENDDO
10079 ENDIF
10080 ENDIF
10081ENDIF
10082IF (PRESENT(nnetwork)) THEN
10083 IF (nnetwork >= 0) THEN
10084 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
10085 ALLOCATE(this%network(nnetwork))
10086 IF (linit) THEN
10087 DO i = 1, nnetwork
10089 ENDDO
10090 ENDIF
10091 ENDIF
10092ENDIF
10093! Dimensioni dei tipi delle variabili
10094CALL vol7d_varvect_alloc(this%anavar, nanavarr, nanavard, &
10095 nanavari, nanavarb, nanavarc, ini)
10096CALL vol7d_varvect_alloc(this%anaattr, nanaattrr, nanaattrd, &
10097 nanaattri, nanaattrb, nanaattrc, ini)
10098CALL vol7d_varvect_alloc(this%anavarattr, nanavarattrr, nanavarattrd, &
10099 nanavarattri, nanavarattrb, nanavarattrc, ini)
10100CALL vol7d_varvect_alloc(this%dativar, ndativarr, ndativard, &
10101 ndativari, ndativarb, ndativarc, ini)
10102CALL vol7d_varvect_alloc(this%datiattr, ndatiattrr, ndatiattrd, &
10103 ndatiattri, ndatiattrb, ndatiattrc, ini)
10104CALL vol7d_varvect_alloc(this%dativarattr, ndativarattrr, ndativarattrd, &
10105 ndativarattri, ndativarattrb, ndativarattrc, ini)
10106
10107END SUBROUTINE vol7d_alloc
10108
10109
10110FUNCTION vol7d_check_alloc_ana(this)
10111TYPE(vol7d),INTENT(in) :: this
10112LOGICAL :: vol7d_check_alloc_ana
10113
10114vol7d_check_alloc_ana = ASSOCIATED(this%ana) .AND. ASSOCIATED(this%network)
10115
10116END FUNCTION vol7d_check_alloc_ana
10117
10118SUBROUTINE vol7d_force_alloc_ana(this, ini)
10119TYPE(vol7d),INTENT(inout) :: this
10120LOGICAL,INTENT(in),OPTIONAL :: ini
10121
10122! Alloco i descrittori minimi per avere un volume di anagrafica
10123IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=1, ini=ini)
10124IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=1, ini=ini)
10125
10126END SUBROUTINE vol7d_force_alloc_ana
10127
10128
10129FUNCTION vol7d_check_alloc_dati(this)
10130TYPE(vol7d),INTENT(in) :: this
10131LOGICAL :: vol7d_check_alloc_dati
10132
10133vol7d_check_alloc_dati = vol7d_check_alloc_ana(this) .AND. &
10134 ASSOCIATED(this%time) .AND. ASSOCIATED(this%level) .AND. &
10135 ASSOCIATED(this%timerange)
10136
10137END FUNCTION vol7d_check_alloc_dati
10138
10139SUBROUTINE vol7d_force_alloc_dati(this, ini)
10140TYPE(vol7d),INTENT(inout) :: this
10141LOGICAL,INTENT(in),OPTIONAL :: ini
10142
10143! Alloco i descrittori minimi per avere un volume di dati
10144CALL vol7d_force_alloc_ana(this, ini)
10145IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=1, ini=ini)
10146IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=1, ini=ini)
10147IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=1, ini=ini)
10148
10149END SUBROUTINE vol7d_force_alloc_dati
10150
10151
10152SUBROUTINE vol7d_force_alloc(this)
10153TYPE(vol7d),INTENT(inout) :: this
10154
10155! If anything really not allocated yet, allocate with size 0
10156IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=0)
10157IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=0)
10158IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=0)
10159IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=0)
10160IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=0)
10161
10162END SUBROUTINE vol7d_force_alloc
10163
10164
10165FUNCTION vol7d_check_vol(this)
10166TYPE(vol7d),INTENT(in) :: this
10167LOGICAL :: vol7d_check_vol
10168
10169vol7d_check_vol = c_e(this)
10170
10171! Anagrafica
10172IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
10173 vol7d_check_vol = .false.
10174ENDIF
10175
10176IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
10177 vol7d_check_vol = .false.
10178ENDIF
10179
10180IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
10181 vol7d_check_vol = .false.
10182ENDIF
10183
10184IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
10185 vol7d_check_vol = .false.
10186ENDIF
10187
10188IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
10189 vol7d_check_vol = .false.
10190ENDIF
10191IF (ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
10192 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
10193 ASSOCIATED(this%anavar%c)) THEN
10194 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_ana(this)
10195ENDIF
10196
10197! Attributi dell'anagrafica
10198IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
10199 .NOT.ASSOCIATED(this%volanaattrr)) THEN
10200 vol7d_check_vol = .false.
10201ENDIF
10202
10203IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
10204 .NOT.ASSOCIATED(this%volanaattrd)) THEN
10205 vol7d_check_vol = .false.
10206ENDIF
10207
10208IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
10209 .NOT.ASSOCIATED(this%volanaattri)) THEN
10210 vol7d_check_vol = .false.
10211ENDIF
10212
10213IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
10214 .NOT.ASSOCIATED(this%volanaattrb)) THEN
10215 vol7d_check_vol = .false.
10216ENDIF
10217
10218IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
10219 .NOT.ASSOCIATED(this%volanaattrc)) THEN
10220 vol7d_check_vol = .false.
10221ENDIF
10222
10223! Dati
10224IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
10225 vol7d_check_vol = .false.
10226ENDIF
10227
10228IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
10229 vol7d_check_vol = .false.
10230ENDIF
10231
10232IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
10233 vol7d_check_vol = .false.
10234ENDIF
10235
10236IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
10237 vol7d_check_vol = .false.
10238ENDIF
10239
10240IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
10241 vol7d_check_vol = .false.
10242ENDIF
10243
10244! Attributi dei dati
10245IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
10246 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
10247 vol7d_check_vol = .false.
10248ENDIF
10249
10250IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
10251 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
10252 vol7d_check_vol = .false.
10253ENDIF
10254
10255IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
10256 .NOT.ASSOCIATED(this%voldatiattri)) THEN
10257 vol7d_check_vol = .false.
10258ENDIF
10259
10260IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
10261 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
10262 vol7d_check_vol = .false.
10263ENDIF
10264
10265IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
10266 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
10267 vol7d_check_vol = .false.
10268ENDIF
10269IF (ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
10270 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
10271 ASSOCIATED(this%dativar%c)) THEN
10272 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_dati(this)
10273ENDIF
10274
10275END FUNCTION vol7d_check_vol
10276
10277
10292SUBROUTINE vol7d_alloc_vol(this, ini, inivol)
10293TYPE(vol7d),INTENT(inout) :: this
10294LOGICAL,INTENT(in),OPTIONAL :: ini
10295LOGICAL,INTENT(in),OPTIONAL :: inivol
10296
10297LOGICAL :: linivol
10298
10299IF (PRESENT(inivol)) THEN
10300 linivol = inivol
10301ELSE
10302 linivol = .true.
10303ENDIF
10304
10305! Anagrafica
10306IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
10307 CALL vol7d_force_alloc_ana(this, ini)
10308 ALLOCATE(this%volanar(SIZE(this%ana), SIZE(this%anavar%r), SIZE(this%network)))
10309 IF (linivol) this%volanar(:,:,:) = rmiss
10310ENDIF
10311
10312IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
10313 CALL vol7d_force_alloc_ana(this, ini)
10314 ALLOCATE(this%volanad(SIZE(this%ana), SIZE(this%anavar%d), SIZE(this%network)))
10315 IF (linivol) this%volanad(:,:,:) = rdmiss
10316ENDIF
10317
10318IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
10319 CALL vol7d_force_alloc_ana(this, ini)
10320 ALLOCATE(this%volanai(SIZE(this%ana), SIZE(this%anavar%i), SIZE(this%network)))
10321 IF (linivol) this%volanai(:,:,:) = imiss
10322ENDIF
10323
10324IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
10325 CALL vol7d_force_alloc_ana(this, ini)
10326 ALLOCATE(this%volanab(SIZE(this%ana), SIZE(this%anavar%b), SIZE(this%network)))
10327 IF (linivol) this%volanab(:,:,:) = ibmiss
10328ENDIF
10329
10330IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
10331 CALL vol7d_force_alloc_ana(this, ini)
10332 ALLOCATE(this%volanac(SIZE(this%ana), SIZE(this%anavar%c), SIZE(this%network)))
10333 IF (linivol) this%volanac(:,:,:) = cmiss
10334ENDIF
10335
10336! Attributi dell'anagrafica
10337IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
10338 .NOT.ASSOCIATED(this%volanaattrr)) THEN
10339 CALL vol7d_force_alloc_ana(this, ini)
10340 ALLOCATE(this%volanaattrr(SIZE(this%ana), SIZE(this%anavarattr%r), &
10341 SIZE(this%network), SIZE(this%anaattr%r)))
10342 IF (linivol) this%volanaattrr(:,:,:,:) = rmiss
10343ENDIF
10344
10345IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
10346 .NOT.ASSOCIATED(this%volanaattrd)) THEN
10347 CALL vol7d_force_alloc_ana(this, ini)
10348 ALLOCATE(this%volanaattrd(SIZE(this%ana), SIZE(this%anavarattr%d), &
10349 SIZE(this%network), SIZE(this%anaattr%d)))
10350 IF (linivol) this%volanaattrd(:,:,:,:) = rdmiss
10351ENDIF
10352
10353IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
10354 .NOT.ASSOCIATED(this%volanaattri)) THEN
10355 CALL vol7d_force_alloc_ana(this, ini)
10356 ALLOCATE(this%volanaattri(SIZE(this%ana), SIZE(this%anavarattr%i), &
10357 SIZE(this%network), SIZE(this%anaattr%i)))
10358 IF (linivol) this%volanaattri(:,:,:,:) = imiss
10359ENDIF
10360
10361IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
10362 .NOT.ASSOCIATED(this%volanaattrb)) THEN
10363 CALL vol7d_force_alloc_ana(this, ini)
10364 ALLOCATE(this%volanaattrb(SIZE(this%ana), SIZE(this%anavarattr%b), &
10365 SIZE(this%network), SIZE(this%anaattr%b)))
10366 IF (linivol) this%volanaattrb(:,:,:,:) = ibmiss
10367ENDIF
10368
10369IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
10370 .NOT.ASSOCIATED(this%volanaattrc)) THEN
10371 CALL vol7d_force_alloc_ana(this, ini)
10372 ALLOCATE(this%volanaattrc(SIZE(this%ana), SIZE(this%anavarattr%c), &
10373 SIZE(this%network), SIZE(this%anaattr%c)))
10374 IF (linivol) this%volanaattrc(:,:,:,:) = cmiss
10375ENDIF
10376
10377! Dati
10378IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
10379 CALL vol7d_force_alloc_dati(this, ini)
10380 ALLOCATE(this%voldatir(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10381 SIZE(this%timerange), SIZE(this%dativar%r), SIZE(this%network)))
10382 IF (linivol) this%voldatir(:,:,:,:,:,:) = rmiss
10383ENDIF
10384
10385IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
10386 CALL vol7d_force_alloc_dati(this, ini)
10387 ALLOCATE(this%voldatid(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10388 SIZE(this%timerange), SIZE(this%dativar%d), SIZE(this%network)))
10389 IF (linivol) this%voldatid(:,:,:,:,:,:) = rdmiss
10390ENDIF
10391
10392IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
10393 CALL vol7d_force_alloc_dati(this, ini)
10394 ALLOCATE(this%voldatii(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10395 SIZE(this%timerange), SIZE(this%dativar%i), SIZE(this%network)))
10396 IF (linivol) this%voldatii(:,:,:,:,:,:) = imiss
10397ENDIF
10398
10399IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
10400 CALL vol7d_force_alloc_dati(this, ini)
10401 ALLOCATE(this%voldatib(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10402 SIZE(this%timerange), SIZE(this%dativar%b), SIZE(this%network)))
10403 IF (linivol) this%voldatib(:,:,:,:,:,:) = ibmiss
10404ENDIF
10405
10406IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
10407 CALL vol7d_force_alloc_dati(this, ini)
10408 ALLOCATE(this%voldatic(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10409 SIZE(this%timerange), SIZE(this%dativar%c), SIZE(this%network)))
10410 IF (linivol) this%voldatic(:,:,:,:,:,:) = cmiss
10411ENDIF
10412
10413! Attributi dei dati
10414IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
10415 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
10416 CALL vol7d_force_alloc_dati(this, ini)
10417 ALLOCATE(this%voldatiattrr(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10418 SIZE(this%timerange), SIZE(this%dativarattr%r), SIZE(this%network), &
10419 SIZE(this%datiattr%r)))
10420 IF (linivol) this%voldatiattrr(:,:,:,:,:,:,:) = rmiss
10421ENDIF
10422
10423IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
10424 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
10425 CALL vol7d_force_alloc_dati(this, ini)
10426 ALLOCATE(this%voldatiattrd(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10427 SIZE(this%timerange), SIZE(this%dativarattr%d), SIZE(this%network), &
10428 SIZE(this%datiattr%d)))
10429 IF (linivol) this%voldatiattrd(:,:,:,:,:,:,:) = rdmiss
10430ENDIF
10431
10432IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
10433 .NOT.ASSOCIATED(this%voldatiattri)) THEN
10434 CALL vol7d_force_alloc_dati(this, ini)
10435 ALLOCATE(this%voldatiattri(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10436 SIZE(this%timerange), SIZE(this%dativarattr%i), SIZE(this%network), &
10437 SIZE(this%datiattr%i)))
10438 IF (linivol) this%voldatiattri(:,:,:,:,:,:,:) = imiss
10439ENDIF
10440
10441IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
10442 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
10443 CALL vol7d_force_alloc_dati(this, ini)
10444 ALLOCATE(this%voldatiattrb(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10445 SIZE(this%timerange), SIZE(this%dativarattr%b), SIZE(this%network), &
10446 SIZE(this%datiattr%b)))
10447 IF (linivol) this%voldatiattrb(:,:,:,:,:,:,:) = ibmiss
10448ENDIF
10449
10450IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
10451 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
10452 CALL vol7d_force_alloc_dati(this, ini)
10453 ALLOCATE(this%voldatiattrc(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10454 SIZE(this%timerange), SIZE(this%dativarattr%c), SIZE(this%network), &
10455 SIZE(this%datiattr%c)))
10456 IF (linivol) this%voldatiattrc(:,:,:,:,:,:,:) = cmiss
10457ENDIF
10458
10459! Catch-all method
10460CALL vol7d_force_alloc(this)
10461
10462! Creo gli indici var-attr
10463
10464#ifdef DEBUG
10465CALL l4f_log(l4f_debug,"calling: vol7d_set_attr_ind")
10466#endif
10467
10468CALL vol7d_set_attr_ind(this)
10469
10470
10471
10472END SUBROUTINE vol7d_alloc_vol
10473
10474
10481SUBROUTINE vol7d_set_attr_ind(this)
10482TYPE(vol7d),INTENT(inout) :: this
10483
10484INTEGER :: i
10485
10486! real
10487IF (ASSOCIATED(this%dativar%r)) THEN
10488 IF (ASSOCIATED(this%dativarattr%r)) THEN
10489 DO i = 1, SIZE(this%dativar%r)
10490 this%dativar%r(i)%r = &
10491 firsttrue(this%dativar%r(i)%btable == this%dativarattr%r(:)%btable)
10492 ENDDO
10493 ENDIF
10494
10495 IF (ASSOCIATED(this%dativarattr%d)) THEN
10496 DO i = 1, SIZE(this%dativar%r)
10497 this%dativar%r(i)%d = &
10498 firsttrue(this%dativar%r(i)%btable == this%dativarattr%d(:)%btable)
10499 ENDDO
10500 ENDIF
10501
10502 IF (ASSOCIATED(this%dativarattr%i)) THEN
10503 DO i = 1, SIZE(this%dativar%r)
10504 this%dativar%r(i)%i = &
10505 firsttrue(this%dativar%r(i)%btable == this%dativarattr%i(:)%btable)
10506 ENDDO
10507 ENDIF
10508
10509 IF (ASSOCIATED(this%dativarattr%b)) THEN
10510 DO i = 1, SIZE(this%dativar%r)
10511 this%dativar%r(i)%b = &
10512 firsttrue(this%dativar%r(i)%btable == this%dativarattr%b(:)%btable)
10513 ENDDO
10514 ENDIF
10515
10516 IF (ASSOCIATED(this%dativarattr%c)) THEN
10517 DO i = 1, SIZE(this%dativar%r)
10518 this%dativar%r(i)%c = &
10519 firsttrue(this%dativar%r(i)%btable == this%dativarattr%c(:)%btable)
10520 ENDDO
10521 ENDIF
10522ENDIF
10523! double
10524IF (ASSOCIATED(this%dativar%d)) THEN
10525 IF (ASSOCIATED(this%dativarattr%r)) THEN
10526 DO i = 1, SIZE(this%dativar%d)
10527 this%dativar%d(i)%r = &
10528 firsttrue(this%dativar%d(i)%btable == this%dativarattr%r(:)%btable)
10529 ENDDO
10530 ENDIF
10531
10532 IF (ASSOCIATED(this%dativarattr%d)) THEN
10533 DO i = 1, SIZE(this%dativar%d)
10534 this%dativar%d(i)%d = &
10535 firsttrue(this%dativar%d(i)%btable == this%dativarattr%d(:)%btable)
10536 ENDDO
10537 ENDIF
10538
10539 IF (ASSOCIATED(this%dativarattr%i)) THEN
10540 DO i = 1, SIZE(this%dativar%d)
10541 this%dativar%d(i)%i = &
10542 firsttrue(this%dativar%d(i)%btable == this%dativarattr%i(:)%btable)
10543 ENDDO
10544 ENDIF
10545
10546 IF (ASSOCIATED(this%dativarattr%b)) THEN
10547 DO i = 1, SIZE(this%dativar%d)
10548 this%dativar%d(i)%b = &
10549 firsttrue(this%dativar%d(i)%btable == this%dativarattr%b(:)%btable)
10550 ENDDO
10551 ENDIF
10552
10553 IF (ASSOCIATED(this%dativarattr%c)) THEN
10554 DO i = 1, SIZE(this%dativar%d)
10555 this%dativar%d(i)%c = &
10556 firsttrue(this%dativar%d(i)%btable == this%dativarattr%c(:)%btable)
10557 ENDDO
10558 ENDIF
10559ENDIF
10560! integer
10561IF (ASSOCIATED(this%dativar%i)) THEN
10562 IF (ASSOCIATED(this%dativarattr%r)) THEN
10563 DO i = 1, SIZE(this%dativar%i)
10564 this%dativar%i(i)%r = &
10565 firsttrue(this%dativar%i(i)%btable == this%dativarattr%r(:)%btable)
10566 ENDDO
10567 ENDIF
10568
10569 IF (ASSOCIATED(this%dativarattr%d)) THEN
10570 DO i = 1, SIZE(this%dativar%i)
10571 this%dativar%i(i)%d = &
10572 firsttrue(this%dativar%i(i)%btable == this%dativarattr%d(:)%btable)
10573 ENDDO
10574 ENDIF
10575
10576 IF (ASSOCIATED(this%dativarattr%i)) THEN
10577 DO i = 1, SIZE(this%dativar%i)
10578 this%dativar%i(i)%i = &
10579 firsttrue(this%dativar%i(i)%btable == this%dativarattr%i(:)%btable)
10580 ENDDO
10581 ENDIF
10582
10583 IF (ASSOCIATED(this%dativarattr%b)) THEN
10584 DO i = 1, SIZE(this%dativar%i)
10585 this%dativar%i(i)%b = &
10586 firsttrue(this%dativar%i(i)%btable == this%dativarattr%b(:)%btable)
10587 ENDDO
10588 ENDIF
10589
10590 IF (ASSOCIATED(this%dativarattr%c)) THEN
10591 DO i = 1, SIZE(this%dativar%i)
10592 this%dativar%i(i)%c = &
10593 firsttrue(this%dativar%i(i)%btable == this%dativarattr%c(:)%btable)
10594 ENDDO
10595 ENDIF
10596ENDIF
10597! byte
10598IF (ASSOCIATED(this%dativar%b)) THEN
10599 IF (ASSOCIATED(this%dativarattr%r)) THEN
10600 DO i = 1, SIZE(this%dativar%b)
10601 this%dativar%b(i)%r = &
10602 firsttrue(this%dativar%b(i)%btable == this%dativarattr%r(:)%btable)
10603 ENDDO
10604 ENDIF
10605
10606 IF (ASSOCIATED(this%dativarattr%d)) THEN
10607 DO i = 1, SIZE(this%dativar%b)
10608 this%dativar%b(i)%d = &
10609 firsttrue(this%dativar%b(i)%btable == this%dativarattr%d(:)%btable)
10610 ENDDO
10611 ENDIF
10612
10613 IF (ASSOCIATED(this%dativarattr%i)) THEN
10614 DO i = 1, SIZE(this%dativar%b)
10615 this%dativar%b(i)%i = &
10616 firsttrue(this%dativar%b(i)%btable == this%dativarattr%i(:)%btable)
10617 ENDDO
10618 ENDIF
10619
10620 IF (ASSOCIATED(this%dativarattr%b)) THEN
10621 DO i = 1, SIZE(this%dativar%b)
10622 this%dativar%b(i)%b = &
10623 firsttrue(this%dativar%b(i)%btable == this%dativarattr%b(:)%btable)
10624 ENDDO
10625 ENDIF
10626
10627 IF (ASSOCIATED(this%dativarattr%c)) THEN
10628 DO i = 1, SIZE(this%dativar%b)
10629 this%dativar%b(i)%c = &
10630 firsttrue(this%dativar%b(i)%btable == this%dativarattr%c(:)%btable)
10631 ENDDO
10632 ENDIF
10633ENDIF
10634! character
10635IF (ASSOCIATED(this%dativar%c)) THEN
10636 IF (ASSOCIATED(this%dativarattr%r)) THEN
10637 DO i = 1, SIZE(this%dativar%c)
10638 this%dativar%c(i)%r = &
10639 firsttrue(this%dativar%c(i)%btable == this%dativarattr%r(:)%btable)
10640 ENDDO
10641 ENDIF
10642
10643 IF (ASSOCIATED(this%dativarattr%d)) THEN
10644 DO i = 1, SIZE(this%dativar%c)
10645 this%dativar%c(i)%d = &
10646 firsttrue(this%dativar%c(i)%btable == this%dativarattr%d(:)%btable)
10647 ENDDO
10648 ENDIF
10649
10650 IF (ASSOCIATED(this%dativarattr%i)) THEN
10651 DO i = 1, SIZE(this%dativar%c)
10652 this%dativar%c(i)%i = &
10653 firsttrue(this%dativar%c(i)%btable == this%dativarattr%i(:)%btable)
10654 ENDDO
10655 ENDIF
10656
10657 IF (ASSOCIATED(this%dativarattr%b)) THEN
10658 DO i = 1, SIZE(this%dativar%c)
10659 this%dativar%c(i)%b = &
10660 firsttrue(this%dativar%c(i)%btable == this%dativarattr%b(:)%btable)
10661 ENDDO
10662 ENDIF
10663
10664 IF (ASSOCIATED(this%dativarattr%c)) THEN
10665 DO i = 1, SIZE(this%dativar%c)
10666 this%dativar%c(i)%c = &
10667 firsttrue(this%dativar%c(i)%btable == this%dativarattr%c(:)%btable)
10668 ENDDO
10669 ENDIF
10670ENDIF
10671
10672END SUBROUTINE vol7d_set_attr_ind
10673
10674
10679SUBROUTINE vol7d_merge(this, that, sort, bestdata, &
10680 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
10681TYPE(vol7d),INTENT(INOUT) :: this
10682TYPE(vol7d),INTENT(INOUT) :: that
10683LOGICAL,INTENT(IN),OPTIONAL :: sort
10684LOGICAL,INTENT(in),OPTIONAL :: bestdata
10685LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple ! experimental, please do not use outside the library now
10686
10687TYPE(vol7d) :: v7d_clean
10688
10689
10691 this = that
10693 that = v7d_clean ! destroy that without deallocating
10694ELSE ! Append that to this and destroy that
10696 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
10698ENDIF
10699
10700END SUBROUTINE vol7d_merge
10701
10702
10731SUBROUTINE vol7d_append(this, that, sort, bestdata, &
10732 ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple)
10733TYPE(vol7d),INTENT(INOUT) :: this
10734TYPE(vol7d),INTENT(IN) :: that
10735LOGICAL,INTENT(IN),OPTIONAL :: sort
10736! experimental, please do not use outside the library now, they force the use
10737! of a simplified mapping algorithm which is valid only whene the dimension
10738! content is the same in both volumes , or when one of them is empty
10739LOGICAL,INTENT(in),OPTIONAL :: bestdata
10740LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple
10741
10742
10743TYPE(vol7d) :: v7dtmp
10744LOGICAL :: lsort, lbestdata
10745INTEGER,POINTER :: remapt1(:), remapt2(:), remaptr1(:), remaptr2(:), &
10746 remapl1(:), remapl2(:), remapa1(:), remapa2(:), remapn1(:), remapn2(:)
10747
10749IF (.NOT.vol7d_check_vol(that)) RETURN ! be safe
10752 RETURN
10753ENDIF
10754
10755IF (this%time_definition /= that%time_definition) THEN
10756 CALL l4f_log(l4f_fatal, &
10757 'in vol7d_append, cannot append volumes with different &
10758 &time definition')
10759 CALL raise_fatal_error()
10760ENDIF
10761
10762! Completo l'allocazione per avere volumi a norma
10763CALL vol7d_alloc_vol(this)
10764
10768
10769! Calcolo le mappature tra volumi vecchi e volume nuovo
10770! I puntatori remap* vengono tutti o allocati o nullificati
10771IF (optio_log(ltimesimple)) THEN
10772 CALL vol7d_remap2simple_datetime(this%time, that%time, v7dtmp%time, &
10773 lsort, remapt1, remapt2)
10774ELSE
10775 CALL vol7d_remap2_datetime(this%time, that%time, v7dtmp%time, &
10776 lsort, remapt1, remapt2)
10777ENDIF
10778IF (optio_log(ltimerangesimple)) THEN
10779 CALL vol7d_remap2simple_vol7d_timerange(this%timerange, that%timerange, &
10780 v7dtmp%timerange, lsort, remaptr1, remaptr2)
10781ELSE
10782 CALL vol7d_remap2_vol7d_timerange(this%timerange, that%timerange, &
10783 v7dtmp%timerange, lsort, remaptr1, remaptr2)
10784ENDIF
10785IF (optio_log(llevelsimple)) THEN
10786 CALL vol7d_remap2simple_vol7d_level(this%level, that%level, v7dtmp%level, &
10787 lsort, remapl1, remapl2)
10788ELSE
10789 CALL vol7d_remap2_vol7d_level(this%level, that%level, v7dtmp%level, &
10790 lsort, remapl1, remapl2)
10791ENDIF
10792IF (optio_log(lanasimple)) THEN
10793 CALL vol7d_remap2simple_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
10794 .false., remapa1, remapa2)
10795ELSE
10796 CALL vol7d_remap2_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
10797 .false., remapa1, remapa2)
10798ENDIF
10799IF (optio_log(lnetworksimple)) THEN
10800 CALL vol7d_remap2simple_vol7d_network(this%network, that%network, v7dtmp%network, &
10801 .false., remapn1, remapn2)
10802ELSE
10803 CALL vol7d_remap2_vol7d_network(this%network, that%network, v7dtmp%network, &
10804 .false., remapn1, remapn2)
10805ENDIF
10806
10807! Faccio la fusione fisica dei volumi
10808CALL vol7d_merge_finalr(this, that, v7dtmp, &
10809 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10810 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10811CALL vol7d_merge_finald(this, that, v7dtmp, &
10812 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10813 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10814CALL vol7d_merge_finali(this, that, v7dtmp, &
10815 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10816 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10817CALL vol7d_merge_finalb(this, that, v7dtmp, &
10818 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10819 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10820CALL vol7d_merge_finalc(this, that, v7dtmp, &
10821 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10822 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10823
10824! Dealloco i vettori di rimappatura
10825IF (ASSOCIATED(remapt1)) DEALLOCATE(remapt1)
10826IF (ASSOCIATED(remapt2)) DEALLOCATE(remapt2)
10827IF (ASSOCIATED(remaptr1)) DEALLOCATE(remaptr1)
10828IF (ASSOCIATED(remaptr2)) DEALLOCATE(remaptr2)
10829IF (ASSOCIATED(remapl1)) DEALLOCATE(remapl1)
10830IF (ASSOCIATED(remapl2)) DEALLOCATE(remapl2)
10831IF (ASSOCIATED(remapa1)) DEALLOCATE(remapa1)
10832IF (ASSOCIATED(remapa2)) DEALLOCATE(remapa2)
10833IF (ASSOCIATED(remapn1)) DEALLOCATE(remapn1)
10834IF (ASSOCIATED(remapn2)) DEALLOCATE(remapn2)
10835
10836! Distruggo il vecchio volume e assegno il nuovo a this
10838this = v7dtmp
10839! Ricreo gli indici var-attr
10840CALL vol7d_set_attr_ind(this)
10841
10842END SUBROUTINE vol7d_append
10843
10844
10877SUBROUTINE vol7d_copy(this, that, sort, unique, miss, &
10878 lsort_time, lsort_timerange, lsort_level, &
10879 ltime, ltimerange, llevel, lana, lnetwork, &
10880 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
10881 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
10882 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
10883 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
10884 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
10885 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
10886TYPE(vol7d),INTENT(IN) :: this
10887TYPE(vol7d),INTENT(INOUT) :: that
10888LOGICAL,INTENT(IN),OPTIONAL :: sort
10889LOGICAL,INTENT(IN),OPTIONAL :: unique
10890LOGICAL,INTENT(IN),OPTIONAL :: miss
10891LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
10892LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
10893LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
10901LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
10903LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
10905LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
10907LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
10909LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
10911LOGICAL,INTENT(in),OPTIONAL :: &
10912 lanavarr(:), lanavard(:), lanavari(:), lanavarb(:), lanavarc(:), &
10913 lanaattrr(:), lanaattrd(:), lanaattri(:), lanaattrb(:), lanaattrc(:), &
10914 lanavarattrr(:), lanavarattrd(:), lanavarattri(:), lanavarattrb(:), lanavarattrc(:), &
10915 ldativarr(:), ldativard(:), ldativari(:), ldativarb(:), ldativarc(:), &
10916 ldatiattrr(:), ldatiattrd(:), ldatiattri(:), ldatiattrb(:), ldatiattrc(:), &
10917 ldativarattrr(:), ldativarattrd(:), ldativarattri(:), ldativarattrb(:), ldativarattrc(:)
10918
10919LOGICAL :: lsort, lunique, lmiss
10920INTEGER,POINTER :: remapt(:), remaptr(:), remapl(:), remapa(:), remapn(:)
10921
10924IF (.NOT.vol7d_check_vol(this)) RETURN ! be safe
10925
10929
10930! Calcolo le mappature tra volume vecchio e volume nuovo
10931! I puntatori remap* vengono tutti o allocati o nullificati
10932CALL vol7d_remap1_datetime(this%time, that%time, &
10933 lsort.OR.optio_log(lsort_time), lunique, lmiss, remapt, ltime)
10934CALL vol7d_remap1_vol7d_timerange(this%timerange, that%timerange, &
10935 lsort.OR.optio_log(lsort_timerange), lunique, lmiss, remaptr, ltimerange)
10936CALL vol7d_remap1_vol7d_level(this%level, that%level, &
10937 lsort.OR.optio_log(lsort_level), lunique, lmiss, remapl, llevel)
10938CALL vol7d_remap1_vol7d_ana(this%ana, that%ana, &
10939 lsort, lunique, lmiss, remapa, lana)
10940CALL vol7d_remap1_vol7d_network(this%network, that%network, &
10941 lsort, lunique, lmiss, remapn, lnetwork)
10942
10943! lanavari, lanavarb, lanavarc, &
10944! lanaattri, lanaattrb, lanaattrc, &
10945! lanavarattri, lanavarattrb, lanavarattrc, &
10946! ldativari, ldativarb, ldativarc, &
10947! ldatiattri, ldatiattrb, ldatiattrc, &
10948! ldativarattri, ldativarattrb, ldativarattrc
10949! Faccio la riforma fisica dei volumi
10950CALL vol7d_reform_finalr(this, that, &
10951 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10952 lanavarr, lanaattrr, lanavarattrr, ldativarr, ldatiattrr, ldativarattrr)
10953CALL vol7d_reform_finald(this, that, &
10954 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10955 lanavard, lanaattrd, lanavarattrd, ldativard, ldatiattrd, ldativarattrd)
10956CALL vol7d_reform_finali(this, that, &
10957 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10958 lanavari, lanaattri, lanavarattri, ldativari, ldatiattri, ldativarattri)
10959CALL vol7d_reform_finalb(this, that, &
10960 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10961 lanavarb, lanaattrb, lanavarattrb, ldativarb, ldatiattrb, ldativarattrb)
10962CALL vol7d_reform_finalc(this, that, &
10963 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10964 lanavarc, lanaattrc, lanavarattrc, ldativarc, ldatiattrc, ldativarattrc)
10965
10966! Dealloco i vettori di rimappatura
10967IF (ASSOCIATED(remapt)) DEALLOCATE(remapt)
10968IF (ASSOCIATED(remaptr)) DEALLOCATE(remaptr)
10969IF (ASSOCIATED(remapl)) DEALLOCATE(remapl)
10970IF (ASSOCIATED(remapa)) DEALLOCATE(remapa)
10971IF (ASSOCIATED(remapn)) DEALLOCATE(remapn)
10972
10973! Ricreo gli indici var-attr
10974CALL vol7d_set_attr_ind(that)
10975that%time_definition = this%time_definition
10976
10977END SUBROUTINE vol7d_copy
10978
10979
10990SUBROUTINE vol7d_reform(this, sort, unique, miss, &
10991 lsort_time, lsort_timerange, lsort_level, &
10992 ltime, ltimerange, llevel, lana, lnetwork, &
10993 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
10994 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
10995 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
10996 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
10997 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
10998 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc&
10999 ,purgeana)
11000TYPE(vol7d),INTENT(INOUT) :: this
11001LOGICAL,INTENT(IN),OPTIONAL :: sort
11002LOGICAL,INTENT(IN),OPTIONAL :: unique
11003LOGICAL,INTENT(IN),OPTIONAL :: miss
11004LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
11005LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
11006LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
11014LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
11015LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
11016LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
11017LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
11018LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
11020LOGICAL,INTENT(in),OPTIONAL :: &
11021 lanavarr(:), lanavard(:), lanavari(:), lanavarb(:), lanavarc(:), &
11022 lanaattrr(:), lanaattrd(:), lanaattri(:), lanaattrb(:), lanaattrc(:), &
11023 lanavarattrr(:), lanavarattrd(:), lanavarattri(:), lanavarattrb(:), lanavarattrc(:), &
11024 ldativarr(:), ldativard(:), ldativari(:), ldativarb(:), ldativarc(:), &
11025 ldatiattrr(:), ldatiattrd(:), ldatiattri(:), ldatiattrb(:), ldatiattrc(:), &
11026 ldativarattrr(:), ldativarattrd(:), ldativarattri(:), ldativarattrb(:), ldativarattrc(:)
11027LOGICAL,INTENT(IN),OPTIONAL :: purgeana
11028
11029TYPE(vol7d) :: v7dtmp
11030logical,allocatable :: llana(:)
11031integer :: i
11032
11034 lsort_time, lsort_timerange, lsort_level, &
11035 ltime, ltimerange, llevel, lana, lnetwork, &
11036 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
11037 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
11038 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
11039 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
11040 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
11041 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
11042
11043! destroy old volume
11045
11046if (optio_log(purgeana)) then
11047 allocate(llana(size(v7dtmp%ana)))
11048 llana =.false.
11049 do i =1,size(v7dtmp%ana)
11050 if (associated(v7dtmp%voldatii)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatii(i,:,:,:,:,:)))
11051 if (associated(v7dtmp%voldatir)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatir(i,:,:,:,:,:)))
11052 if (associated(v7dtmp%voldatid)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatid(i,:,:,:,:,:)))
11053 if (associated(v7dtmp%voldatib)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatib(i,:,:,:,:,:)))
11054 if (associated(v7dtmp%voldatic)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatic(i,:,:,:,:,:)))
11055 end do
11056 CALL vol7d_copy(v7dtmp, this,lana=llana)
11058 deallocate(llana)
11059else
11060 this=v7dtmp
11061end if
11062
11063END SUBROUTINE vol7d_reform
11064
11065
11073SUBROUTINE vol7d_smart_sort(this, lsort_time, lsort_timerange, lsort_level)
11074TYPE(vol7d),INTENT(INOUT) :: this
11075LOGICAL,OPTIONAL,INTENT(in) :: lsort_time
11076LOGICAL,OPTIONAL,INTENT(in) :: lsort_timerange
11077LOGICAL,OPTIONAL,INTENT(in) :: lsort_level
11078
11079INTEGER :: i
11080LOGICAL :: to_be_sorted
11081
11082to_be_sorted = .false.
11083CALL vol7d_alloc_vol(this) ! usual safety check
11084
11085IF (optio_log(lsort_time)) THEN
11086 DO i = 2, SIZE(this%time)
11087 IF (this%time(i) < this%time(i-1)) THEN
11088 to_be_sorted = .true.
11089 EXIT
11090 ENDIF
11091 ENDDO
11092ENDIF
11093IF (optio_log(lsort_timerange)) THEN
11094 DO i = 2, SIZE(this%timerange)
11095 IF (this%timerange(i) < this%timerange(i-1)) THEN
11096 to_be_sorted = .true.
11097 EXIT
11098 ENDIF
11099 ENDDO
11100ENDIF
11101IF (optio_log(lsort_level)) THEN
11102 DO i = 2, SIZE(this%level)
11103 IF (this%level(i) < this%level(i-1)) THEN
11104 to_be_sorted = .true.
11105 EXIT
11106 ENDIF
11107 ENDDO
11108ENDIF
11109
11110IF (to_be_sorted) CALL vol7d_reform(this, &
11111 lsort_time=lsort_time, lsort_timerange=lsort_timerange, lsort_level=lsort_level )
11112
11113END SUBROUTINE vol7d_smart_sort
11114
11122SUBROUTINE vol7d_filter(this, avl, vl, nl, s_d, e_d)
11123TYPE(vol7d),INTENT(inout) :: this
11124CHARACTER(len=*),INTENT(in),OPTIONAL :: avl(:)
11125CHARACTER(len=*),INTENT(in),OPTIONAL :: vl(:)
11126TYPE(vol7d_network),OPTIONAL :: nl(:)
11127TYPE(datetime),INTENT(in),OPTIONAL :: s_d
11128TYPE(datetime),INTENT(in),OPTIONAL :: e_d
11129
11130INTEGER :: i
11131
11132IF (PRESENT(avl)) THEN
11133 IF (SIZE(avl) > 0) THEN
11134
11135 IF (ASSOCIATED(this%anavar%r)) THEN
11136 DO i = 1, SIZE(this%anavar%r)
11137 IF (all(this%anavar%r(i)%btable /= avl)) this%anavar%r(i) = vol7d_var_miss
11138 ENDDO
11139 ENDIF
11140
11141 IF (ASSOCIATED(this%anavar%i)) THEN
11142 DO i = 1, SIZE(this%anavar%i)
11143 IF (all(this%anavar%i(i)%btable /= avl)) this%anavar%i(i) = vol7d_var_miss
11144 ENDDO
11145 ENDIF
11146
11147 IF (ASSOCIATED(this%anavar%b)) THEN
11148 DO i = 1, SIZE(this%anavar%b)
11149 IF (all(this%anavar%b(i)%btable /= avl)) this%anavar%b(i) = vol7d_var_miss
11150 ENDDO
11151 ENDIF
11152
11153 IF (ASSOCIATED(this%anavar%d)) THEN
11154 DO i = 1, SIZE(this%anavar%d)
11155 IF (all(this%anavar%d(i)%btable /= avl)) this%anavar%d(i) = vol7d_var_miss
11156 ENDDO
11157 ENDIF
11158
11159 IF (ASSOCIATED(this%anavar%c)) THEN
11160 DO i = 1, SIZE(this%anavar%c)
11161 IF (all(this%anavar%c(i)%btable /= avl)) this%anavar%c(i) = vol7d_var_miss
11162 ENDDO
11163 ENDIF
11164
11165 ENDIF
11166ENDIF
11167
11168
11169IF (PRESENT(vl)) THEN
11170 IF (size(vl) > 0) THEN
11171 IF (ASSOCIATED(this%dativar%r)) THEN
11172 DO i = 1, SIZE(this%dativar%r)
11173 IF (all(this%dativar%r(i)%btable /= vl)) this%dativar%r(i) = vol7d_var_miss
11174 ENDDO
11175 ENDIF
11176
11177 IF (ASSOCIATED(this%dativar%i)) THEN
11178 DO i = 1, SIZE(this%dativar%i)
11179 IF (all(this%dativar%i(i)%btable /= vl)) this%dativar%i(i) = vol7d_var_miss
11180 ENDDO
11181 ENDIF
11182
11183 IF (ASSOCIATED(this%dativar%b)) THEN
11184 DO i = 1, SIZE(this%dativar%b)
11185 IF (all(this%dativar%b(i)%btable /= vl)) this%dativar%b(i) = vol7d_var_miss
11186 ENDDO
11187 ENDIF
11188
11189 IF (ASSOCIATED(this%dativar%d)) THEN
11190 DO i = 1, SIZE(this%dativar%d)
11191 IF (all(this%dativar%d(i)%btable /= vl)) this%dativar%d(i) = vol7d_var_miss
11192 ENDDO
11193 ENDIF
11194
11195 IF (ASSOCIATED(this%dativar%c)) THEN
11196 DO i = 1, SIZE(this%dativar%c)
11197 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
11198 ENDDO
11199 ENDIF
11200
11201 IF (ASSOCIATED(this%dativar%c)) THEN
11202 DO i = 1, SIZE(this%dativar%c)
11203 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
11204 ENDDO
11205 ENDIF
11206
11207 ENDIF
11208ENDIF
11209
11210IF (PRESENT(nl)) THEN
11211 IF (SIZE(nl) > 0) THEN
11212 DO i = 1, SIZE(this%network)
11213 IF (all(this%network(i) /= nl)) this%network(i) = vol7d_network_miss
11214 ENDDO
11215 ENDIF
11216ENDIF
11217
11218IF (PRESENT(s_d)) THEN
11220 WHERE (this%time < s_d)
11221 this%time = datetime_miss
11222 END WHERE
11223 ENDIF
11224ENDIF
11225
11226IF (PRESENT(e_d)) THEN
11228 WHERE (this%time > e_d)
11229 this%time = datetime_miss
11230 END WHERE
11231 ENDIF
11232ENDIF
11233
11234CALL vol7d_reform(this, miss=.true.)
11235
11236END SUBROUTINE vol7d_filter
11237
11238
11245SUBROUTINE vol7d_convr(this, that, anaconv)
11246TYPE(vol7d),INTENT(IN) :: this
11247TYPE(vol7d),INTENT(INOUT) :: that
11248LOGICAL,OPTIONAL,INTENT(in) :: anaconv
11249INTEGER :: i
11250LOGICAL :: fv(1)=(/.false./), tv(1)=(/.true./), acp(1), acn(1)
11251TYPE(vol7d) :: v7d_tmp
11252
11253IF (optio_log(anaconv)) THEN
11254 acp=fv
11255 acn=tv
11256ELSE
11257 acp=tv
11258 acn=fv
11259ENDIF
11260
11261! Volume con solo i dati reali e tutti gli attributi
11262! l'anagrafica e` copiata interamente se necessario
11263CALL vol7d_copy(this, that, &
11264 lanavarr=tv, lanavard=acp, lanavari=acp, lanavarb=acp, lanavarc=acp, &
11265 ldativarr=tv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=fv)
11266
11267! Volume solo di dati double
11268CALL vol7d_copy(this, v7d_tmp, &
11269 lanavarr=fv, lanavard=acn, lanavari=fv, lanavarb=fv, lanavarc=fv, &
11270 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11271 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11272 ldativarr=fv, ldativard=tv, ldativari=fv, ldativarb=fv, ldativarc=fv, &
11273 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11274 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11275
11276! converto a dati reali
11277IF (ASSOCIATED(v7d_tmp%anavar%d) .OR. ASSOCIATED(v7d_tmp%dativar%d)) THEN
11278
11279 IF (ASSOCIATED(v7d_tmp%anavar%d)) THEN
11280! alloco i dati reali e vi trasferisco i double
11281 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanad, 1), SIZE(v7d_tmp%volanad, 2), &
11282 SIZE(v7d_tmp%volanad, 3)))
11283 DO i = 1, SIZE(v7d_tmp%anavar%d)
11284 v7d_tmp%volanar(:,i,:) = &
11285 realdat(v7d_tmp%volanad(:,i,:), v7d_tmp%anavar%d(i))
11286 ENDDO
11287 DEALLOCATE(v7d_tmp%volanad)
11288! trasferisco le variabili
11289 v7d_tmp%anavar%r => v7d_tmp%anavar%d
11290 NULLIFY(v7d_tmp%anavar%d)
11291 ENDIF
11292
11293 IF (ASSOCIATED(v7d_tmp%dativar%d)) THEN
11294! alloco i dati reali e vi trasferisco i double
11295 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatid, 1), SIZE(v7d_tmp%voldatid, 2), &
11296 SIZE(v7d_tmp%voldatid, 3), SIZE(v7d_tmp%voldatid, 4), SIZE(v7d_tmp%voldatid, 5), &
11297 SIZE(v7d_tmp%voldatid, 6)))
11298 DO i = 1, SIZE(v7d_tmp%dativar%d)
11299 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11300 realdat(v7d_tmp%voldatid(:,:,:,:,i,:), v7d_tmp%dativar%d(i))
11301 ENDDO
11302 DEALLOCATE(v7d_tmp%voldatid)
11303! trasferisco le variabili
11304 v7d_tmp%dativar%r => v7d_tmp%dativar%d
11305 NULLIFY(v7d_tmp%dativar%d)
11306 ENDIF
11307
11308! fondo con il volume definitivo
11309 CALL vol7d_merge(that, v7d_tmp)
11310ELSE
11312ENDIF
11313
11314
11315! Volume solo di dati interi
11316CALL vol7d_copy(this, v7d_tmp, &
11317 lanavarr=fv, lanavard=fv, lanavari=acn, lanavarb=fv, lanavarc=fv, &
11318 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11319 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11320 ldativarr=fv, ldativard=fv, ldativari=tv, ldativarb=fv, ldativarc=fv, &
11321 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11322 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11323
11324! converto a dati reali
11325IF (ASSOCIATED(v7d_tmp%anavar%i) .OR. ASSOCIATED(v7d_tmp%dativar%i)) THEN
11326
11327 IF (ASSOCIATED(v7d_tmp%anavar%i)) THEN
11328! alloco i dati reali e vi trasferisco gli interi
11329 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanai, 1), SIZE(v7d_tmp%volanai, 2), &
11330 SIZE(v7d_tmp%volanai, 3)))
11331 DO i = 1, SIZE(v7d_tmp%anavar%i)
11332 v7d_tmp%volanar(:,i,:) = &
11333 realdat(v7d_tmp%volanai(:,i,:), v7d_tmp%anavar%i(i))
11334 ENDDO
11335 DEALLOCATE(v7d_tmp%volanai)
11336! trasferisco le variabili
11337 v7d_tmp%anavar%r => v7d_tmp%anavar%i
11338 NULLIFY(v7d_tmp%anavar%i)
11339 ENDIF
11340
11341 IF (ASSOCIATED(v7d_tmp%dativar%i)) THEN
11342! alloco i dati reali e vi trasferisco gli interi
11343 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatii, 1), SIZE(v7d_tmp%voldatii, 2), &
11344 SIZE(v7d_tmp%voldatii, 3), SIZE(v7d_tmp%voldatii, 4), SIZE(v7d_tmp%voldatii, 5), &
11345 SIZE(v7d_tmp%voldatii, 6)))
11346 DO i = 1, SIZE(v7d_tmp%dativar%i)
11347 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11348 realdat(v7d_tmp%voldatii(:,:,:,:,i,:), v7d_tmp%dativar%i(i))
11349 ENDDO
11350 DEALLOCATE(v7d_tmp%voldatii)
11351! trasferisco le variabili
11352 v7d_tmp%dativar%r => v7d_tmp%dativar%i
11353 NULLIFY(v7d_tmp%dativar%i)
11354 ENDIF
11355
11356! fondo con il volume definitivo
11357 CALL vol7d_merge(that, v7d_tmp)
11358ELSE
11360ENDIF
11361
11362
11363! Volume solo di dati byte
11364CALL vol7d_copy(this, v7d_tmp, &
11365 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=acn, lanavarc=fv, &
11366 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11367 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11368 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=tv, ldativarc=fv, &
11369 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11370 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11371
11372! converto a dati reali
11373IF (ASSOCIATED(v7d_tmp%anavar%b) .OR. ASSOCIATED(v7d_tmp%dativar%b)) THEN
11374
11375 IF (ASSOCIATED(v7d_tmp%anavar%b)) THEN
11376! alloco i dati reali e vi trasferisco i byte
11377 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanab, 1), SIZE(v7d_tmp%volanab, 2), &
11378 SIZE(v7d_tmp%volanab, 3)))
11379 DO i = 1, SIZE(v7d_tmp%anavar%b)
11380 v7d_tmp%volanar(:,i,:) = &
11381 realdat(v7d_tmp%volanab(:,i,:), v7d_tmp%anavar%b(i))
11382 ENDDO
11383 DEALLOCATE(v7d_tmp%volanab)
11384! trasferisco le variabili
11385 v7d_tmp%anavar%r => v7d_tmp%anavar%b
11386 NULLIFY(v7d_tmp%anavar%b)
11387 ENDIF
11388
11389 IF (ASSOCIATED(v7d_tmp%dativar%b)) THEN
11390! alloco i dati reali e vi trasferisco i byte
11391 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatib, 1), SIZE(v7d_tmp%voldatib, 2), &
11392 SIZE(v7d_tmp%voldatib, 3), SIZE(v7d_tmp%voldatib, 4), SIZE(v7d_tmp%voldatib, 5), &
11393 SIZE(v7d_tmp%voldatib, 6)))
11394 DO i = 1, SIZE(v7d_tmp%dativar%b)
11395 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11396 realdat(v7d_tmp%voldatib(:,:,:,:,i,:), v7d_tmp%dativar%b(i))
11397 ENDDO
11398 DEALLOCATE(v7d_tmp%voldatib)
11399! trasferisco le variabili
11400 v7d_tmp%dativar%r => v7d_tmp%dativar%b
11401 NULLIFY(v7d_tmp%dativar%b)
11402 ENDIF
11403
11404! fondo con il volume definitivo
11405 CALL vol7d_merge(that, v7d_tmp)
11406ELSE
11408ENDIF
11409
11410
11411! Volume solo di dati character
11412CALL vol7d_copy(this, v7d_tmp, &
11413 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=fv, lanavarc=acn, &
11414 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11415 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11416 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=tv, &
11417 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11418 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11419
11420! converto a dati reali
11421IF (ASSOCIATED(v7d_tmp%anavar%c) .OR. ASSOCIATED(v7d_tmp%dativar%c)) THEN
11422
11423 IF (ASSOCIATED(v7d_tmp%anavar%c)) THEN
11424! alloco i dati reali e vi trasferisco i character
11425 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanac, 1), SIZE(v7d_tmp%volanac, 2), &
11426 SIZE(v7d_tmp%volanac, 3)))
11427 DO i = 1, SIZE(v7d_tmp%anavar%c)
11428 v7d_tmp%volanar(:,i,:) = &
11429 realdat(v7d_tmp%volanac(:,i,:), v7d_tmp%anavar%c(i))
11430 ENDDO
11431 DEALLOCATE(v7d_tmp%volanac)
11432! trasferisco le variabili
11433 v7d_tmp%anavar%r => v7d_tmp%anavar%c
11434 NULLIFY(v7d_tmp%anavar%c)
11435 ENDIF
11436
11437 IF (ASSOCIATED(v7d_tmp%dativar%c)) THEN
11438! alloco i dati reali e vi trasferisco i character
11439 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatic, 1), SIZE(v7d_tmp%voldatic, 2), &
11440 SIZE(v7d_tmp%voldatic, 3), SIZE(v7d_tmp%voldatic, 4), SIZE(v7d_tmp%voldatic, 5), &
11441 SIZE(v7d_tmp%voldatic, 6)))
11442 DO i = 1, SIZE(v7d_tmp%dativar%c)
11443 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11444 realdat(v7d_tmp%voldatic(:,:,:,:,i,:), v7d_tmp%dativar%c(i))
11445 ENDDO
11446 DEALLOCATE(v7d_tmp%voldatic)
11447! trasferisco le variabili
11448 v7d_tmp%dativar%r => v7d_tmp%dativar%c
11449 NULLIFY(v7d_tmp%dativar%c)
11450 ENDIF
11451
11452! fondo con il volume definitivo
11453 CALL vol7d_merge(that, v7d_tmp)
11454ELSE
11456ENDIF
11457
11458END SUBROUTINE vol7d_convr
11459
11460
11464SUBROUTINE vol7d_diff_only (this, that, data_only,ana)
11465TYPE(vol7d),INTENT(IN) :: this
11466TYPE(vol7d),INTENT(OUT) :: that
11467logical , optional, intent(in) :: data_only
11468logical , optional, intent(in) :: ana
11469logical :: ldata_only,lana
11470
11471IF (PRESENT(data_only)) THEN
11472 ldata_only = data_only
11473ELSE
11474 ldata_only = .false.
11475ENDIF
11476
11477IF (PRESENT(ana)) THEN
11478 lana = ana
11479ELSE
11480 lana = .false.
11481ENDIF
11482
11483
11484#undef VOL7D_POLY_ARRAY
11485#define VOL7D_POLY_ARRAY voldati
11486#include "vol7d_class_diff.F90"
11487#undef VOL7D_POLY_ARRAY
11488#define VOL7D_POLY_ARRAY voldatiattr
11489#include "vol7d_class_diff.F90"
11490#undef VOL7D_POLY_ARRAY
11491
11492if ( .not. ldata_only) then
11493
11494#define VOL7D_POLY_ARRAY volana
11495#include "vol7d_class_diff.F90"
11496#undef VOL7D_POLY_ARRAY
11497#define VOL7D_POLY_ARRAY volanaattr
11498#include "vol7d_class_diff.F90"
11499#undef VOL7D_POLY_ARRAY
11500
11501 if(lana)then
11502 where ( this%ana == that%ana )
11503 that%ana = vol7d_ana_miss
11504 end where
11505 end if
11506
11507end if
11508
11509
11510
11511END SUBROUTINE vol7d_diff_only
11512
11513
11514
11515! Creo le routine da ripetere per i vari tipi di dati di v7d
11516! tramite un template e il preprocessore
11517#undef VOL7D_POLY_TYPE
11518#undef VOL7D_POLY_TYPES
11519#define VOL7D_POLY_TYPE REAL
11520#define VOL7D_POLY_TYPES r
11521#include "vol7d_class_type_templ.F90"
11522#undef VOL7D_POLY_TYPE
11523#undef VOL7D_POLY_TYPES
11524#define VOL7D_POLY_TYPE DOUBLE PRECISION
11525#define VOL7D_POLY_TYPES d
11526#include "vol7d_class_type_templ.F90"
11527#undef VOL7D_POLY_TYPE
11528#undef VOL7D_POLY_TYPES
11529#define VOL7D_POLY_TYPE INTEGER
11530#define VOL7D_POLY_TYPES i
11531#include "vol7d_class_type_templ.F90"
11532#undef VOL7D_POLY_TYPE
11533#undef VOL7D_POLY_TYPES
11534#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
11535#define VOL7D_POLY_TYPES b
11536#include "vol7d_class_type_templ.F90"
11537#undef VOL7D_POLY_TYPE
11538#undef VOL7D_POLY_TYPES
11539#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
11540#define VOL7D_POLY_TYPES c
11541#include "vol7d_class_type_templ.F90"
11542
11543! Creo le routine da ripetere per i vari descrittori di dimensioni di v7d
11544! tramite un template e il preprocessore
11545#define VOL7D_SORT
11546#undef VOL7D_NO_ZERO_ALLOC
11547#undef VOL7D_POLY_TYPE
11548#define VOL7D_POLY_TYPE datetime
11549#include "vol7d_class_desc_templ.F90"
11550#undef VOL7D_POLY_TYPE
11551#define VOL7D_POLY_TYPE vol7d_timerange
11552#include "vol7d_class_desc_templ.F90"
11553#undef VOL7D_POLY_TYPE
11554#define VOL7D_POLY_TYPE vol7d_level
11555#include "vol7d_class_desc_templ.F90"
11556#undef VOL7D_SORT
11557#undef VOL7D_POLY_TYPE
11558#define VOL7D_POLY_TYPE vol7d_network
11559#include "vol7d_class_desc_templ.F90"
11560#undef VOL7D_POLY_TYPE
11561#define VOL7D_POLY_TYPE vol7d_ana
11562#include "vol7d_class_desc_templ.F90"
11563#define VOL7D_NO_ZERO_ALLOC
11564#undef VOL7D_POLY_TYPE
11565#define VOL7D_POLY_TYPE vol7d_var
11566#include "vol7d_class_desc_templ.F90"
11567
11577subroutine vol7d_write_on_file (this,unit,description,filename,filename_auto)
11578
11579TYPE(vol7d),INTENT(IN) :: this
11580integer,optional,intent(inout) :: unit
11581character(len=*),intent(in),optional :: filename
11582character(len=*),intent(out),optional :: filename_auto
11583character(len=*),INTENT(IN),optional :: description
11584
11585integer :: lunit
11586character(len=254) :: ldescription,arg,lfilename
11587integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
11588 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11589 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11590 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11591 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11592 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11593 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
11594!integer :: im,id,iy
11595integer :: tarray(8)
11596logical :: opened,exist
11597
11598 nana=0
11599 ntime=0
11600 ntimerange=0
11601 nlevel=0
11602 nnetwork=0
11603 ndativarr=0
11604 ndativari=0
11605 ndativarb=0
11606 ndativard=0
11607 ndativarc=0
11608 ndatiattrr=0
11609 ndatiattri=0
11610 ndatiattrb=0
11611 ndatiattrd=0
11612 ndatiattrc=0
11613 ndativarattrr=0
11614 ndativarattri=0
11615 ndativarattrb=0
11616 ndativarattrd=0
11617 ndativarattrc=0
11618 nanavarr=0
11619 nanavari=0
11620 nanavarb=0
11621 nanavard=0
11622 nanavarc=0
11623 nanaattrr=0
11624 nanaattri=0
11625 nanaattrb=0
11626 nanaattrd=0
11627 nanaattrc=0
11628 nanavarattrr=0
11629 nanavarattri=0
11630 nanavarattrb=0
11631 nanavarattrd=0
11632 nanavarattrc=0
11633
11634
11635!call idate(im,id,iy)
11636call date_and_time(values=tarray)
11637call getarg(0,arg)
11638
11639if (present(description))then
11640 ldescription=description
11641else
11642 ldescription="Vol7d generated by: "//trim(arg)
11643end if
11644
11645if (.not. present(unit))then
11646 lunit=getunit()
11647else
11648 if (unit==0)then
11649 lunit=getunit()
11650 unit=lunit
11651 else
11652 lunit=unit
11653 end if
11654end if
11655
11656lfilename=trim(arg)//".v7d"
11658
11659if (present(filename))then
11660 if (filename /= "")then
11661 lfilename=filename
11662 end if
11663end if
11664
11665if (present(filename_auto))filename_auto=lfilename
11666
11667
11668inquire(unit=lunit,opened=opened)
11669if (.not. opened) then
11670! inquire(file=lfilename, EXIST=exist)
11671! IF (exist) THEN
11672! CALL l4f_log(L4F_FATAL, &
11673! 'in vol7d_write_on_file, file exists, cannot open file '//TRIM(lfilename))
11674! CALL raise_fatal_error()
11675! ENDIF
11676 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM')
11677 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
11678end if
11679
11680if (associated(this%ana)) nana=size(this%ana)
11681if (associated(this%time)) ntime=size(this%time)
11682if (associated(this%timerange)) ntimerange=size(this%timerange)
11683if (associated(this%level)) nlevel=size(this%level)
11684if (associated(this%network)) nnetwork=size(this%network)
11685
11686if (associated(this%dativar%r)) ndativarr=size(this%dativar%r)
11687if (associated(this%dativar%i)) ndativari=size(this%dativar%i)
11688if (associated(this%dativar%b)) ndativarb=size(this%dativar%b)
11689if (associated(this%dativar%d)) ndativard=size(this%dativar%d)
11690if (associated(this%dativar%c)) ndativarc=size(this%dativar%c)
11691
11692if (associated(this%datiattr%r)) ndatiattrr=size(this%datiattr%r)
11693if (associated(this%datiattr%i)) ndatiattri=size(this%datiattr%i)
11694if (associated(this%datiattr%b)) ndatiattrb=size(this%datiattr%b)
11695if (associated(this%datiattr%d)) ndatiattrd=size(this%datiattr%d)
11696if (associated(this%datiattr%c)) ndatiattrc=size(this%datiattr%c)
11697
11698if (associated(this%dativarattr%r)) ndativarattrr=size(this%dativarattr%r)
11699if (associated(this%dativarattr%i)) ndativarattri=size(this%dativarattr%i)
11700if (associated(this%dativarattr%b)) ndativarattrb=size(this%dativarattr%b)
11701if (associated(this%dativarattr%d)) ndativarattrd=size(this%dativarattr%d)
11702if (associated(this%dativarattr%c)) ndativarattrc=size(this%dativarattr%c)
11703
11704if (associated(this%anavar%r)) nanavarr=size(this%anavar%r)
11705if (associated(this%anavar%i)) nanavari=size(this%anavar%i)
11706if (associated(this%anavar%b)) nanavarb=size(this%anavar%b)
11707if (associated(this%anavar%d)) nanavard=size(this%anavar%d)
11708if (associated(this%anavar%c)) nanavarc=size(this%anavar%c)
11709
11710if (associated(this%anaattr%r)) nanaattrr=size(this%anaattr%r)
11711if (associated(this%anaattr%i)) nanaattri=size(this%anaattr%i)
11712if (associated(this%anaattr%b)) nanaattrb=size(this%anaattr%b)
11713if (associated(this%anaattr%d)) nanaattrd=size(this%anaattr%d)
11714if (associated(this%anaattr%c)) nanaattrc=size(this%anaattr%c)
11715
11716if (associated(this%anavarattr%r)) nanavarattrr=size(this%anavarattr%r)
11717if (associated(this%anavarattr%i)) nanavarattri=size(this%anavarattr%i)
11718if (associated(this%anavarattr%b)) nanavarattrb=size(this%anavarattr%b)
11719if (associated(this%anavarattr%d)) nanavarattrd=size(this%anavarattr%d)
11720if (associated(this%anavarattr%c)) nanavarattrc=size(this%anavarattr%c)
11721
11722write(unit=lunit)ldescription
11723write(unit=lunit)tarray
11724
11725write(unit=lunit)&
11726 nana, ntime, ntimerange, nlevel, nnetwork, &
11727 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11728 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11729 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11730 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11731 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11732 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
11733 this%time_definition
11734
11735
11736!write(unit=lunit)this
11737
11738
11739!! prime 5 dimensioni
11742if (associated(this%level)) write(unit=lunit)this%level
11743if (associated(this%timerange)) write(unit=lunit)this%timerange
11744if (associated(this%network)) write(unit=lunit)this%network
11745
11746 !! 6a dimensione: variabile dell'anagrafica e dei dati
11747 !! con relativi attributi e in 5 tipi diversi
11748
11749if (associated(this%anavar%r)) write(unit=lunit)this%anavar%r
11750if (associated(this%anavar%i)) write(unit=lunit)this%anavar%i
11751if (associated(this%anavar%b)) write(unit=lunit)this%anavar%b
11752if (associated(this%anavar%d)) write(unit=lunit)this%anavar%d
11753if (associated(this%anavar%c)) write(unit=lunit)this%anavar%c
11754
11755if (associated(this%anaattr%r)) write(unit=lunit)this%anaattr%r
11756if (associated(this%anaattr%i)) write(unit=lunit)this%anaattr%i
11757if (associated(this%anaattr%b)) write(unit=lunit)this%anaattr%b
11758if (associated(this%anaattr%d)) write(unit=lunit)this%anaattr%d
11759if (associated(this%anaattr%c)) write(unit=lunit)this%anaattr%c
11760
11761if (associated(this%anavarattr%r)) write(unit=lunit)this%anavarattr%r
11762if (associated(this%anavarattr%i)) write(unit=lunit)this%anavarattr%i
11763if (associated(this%anavarattr%b)) write(unit=lunit)this%anavarattr%b
11764if (associated(this%anavarattr%d)) write(unit=lunit)this%anavarattr%d
11765if (associated(this%anavarattr%c)) write(unit=lunit)this%anavarattr%c
11766
11767if (associated(this%dativar%r)) write(unit=lunit)this%dativar%r
11768if (associated(this%dativar%i)) write(unit=lunit)this%dativar%i
11769if (associated(this%dativar%b)) write(unit=lunit)this%dativar%b
11770if (associated(this%dativar%d)) write(unit=lunit)this%dativar%d
11771if (associated(this%dativar%c)) write(unit=lunit)this%dativar%c
11772
11773if (associated(this%datiattr%r)) write(unit=lunit)this%datiattr%r
11774if (associated(this%datiattr%i)) write(unit=lunit)this%datiattr%i
11775if (associated(this%datiattr%b)) write(unit=lunit)this%datiattr%b
11776if (associated(this%datiattr%d)) write(unit=lunit)this%datiattr%d
11777if (associated(this%datiattr%c)) write(unit=lunit)this%datiattr%c
11778
11779if (associated(this%dativarattr%r)) write(unit=lunit)this%dativarattr%r
11780if (associated(this%dativarattr%i)) write(unit=lunit)this%dativarattr%i
11781if (associated(this%dativarattr%b)) write(unit=lunit)this%dativarattr%b
11782if (associated(this%dativarattr%d)) write(unit=lunit)this%dativarattr%d
11783if (associated(this%dativarattr%c)) write(unit=lunit)this%dativarattr%c
11784
11785!! Volumi di valori e attributi per anagrafica e dati
11786
11787if (associated(this%volanar)) write(unit=lunit)this%volanar
11788if (associated(this%volanaattrr)) write(unit=lunit)this%volanaattrr
11789if (associated(this%voldatir)) write(unit=lunit)this%voldatir
11790if (associated(this%voldatiattrr)) write(unit=lunit)this%voldatiattrr
11791
11792if (associated(this%volanai)) write(unit=lunit)this%volanai
11793if (associated(this%volanaattri)) write(unit=lunit)this%volanaattri
11794if (associated(this%voldatii)) write(unit=lunit)this%voldatii
11795if (associated(this%voldatiattri)) write(unit=lunit)this%voldatiattri
11796
11797if (associated(this%volanab)) write(unit=lunit)this%volanab
11798if (associated(this%volanaattrb)) write(unit=lunit)this%volanaattrb
11799if (associated(this%voldatib)) write(unit=lunit)this%voldatib
11800if (associated(this%voldatiattrb)) write(unit=lunit)this%voldatiattrb
11801
11802if (associated(this%volanad)) write(unit=lunit)this%volanad
11803if (associated(this%volanaattrd)) write(unit=lunit)this%volanaattrd
11804if (associated(this%voldatid)) write(unit=lunit)this%voldatid
11805if (associated(this%voldatiattrd)) write(unit=lunit)this%voldatiattrd
11806
11807if (associated(this%volanac)) write(unit=lunit)this%volanac
11808if (associated(this%volanaattrc)) write(unit=lunit)this%volanaattrc
11809if (associated(this%voldatic)) write(unit=lunit)this%voldatic
11810if (associated(this%voldatiattrc)) write(unit=lunit)this%voldatiattrc
11811
11812if (.not. present(unit)) close(unit=lunit)
11813
11814end subroutine vol7d_write_on_file
11815
11816
11823
11824
11825subroutine vol7d_read_from_file (this,unit,filename,description,tarray,filename_auto)
11826
11827TYPE(vol7d),INTENT(OUT) :: this
11828integer,intent(inout),optional :: unit
11829character(len=*),INTENT(in),optional :: filename
11830character(len=*),intent(out),optional :: filename_auto
11831character(len=*),INTENT(out),optional :: description
11832integer,intent(out),optional :: tarray(8)
11833
11834
11835integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
11836 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11837 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11838 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11839 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11840 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11841 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
11842
11843character(len=254) :: ldescription,lfilename,arg
11844integer :: ltarray(8),lunit,ios
11845logical :: opened,exist
11846
11847
11848call getarg(0,arg)
11849
11850if (.not. present(unit))then
11851 lunit=getunit()
11852else
11853 if (unit==0)then
11854 lunit=getunit()
11855 unit=lunit
11856 else
11857 lunit=unit
11858 end if
11859end if
11860
11861lfilename=trim(arg)//".v7d"
11863
11864if (present(filename))then
11865 if (filename /= "")then
11866 lfilename=filename
11867 end if
11868end if
11869
11870if (present(filename_auto))filename_auto=lfilename
11871
11872
11873inquire(unit=lunit,opened=opened)
11874IF (.NOT. opened) THEN
11875 inquire(file=lfilename,exist=exist)
11876 IF (.NOT.exist) THEN
11877 CALL l4f_log(l4f_fatal, &
11878 'in vol7d_read_from_file, file does not exists, cannot open')
11879 CALL raise_fatal_error()
11880 ENDIF
11881 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM', &
11882 status='OLD', action='READ')
11883 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
11884end if
11885
11886
11888read(unit=lunit,iostat=ios)ldescription
11889
11890if (ios < 0) then ! A negative value indicates that the End of File or End of Record
11891 call vol7d_alloc (this)
11892 call vol7d_alloc_vol (this)
11893 if (present(description))description=ldescription
11894 if (present(tarray))tarray=ltarray
11895 if (.not. present(unit)) close(unit=lunit)
11896end if
11897
11898read(unit=lunit)ltarray
11899
11900CALL l4f_log(l4f_info, 'Reading vol7d from file')
11901CALL l4f_log(l4f_info, 'description: '//trim(ldescription))
11904
11905if (present(description))description=ldescription
11906if (present(tarray))tarray=ltarray
11907
11908read(unit=lunit)&
11909 nana, ntime, ntimerange, nlevel, nnetwork, &
11910 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11911 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11912 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11913 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11914 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11915 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
11916 this%time_definition
11917
11918call vol7d_alloc (this, &
11919 nana=nana, ntime=ntime, ntimerange=ntimerange, nlevel=nlevel, nnetwork=nnetwork,&
11920 ndativarr=ndativarr, ndativari=ndativari, ndativarb=ndativarb,&
11921 ndativard=ndativard, ndativarc=ndativarc,&
11922 ndatiattrr=ndatiattrr, ndatiattri=ndatiattri, ndatiattrb=ndatiattrb,&
11923 ndatiattrd=ndatiattrd, ndatiattrc=ndatiattrc,&
11924 ndativarattrr=ndativarattrr, ndativarattri=ndativarattri, ndativarattrb=ndativarattrb, &
11925 ndativarattrd=ndativarattrd, ndativarattrc=ndativarattrc,&
11926 nanavarr=nanavarr, nanavari=nanavari, nanavarb=nanavarb, &
11927 nanavard=nanavard, nanavarc=nanavarc,&
11928 nanaattrr=nanaattrr, nanaattri=nanaattri, nanaattrb=nanaattrb,&
11929 nanaattrd=nanaattrd, nanaattrc=nanaattrc,&
11930 nanavarattrr=nanavarattrr, nanavarattri=nanavarattri, nanavarattrb=nanavarattrb, &
11931 nanavarattrd=nanavarattrd, nanavarattrc=nanavarattrc)
11932
11933
11936if (associated(this%level)) read(unit=lunit)this%level
11937if (associated(this%timerange)) read(unit=lunit)this%timerange
11938if (associated(this%network)) read(unit=lunit)this%network
11939
11940if (associated(this%anavar%r)) read(unit=lunit)this%anavar%r
11941if (associated(this%anavar%i)) read(unit=lunit)this%anavar%i
11942if (associated(this%anavar%b)) read(unit=lunit)this%anavar%b
11943if (associated(this%anavar%d)) read(unit=lunit)this%anavar%d
11944if (associated(this%anavar%c)) read(unit=lunit)this%anavar%c
11945
11946if (associated(this%anaattr%r)) read(unit=lunit)this%anaattr%r
11947if (associated(this%anaattr%i)) read(unit=lunit)this%anaattr%i
11948if (associated(this%anaattr%b)) read(unit=lunit)this%anaattr%b
11949if (associated(this%anaattr%d)) read(unit=lunit)this%anaattr%d
11950if (associated(this%anaattr%c)) read(unit=lunit)this%anaattr%c
11951
11952if (associated(this%anavarattr%r)) read(unit=lunit)this%anavarattr%r
11953if (associated(this%anavarattr%i)) read(unit=lunit)this%anavarattr%i
11954if (associated(this%anavarattr%b)) read(unit=lunit)this%anavarattr%b
11955if (associated(this%anavarattr%d)) read(unit=lunit)this%anavarattr%d
11956if (associated(this%anavarattr%c)) read(unit=lunit)this%anavarattr%c
11957
11958if (associated(this%dativar%r)) read(unit=lunit)this%dativar%r
11959if (associated(this%dativar%i)) read(unit=lunit)this%dativar%i
11960if (associated(this%dativar%b)) read(unit=lunit)this%dativar%b
11961if (associated(this%dativar%d)) read(unit=lunit)this%dativar%d
11962if (associated(this%dativar%c)) read(unit=lunit)this%dativar%c
11963
11964if (associated(this%datiattr%r)) read(unit=lunit)this%datiattr%r
11965if (associated(this%datiattr%i)) read(unit=lunit)this%datiattr%i
11966if (associated(this%datiattr%b)) read(unit=lunit)this%datiattr%b
11967if (associated(this%datiattr%d)) read(unit=lunit)this%datiattr%d
11968if (associated(this%datiattr%c)) read(unit=lunit)this%datiattr%c
11969
11970if (associated(this%dativarattr%r)) read(unit=lunit)this%dativarattr%r
11971if (associated(this%dativarattr%i)) read(unit=lunit)this%dativarattr%i
11972if (associated(this%dativarattr%b)) read(unit=lunit)this%dativarattr%b
11973if (associated(this%dativarattr%d)) read(unit=lunit)this%dativarattr%d
11974if (associated(this%dativarattr%c)) read(unit=lunit)this%dativarattr%c
11975
11976call vol7d_alloc_vol (this)
11977
11978!! Volumi di valori e attributi per anagrafica e dati
11979
11980if (associated(this%volanar)) read(unit=lunit)this%volanar
11981if (associated(this%volanaattrr)) read(unit=lunit)this%volanaattrr
11982if (associated(this%voldatir)) read(unit=lunit)this%voldatir
11983if (associated(this%voldatiattrr)) read(unit=lunit)this%voldatiattrr
11984
11985if (associated(this%volanai)) read(unit=lunit)this%volanai
11986if (associated(this%volanaattri)) read(unit=lunit)this%volanaattri
11987if (associated(this%voldatii)) read(unit=lunit)this%voldatii
11988if (associated(this%voldatiattri)) read(unit=lunit)this%voldatiattri
11989
11990if (associated(this%volanab)) read(unit=lunit)this%volanab
11991if (associated(this%volanaattrb)) read(unit=lunit)this%volanaattrb
11992if (associated(this%voldatib)) read(unit=lunit)this%voldatib
11993if (associated(this%voldatiattrb)) read(unit=lunit)this%voldatiattrb
11994
11995if (associated(this%volanad)) read(unit=lunit)this%volanad
11996if (associated(this%volanaattrd)) read(unit=lunit)this%volanaattrd
11997if (associated(this%voldatid)) read(unit=lunit)this%voldatid
11998if (associated(this%voldatiattrd)) read(unit=lunit)this%voldatiattrd
11999
12000if (associated(this%volanac)) read(unit=lunit)this%volanac
12001if (associated(this%volanaattrc)) read(unit=lunit)this%volanaattrc
12002if (associated(this%voldatic)) read(unit=lunit)this%voldatic
12003if (associated(this%voldatiattrc)) read(unit=lunit)this%voldatiattrc
12004
12005if (.not. present(unit)) close(unit=lunit)
12006
12007end subroutine vol7d_read_from_file
12008
12009
12010! to double precision
12011elemental doubleprecision function doubledatd(voldat,var)
12012doubleprecision,intent(in) :: voldat
12013type(vol7d_var),intent(in) :: var
12014
12015doubledatd=voldat
12016
12017end function doubledatd
12018
12019
12020elemental doubleprecision function doubledatr(voldat,var)
12021real,intent(in) :: voldat
12022type(vol7d_var),intent(in) :: var
12023
12025 doubledatr=dble(voldat)
12026else
12027 doubledatr=dmiss
12028end if
12029
12030end function doubledatr
12031
12032
12033elemental doubleprecision function doubledati(voldat,var)
12034integer,intent(in) :: voldat
12035type(vol7d_var),intent(in) :: var
12036
12039 doubledati=dble(voldat)/10.d0**var%scalefactor
12040 else
12041 doubledati=dble(voldat)
12042 endif
12043else
12044 doubledati=dmiss
12045end if
12046
12047end function doubledati
12048
12049
12050elemental doubleprecision function doubledatb(voldat,var)
12051integer(kind=int_b),intent(in) :: voldat
12052type(vol7d_var),intent(in) :: var
12053
12056 doubledatb=dble(voldat)/10.d0**var%scalefactor
12057 else
12058 doubledatb=dble(voldat)
12059 endif
12060else
12061 doubledatb=dmiss
12062end if
12063
12064end function doubledatb
12065
12066
12067elemental doubleprecision function doubledatc(voldat,var)
12068CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
12069type(vol7d_var),intent(in) :: var
12070
12071doubledatc = c2d(voldat)
12073 doubledatc=doubledatc/10.d0**var%scalefactor
12074end if
12075
12076end function doubledatc
12077
12078
12079! to integer
12080elemental integer function integerdatd(voldat,var)
12081doubleprecision,intent(in) :: voldat
12082type(vol7d_var),intent(in) :: var
12083
12086 integerdatd=nint(voldat*10d0**var%scalefactor)
12087 else
12088 integerdatd=nint(voldat)
12089 endif
12090else
12091 integerdatd=imiss
12092end if
12093
12094end function integerdatd
12095
12096
12097elemental integer function integerdatr(voldat,var)
12098real,intent(in) :: voldat
12099type(vol7d_var),intent(in) :: var
12100
12103 integerdatr=nint(voldat*10d0**var%scalefactor)
12104 else
12105 integerdatr=nint(voldat)
12106 endif
12107else
12108 integerdatr=imiss
12109end if
12110
12111end function integerdatr
12112
12113
12114elemental integer function integerdati(voldat,var)
12115integer,intent(in) :: voldat
12116type(vol7d_var),intent(in) :: var
12117
12118integerdati=voldat
12119
12120end function integerdati
12121
12122
12123elemental integer function integerdatb(voldat,var)
12124integer(kind=int_b),intent(in) :: voldat
12125type(vol7d_var),intent(in) :: var
12126
12128 integerdatb=voldat
12129else
12130 integerdatb=imiss
12131end if
12132
12133end function integerdatb
12134
12135
12136elemental integer function integerdatc(voldat,var)
12137CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
12138type(vol7d_var),intent(in) :: var
12139
12140integerdatc=c2i(voldat)
12141
12142end function integerdatc
12143
12144
12145! to real
12146elemental real function realdatd(voldat,var)
12147doubleprecision,intent(in) :: voldat
12148type(vol7d_var),intent(in) :: var
12149
12151 realdatd=real(voldat)
12152else
12153 realdatd=rmiss
12154end if
12155
12156end function realdatd
12157
12158
12159elemental real function realdatr(voldat,var)
12160real,intent(in) :: voldat
12161type(vol7d_var),intent(in) :: var
12162
12163realdatr=voldat
12164
12165end function realdatr
12166
12167
12168elemental real function realdati(voldat,var)
12169integer,intent(in) :: voldat
12170type(vol7d_var),intent(in) :: var
12171
12174 realdati=float(voldat)/10.**var%scalefactor
12175 else
12176 realdati=float(voldat)
12177 endif
12178else
12179 realdati=rmiss
12180end if
12181
12182end function realdati
12183
12184
12185elemental real function realdatb(voldat,var)
12186integer(kind=int_b),intent(in) :: voldat
12187type(vol7d_var),intent(in) :: var
12188
12191 realdatb=float(voldat)/10**var%scalefactor
12192 else
12193 realdatb=float(voldat)
12194 endif
12195else
12196 realdatb=rmiss
12197end if
12198
12199end function realdatb
12200
12201
12202elemental real function realdatc(voldat,var)
12203CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
12204type(vol7d_var),intent(in) :: var
12205
12206realdatc=c2r(voldat)
12208 realdatc=realdatc/10.**var%scalefactor
12209end if
12210
12211end function realdatc
12212
12213
12219FUNCTION realanavol(this, var) RESULT(vol)
12220TYPE(vol7d),INTENT(in) :: this
12221TYPE(vol7d_var),INTENT(in) :: var
12222REAL :: vol(SIZE(this%ana),size(this%network))
12223
12224CHARACTER(len=1) :: dtype
12225INTEGER :: indvar
12226
12227dtype = cmiss
12228indvar = index(this%anavar, var, type=dtype)
12229
12230IF (indvar > 0) THEN
12231 SELECT CASE (dtype)
12232 CASE("d")
12233 vol = realdat(this%volanad(:,indvar,:), var)
12234 CASE("r")
12235 vol = this%volanar(:,indvar,:)
12236 CASE("i")
12237 vol = realdat(this%volanai(:,indvar,:), var)
12238 CASE("b")
12239 vol = realdat(this%volanab(:,indvar,:), var)
12240 CASE("c")
12241 vol = realdat(this%volanac(:,indvar,:), var)
12242 CASE default
12243 vol = rmiss
12244 END SELECT
12245ELSE
12246 vol = rmiss
12247ENDIF
12248
12249END FUNCTION realanavol
12250
12251
12257FUNCTION integeranavol(this, var) RESULT(vol)
12258TYPE(vol7d),INTENT(in) :: this
12259TYPE(vol7d_var),INTENT(in) :: var
12260INTEGER :: vol(SIZE(this%ana),size(this%network))
12261
12262CHARACTER(len=1) :: dtype
12263INTEGER :: indvar
12264
12265dtype = cmiss
12266indvar = index(this%anavar, var, type=dtype)
12267
12268IF (indvar > 0) THEN
12269 SELECT CASE (dtype)
12270 CASE("d")
12271 vol = integerdat(this%volanad(:,indvar,:), var)
12272 CASE("r")
12273 vol = integerdat(this%volanar(:,indvar,:), var)
12274 CASE("i")
12275 vol = this%volanai(:,indvar,:)
12276 CASE("b")
12277 vol = integerdat(this%volanab(:,indvar,:), var)
12278 CASE("c")
12279 vol = integerdat(this%volanac(:,indvar,:), var)
12280 CASE default
12281 vol = imiss
12282 END SELECT
12283ELSE
12284 vol = imiss
12285ENDIF
12286
12287END FUNCTION integeranavol
12288
12289
12295subroutine move_datac (v7d,&
12296 indana,indtime,indlevel,indtimerange,indnetwork,&
12297 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
12298
12299TYPE(vol7d),intent(inout) :: v7d
12300
12301integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
12302integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
12303integer :: inddativar,inddativarattr
12304
12305
12306do inddativar=1,size(v7d%dativar%c)
12307
12309 .not. c_e(v7d%voldatic(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
12310 ) then
12311
12312 ! dati
12313 v7d%voldatic &
12314 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
12315 v7d%voldatic &
12316 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
12317
12318
12319 ! attributi
12320 if (associated (v7d%dativarattr%i)) then
12321 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%c(inddativar))
12322 if (inddativarattr > 0 ) then
12323 v7d%voldatiattri &
12324 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12325 v7d%voldatiattri &
12326 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12327 end if
12328 end if
12329
12330 if (associated (v7d%dativarattr%r)) then
12331 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%c(inddativar))
12332 if (inddativarattr > 0 ) then
12333 v7d%voldatiattrr &
12334 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12335 v7d%voldatiattrr &
12336 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12337 end if
12338 end if
12339
12340 if (associated (v7d%dativarattr%d)) then
12341 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%c(inddativar))
12342 if (inddativarattr > 0 ) then
12343 v7d%voldatiattrd &
12344 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12345 v7d%voldatiattrd &
12346 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12347 end if
12348 end if
12349
12350 if (associated (v7d%dativarattr%b)) then
12351 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%c(inddativar))
12352 if (inddativarattr > 0 ) then
12353 v7d%voldatiattrb &
12354 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12355 v7d%voldatiattrb &
12356 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12357 end if
12358 end if
12359
12360 if (associated (v7d%dativarattr%c)) then
12361 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%c(inddativar))
12362 if (inddativarattr > 0 ) then
12363 v7d%voldatiattrc &
12364 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12365 v7d%voldatiattrc &
12366 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12367 end if
12368 end if
12369
12370 end if
12371
12372end do
12373
12374end subroutine move_datac
12375
12381subroutine move_datar (v7d,&
12382 indana,indtime,indlevel,indtimerange,indnetwork,&
12383 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
12384
12385TYPE(vol7d),intent(inout) :: v7d
12386
12387integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
12388integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
12389integer :: inddativar,inddativarattr
12390
12391
12392do inddativar=1,size(v7d%dativar%r)
12393
12395 .not. c_e(v7d%voldatir(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
12396 ) then
12397
12398 ! dati
12399 v7d%voldatir &
12400 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
12401 v7d%voldatir &
12402 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
12403
12404
12405 ! attributi
12406 if (associated (v7d%dativarattr%i)) then
12407 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%r(inddativar))
12408 if (inddativarattr > 0 ) then
12409 v7d%voldatiattri &
12410 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12411 v7d%voldatiattri &
12412 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12413 end if
12414 end if
12415
12416 if (associated (v7d%dativarattr%r)) then
12417 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%r(inddativar))
12418 if (inddativarattr > 0 ) then
12419 v7d%voldatiattrr &
12420 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12421 v7d%voldatiattrr &
12422 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12423 end if
12424 end if
12425
12426 if (associated (v7d%dativarattr%d)) then
12427 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%r(inddativar))
12428 if (inddativarattr > 0 ) then
12429 v7d%voldatiattrd &
12430 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12431 v7d%voldatiattrd &
12432 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12433 end if
12434 end if
12435
12436 if (associated (v7d%dativarattr%b)) then
12437 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%r(inddativar))
12438 if (inddativarattr > 0 ) then
12439 v7d%voldatiattrb &
12440 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12441 v7d%voldatiattrb &
12442 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12443 end if
12444 end if
12445
12446 if (associated (v7d%dativarattr%c)) then
12447 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%r(inddativar))
12448 if (inddativarattr > 0 ) then
12449 v7d%voldatiattrc &
12450 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12451 v7d%voldatiattrc &
12452 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12453 end if
12454 end if
12455
12456 end if
12457
12458end do
12459
12460end subroutine move_datar
12461
12462
12476subroutine v7d_rounding(v7din,v7dout,level,timerange,nostatproc)
12477type(vol7d),intent(inout) :: v7din
12478type(vol7d),intent(out) :: v7dout
12479type(vol7d_level),intent(in),optional :: level(:)
12480type(vol7d_timerange),intent(in),optional :: timerange(:)
12481!logical,intent(in),optional :: merge !< if there are data on more then one almost equal levels or timeranges
12482!! will be merged POINT BY POINT with priority for the fird data found ordered by icreasing var index
12483logical,intent(in),optional :: nostatproc
12484
12485integer :: nana,nlevel,ntime,ntimerange,nnetwork,nbin
12486integer :: iana,ilevel,itimerange,indl,indt,itime,inetwork
12487type(vol7d_level) :: roundlevel(size(v7din%level))
12488type(vol7d_timerange) :: roundtimerange(size(v7din%timerange))
12489type(vol7d) :: v7d_tmp
12490
12491
12492nbin=0
12493
12494if (associated(v7din%dativar%r)) nbin = nbin + size(v7din%dativar%r)
12495if (associated(v7din%dativar%i)) nbin = nbin + size(v7din%dativar%i)
12496if (associated(v7din%dativar%d)) nbin = nbin + size(v7din%dativar%d)
12497if (associated(v7din%dativar%b)) nbin = nbin + size(v7din%dativar%b)
12498
12500
12501roundlevel=v7din%level
12502
12503if (present(level))then
12504 do ilevel = 1, size(v7din%level)
12505 if ((any(v7din%level(ilevel) .almosteq. level))) then
12506 roundlevel(ilevel)=level(1)
12507 end if
12508 end do
12509end if
12510
12511roundtimerange=v7din%timerange
12512
12513if (present(timerange))then
12514 do itimerange = 1, size(v7din%timerange)
12515 if ((any(v7din%timerange(itimerange) .almosteq. timerange))) then
12516 roundtimerange(itimerange)=timerange(1)
12517 end if
12518 end do
12519end if
12520
12521!set istantaneous values everywere
12522!preserve p1 for forecast time
12523if (optio_log(nostatproc)) then
12524 roundtimerange(:)%timerange=254
12525 roundtimerange(:)%p2=0
12526end if
12527
12528
12529nana=size(v7din%ana)
12530nlevel=count_distinct(roundlevel,back=.true.)
12531ntime=size(v7din%time)
12532ntimerange=count_distinct(roundtimerange,back=.true.)
12533nnetwork=size(v7din%network)
12534
12536
12537if (nbin == 0) then
12539else
12540 call vol7d_convr(v7din,v7d_tmp)
12541end if
12542
12543v7d_tmp%level=roundlevel
12544v7d_tmp%timerange=roundtimerange
12545
12546do ilevel=1, size(v7d_tmp%level)
12547 indl=index(v7d_tmp%level,roundlevel(ilevel))
12548 do itimerange=1,size(v7d_tmp%timerange)
12549 indt=index(v7d_tmp%timerange,roundtimerange(itimerange))
12550
12551 if (indl /= ilevel .or. indt /= itimerange) then
12552
12553 do iana=1, nana
12554 do itime=1,ntime
12555 do inetwork=1,nnetwork
12556
12557 if (nbin > 0) then
12558 call move_datar (v7d_tmp,&
12559 iana,itime,ilevel,itimerange,inetwork,&
12560 iana,itime,indl,indt,inetwork)
12561 else
12562 call move_datac (v7d_tmp,&
12563 iana,itime,ilevel,itimerange,inetwork,&
12564 iana,itime,indl,indt,inetwork)
12565 end if
12566
12567 end do
12568 end do
12569 end do
12570
12571 end if
12572
12573 end do
12574end do
12575
12576! set to missing level and time > nlevel
12577do ilevel=nlevel+1,size(v7d_tmp%level)
12579end do
12580
12581do itimerange=ntimerange+1,size(v7d_tmp%timerange)
12583end do
12584
12585!copy with remove
12588
12589!call display(v7dout)
12590
12591end subroutine v7d_rounding
12592
12593
12595
12601
12602
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 |