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