libsim Versione 7.1.11
|
◆ realanavol()
Return an ana volume of a requested variable as real data. It returns a 2-d array of the proper shape (ana x network) for the ana variable requested, converted to real type. If the conversion fails or if the variable is not contained in the ana volume, missing data are returned.
Definizione alla linea 9116 del file vol7d_class.F90. 9117! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
9118! authors:
9119! Davide Cesari <dcesari@arpa.emr.it>
9120! Paolo Patruno <ppatruno@arpa.emr.it>
9121
9122! This program is free software; you can redistribute it and/or
9123! modify it under the terms of the GNU General Public License as
9124! published by the Free Software Foundation; either version 2 of
9125! the License, or (at your option) any later version.
9126
9127! This program is distributed in the hope that it will be useful,
9128! but WITHOUT ANY WARRANTY; without even the implied warranty of
9129! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
9130! GNU General Public License for more details.
9131
9132! You should have received a copy of the GNU General Public License
9133! along with this program. If not, see <http://www.gnu.org/licenses/>.
9134#include "config.h"
9135
9147
9215IMPLICIT NONE
9216
9217
9218INTEGER, PARAMETER :: vol7d_maxdim_a = 3, vol7d_maxdim_aa = 4, &
9219 vol7d_maxdim_d = 6, vol7d_maxdim_ad = 7
9220
9221INTEGER, PARAMETER :: vol7d_ana_a=1
9222INTEGER, PARAMETER :: vol7d_var_a=2
9223INTEGER, PARAMETER :: vol7d_network_a=3
9224INTEGER, PARAMETER :: vol7d_attr_a=4
9225INTEGER, PARAMETER :: vol7d_ana_d=1
9226INTEGER, PARAMETER :: vol7d_time_d=2
9227INTEGER, PARAMETER :: vol7d_level_d=3
9228INTEGER, PARAMETER :: vol7d_timerange_d=4
9229INTEGER, PARAMETER :: vol7d_var_d=5
9230INTEGER, PARAMETER :: vol7d_network_d=6
9231INTEGER, PARAMETER :: vol7d_attr_d=7
9232INTEGER, PARAMETER :: vol7d_cdatalen=32
9233
9234TYPE vol7d_varmap
9235 INTEGER :: r, d, i, b, c
9236END TYPE vol7d_varmap
9237
9242 TYPE(vol7d_ana),POINTER :: ana(:)
9244 TYPE(datetime),POINTER :: time(:)
9246 TYPE(vol7d_level),POINTER :: level(:)
9248 TYPE(vol7d_timerange),POINTER :: timerange(:)
9250 TYPE(vol7d_network),POINTER :: network(:)
9252 TYPE(vol7d_varvect) :: anavar
9254 TYPE(vol7d_varvect) :: anaattr
9256 TYPE(vol7d_varvect) :: anavarattr
9258 TYPE(vol7d_varvect) :: dativar
9260 TYPE(vol7d_varvect) :: datiattr
9262 TYPE(vol7d_varvect) :: dativarattr
9263
9265 REAL,POINTER :: volanar(:,:,:)
9267 DOUBLE PRECISION,POINTER :: volanad(:,:,:)
9269 INTEGER,POINTER :: volanai(:,:,:)
9271 INTEGER(kind=int_b),POINTER :: volanab(:,:,:)
9273 CHARACTER(len=vol7d_cdatalen),POINTER :: volanac(:,:,:)
9274
9276 REAL,POINTER :: volanaattrr(:,:,:,:)
9278 DOUBLE PRECISION,POINTER :: volanaattrd(:,:,:,:)
9280 INTEGER,POINTER :: volanaattri(:,:,:,:)
9282 INTEGER(kind=int_b),POINTER :: volanaattrb(:,:,:,:)
9284 CHARACTER(len=vol7d_cdatalen),POINTER :: volanaattrc(:,:,:,:)
9285
9287 REAL,POINTER :: voldatir(:,:,:,:,:,:) ! sono i dati
9289 DOUBLE PRECISION,POINTER :: voldatid(:,:,:,:,:,:)
9291 INTEGER,POINTER :: voldatii(:,:,:,:,:,:)
9293 INTEGER(kind=int_b),POINTER :: voldatib(:,:,:,:,:,:)
9295 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatic(:,:,:,:,:,:)
9296
9298 REAL,POINTER :: voldatiattrr(:,:,:,:,:,:,:)
9300 DOUBLE PRECISION,POINTER :: voldatiattrd(:,:,:,:,:,:,:)
9302 INTEGER,POINTER :: voldatiattri(:,:,:,:,:,:,:)
9304 INTEGER(kind=int_b),POINTER :: voldatiattrb(:,:,:,:,:,:,:)
9306 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatiattrc(:,:,:,:,:,:,:)
9307
9309 integer :: time_definition
9310
9312
9317 MODULE PROCEDURE vol7d_init
9318END INTERFACE
9319
9322 MODULE PROCEDURE vol7d_delete
9323END INTERFACE
9324
9327 MODULE PROCEDURE vol7d_write_on_file
9328END INTERFACE
9329
9331INTERFACE import
9332 MODULE PROCEDURE vol7d_read_from_file
9333END INTERFACE
9334
9337 MODULE PROCEDURE vol7d_display, dat_display, dat_vect_display
9338END INTERFACE
9339
9342 MODULE PROCEDURE to_char_dat
9343END INTERFACE
9344
9347 MODULE PROCEDURE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
9348END INTERFACE
9349
9352 MODULE PROCEDURE realdatd,realdatr,realdati,realdatb,realdatc
9353END INTERFACE
9354
9357 MODULE PROCEDURE integerdatd,integerdatr,integerdati,integerdatb,integerdatc
9358END INTERFACE
9359
9362 MODULE PROCEDURE vol7d_copy
9363END INTERFACE
9364
9367 MODULE PROCEDURE vol7d_c_e
9368END INTERFACE
9369
9374 MODULE PROCEDURE vol7d_check
9375END INTERFACE
9376
9391 MODULE PROCEDURE v7d_rounding
9392END INTERFACE
9393
9394!!$INTERFACE get_volana
9395!!$ MODULE PROCEDURE vol7d_get_volanar, vol7d_get_volanad, vol7d_get_volanai, &
9396!!$ vol7d_get_volanab, vol7d_get_volanac
9397!!$END INTERFACE
9398!!$
9399!!$INTERFACE get_voldati
9400!!$ MODULE PROCEDURE vol7d_get_voldatir, vol7d_get_voldatid, vol7d_get_voldatii, &
9401!!$ vol7d_get_voldatib, vol7d_get_voldatic
9402!!$END INTERFACE
9403!!$
9404!!$INTERFACE get_volanaattr
9405!!$ MODULE PROCEDURE vol7d_get_volanaattrr, vol7d_get_volanaattrd, &
9406!!$ vol7d_get_volanaattri, vol7d_get_volanaattrb, vol7d_get_volanaattrc
9407!!$END INTERFACE
9408!!$
9409!!$INTERFACE get_voldatiattr
9410!!$ MODULE PROCEDURE vol7d_get_voldatiattrr, vol7d_get_voldatiattrd, &
9411!!$ vol7d_get_voldatiattri, vol7d_get_voldatiattrb, vol7d_get_voldatiattrc
9412!!$END INTERFACE
9413
9414PRIVATE vol7d_get_volr, vol7d_get_vold, vol7d_get_voli, vol7d_get_volb, &
9415 vol7d_get_volc, &
9416 volptr1dr, volptr2dr, volptr3dr, volptr4dr, volptr5dr, volptr6dr, volptr7dr, &
9417 volptr1dd, volptr2dd, volptr3dd, volptr4dd, volptr5dd, volptr6dd, volptr7dd, &
9418 volptr1di, volptr2di, volptr3di, volptr4di, volptr5di, volptr6di, volptr7di, &
9419 volptr1db, volptr2db, volptr3db, volptr4db, volptr5db, volptr6db, volptr7db, &
9420 volptr1dc, volptr2dc, volptr3dc, volptr4dc, volptr5dc, volptr6dc, volptr7dc, &
9421 vol7d_nullifyr, vol7d_nullifyd, vol7d_nullifyi, vol7d_nullifyb, vol7d_nullifyc, &
9422 vol7d_init, vol7d_delete, vol7d_write_on_file, vol7d_read_from_file, &
9423 vol7d_check_alloc_ana, vol7d_force_alloc_ana, &
9424 vol7d_check_alloc_dati, vol7d_force_alloc_dati, vol7d_force_alloc, &
9425 vol7d_display, dat_display, dat_vect_display, &
9426 to_char_dat, vol7d_check
9427
9428PRIVATE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
9429
9430PRIVATE vol7d_c_e
9431
9432CONTAINS
9433
9434
9439SUBROUTINE vol7d_init(this,time_definition)
9440TYPE(vol7d),intent(out) :: this
9441integer,INTENT(IN),OPTIONAL :: time_definition
9442
9449CALL vol7d_var_features_init() ! initialise var features table once
9450
9451NULLIFY(this%ana, this%time, this%level, this%timerange, this%network)
9452
9453NULLIFY(this%volanar, this%volanaattrr, this%voldatir, this%voldatiattrr)
9454NULLIFY(this%volanad, this%volanaattrd, this%voldatid, this%voldatiattrd)
9455NULLIFY(this%volanai, this%volanaattri, this%voldatii, this%voldatiattri)
9456NULLIFY(this%volanab, this%volanaattrb, this%voldatib, this%voldatiattrb)
9457NULLIFY(this%volanac, this%volanaattrc, this%voldatic, this%voldatiattrc)
9458
9459if(present(time_definition)) then
9460 this%time_definition=time_definition
9461else
9462 this%time_definition=1 !default to validity time
9463end if
9464
9465END SUBROUTINE vol7d_init
9466
9467
9471ELEMENTAL SUBROUTINE vol7d_delete(this, dataonly)
9472TYPE(vol7d),intent(inout) :: this
9473LOGICAL, INTENT(in), OPTIONAL :: dataonly
9474
9475
9476IF (.NOT. optio_log(dataonly)) THEN
9477 IF (ASSOCIATED(this%volanar)) DEALLOCATE(this%volanar)
9478 IF (ASSOCIATED(this%volanad)) DEALLOCATE(this%volanad)
9479 IF (ASSOCIATED(this%volanai)) DEALLOCATE(this%volanai)
9480 IF (ASSOCIATED(this%volanab)) DEALLOCATE(this%volanab)
9481 IF (ASSOCIATED(this%volanac)) DEALLOCATE(this%volanac)
9482 IF (ASSOCIATED(this%volanaattrr)) DEALLOCATE(this%volanaattrr)
9483 IF (ASSOCIATED(this%volanaattrd)) DEALLOCATE(this%volanaattrd)
9484 IF (ASSOCIATED(this%volanaattri)) DEALLOCATE(this%volanaattri)
9485 IF (ASSOCIATED(this%volanaattrb)) DEALLOCATE(this%volanaattrb)
9486 IF (ASSOCIATED(this%volanaattrc)) DEALLOCATE(this%volanaattrc)
9487ENDIF
9488IF (ASSOCIATED(this%voldatir)) DEALLOCATE(this%voldatir)
9489IF (ASSOCIATED(this%voldatid)) DEALLOCATE(this%voldatid)
9490IF (ASSOCIATED(this%voldatii)) DEALLOCATE(this%voldatii)
9491IF (ASSOCIATED(this%voldatib)) DEALLOCATE(this%voldatib)
9492IF (ASSOCIATED(this%voldatic)) DEALLOCATE(this%voldatic)
9493IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
9494IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
9495IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
9496IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
9497IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
9498
9499IF (.NOT. optio_log(dataonly)) THEN
9500 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
9501 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
9502ENDIF
9503IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
9504IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
9505IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
9506
9507IF (.NOT. optio_log(dataonly)) THEN
9511ENDIF
9515
9516END SUBROUTINE vol7d_delete
9517
9518
9519
9520integer function vol7d_check(this)
9521TYPE(vol7d),intent(in) :: this
9522integer :: i,j,k,l,m,n
9523
9524vol7d_check=0
9525
9526if (associated(this%voldatii)) then
9527do i = 1,size(this%voldatii,1)
9528 do j = 1,size(this%voldatii,2)
9529 do k = 1,size(this%voldatii,3)
9530 do l = 1,size(this%voldatii,4)
9531 do m = 1,size(this%voldatii,5)
9532 do n = 1,size(this%voldatii,6)
9533 if (this%voldatii(i,j,k,l,m,n) /= this%voldatii(i,j,k,l,m,n) ) then
9534 CALL l4f_log(l4f_warn,"check: abnormal value at voldatii("&
9536 vol7d_check=1
9537 end if
9538 end do
9539 end do
9540 end do
9541 end do
9542 end do
9543end do
9544end if
9545
9546
9547if (associated(this%voldatir)) then
9548do i = 1,size(this%voldatir,1)
9549 do j = 1,size(this%voldatir,2)
9550 do k = 1,size(this%voldatir,3)
9551 do l = 1,size(this%voldatir,4)
9552 do m = 1,size(this%voldatir,5)
9553 do n = 1,size(this%voldatir,6)
9554 if (this%voldatir(i,j,k,l,m,n) /= this%voldatir(i,j,k,l,m,n) ) then
9555 CALL l4f_log(l4f_warn,"check: abnormal value at voldatir("&
9557 vol7d_check=2
9558 end if
9559 end do
9560 end do
9561 end do
9562 end do
9563 end do
9564end do
9565end if
9566
9567if (associated(this%voldatid)) then
9568do i = 1,size(this%voldatid,1)
9569 do j = 1,size(this%voldatid,2)
9570 do k = 1,size(this%voldatid,3)
9571 do l = 1,size(this%voldatid,4)
9572 do m = 1,size(this%voldatid,5)
9573 do n = 1,size(this%voldatid,6)
9574 if (this%voldatid(i,j,k,l,m,n) /= this%voldatid(i,j,k,l,m,n) ) then
9575 CALL l4f_log(l4f_warn,"check: abnormal value at voldatid("&
9577 vol7d_check=3
9578 end if
9579 end do
9580 end do
9581 end do
9582 end do
9583 end do
9584end do
9585end if
9586
9587if (associated(this%voldatib)) then
9588do i = 1,size(this%voldatib,1)
9589 do j = 1,size(this%voldatib,2)
9590 do k = 1,size(this%voldatib,3)
9591 do l = 1,size(this%voldatib,4)
9592 do m = 1,size(this%voldatib,5)
9593 do n = 1,size(this%voldatib,6)
9594 if (this%voldatib(i,j,k,l,m,n) /= this%voldatib(i,j,k,l,m,n) ) then
9595 CALL l4f_log(l4f_warn,"check: abnormal value at voldatib("&
9597 vol7d_check=4
9598 end if
9599 end do
9600 end do
9601 end do
9602 end do
9603 end do
9604end do
9605end if
9606
9607end function vol7d_check
9608
9609
9610
9611!TODO da completare ! aborta se i volumi sono allocati a dimensione 0
9613SUBROUTINE vol7d_display(this)
9614TYPE(vol7d),intent(in) :: this
9615integer :: i
9616
9617REAL :: rdat
9618DOUBLE PRECISION :: ddat
9619INTEGER :: idat
9620INTEGER(kind=int_b) :: bdat
9621CHARACTER(len=vol7d_cdatalen) :: cdat
9622
9623
9624print*,"<<<<<<<<<<<<<<<<<<< vol7d object >>>>>>>>>>>>>>>>>>>>"
9625if (this%time_definition == 0) then
9626 print*,"TIME DEFINITION: time is reference time"
9627else if (this%time_definition == 1) then
9628 print*,"TIME DEFINITION: time is validity time"
9629else
9630 print*,"Time definition have a wrong walue:", this%time_definition
9631end if
9632
9633IF (ASSOCIATED(this%network))then
9634 print*,"---- network vector ----"
9635 print*,"elements=",size(this%network)
9636 do i=1, size(this%network)
9638 end do
9639end IF
9640
9641IF (ASSOCIATED(this%ana))then
9642 print*,"---- ana vector ----"
9643 print*,"elements=",size(this%ana)
9644 do i=1, size(this%ana)
9646 end do
9647end IF
9648
9649IF (ASSOCIATED(this%time))then
9650 print*,"---- time vector ----"
9651 print*,"elements=",size(this%time)
9652 do i=1, size(this%time)
9654 end do
9655end if
9656
9657IF (ASSOCIATED(this%level)) then
9658 print*,"---- level vector ----"
9659 print*,"elements=",size(this%level)
9660 do i =1,size(this%level)
9662 end do
9663end if
9664
9665IF (ASSOCIATED(this%timerange))then
9666 print*,"---- timerange vector ----"
9667 print*,"elements=",size(this%timerange)
9668 do i =1,size(this%timerange)
9670 end do
9671end if
9672
9673
9674print*,"---- ana vector ----"
9675print*,""
9676print*,"->>>>>>>>> anavar -"
9678print*,""
9679print*,"->>>>>>>>> anaattr -"
9681print*,""
9682print*,"->>>>>>>>> anavarattr -"
9684
9685print*,"-- ana data section (first point) --"
9686
9687idat=imiss
9688rdat=rmiss
9689ddat=dmiss
9690bdat=ibmiss
9691cdat=cmiss
9692
9693!ntime = MIN(SIZE(this%time),nprint)
9694!ntimerange = MIN(SIZE(this%timerange),nprint)
9695!nlevel = MIN(SIZE(this%level),nprint)
9696!nnetwork = MIN(SIZE(this%network),nprint)
9697!nana = MIN(SIZE(this%ana),nprint)
9698
9699IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0) THEN
9700if (associated(this%volanai)) then
9701 do i=1,size(this%anavar%i)
9702 idat=this%volanai(1,i,1)
9704 end do
9705end if
9706idat=imiss
9707
9708if (associated(this%volanar)) then
9709 do i=1,size(this%anavar%r)
9710 rdat=this%volanar(1,i,1)
9712 end do
9713end if
9714rdat=rmiss
9715
9716if (associated(this%volanad)) then
9717 do i=1,size(this%anavar%d)
9718 ddat=this%volanad(1,i,1)
9720 end do
9721end if
9722ddat=dmiss
9723
9724if (associated(this%volanab)) then
9725 do i=1,size(this%anavar%b)
9726 bdat=this%volanab(1,i,1)
9728 end do
9729end if
9730bdat=ibmiss
9731
9732if (associated(this%volanac)) then
9733 do i=1,size(this%anavar%c)
9734 cdat=this%volanac(1,i,1)
9736 end do
9737end if
9738cdat=cmiss
9739ENDIF
9740
9741print*,"---- data vector ----"
9742print*,""
9743print*,"->>>>>>>>> dativar -"
9745print*,""
9746print*,"->>>>>>>>> datiattr -"
9748print*,""
9749print*,"->>>>>>>>> dativarattr -"
9751
9752print*,"-- data data section (first point) --"
9753
9754idat=imiss
9755rdat=rmiss
9756ddat=dmiss
9757bdat=ibmiss
9758cdat=cmiss
9759
9760IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0 .AND. size(this%time) > 0 &
9761 .AND. size(this%level) > 0 .AND. size(this%timerange) > 0) THEN
9762if (associated(this%voldatii)) then
9763 do i=1,size(this%dativar%i)
9764 idat=this%voldatii(1,1,1,1,i,1)
9766 end do
9767end if
9768idat=imiss
9769
9770if (associated(this%voldatir)) then
9771 do i=1,size(this%dativar%r)
9772 rdat=this%voldatir(1,1,1,1,i,1)
9774 end do
9775end if
9776rdat=rmiss
9777
9778if (associated(this%voldatid)) then
9779 do i=1,size(this%dativar%d)
9780 ddat=this%voldatid(1,1,1,1,i,1)
9782 end do
9783end if
9784ddat=dmiss
9785
9786if (associated(this%voldatib)) then
9787 do i=1,size(this%dativar%b)
9788 bdat=this%voldatib(1,1,1,1,i,1)
9790 end do
9791end if
9792bdat=ibmiss
9793
9794if (associated(this%voldatic)) then
9795 do i=1,size(this%dativar%c)
9796 cdat=this%voldatic(1,1,1,1,i,1)
9798 end do
9799end if
9800cdat=cmiss
9801ENDIF
9802
9803print*,"<<<<<<<<<<<<<<<<<<< END vol7d object >>>>>>>>>>>>>>>>>>>>"
9804
9805END SUBROUTINE vol7d_display
9806
9807
9809SUBROUTINE dat_display(this,idat,rdat,ddat,bdat,cdat)
9810TYPE(vol7d_var),intent(in) :: this
9812REAL :: rdat
9814DOUBLE PRECISION :: ddat
9816INTEGER :: idat
9818INTEGER(kind=int_b) :: bdat
9820CHARACTER(len=*) :: cdat
9821
9822print *, to_char_dat(this,idat,rdat,ddat,bdat,cdat)
9823
9824end SUBROUTINE dat_display
9825
9827SUBROUTINE dat_vect_display(this,idat,rdat,ddat,bdat,cdat)
9828
9829TYPE(vol7d_var),intent(in) :: this(:)
9831REAL :: rdat(:)
9833DOUBLE PRECISION :: ddat(:)
9835INTEGER :: idat(:)
9837INTEGER(kind=int_b) :: bdat(:)
9839CHARACTER(len=*):: cdat(:)
9840
9841integer :: i
9842
9843do i =1,size(this)
9845end do
9846
9847end SUBROUTINE dat_vect_display
9848
9849
9850FUNCTION to_char_dat(this,idat,rdat,ddat,bdat,cdat)
9851#ifdef HAVE_DBALLE
9852USE dballef
9853#endif
9854TYPE(vol7d_var),INTENT(in) :: this
9856REAL :: rdat
9858DOUBLE PRECISION :: ddat
9860INTEGER :: idat
9862INTEGER(kind=int_b) :: bdat
9864CHARACTER(len=*) :: cdat
9865CHARACTER(len=80) :: to_char_dat
9866
9867CHARACTER(len=LEN(to_char_dat)) :: to_char_tmp
9868
9869
9870#ifdef HAVE_DBALLE
9871INTEGER :: handle, ier
9872
9873handle = 0
9874to_char_dat="VALUE: "
9875
9880
9882 ier = idba_messaggi(handle,"/dev/null", "w", "BUFR")
9883 ier = idba_spiegab(handle,this%btable,cdat,to_char_tmp)
9884 ier = idba_fatto(handle)
9885 to_char_dat=trim(to_char_dat)//" ;char> "//trim(to_char_tmp)
9886endif
9887
9888#else
9889
9890to_char_dat="VALUE: "
9896
9897#endif
9898
9899END FUNCTION to_char_dat
9900
9901
9904FUNCTION vol7d_c_e(this) RESULT(c_e)
9905TYPE(vol7d), INTENT(in) :: this
9906
9907LOGICAL :: c_e
9908
9910 ASSOCIATED(this%level) .OR. ASSOCIATED(this%timerange) .OR. &
9911 ASSOCIATED(this%network) .OR. &
9912 ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
9913 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
9914 ASSOCIATED(this%anavar%c) .OR. &
9915 ASSOCIATED(this%anaattr%r) .OR. ASSOCIATED(this%anaattr%d) .OR. &
9916 ASSOCIATED(this%anaattr%i) .OR. ASSOCIATED(this%anaattr%b) .OR. &
9917 ASSOCIATED(this%anaattr%c) .OR. &
9918 ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
9919 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
9920 ASSOCIATED(this%dativar%c) .OR. &
9921 ASSOCIATED(this%datiattr%r) .OR. ASSOCIATED(this%datiattr%d) .OR. &
9922 ASSOCIATED(this%datiattr%i) .OR. ASSOCIATED(this%datiattr%b) .OR. &
9923 ASSOCIATED(this%datiattr%c)
9924
9925END FUNCTION vol7d_c_e
9926
9927
9966SUBROUTINE vol7d_alloc(this, nana, ntime, nlevel, ntimerange, nnetwork, &
9967 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
9968 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
9969 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
9970 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
9971 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
9972 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc, &
9973 ini)
9974TYPE(vol7d),INTENT(inout) :: this
9975INTEGER,INTENT(in),OPTIONAL :: nana
9976INTEGER,INTENT(in),OPTIONAL :: ntime
9977INTEGER,INTENT(in),OPTIONAL :: nlevel
9978INTEGER,INTENT(in),OPTIONAL :: ntimerange
9979INTEGER,INTENT(in),OPTIONAL :: nnetwork
9981INTEGER,INTENT(in),OPTIONAL :: &
9982 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
9983 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
9984 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
9985 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
9986 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
9987 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc
9988LOGICAL,INTENT(in),OPTIONAL :: ini
9989
9990INTEGER :: i
9991LOGICAL :: linit
9992
9993IF (PRESENT(ini)) THEN
9994 linit = ini
9995ELSE
9996 linit = .false.
9997ENDIF
9998
9999! Dimensioni principali
10000IF (PRESENT(nana)) THEN
10001 IF (nana >= 0) THEN
10002 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
10003 ALLOCATE(this%ana(nana))
10004 IF (linit) THEN
10005 DO i = 1, nana
10007 ENDDO
10008 ENDIF
10009 ENDIF
10010ENDIF
10011IF (PRESENT(ntime)) THEN
10012 IF (ntime >= 0) THEN
10013 IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
10014 ALLOCATE(this%time(ntime))
10015 IF (linit) THEN
10016 DO i = 1, ntime
10018 ENDDO
10019 ENDIF
10020 ENDIF
10021ENDIF
10022IF (PRESENT(nlevel)) THEN
10023 IF (nlevel >= 0) THEN
10024 IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
10025 ALLOCATE(this%level(nlevel))
10026 IF (linit) THEN
10027 DO i = 1, nlevel
10029 ENDDO
10030 ENDIF
10031 ENDIF
10032ENDIF
10033IF (PRESENT(ntimerange)) THEN
10034 IF (ntimerange >= 0) THEN
10035 IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
10036 ALLOCATE(this%timerange(ntimerange))
10037 IF (linit) THEN
10038 DO i = 1, ntimerange
10040 ENDDO
10041 ENDIF
10042 ENDIF
10043ENDIF
10044IF (PRESENT(nnetwork)) THEN
10045 IF (nnetwork >= 0) THEN
10046 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
10047 ALLOCATE(this%network(nnetwork))
10048 IF (linit) THEN
10049 DO i = 1, nnetwork
10051 ENDDO
10052 ENDIF
10053 ENDIF
10054ENDIF
10055! Dimensioni dei tipi delle variabili
10056CALL vol7d_varvect_alloc(this%anavar, nanavarr, nanavard, &
10057 nanavari, nanavarb, nanavarc, ini)
10058CALL vol7d_varvect_alloc(this%anaattr, nanaattrr, nanaattrd, &
10059 nanaattri, nanaattrb, nanaattrc, ini)
10060CALL vol7d_varvect_alloc(this%anavarattr, nanavarattrr, nanavarattrd, &
10061 nanavarattri, nanavarattrb, nanavarattrc, ini)
10062CALL vol7d_varvect_alloc(this%dativar, ndativarr, ndativard, &
10063 ndativari, ndativarb, ndativarc, ini)
10064CALL vol7d_varvect_alloc(this%datiattr, ndatiattrr, ndatiattrd, &
10065 ndatiattri, ndatiattrb, ndatiattrc, ini)
10066CALL vol7d_varvect_alloc(this%dativarattr, ndativarattrr, ndativarattrd, &
10067 ndativarattri, ndativarattrb, ndativarattrc, ini)
10068
10069END SUBROUTINE vol7d_alloc
10070
10071
10072FUNCTION vol7d_check_alloc_ana(this)
10073TYPE(vol7d),INTENT(in) :: this
10074LOGICAL :: vol7d_check_alloc_ana
10075
10076vol7d_check_alloc_ana = ASSOCIATED(this%ana) .AND. ASSOCIATED(this%network)
10077
10078END FUNCTION vol7d_check_alloc_ana
10079
10080SUBROUTINE vol7d_force_alloc_ana(this, ini)
10081TYPE(vol7d),INTENT(inout) :: this
10082LOGICAL,INTENT(in),OPTIONAL :: ini
10083
10084! Alloco i descrittori minimi per avere un volume di anagrafica
10085IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=1, ini=ini)
10086IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=1, ini=ini)
10087
10088END SUBROUTINE vol7d_force_alloc_ana
10089
10090
10091FUNCTION vol7d_check_alloc_dati(this)
10092TYPE(vol7d),INTENT(in) :: this
10093LOGICAL :: vol7d_check_alloc_dati
10094
10095vol7d_check_alloc_dati = vol7d_check_alloc_ana(this) .AND. &
10096 ASSOCIATED(this%time) .AND. ASSOCIATED(this%level) .AND. &
10097 ASSOCIATED(this%timerange)
10098
10099END FUNCTION vol7d_check_alloc_dati
10100
10101SUBROUTINE vol7d_force_alloc_dati(this, ini)
10102TYPE(vol7d),INTENT(inout) :: this
10103LOGICAL,INTENT(in),OPTIONAL :: ini
10104
10105! Alloco i descrittori minimi per avere un volume di dati
10106CALL vol7d_force_alloc_ana(this, ini)
10107IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=1, ini=ini)
10108IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=1, ini=ini)
10109IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=1, ini=ini)
10110
10111END SUBROUTINE vol7d_force_alloc_dati
10112
10113
10114SUBROUTINE vol7d_force_alloc(this)
10115TYPE(vol7d),INTENT(inout) :: this
10116
10117! If anything really not allocated yet, allocate with size 0
10118IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=0)
10119IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=0)
10120IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=0)
10121IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=0)
10122IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=0)
10123
10124END SUBROUTINE vol7d_force_alloc
10125
10126
10127FUNCTION vol7d_check_vol(this)
10128TYPE(vol7d),INTENT(in) :: this
10129LOGICAL :: vol7d_check_vol
10130
10131vol7d_check_vol = c_e(this)
10132
10133! Anagrafica
10134IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
10135 vol7d_check_vol = .false.
10136ENDIF
10137
10138IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
10139 vol7d_check_vol = .false.
10140ENDIF
10141
10142IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
10143 vol7d_check_vol = .false.
10144ENDIF
10145
10146IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
10147 vol7d_check_vol = .false.
10148ENDIF
10149
10150IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
10151 vol7d_check_vol = .false.
10152ENDIF
10153IF (ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
10154 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
10155 ASSOCIATED(this%anavar%c)) THEN
10156 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_ana(this)
10157ENDIF
10158
10159! Attributi dell'anagrafica
10160IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
10161 .NOT.ASSOCIATED(this%volanaattrr)) THEN
10162 vol7d_check_vol = .false.
10163ENDIF
10164
10165IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
10166 .NOT.ASSOCIATED(this%volanaattrd)) THEN
10167 vol7d_check_vol = .false.
10168ENDIF
10169
10170IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
10171 .NOT.ASSOCIATED(this%volanaattri)) THEN
10172 vol7d_check_vol = .false.
10173ENDIF
10174
10175IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
10176 .NOT.ASSOCIATED(this%volanaattrb)) THEN
10177 vol7d_check_vol = .false.
10178ENDIF
10179
10180IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
10181 .NOT.ASSOCIATED(this%volanaattrc)) THEN
10182 vol7d_check_vol = .false.
10183ENDIF
10184
10185! Dati
10186IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
10187 vol7d_check_vol = .false.
10188ENDIF
10189
10190IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
10191 vol7d_check_vol = .false.
10192ENDIF
10193
10194IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
10195 vol7d_check_vol = .false.
10196ENDIF
10197
10198IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
10199 vol7d_check_vol = .false.
10200ENDIF
10201
10202IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
10203 vol7d_check_vol = .false.
10204ENDIF
10205
10206! Attributi dei dati
10207IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
10208 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
10209 vol7d_check_vol = .false.
10210ENDIF
10211
10212IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
10213 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
10214 vol7d_check_vol = .false.
10215ENDIF
10216
10217IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
10218 .NOT.ASSOCIATED(this%voldatiattri)) THEN
10219 vol7d_check_vol = .false.
10220ENDIF
10221
10222IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
10223 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
10224 vol7d_check_vol = .false.
10225ENDIF
10226
10227IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
10228 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
10229 vol7d_check_vol = .false.
10230ENDIF
10231IF (ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
10232 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
10233 ASSOCIATED(this%dativar%c)) THEN
10234 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_dati(this)
10235ENDIF
10236
10237END FUNCTION vol7d_check_vol
10238
10239
10254SUBROUTINE vol7d_alloc_vol(this, ini, inivol)
10255TYPE(vol7d),INTENT(inout) :: this
10256LOGICAL,INTENT(in),OPTIONAL :: ini
10257LOGICAL,INTENT(in),OPTIONAL :: inivol
10258
10259LOGICAL :: linivol
10260
10261IF (PRESENT(inivol)) THEN
10262 linivol = inivol
10263ELSE
10264 linivol = .true.
10265ENDIF
10266
10267! Anagrafica
10268IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
10269 CALL vol7d_force_alloc_ana(this, ini)
10270 ALLOCATE(this%volanar(SIZE(this%ana), SIZE(this%anavar%r), SIZE(this%network)))
10271 IF (linivol) this%volanar(:,:,:) = rmiss
10272ENDIF
10273
10274IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
10275 CALL vol7d_force_alloc_ana(this, ini)
10276 ALLOCATE(this%volanad(SIZE(this%ana), SIZE(this%anavar%d), SIZE(this%network)))
10277 IF (linivol) this%volanad(:,:,:) = rdmiss
10278ENDIF
10279
10280IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
10281 CALL vol7d_force_alloc_ana(this, ini)
10282 ALLOCATE(this%volanai(SIZE(this%ana), SIZE(this%anavar%i), SIZE(this%network)))
10283 IF (linivol) this%volanai(:,:,:) = imiss
10284ENDIF
10285
10286IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
10287 CALL vol7d_force_alloc_ana(this, ini)
10288 ALLOCATE(this%volanab(SIZE(this%ana), SIZE(this%anavar%b), SIZE(this%network)))
10289 IF (linivol) this%volanab(:,:,:) = ibmiss
10290ENDIF
10291
10292IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
10293 CALL vol7d_force_alloc_ana(this, ini)
10294 ALLOCATE(this%volanac(SIZE(this%ana), SIZE(this%anavar%c), SIZE(this%network)))
10295 IF (linivol) this%volanac(:,:,:) = cmiss
10296ENDIF
10297
10298! Attributi dell'anagrafica
10299IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
10300 .NOT.ASSOCIATED(this%volanaattrr)) THEN
10301 CALL vol7d_force_alloc_ana(this, ini)
10302 ALLOCATE(this%volanaattrr(SIZE(this%ana), SIZE(this%anavarattr%r), &
10303 SIZE(this%network), SIZE(this%anaattr%r)))
10304 IF (linivol) this%volanaattrr(:,:,:,:) = rmiss
10305ENDIF
10306
10307IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
10308 .NOT.ASSOCIATED(this%volanaattrd)) THEN
10309 CALL vol7d_force_alloc_ana(this, ini)
10310 ALLOCATE(this%volanaattrd(SIZE(this%ana), SIZE(this%anavarattr%d), &
10311 SIZE(this%network), SIZE(this%anaattr%d)))
10312 IF (linivol) this%volanaattrd(:,:,:,:) = rdmiss
10313ENDIF
10314
10315IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
10316 .NOT.ASSOCIATED(this%volanaattri)) THEN
10317 CALL vol7d_force_alloc_ana(this, ini)
10318 ALLOCATE(this%volanaattri(SIZE(this%ana), SIZE(this%anavarattr%i), &
10319 SIZE(this%network), SIZE(this%anaattr%i)))
10320 IF (linivol) this%volanaattri(:,:,:,:) = imiss
10321ENDIF
10322
10323IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
10324 .NOT.ASSOCIATED(this%volanaattrb)) THEN
10325 CALL vol7d_force_alloc_ana(this, ini)
10326 ALLOCATE(this%volanaattrb(SIZE(this%ana), SIZE(this%anavarattr%b), &
10327 SIZE(this%network), SIZE(this%anaattr%b)))
10328 IF (linivol) this%volanaattrb(:,:,:,:) = ibmiss
10329ENDIF
10330
10331IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
10332 .NOT.ASSOCIATED(this%volanaattrc)) THEN
10333 CALL vol7d_force_alloc_ana(this, ini)
10334 ALLOCATE(this%volanaattrc(SIZE(this%ana), SIZE(this%anavarattr%c), &
10335 SIZE(this%network), SIZE(this%anaattr%c)))
10336 IF (linivol) this%volanaattrc(:,:,:,:) = cmiss
10337ENDIF
10338
10339! Dati
10340IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
10341 CALL vol7d_force_alloc_dati(this, ini)
10342 ALLOCATE(this%voldatir(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10343 SIZE(this%timerange), SIZE(this%dativar%r), SIZE(this%network)))
10344 IF (linivol) this%voldatir(:,:,:,:,:,:) = rmiss
10345ENDIF
10346
10347IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
10348 CALL vol7d_force_alloc_dati(this, ini)
10349 ALLOCATE(this%voldatid(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10350 SIZE(this%timerange), SIZE(this%dativar%d), SIZE(this%network)))
10351 IF (linivol) this%voldatid(:,:,:,:,:,:) = rdmiss
10352ENDIF
10353
10354IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
10355 CALL vol7d_force_alloc_dati(this, ini)
10356 ALLOCATE(this%voldatii(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10357 SIZE(this%timerange), SIZE(this%dativar%i), SIZE(this%network)))
10358 IF (linivol) this%voldatii(:,:,:,:,:,:) = imiss
10359ENDIF
10360
10361IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
10362 CALL vol7d_force_alloc_dati(this, ini)
10363 ALLOCATE(this%voldatib(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10364 SIZE(this%timerange), SIZE(this%dativar%b), SIZE(this%network)))
10365 IF (linivol) this%voldatib(:,:,:,:,:,:) = ibmiss
10366ENDIF
10367
10368IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
10369 CALL vol7d_force_alloc_dati(this, ini)
10370 ALLOCATE(this%voldatic(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10371 SIZE(this%timerange), SIZE(this%dativar%c), SIZE(this%network)))
10372 IF (linivol) this%voldatic(:,:,:,:,:,:) = cmiss
10373ENDIF
10374
10375! Attributi dei dati
10376IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
10377 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
10378 CALL vol7d_force_alloc_dati(this, ini)
10379 ALLOCATE(this%voldatiattrr(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10380 SIZE(this%timerange), SIZE(this%dativarattr%r), SIZE(this%network), &
10381 SIZE(this%datiattr%r)))
10382 IF (linivol) this%voldatiattrr(:,:,:,:,:,:,:) = rmiss
10383ENDIF
10384
10385IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
10386 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
10387 CALL vol7d_force_alloc_dati(this, ini)
10388 ALLOCATE(this%voldatiattrd(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10389 SIZE(this%timerange), SIZE(this%dativarattr%d), SIZE(this%network), &
10390 SIZE(this%datiattr%d)))
10391 IF (linivol) this%voldatiattrd(:,:,:,:,:,:,:) = rdmiss
10392ENDIF
10393
10394IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
10395 .NOT.ASSOCIATED(this%voldatiattri)) THEN
10396 CALL vol7d_force_alloc_dati(this, ini)
10397 ALLOCATE(this%voldatiattri(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10398 SIZE(this%timerange), SIZE(this%dativarattr%i), SIZE(this%network), &
10399 SIZE(this%datiattr%i)))
10400 IF (linivol) this%voldatiattri(:,:,:,:,:,:,:) = imiss
10401ENDIF
10402
10403IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
10404 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
10405 CALL vol7d_force_alloc_dati(this, ini)
10406 ALLOCATE(this%voldatiattrb(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10407 SIZE(this%timerange), SIZE(this%dativarattr%b), SIZE(this%network), &
10408 SIZE(this%datiattr%b)))
10409 IF (linivol) this%voldatiattrb(:,:,:,:,:,:,:) = ibmiss
10410ENDIF
10411
10412IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
10413 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
10414 CALL vol7d_force_alloc_dati(this, ini)
10415 ALLOCATE(this%voldatiattrc(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10416 SIZE(this%timerange), SIZE(this%dativarattr%c), SIZE(this%network), &
10417 SIZE(this%datiattr%c)))
10418 IF (linivol) this%voldatiattrc(:,:,:,:,:,:,:) = cmiss
10419ENDIF
10420
10421! Catch-all method
10422CALL vol7d_force_alloc(this)
10423
10424! Creo gli indici var-attr
10425
10426#ifdef DEBUG
10427CALL l4f_log(l4f_debug,"calling: vol7d_set_attr_ind")
10428#endif
10429
10430CALL vol7d_set_attr_ind(this)
10431
10432
10433
10434END SUBROUTINE vol7d_alloc_vol
10435
10436
10443SUBROUTINE vol7d_set_attr_ind(this)
10444TYPE(vol7d),INTENT(inout) :: this
10445
10446INTEGER :: i
10447
10448! real
10449IF (ASSOCIATED(this%dativar%r)) THEN
10450 IF (ASSOCIATED(this%dativarattr%r)) THEN
10451 DO i = 1, SIZE(this%dativar%r)
10452 this%dativar%r(i)%r = &
10453 firsttrue(this%dativar%r(i)%btable == this%dativarattr%r(:)%btable)
10454 ENDDO
10455 ENDIF
10456
10457 IF (ASSOCIATED(this%dativarattr%d)) THEN
10458 DO i = 1, SIZE(this%dativar%r)
10459 this%dativar%r(i)%d = &
10460 firsttrue(this%dativar%r(i)%btable == this%dativarattr%d(:)%btable)
10461 ENDDO
10462 ENDIF
10463
10464 IF (ASSOCIATED(this%dativarattr%i)) THEN
10465 DO i = 1, SIZE(this%dativar%r)
10466 this%dativar%r(i)%i = &
10467 firsttrue(this%dativar%r(i)%btable == this%dativarattr%i(:)%btable)
10468 ENDDO
10469 ENDIF
10470
10471 IF (ASSOCIATED(this%dativarattr%b)) THEN
10472 DO i = 1, SIZE(this%dativar%r)
10473 this%dativar%r(i)%b = &
10474 firsttrue(this%dativar%r(i)%btable == this%dativarattr%b(:)%btable)
10475 ENDDO
10476 ENDIF
10477
10478 IF (ASSOCIATED(this%dativarattr%c)) THEN
10479 DO i = 1, SIZE(this%dativar%r)
10480 this%dativar%r(i)%c = &
10481 firsttrue(this%dativar%r(i)%btable == this%dativarattr%c(:)%btable)
10482 ENDDO
10483 ENDIF
10484ENDIF
10485! double
10486IF (ASSOCIATED(this%dativar%d)) THEN
10487 IF (ASSOCIATED(this%dativarattr%r)) THEN
10488 DO i = 1, SIZE(this%dativar%d)
10489 this%dativar%d(i)%r = &
10490 firsttrue(this%dativar%d(i)%btable == this%dativarattr%r(:)%btable)
10491 ENDDO
10492 ENDIF
10493
10494 IF (ASSOCIATED(this%dativarattr%d)) THEN
10495 DO i = 1, SIZE(this%dativar%d)
10496 this%dativar%d(i)%d = &
10497 firsttrue(this%dativar%d(i)%btable == this%dativarattr%d(:)%btable)
10498 ENDDO
10499 ENDIF
10500
10501 IF (ASSOCIATED(this%dativarattr%i)) THEN
10502 DO i = 1, SIZE(this%dativar%d)
10503 this%dativar%d(i)%i = &
10504 firsttrue(this%dativar%d(i)%btable == this%dativarattr%i(:)%btable)
10505 ENDDO
10506 ENDIF
10507
10508 IF (ASSOCIATED(this%dativarattr%b)) THEN
10509 DO i = 1, SIZE(this%dativar%d)
10510 this%dativar%d(i)%b = &
10511 firsttrue(this%dativar%d(i)%btable == this%dativarattr%b(:)%btable)
10512 ENDDO
10513 ENDIF
10514
10515 IF (ASSOCIATED(this%dativarattr%c)) THEN
10516 DO i = 1, SIZE(this%dativar%d)
10517 this%dativar%d(i)%c = &
10518 firsttrue(this%dativar%d(i)%btable == this%dativarattr%c(:)%btable)
10519 ENDDO
10520 ENDIF
10521ENDIF
10522! integer
10523IF (ASSOCIATED(this%dativar%i)) THEN
10524 IF (ASSOCIATED(this%dativarattr%r)) THEN
10525 DO i = 1, SIZE(this%dativar%i)
10526 this%dativar%i(i)%r = &
10527 firsttrue(this%dativar%i(i)%btable == this%dativarattr%r(:)%btable)
10528 ENDDO
10529 ENDIF
10530
10531 IF (ASSOCIATED(this%dativarattr%d)) THEN
10532 DO i = 1, SIZE(this%dativar%i)
10533 this%dativar%i(i)%d = &
10534 firsttrue(this%dativar%i(i)%btable == this%dativarattr%d(:)%btable)
10535 ENDDO
10536 ENDIF
10537
10538 IF (ASSOCIATED(this%dativarattr%i)) THEN
10539 DO i = 1, SIZE(this%dativar%i)
10540 this%dativar%i(i)%i = &
10541 firsttrue(this%dativar%i(i)%btable == this%dativarattr%i(:)%btable)
10542 ENDDO
10543 ENDIF
10544
10545 IF (ASSOCIATED(this%dativarattr%b)) THEN
10546 DO i = 1, SIZE(this%dativar%i)
10547 this%dativar%i(i)%b = &
10548 firsttrue(this%dativar%i(i)%btable == this%dativarattr%b(:)%btable)
10549 ENDDO
10550 ENDIF
10551
10552 IF (ASSOCIATED(this%dativarattr%c)) THEN
10553 DO i = 1, SIZE(this%dativar%i)
10554 this%dativar%i(i)%c = &
10555 firsttrue(this%dativar%i(i)%btable == this%dativarattr%c(:)%btable)
10556 ENDDO
10557 ENDIF
10558ENDIF
10559! byte
10560IF (ASSOCIATED(this%dativar%b)) THEN
10561 IF (ASSOCIATED(this%dativarattr%r)) THEN
10562 DO i = 1, SIZE(this%dativar%b)
10563 this%dativar%b(i)%r = &
10564 firsttrue(this%dativar%b(i)%btable == this%dativarattr%r(:)%btable)
10565 ENDDO
10566 ENDIF
10567
10568 IF (ASSOCIATED(this%dativarattr%d)) THEN
10569 DO i = 1, SIZE(this%dativar%b)
10570 this%dativar%b(i)%d = &
10571 firsttrue(this%dativar%b(i)%btable == this%dativarattr%d(:)%btable)
10572 ENDDO
10573 ENDIF
10574
10575 IF (ASSOCIATED(this%dativarattr%i)) THEN
10576 DO i = 1, SIZE(this%dativar%b)
10577 this%dativar%b(i)%i = &
10578 firsttrue(this%dativar%b(i)%btable == this%dativarattr%i(:)%btable)
10579 ENDDO
10580 ENDIF
10581
10582 IF (ASSOCIATED(this%dativarattr%b)) THEN
10583 DO i = 1, SIZE(this%dativar%b)
10584 this%dativar%b(i)%b = &
10585 firsttrue(this%dativar%b(i)%btable == this%dativarattr%b(:)%btable)
10586 ENDDO
10587 ENDIF
10588
10589 IF (ASSOCIATED(this%dativarattr%c)) THEN
10590 DO i = 1, SIZE(this%dativar%b)
10591 this%dativar%b(i)%c = &
10592 firsttrue(this%dativar%b(i)%btable == this%dativarattr%c(:)%btable)
10593 ENDDO
10594 ENDIF
10595ENDIF
10596! character
10597IF (ASSOCIATED(this%dativar%c)) THEN
10598 IF (ASSOCIATED(this%dativarattr%r)) THEN
10599 DO i = 1, SIZE(this%dativar%c)
10600 this%dativar%c(i)%r = &
10601 firsttrue(this%dativar%c(i)%btable == this%dativarattr%r(:)%btable)
10602 ENDDO
10603 ENDIF
10604
10605 IF (ASSOCIATED(this%dativarattr%d)) THEN
10606 DO i = 1, SIZE(this%dativar%c)
10607 this%dativar%c(i)%d = &
10608 firsttrue(this%dativar%c(i)%btable == this%dativarattr%d(:)%btable)
10609 ENDDO
10610 ENDIF
10611
10612 IF (ASSOCIATED(this%dativarattr%i)) THEN
10613 DO i = 1, SIZE(this%dativar%c)
10614 this%dativar%c(i)%i = &
10615 firsttrue(this%dativar%c(i)%btable == this%dativarattr%i(:)%btable)
10616 ENDDO
10617 ENDIF
10618
10619 IF (ASSOCIATED(this%dativarattr%b)) THEN
10620 DO i = 1, SIZE(this%dativar%c)
10621 this%dativar%c(i)%b = &
10622 firsttrue(this%dativar%c(i)%btable == this%dativarattr%b(:)%btable)
10623 ENDDO
10624 ENDIF
10625
10626 IF (ASSOCIATED(this%dativarattr%c)) THEN
10627 DO i = 1, SIZE(this%dativar%c)
10628 this%dativar%c(i)%c = &
10629 firsttrue(this%dativar%c(i)%btable == this%dativarattr%c(:)%btable)
10630 ENDDO
10631 ENDIF
10632ENDIF
10633
10634END SUBROUTINE vol7d_set_attr_ind
10635
10636
10641SUBROUTINE vol7d_merge(this, that, sort, bestdata, &
10642 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
10643TYPE(vol7d),INTENT(INOUT) :: this
10644TYPE(vol7d),INTENT(INOUT) :: that
10645LOGICAL,INTENT(IN),OPTIONAL :: sort
10646LOGICAL,INTENT(in),OPTIONAL :: bestdata
10647LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple ! experimental, please do not use outside the library now
10648
10649TYPE(vol7d) :: v7d_clean
10650
10651
10653 this = that
10655 that = v7d_clean ! destroy that without deallocating
10656ELSE ! Append that to this and destroy that
10658 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
10660ENDIF
10661
10662END SUBROUTINE vol7d_merge
10663
10664
10693SUBROUTINE vol7d_append(this, that, sort, bestdata, &
10694 ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple)
10695TYPE(vol7d),INTENT(INOUT) :: this
10696TYPE(vol7d),INTENT(IN) :: that
10697LOGICAL,INTENT(IN),OPTIONAL :: sort
10698! experimental, please do not use outside the library now, they force the use
10699! of a simplified mapping algorithm which is valid only whene the dimension
10700! content is the same in both volumes , or when one of them is empty
10701LOGICAL,INTENT(in),OPTIONAL :: bestdata
10702LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple
10703
10704
10705TYPE(vol7d) :: v7dtmp
10706LOGICAL :: lsort, lbestdata
10707INTEGER,POINTER :: remapt1(:), remapt2(:), remaptr1(:), remaptr2(:), &
10708 remapl1(:), remapl2(:), remapa1(:), remapa2(:), remapn1(:), remapn2(:)
10709
10711IF (.NOT.vol7d_check_vol(that)) RETURN ! be safe
10714 RETURN
10715ENDIF
10716
10717IF (this%time_definition /= that%time_definition) THEN
10718 CALL l4f_log(l4f_fatal, &
10719 'in vol7d_append, cannot append volumes with different &
10720 &time definition')
10721 CALL raise_fatal_error()
10722ENDIF
10723
10724! Completo l'allocazione per avere volumi a norma
10725CALL vol7d_alloc_vol(this)
10726
10730
10731! Calcolo le mappature tra volumi vecchi e volume nuovo
10732! I puntatori remap* vengono tutti o allocati o nullificati
10733IF (optio_log(ltimesimple)) THEN
10734 CALL vol7d_remap2simple_datetime(this%time, that%time, v7dtmp%time, &
10735 lsort, remapt1, remapt2)
10736ELSE
10737 CALL vol7d_remap2_datetime(this%time, that%time, v7dtmp%time, &
10738 lsort, remapt1, remapt2)
10739ENDIF
10740IF (optio_log(ltimerangesimple)) THEN
10741 CALL vol7d_remap2simple_vol7d_timerange(this%timerange, that%timerange, &
10742 v7dtmp%timerange, lsort, remaptr1, remaptr2)
10743ELSE
10744 CALL vol7d_remap2_vol7d_timerange(this%timerange, that%timerange, &
10745 v7dtmp%timerange, lsort, remaptr1, remaptr2)
10746ENDIF
10747IF (optio_log(llevelsimple)) THEN
10748 CALL vol7d_remap2simple_vol7d_level(this%level, that%level, v7dtmp%level, &
10749 lsort, remapl1, remapl2)
10750ELSE
10751 CALL vol7d_remap2_vol7d_level(this%level, that%level, v7dtmp%level, &
10752 lsort, remapl1, remapl2)
10753ENDIF
10754IF (optio_log(lanasimple)) THEN
10755 CALL vol7d_remap2simple_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
10756 .false., remapa1, remapa2)
10757ELSE
10758 CALL vol7d_remap2_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
10759 .false., remapa1, remapa2)
10760ENDIF
10761IF (optio_log(lnetworksimple)) THEN
10762 CALL vol7d_remap2simple_vol7d_network(this%network, that%network, v7dtmp%network, &
10763 .false., remapn1, remapn2)
10764ELSE
10765 CALL vol7d_remap2_vol7d_network(this%network, that%network, v7dtmp%network, &
10766 .false., remapn1, remapn2)
10767ENDIF
10768
10769! Faccio la fusione fisica dei volumi
10770CALL vol7d_merge_finalr(this, that, v7dtmp, &
10771 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10772 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10773CALL vol7d_merge_finald(this, that, v7dtmp, &
10774 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10775 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10776CALL vol7d_merge_finali(this, that, v7dtmp, &
10777 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10778 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10779CALL vol7d_merge_finalb(this, that, v7dtmp, &
10780 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10781 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10782CALL vol7d_merge_finalc(this, that, v7dtmp, &
10783 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10784 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10785
10786! Dealloco i vettori di rimappatura
10787IF (ASSOCIATED(remapt1)) DEALLOCATE(remapt1)
10788IF (ASSOCIATED(remapt2)) DEALLOCATE(remapt2)
10789IF (ASSOCIATED(remaptr1)) DEALLOCATE(remaptr1)
10790IF (ASSOCIATED(remaptr2)) DEALLOCATE(remaptr2)
10791IF (ASSOCIATED(remapl1)) DEALLOCATE(remapl1)
10792IF (ASSOCIATED(remapl2)) DEALLOCATE(remapl2)
10793IF (ASSOCIATED(remapa1)) DEALLOCATE(remapa1)
10794IF (ASSOCIATED(remapa2)) DEALLOCATE(remapa2)
10795IF (ASSOCIATED(remapn1)) DEALLOCATE(remapn1)
10796IF (ASSOCIATED(remapn2)) DEALLOCATE(remapn2)
10797
10798! Distruggo il vecchio volume e assegno il nuovo a this
10800this = v7dtmp
10801! Ricreo gli indici var-attr
10802CALL vol7d_set_attr_ind(this)
10803
10804END SUBROUTINE vol7d_append
10805
10806
10839SUBROUTINE vol7d_copy(this, that, sort, unique, miss, &
10840 lsort_time, lsort_timerange, lsort_level, &
10841 ltime, ltimerange, llevel, lana, lnetwork, &
10842 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
10843 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
10844 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
10845 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
10846 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
10847 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
10848TYPE(vol7d),INTENT(IN) :: this
10849TYPE(vol7d),INTENT(INOUT) :: that
10850LOGICAL,INTENT(IN),OPTIONAL :: sort
10851LOGICAL,INTENT(IN),OPTIONAL :: unique
10852LOGICAL,INTENT(IN),OPTIONAL :: miss
10853LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
10854LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
10855LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
10863LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
10865LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
10867LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
10869LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
10871LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
10873LOGICAL,INTENT(in),OPTIONAL :: &
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(:)
10880
10881LOGICAL :: lsort, lunique, lmiss
10882INTEGER,POINTER :: remapt(:), remaptr(:), remapl(:), remapa(:), remapn(:)
10883
10886IF (.NOT.vol7d_check_vol(this)) RETURN ! be safe
10887
10891
10892! Calcolo le mappature tra volume vecchio e volume nuovo
10893! I puntatori remap* vengono tutti o allocati o nullificati
10894CALL vol7d_remap1_datetime(this%time, that%time, &
10895 lsort.OR.optio_log(lsort_time), lunique, lmiss, remapt, ltime)
10896CALL vol7d_remap1_vol7d_timerange(this%timerange, that%timerange, &
10897 lsort.OR.optio_log(lsort_timerange), lunique, lmiss, remaptr, ltimerange)
10898CALL vol7d_remap1_vol7d_level(this%level, that%level, &
10899 lsort.OR.optio_log(lsort_level), lunique, lmiss, remapl, llevel)
10900CALL vol7d_remap1_vol7d_ana(this%ana, that%ana, &
10901 lsort, lunique, lmiss, remapa, lana)
10902CALL vol7d_remap1_vol7d_network(this%network, that%network, &
10903 lsort, lunique, lmiss, remapn, lnetwork)
10904
10905! lanavari, lanavarb, lanavarc, &
10906! lanaattri, lanaattrb, lanaattrc, &
10907! lanavarattri, lanavarattrb, lanavarattrc, &
10908! ldativari, ldativarb, ldativarc, &
10909! ldatiattri, ldatiattrb, ldatiattrc, &
10910! ldativarattri, ldativarattrb, ldativarattrc
10911! Faccio la riforma fisica dei volumi
10912CALL vol7d_reform_finalr(this, that, &
10913 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10914 lanavarr, lanaattrr, lanavarattrr, ldativarr, ldatiattrr, ldativarattrr)
10915CALL vol7d_reform_finald(this, that, &
10916 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10917 lanavard, lanaattrd, lanavarattrd, ldativard, ldatiattrd, ldativarattrd)
10918CALL vol7d_reform_finali(this, that, &
10919 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10920 lanavari, lanaattri, lanavarattri, ldativari, ldatiattri, ldativarattri)
10921CALL vol7d_reform_finalb(this, that, &
10922 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10923 lanavarb, lanaattrb, lanavarattrb, ldativarb, ldatiattrb, ldativarattrb)
10924CALL vol7d_reform_finalc(this, that, &
10925 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10926 lanavarc, lanaattrc, lanavarattrc, ldativarc, ldatiattrc, ldativarattrc)
10927
10928! Dealloco i vettori di rimappatura
10929IF (ASSOCIATED(remapt)) DEALLOCATE(remapt)
10930IF (ASSOCIATED(remaptr)) DEALLOCATE(remaptr)
10931IF (ASSOCIATED(remapl)) DEALLOCATE(remapl)
10932IF (ASSOCIATED(remapa)) DEALLOCATE(remapa)
10933IF (ASSOCIATED(remapn)) DEALLOCATE(remapn)
10934
10935! Ricreo gli indici var-attr
10936CALL vol7d_set_attr_ind(that)
10937that%time_definition = this%time_definition
10938
10939END SUBROUTINE vol7d_copy
10940
10941
10952SUBROUTINE vol7d_reform(this, sort, unique, miss, &
10953 lsort_time, lsort_timerange, lsort_level, &
10954 ltime, ltimerange, llevel, lana, lnetwork, &
10955 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
10956 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
10957 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
10958 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
10959 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
10960 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc&
10961 ,purgeana)
10962TYPE(vol7d),INTENT(INOUT) :: this
10963LOGICAL,INTENT(IN),OPTIONAL :: sort
10964LOGICAL,INTENT(IN),OPTIONAL :: unique
10965LOGICAL,INTENT(IN),OPTIONAL :: miss
10966LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
10967LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
10968LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
10976LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
10977LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
10978LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
10979LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
10980LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
10982LOGICAL,INTENT(in),OPTIONAL :: &
10983 lanavarr(:), lanavard(:), lanavari(:), lanavarb(:), lanavarc(:), &
10984 lanaattrr(:), lanaattrd(:), lanaattri(:), lanaattrb(:), lanaattrc(:), &
10985 lanavarattrr(:), lanavarattrd(:), lanavarattri(:), lanavarattrb(:), lanavarattrc(:), &
10986 ldativarr(:), ldativard(:), ldativari(:), ldativarb(:), ldativarc(:), &
10987 ldatiattrr(:), ldatiattrd(:), ldatiattri(:), ldatiattrb(:), ldatiattrc(:), &
10988 ldativarattrr(:), ldativarattrd(:), ldativarattri(:), ldativarattrb(:), ldativarattrc(:)
10989LOGICAL,INTENT(IN),OPTIONAL :: purgeana
10990
10991TYPE(vol7d) :: v7dtmp
10992logical,allocatable :: llana(:)
10993integer :: i
10994
10996 lsort_time, lsort_timerange, lsort_level, &
10997 ltime, ltimerange, llevel, lana, lnetwork, &
10998 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
10999 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
11000 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
11001 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
11002 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
11003 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
11004
11005! destroy old volume
11007
11008if (optio_log(purgeana)) then
11009 allocate(llana(size(v7dtmp%ana)))
11010 llana =.false.
11011 do i =1,size(v7dtmp%ana)
11012 if (associated(v7dtmp%voldatii)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatii(i,:,:,:,:,:)))
11013 if (associated(v7dtmp%voldatir)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatir(i,:,:,:,:,:)))
11014 if (associated(v7dtmp%voldatid)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatid(i,:,:,:,:,:)))
11015 if (associated(v7dtmp%voldatib)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatib(i,:,:,:,:,:)))
11016 if (associated(v7dtmp%voldatic)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatic(i,:,:,:,:,:)))
11017 end do
11018 CALL vol7d_copy(v7dtmp, this,lana=llana)
11020 deallocate(llana)
11021else
11022 this=v7dtmp
11023end if
11024
11025END SUBROUTINE vol7d_reform
11026
11027
11035SUBROUTINE vol7d_smart_sort(this, lsort_time, lsort_timerange, lsort_level)
11036TYPE(vol7d),INTENT(INOUT) :: this
11037LOGICAL,OPTIONAL,INTENT(in) :: lsort_time
11038LOGICAL,OPTIONAL,INTENT(in) :: lsort_timerange
11039LOGICAL,OPTIONAL,INTENT(in) :: lsort_level
11040
11041INTEGER :: i
11042LOGICAL :: to_be_sorted
11043
11044to_be_sorted = .false.
11045CALL vol7d_alloc_vol(this) ! usual safety check
11046
11047IF (optio_log(lsort_time)) THEN
11048 DO i = 2, SIZE(this%time)
11049 IF (this%time(i) < this%time(i-1)) THEN
11050 to_be_sorted = .true.
11051 EXIT
11052 ENDIF
11053 ENDDO
11054ENDIF
11055IF (optio_log(lsort_timerange)) THEN
11056 DO i = 2, SIZE(this%timerange)
11057 IF (this%timerange(i) < this%timerange(i-1)) THEN
11058 to_be_sorted = .true.
11059 EXIT
11060 ENDIF
11061 ENDDO
11062ENDIF
11063IF (optio_log(lsort_level)) THEN
11064 DO i = 2, SIZE(this%level)
11065 IF (this%level(i) < this%level(i-1)) THEN
11066 to_be_sorted = .true.
11067 EXIT
11068 ENDIF
11069 ENDDO
11070ENDIF
11071
11072IF (to_be_sorted) CALL vol7d_reform(this, &
11073 lsort_time=lsort_time, lsort_timerange=lsort_timerange, lsort_level=lsort_level )
11074
11075END SUBROUTINE vol7d_smart_sort
11076
11084SUBROUTINE vol7d_filter(this, avl, vl, nl, s_d, e_d)
11085TYPE(vol7d),INTENT(inout) :: this
11086CHARACTER(len=*),INTENT(in),OPTIONAL :: avl(:)
11087CHARACTER(len=*),INTENT(in),OPTIONAL :: vl(:)
11088TYPE(vol7d_network),OPTIONAL :: nl(:)
11089TYPE(datetime),INTENT(in),OPTIONAL :: s_d
11090TYPE(datetime),INTENT(in),OPTIONAL :: e_d
11091
11092INTEGER :: i
11093
11094IF (PRESENT(avl)) THEN
11095 IF (SIZE(avl) > 0) THEN
11096
11097 IF (ASSOCIATED(this%anavar%r)) THEN
11098 DO i = 1, SIZE(this%anavar%r)
11099 IF (all(this%anavar%r(i)%btable /= avl)) this%anavar%r(i) = vol7d_var_miss
11100 ENDDO
11101 ENDIF
11102
11103 IF (ASSOCIATED(this%anavar%i)) THEN
11104 DO i = 1, SIZE(this%anavar%i)
11105 IF (all(this%anavar%i(i)%btable /= avl)) this%anavar%i(i) = vol7d_var_miss
11106 ENDDO
11107 ENDIF
11108
11109 IF (ASSOCIATED(this%anavar%b)) THEN
11110 DO i = 1, SIZE(this%anavar%b)
11111 IF (all(this%anavar%b(i)%btable /= avl)) this%anavar%b(i) = vol7d_var_miss
11112 ENDDO
11113 ENDIF
11114
11115 IF (ASSOCIATED(this%anavar%d)) THEN
11116 DO i = 1, SIZE(this%anavar%d)
11117 IF (all(this%anavar%d(i)%btable /= avl)) this%anavar%d(i) = vol7d_var_miss
11118 ENDDO
11119 ENDIF
11120
11121 IF (ASSOCIATED(this%anavar%c)) THEN
11122 DO i = 1, SIZE(this%anavar%c)
11123 IF (all(this%anavar%c(i)%btable /= avl)) this%anavar%c(i) = vol7d_var_miss
11124 ENDDO
11125 ENDIF
11126
11127 ENDIF
11128ENDIF
11129
11130
11131IF (PRESENT(vl)) THEN
11132 IF (size(vl) > 0) THEN
11133 IF (ASSOCIATED(this%dativar%r)) THEN
11134 DO i = 1, SIZE(this%dativar%r)
11135 IF (all(this%dativar%r(i)%btable /= vl)) this%dativar%r(i) = vol7d_var_miss
11136 ENDDO
11137 ENDIF
11138
11139 IF (ASSOCIATED(this%dativar%i)) THEN
11140 DO i = 1, SIZE(this%dativar%i)
11141 IF (all(this%dativar%i(i)%btable /= vl)) this%dativar%i(i) = vol7d_var_miss
11142 ENDDO
11143 ENDIF
11144
11145 IF (ASSOCIATED(this%dativar%b)) THEN
11146 DO i = 1, SIZE(this%dativar%b)
11147 IF (all(this%dativar%b(i)%btable /= vl)) this%dativar%b(i) = vol7d_var_miss
11148 ENDDO
11149 ENDIF
11150
11151 IF (ASSOCIATED(this%dativar%d)) THEN
11152 DO i = 1, SIZE(this%dativar%d)
11153 IF (all(this%dativar%d(i)%btable /= vl)) this%dativar%d(i) = vol7d_var_miss
11154 ENDDO
11155 ENDIF
11156
11157 IF (ASSOCIATED(this%dativar%c)) THEN
11158 DO i = 1, SIZE(this%dativar%c)
11159 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
11160 ENDDO
11161 ENDIF
11162
11163 IF (ASSOCIATED(this%dativar%c)) THEN
11164 DO i = 1, SIZE(this%dativar%c)
11165 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
11166 ENDDO
11167 ENDIF
11168
11169 ENDIF
11170ENDIF
11171
11172IF (PRESENT(nl)) THEN
11173 IF (SIZE(nl) > 0) THEN
11174 DO i = 1, SIZE(this%network)
11175 IF (all(this%network(i) /= nl)) this%network(i) = vol7d_network_miss
11176 ENDDO
11177 ENDIF
11178ENDIF
11179
11180IF (PRESENT(s_d)) THEN
11182 WHERE (this%time < s_d)
11183 this%time = datetime_miss
11184 END WHERE
11185 ENDIF
11186ENDIF
11187
11188IF (PRESENT(e_d)) THEN
11190 WHERE (this%time > e_d)
11191 this%time = datetime_miss
11192 END WHERE
11193 ENDIF
11194ENDIF
11195
11196CALL vol7d_reform(this, miss=.true.)
11197
11198END SUBROUTINE vol7d_filter
11199
11200
11207SUBROUTINE vol7d_convr(this, that, anaconv)
11208TYPE(vol7d),INTENT(IN) :: this
11209TYPE(vol7d),INTENT(INOUT) :: that
11210LOGICAL,OPTIONAL,INTENT(in) :: anaconv
11211INTEGER :: i
11212LOGICAL :: fv(1)=(/.false./), tv(1)=(/.true./), acp(1), acn(1)
11213TYPE(vol7d) :: v7d_tmp
11214
11215IF (optio_log(anaconv)) THEN
11216 acp=fv
11217 acn=tv
11218ELSE
11219 acp=tv
11220 acn=fv
11221ENDIF
11222
11223! Volume con solo i dati reali e tutti gli attributi
11224! l'anagrafica e` copiata interamente se necessario
11225CALL vol7d_copy(this, that, &
11226 lanavarr=tv, lanavard=acp, lanavari=acp, lanavarb=acp, lanavarc=acp, &
11227 ldativarr=tv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=fv)
11228
11229! Volume solo di dati double
11230CALL vol7d_copy(this, v7d_tmp, &
11231 lanavarr=fv, lanavard=acn, lanavari=fv, lanavarb=fv, lanavarc=fv, &
11232 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11233 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11234 ldativarr=fv, ldativard=tv, ldativari=fv, ldativarb=fv, ldativarc=fv, &
11235 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11236 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11237
11238! converto a dati reali
11239IF (ASSOCIATED(v7d_tmp%anavar%d) .OR. ASSOCIATED(v7d_tmp%dativar%d)) THEN
11240
11241 IF (ASSOCIATED(v7d_tmp%anavar%d)) THEN
11242! alloco i dati reali e vi trasferisco i double
11243 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanad, 1), SIZE(v7d_tmp%volanad, 2), &
11244 SIZE(v7d_tmp%volanad, 3)))
11245 DO i = 1, SIZE(v7d_tmp%anavar%d)
11246 v7d_tmp%volanar(:,i,:) = &
11247 realdat(v7d_tmp%volanad(:,i,:), v7d_tmp%anavar%d(i))
11248 ENDDO
11249 DEALLOCATE(v7d_tmp%volanad)
11250! trasferisco le variabili
11251 v7d_tmp%anavar%r => v7d_tmp%anavar%d
11252 NULLIFY(v7d_tmp%anavar%d)
11253 ENDIF
11254
11255 IF (ASSOCIATED(v7d_tmp%dativar%d)) THEN
11256! alloco i dati reali e vi trasferisco i double
11257 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatid, 1), SIZE(v7d_tmp%voldatid, 2), &
11258 SIZE(v7d_tmp%voldatid, 3), SIZE(v7d_tmp%voldatid, 4), SIZE(v7d_tmp%voldatid, 5), &
11259 SIZE(v7d_tmp%voldatid, 6)))
11260 DO i = 1, SIZE(v7d_tmp%dativar%d)
11261 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11262 realdat(v7d_tmp%voldatid(:,:,:,:,i,:), v7d_tmp%dativar%d(i))
11263 ENDDO
11264 DEALLOCATE(v7d_tmp%voldatid)
11265! trasferisco le variabili
11266 v7d_tmp%dativar%r => v7d_tmp%dativar%d
11267 NULLIFY(v7d_tmp%dativar%d)
11268 ENDIF
11269
11270! fondo con il volume definitivo
11271 CALL vol7d_merge(that, v7d_tmp)
11272ELSE
11274ENDIF
11275
11276
11277! Volume solo di dati interi
11278CALL vol7d_copy(this, v7d_tmp, &
11279 lanavarr=fv, lanavard=fv, lanavari=acn, lanavarb=fv, lanavarc=fv, &
11280 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11281 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11282 ldativarr=fv, ldativard=fv, ldativari=tv, ldativarb=fv, ldativarc=fv, &
11283 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11284 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11285
11286! converto a dati reali
11287IF (ASSOCIATED(v7d_tmp%anavar%i) .OR. ASSOCIATED(v7d_tmp%dativar%i)) THEN
11288
11289 IF (ASSOCIATED(v7d_tmp%anavar%i)) THEN
11290! alloco i dati reali e vi trasferisco gli interi
11291 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanai, 1), SIZE(v7d_tmp%volanai, 2), &
11292 SIZE(v7d_tmp%volanai, 3)))
11293 DO i = 1, SIZE(v7d_tmp%anavar%i)
11294 v7d_tmp%volanar(:,i,:) = &
11295 realdat(v7d_tmp%volanai(:,i,:), v7d_tmp%anavar%i(i))
11296 ENDDO
11297 DEALLOCATE(v7d_tmp%volanai)
11298! trasferisco le variabili
11299 v7d_tmp%anavar%r => v7d_tmp%anavar%i
11300 NULLIFY(v7d_tmp%anavar%i)
11301 ENDIF
11302
11303 IF (ASSOCIATED(v7d_tmp%dativar%i)) THEN
11304! alloco i dati reali e vi trasferisco gli interi
11305 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatii, 1), SIZE(v7d_tmp%voldatii, 2), &
11306 SIZE(v7d_tmp%voldatii, 3), SIZE(v7d_tmp%voldatii, 4), SIZE(v7d_tmp%voldatii, 5), &
11307 SIZE(v7d_tmp%voldatii, 6)))
11308 DO i = 1, SIZE(v7d_tmp%dativar%i)
11309 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11310 realdat(v7d_tmp%voldatii(:,:,:,:,i,:), v7d_tmp%dativar%i(i))
11311 ENDDO
11312 DEALLOCATE(v7d_tmp%voldatii)
11313! trasferisco le variabili
11314 v7d_tmp%dativar%r => v7d_tmp%dativar%i
11315 NULLIFY(v7d_tmp%dativar%i)
11316 ENDIF
11317
11318! fondo con il volume definitivo
11319 CALL vol7d_merge(that, v7d_tmp)
11320ELSE
11322ENDIF
11323
11324
11325! Volume solo di dati byte
11326CALL vol7d_copy(this, v7d_tmp, &
11327 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=acn, lanavarc=fv, &
11328 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11329 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11330 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=tv, ldativarc=fv, &
11331 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11332 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11333
11334! converto a dati reali
11335IF (ASSOCIATED(v7d_tmp%anavar%b) .OR. ASSOCIATED(v7d_tmp%dativar%b)) THEN
11336
11337 IF (ASSOCIATED(v7d_tmp%anavar%b)) THEN
11338! alloco i dati reali e vi trasferisco i byte
11339 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanab, 1), SIZE(v7d_tmp%volanab, 2), &
11340 SIZE(v7d_tmp%volanab, 3)))
11341 DO i = 1, SIZE(v7d_tmp%anavar%b)
11342 v7d_tmp%volanar(:,i,:) = &
11343 realdat(v7d_tmp%volanab(:,i,:), v7d_tmp%anavar%b(i))
11344 ENDDO
11345 DEALLOCATE(v7d_tmp%volanab)
11346! trasferisco le variabili
11347 v7d_tmp%anavar%r => v7d_tmp%anavar%b
11348 NULLIFY(v7d_tmp%anavar%b)
11349 ENDIF
11350
11351 IF (ASSOCIATED(v7d_tmp%dativar%b)) THEN
11352! alloco i dati reali e vi trasferisco i byte
11353 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatib, 1), SIZE(v7d_tmp%voldatib, 2), &
11354 SIZE(v7d_tmp%voldatib, 3), SIZE(v7d_tmp%voldatib, 4), SIZE(v7d_tmp%voldatib, 5), &
11355 SIZE(v7d_tmp%voldatib, 6)))
11356 DO i = 1, SIZE(v7d_tmp%dativar%b)
11357 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11358 realdat(v7d_tmp%voldatib(:,:,:,:,i,:), v7d_tmp%dativar%b(i))
11359 ENDDO
11360 DEALLOCATE(v7d_tmp%voldatib)
11361! trasferisco le variabili
11362 v7d_tmp%dativar%r => v7d_tmp%dativar%b
11363 NULLIFY(v7d_tmp%dativar%b)
11364 ENDIF
11365
11366! fondo con il volume definitivo
11367 CALL vol7d_merge(that, v7d_tmp)
11368ELSE
11370ENDIF
11371
11372
11373! Volume solo di dati character
11374CALL vol7d_copy(this, v7d_tmp, &
11375 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=fv, lanavarc=acn, &
11376 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11377 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11378 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=tv, &
11379 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11380 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11381
11382! converto a dati reali
11383IF (ASSOCIATED(v7d_tmp%anavar%c) .OR. ASSOCIATED(v7d_tmp%dativar%c)) THEN
11384
11385 IF (ASSOCIATED(v7d_tmp%anavar%c)) THEN
11386! alloco i dati reali e vi trasferisco i character
11387 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanac, 1), SIZE(v7d_tmp%volanac, 2), &
11388 SIZE(v7d_tmp%volanac, 3)))
11389 DO i = 1, SIZE(v7d_tmp%anavar%c)
11390 v7d_tmp%volanar(:,i,:) = &
11391 realdat(v7d_tmp%volanac(:,i,:), v7d_tmp%anavar%c(i))
11392 ENDDO
11393 DEALLOCATE(v7d_tmp%volanac)
11394! trasferisco le variabili
11395 v7d_tmp%anavar%r => v7d_tmp%anavar%c
11396 NULLIFY(v7d_tmp%anavar%c)
11397 ENDIF
11398
11399 IF (ASSOCIATED(v7d_tmp%dativar%c)) THEN
11400! alloco i dati reali e vi trasferisco i character
11401 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatic, 1), SIZE(v7d_tmp%voldatic, 2), &
11402 SIZE(v7d_tmp%voldatic, 3), SIZE(v7d_tmp%voldatic, 4), SIZE(v7d_tmp%voldatic, 5), &
11403 SIZE(v7d_tmp%voldatic, 6)))
11404 DO i = 1, SIZE(v7d_tmp%dativar%c)
11405 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11406 realdat(v7d_tmp%voldatic(:,:,:,:,i,:), v7d_tmp%dativar%c(i))
11407 ENDDO
11408 DEALLOCATE(v7d_tmp%voldatic)
11409! trasferisco le variabili
11410 v7d_tmp%dativar%r => v7d_tmp%dativar%c
11411 NULLIFY(v7d_tmp%dativar%c)
11412 ENDIF
11413
11414! fondo con il volume definitivo
11415 CALL vol7d_merge(that, v7d_tmp)
11416ELSE
11418ENDIF
11419
11420END SUBROUTINE vol7d_convr
11421
11422
11426SUBROUTINE vol7d_diff_only (this, that, data_only,ana)
11427TYPE(vol7d),INTENT(IN) :: this
11428TYPE(vol7d),INTENT(OUT) :: that
11429logical , optional, intent(in) :: data_only
11430logical , optional, intent(in) :: ana
11431logical :: ldata_only,lana
11432
11433IF (PRESENT(data_only)) THEN
11434 ldata_only = data_only
11435ELSE
11436 ldata_only = .false.
11437ENDIF
11438
11439IF (PRESENT(ana)) THEN
11440 lana = ana
11441ELSE
11442 lana = .false.
11443ENDIF
11444
11445
11446#undef VOL7D_POLY_ARRAY
11447#define VOL7D_POLY_ARRAY voldati
11448#include "vol7d_class_diff.F90"
11449#undef VOL7D_POLY_ARRAY
11450#define VOL7D_POLY_ARRAY voldatiattr
11451#include "vol7d_class_diff.F90"
11452#undef VOL7D_POLY_ARRAY
11453
11454if ( .not. ldata_only) then
11455
11456#define VOL7D_POLY_ARRAY volana
11457#include "vol7d_class_diff.F90"
11458#undef VOL7D_POLY_ARRAY
11459#define VOL7D_POLY_ARRAY volanaattr
11460#include "vol7d_class_diff.F90"
11461#undef VOL7D_POLY_ARRAY
11462
11463 if(lana)then
11464 where ( this%ana == that%ana )
11465 that%ana = vol7d_ana_miss
11466 end where
11467 end if
11468
11469end if
11470
11471
11472
11473END SUBROUTINE vol7d_diff_only
11474
11475
11476
11477! Creo le routine da ripetere per i vari tipi di dati di v7d
11478! tramite un template e il preprocessore
11479#undef VOL7D_POLY_TYPE
11480#undef VOL7D_POLY_TYPES
11481#define VOL7D_POLY_TYPE REAL
11482#define VOL7D_POLY_TYPES r
11483#include "vol7d_class_type_templ.F90"
11484#undef VOL7D_POLY_TYPE
11485#undef VOL7D_POLY_TYPES
11486#define VOL7D_POLY_TYPE DOUBLE PRECISION
11487#define VOL7D_POLY_TYPES d
11488#include "vol7d_class_type_templ.F90"
11489#undef VOL7D_POLY_TYPE
11490#undef VOL7D_POLY_TYPES
11491#define VOL7D_POLY_TYPE INTEGER
11492#define VOL7D_POLY_TYPES i
11493#include "vol7d_class_type_templ.F90"
11494#undef VOL7D_POLY_TYPE
11495#undef VOL7D_POLY_TYPES
11496#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
11497#define VOL7D_POLY_TYPES b
11498#include "vol7d_class_type_templ.F90"
11499#undef VOL7D_POLY_TYPE
11500#undef VOL7D_POLY_TYPES
11501#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
11502#define VOL7D_POLY_TYPES c
11503#include "vol7d_class_type_templ.F90"
11504
11505! Creo le routine da ripetere per i vari descrittori di dimensioni di v7d
11506! tramite un template e il preprocessore
11507#define VOL7D_SORT
11508#undef VOL7D_NO_ZERO_ALLOC
11509#undef VOL7D_POLY_TYPE
11510#define VOL7D_POLY_TYPE datetime
11511#include "vol7d_class_desc_templ.F90"
11512#undef VOL7D_POLY_TYPE
11513#define VOL7D_POLY_TYPE vol7d_timerange
11514#include "vol7d_class_desc_templ.F90"
11515#undef VOL7D_POLY_TYPE
11516#define VOL7D_POLY_TYPE vol7d_level
11517#include "vol7d_class_desc_templ.F90"
11518#undef VOL7D_SORT
11519#undef VOL7D_POLY_TYPE
11520#define VOL7D_POLY_TYPE vol7d_network
11521#include "vol7d_class_desc_templ.F90"
11522#undef VOL7D_POLY_TYPE
11523#define VOL7D_POLY_TYPE vol7d_ana
11524#include "vol7d_class_desc_templ.F90"
11525#define VOL7D_NO_ZERO_ALLOC
11526#undef VOL7D_POLY_TYPE
11527#define VOL7D_POLY_TYPE vol7d_var
11528#include "vol7d_class_desc_templ.F90"
11529
11539subroutine vol7d_write_on_file (this,unit,description,filename,filename_auto)
11540
11541TYPE(vol7d),INTENT(IN) :: this
11542integer,optional,intent(inout) :: unit
11543character(len=*),intent(in),optional :: filename
11544character(len=*),intent(out),optional :: filename_auto
11545character(len=*),INTENT(IN),optional :: description
11546
11547integer :: lunit
11548character(len=254) :: ldescription,arg,lfilename
11549integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
11550 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11551 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11552 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11553 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11554 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11555 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
11556!integer :: im,id,iy
11557integer :: tarray(8)
11558logical :: opened,exist
11559
11560 nana=0
11561 ntime=0
11562 ntimerange=0
11563 nlevel=0
11564 nnetwork=0
11565 ndativarr=0
11566 ndativari=0
11567 ndativarb=0
11568 ndativard=0
11569 ndativarc=0
11570 ndatiattrr=0
11571 ndatiattri=0
11572 ndatiattrb=0
11573 ndatiattrd=0
11574 ndatiattrc=0
11575 ndativarattrr=0
11576 ndativarattri=0
11577 ndativarattrb=0
11578 ndativarattrd=0
11579 ndativarattrc=0
11580 nanavarr=0
11581 nanavari=0
11582 nanavarb=0
11583 nanavard=0
11584 nanavarc=0
11585 nanaattrr=0
11586 nanaattri=0
11587 nanaattrb=0
11588 nanaattrd=0
11589 nanaattrc=0
11590 nanavarattrr=0
11591 nanavarattri=0
11592 nanavarattrb=0
11593 nanavarattrd=0
11594 nanavarattrc=0
11595
11596
11597!call idate(im,id,iy)
11598call date_and_time(values=tarray)
11599call getarg(0,arg)
11600
11601if (present(description))then
11602 ldescription=description
11603else
11604 ldescription="Vol7d generated by: "//trim(arg)
11605end if
11606
11607if (.not. present(unit))then
11608 lunit=getunit()
11609else
11610 if (unit==0)then
11611 lunit=getunit()
11612 unit=lunit
11613 else
11614 lunit=unit
11615 end if
11616end if
11617
11618lfilename=trim(arg)//".v7d"
11620
11621if (present(filename))then
11622 if (filename /= "")then
11623 lfilename=filename
11624 end if
11625end if
11626
11627if (present(filename_auto))filename_auto=lfilename
11628
11629
11630inquire(unit=lunit,opened=opened)
11631if (.not. opened) then
11632! inquire(file=lfilename, EXIST=exist)
11633! IF (exist) THEN
11634! CALL l4f_log(L4F_FATAL, &
11635! 'in vol7d_write_on_file, file exists, cannot open file '//TRIM(lfilename))
11636! CALL raise_fatal_error()
11637! ENDIF
11638 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM')
11639 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
11640end if
11641
11642if (associated(this%ana)) nana=size(this%ana)
11643if (associated(this%time)) ntime=size(this%time)
11644if (associated(this%timerange)) ntimerange=size(this%timerange)
11645if (associated(this%level)) nlevel=size(this%level)
11646if (associated(this%network)) nnetwork=size(this%network)
11647
11648if (associated(this%dativar%r)) ndativarr=size(this%dativar%r)
11649if (associated(this%dativar%i)) ndativari=size(this%dativar%i)
11650if (associated(this%dativar%b)) ndativarb=size(this%dativar%b)
11651if (associated(this%dativar%d)) ndativard=size(this%dativar%d)
11652if (associated(this%dativar%c)) ndativarc=size(this%dativar%c)
11653
11654if (associated(this%datiattr%r)) ndatiattrr=size(this%datiattr%r)
11655if (associated(this%datiattr%i)) ndatiattri=size(this%datiattr%i)
11656if (associated(this%datiattr%b)) ndatiattrb=size(this%datiattr%b)
11657if (associated(this%datiattr%d)) ndatiattrd=size(this%datiattr%d)
11658if (associated(this%datiattr%c)) ndatiattrc=size(this%datiattr%c)
11659
11660if (associated(this%dativarattr%r)) ndativarattrr=size(this%dativarattr%r)
11661if (associated(this%dativarattr%i)) ndativarattri=size(this%dativarattr%i)
11662if (associated(this%dativarattr%b)) ndativarattrb=size(this%dativarattr%b)
11663if (associated(this%dativarattr%d)) ndativarattrd=size(this%dativarattr%d)
11664if (associated(this%dativarattr%c)) ndativarattrc=size(this%dativarattr%c)
11665
11666if (associated(this%anavar%r)) nanavarr=size(this%anavar%r)
11667if (associated(this%anavar%i)) nanavari=size(this%anavar%i)
11668if (associated(this%anavar%b)) nanavarb=size(this%anavar%b)
11669if (associated(this%anavar%d)) nanavard=size(this%anavar%d)
11670if (associated(this%anavar%c)) nanavarc=size(this%anavar%c)
11671
11672if (associated(this%anaattr%r)) nanaattrr=size(this%anaattr%r)
11673if (associated(this%anaattr%i)) nanaattri=size(this%anaattr%i)
11674if (associated(this%anaattr%b)) nanaattrb=size(this%anaattr%b)
11675if (associated(this%anaattr%d)) nanaattrd=size(this%anaattr%d)
11676if (associated(this%anaattr%c)) nanaattrc=size(this%anaattr%c)
11677
11678if (associated(this%anavarattr%r)) nanavarattrr=size(this%anavarattr%r)
11679if (associated(this%anavarattr%i)) nanavarattri=size(this%anavarattr%i)
11680if (associated(this%anavarattr%b)) nanavarattrb=size(this%anavarattr%b)
11681if (associated(this%anavarattr%d)) nanavarattrd=size(this%anavarattr%d)
11682if (associated(this%anavarattr%c)) nanavarattrc=size(this%anavarattr%c)
11683
11684write(unit=lunit)ldescription
11685write(unit=lunit)tarray
11686
11687write(unit=lunit)&
11688 nana, ntime, ntimerange, nlevel, nnetwork, &
11689 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11690 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11691 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11692 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11693 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11694 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
11695 this%time_definition
11696
11697
11698!write(unit=lunit)this
11699
11700
11701!! prime 5 dimensioni
11704if (associated(this%level)) write(unit=lunit)this%level
11705if (associated(this%timerange)) write(unit=lunit)this%timerange
11706if (associated(this%network)) write(unit=lunit)this%network
11707
11708 !! 6a dimensione: variabile dell'anagrafica e dei dati
11709 !! con relativi attributi e in 5 tipi diversi
11710
11711if (associated(this%anavar%r)) write(unit=lunit)this%anavar%r
11712if (associated(this%anavar%i)) write(unit=lunit)this%anavar%i
11713if (associated(this%anavar%b)) write(unit=lunit)this%anavar%b
11714if (associated(this%anavar%d)) write(unit=lunit)this%anavar%d
11715if (associated(this%anavar%c)) write(unit=lunit)this%anavar%c
11716
11717if (associated(this%anaattr%r)) write(unit=lunit)this%anaattr%r
11718if (associated(this%anaattr%i)) write(unit=lunit)this%anaattr%i
11719if (associated(this%anaattr%b)) write(unit=lunit)this%anaattr%b
11720if (associated(this%anaattr%d)) write(unit=lunit)this%anaattr%d
11721if (associated(this%anaattr%c)) write(unit=lunit)this%anaattr%c
11722
11723if (associated(this%anavarattr%r)) write(unit=lunit)this%anavarattr%r
11724if (associated(this%anavarattr%i)) write(unit=lunit)this%anavarattr%i
11725if (associated(this%anavarattr%b)) write(unit=lunit)this%anavarattr%b
11726if (associated(this%anavarattr%d)) write(unit=lunit)this%anavarattr%d
11727if (associated(this%anavarattr%c)) write(unit=lunit)this%anavarattr%c
11728
11729if (associated(this%dativar%r)) write(unit=lunit)this%dativar%r
11730if (associated(this%dativar%i)) write(unit=lunit)this%dativar%i
11731if (associated(this%dativar%b)) write(unit=lunit)this%dativar%b
11732if (associated(this%dativar%d)) write(unit=lunit)this%dativar%d
11733if (associated(this%dativar%c)) write(unit=lunit)this%dativar%c
11734
11735if (associated(this%datiattr%r)) write(unit=lunit)this%datiattr%r
11736if (associated(this%datiattr%i)) write(unit=lunit)this%datiattr%i
11737if (associated(this%datiattr%b)) write(unit=lunit)this%datiattr%b
11738if (associated(this%datiattr%d)) write(unit=lunit)this%datiattr%d
11739if (associated(this%datiattr%c)) write(unit=lunit)this%datiattr%c
11740
11741if (associated(this%dativarattr%r)) write(unit=lunit)this%dativarattr%r
11742if (associated(this%dativarattr%i)) write(unit=lunit)this%dativarattr%i
11743if (associated(this%dativarattr%b)) write(unit=lunit)this%dativarattr%b
11744if (associated(this%dativarattr%d)) write(unit=lunit)this%dativarattr%d
11745if (associated(this%dativarattr%c)) write(unit=lunit)this%dativarattr%c
11746
11747!! Volumi di valori e attributi per anagrafica e dati
11748
11749if (associated(this%volanar)) write(unit=lunit)this%volanar
11750if (associated(this%volanaattrr)) write(unit=lunit)this%volanaattrr
11751if (associated(this%voldatir)) write(unit=lunit)this%voldatir
11752if (associated(this%voldatiattrr)) write(unit=lunit)this%voldatiattrr
11753
11754if (associated(this%volanai)) write(unit=lunit)this%volanai
11755if (associated(this%volanaattri)) write(unit=lunit)this%volanaattri
11756if (associated(this%voldatii)) write(unit=lunit)this%voldatii
11757if (associated(this%voldatiattri)) write(unit=lunit)this%voldatiattri
11758
11759if (associated(this%volanab)) write(unit=lunit)this%volanab
11760if (associated(this%volanaattrb)) write(unit=lunit)this%volanaattrb
11761if (associated(this%voldatib)) write(unit=lunit)this%voldatib
11762if (associated(this%voldatiattrb)) write(unit=lunit)this%voldatiattrb
11763
11764if (associated(this%volanad)) write(unit=lunit)this%volanad
11765if (associated(this%volanaattrd)) write(unit=lunit)this%volanaattrd
11766if (associated(this%voldatid)) write(unit=lunit)this%voldatid
11767if (associated(this%voldatiattrd)) write(unit=lunit)this%voldatiattrd
11768
11769if (associated(this%volanac)) write(unit=lunit)this%volanac
11770if (associated(this%volanaattrc)) write(unit=lunit)this%volanaattrc
11771if (associated(this%voldatic)) write(unit=lunit)this%voldatic
11772if (associated(this%voldatiattrc)) write(unit=lunit)this%voldatiattrc
11773
11774if (.not. present(unit)) close(unit=lunit)
11775
11776end subroutine vol7d_write_on_file
11777
11778
11785
11786
11787subroutine vol7d_read_from_file (this,unit,filename,description,tarray,filename_auto)
11788
11789TYPE(vol7d),INTENT(OUT) :: this
11790integer,intent(inout),optional :: unit
11791character(len=*),INTENT(in),optional :: filename
11792character(len=*),intent(out),optional :: filename_auto
11793character(len=*),INTENT(out),optional :: description
11794integer,intent(out),optional :: tarray(8)
11795
11796
11797integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
11798 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11799 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11800 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11801 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11802 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11803 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
11804
11805character(len=254) :: ldescription,lfilename,arg
11806integer :: ltarray(8),lunit,ios
11807logical :: opened,exist
11808
11809
11810call getarg(0,arg)
11811
11812if (.not. present(unit))then
11813 lunit=getunit()
11814else
11815 if (unit==0)then
11816 lunit=getunit()
11817 unit=lunit
11818 else
11819 lunit=unit
11820 end if
11821end if
11822
11823lfilename=trim(arg)//".v7d"
11825
11826if (present(filename))then
11827 if (filename /= "")then
11828 lfilename=filename
11829 end if
11830end if
11831
11832if (present(filename_auto))filename_auto=lfilename
11833
11834
11835inquire(unit=lunit,opened=opened)
11836IF (.NOT. opened) THEN
11837 inquire(file=lfilename,exist=exist)
11838 IF (.NOT.exist) THEN
11839 CALL l4f_log(l4f_fatal, &
11840 'in vol7d_read_from_file, file does not exists, cannot open')
11841 CALL raise_fatal_error()
11842 ENDIF
11843 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM', &
11844 status='OLD', action='READ')
11845 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
11846end if
11847
11848
11850read(unit=lunit,iostat=ios)ldescription
11851
11852if (ios < 0) then ! A negative value indicates that the End of File or End of Record
11853 call vol7d_alloc (this)
11854 call vol7d_alloc_vol (this)
11855 if (present(description))description=ldescription
11856 if (present(tarray))tarray=ltarray
11857 if (.not. present(unit)) close(unit=lunit)
11858end if
11859
11860read(unit=lunit)ltarray
11861
11862CALL l4f_log(l4f_info, 'Reading vol7d from file')
11863CALL l4f_log(l4f_info, 'description: '//trim(ldescription))
11866
11867if (present(description))description=ldescription
11868if (present(tarray))tarray=ltarray
11869
11870read(unit=lunit)&
11871 nana, ntime, ntimerange, nlevel, nnetwork, &
11872 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11873 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11874 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11875 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11876 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11877 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
11878 this%time_definition
11879
11880call vol7d_alloc (this, &
11881 nana=nana, ntime=ntime, ntimerange=ntimerange, nlevel=nlevel, nnetwork=nnetwork,&
11882 ndativarr=ndativarr, ndativari=ndativari, ndativarb=ndativarb,&
11883 ndativard=ndativard, ndativarc=ndativarc,&
11884 ndatiattrr=ndatiattrr, ndatiattri=ndatiattri, ndatiattrb=ndatiattrb,&
11885 ndatiattrd=ndatiattrd, ndatiattrc=ndatiattrc,&
11886 ndativarattrr=ndativarattrr, ndativarattri=ndativarattri, ndativarattrb=ndativarattrb, &
11887 ndativarattrd=ndativarattrd, ndativarattrc=ndativarattrc,&
11888 nanavarr=nanavarr, nanavari=nanavari, nanavarb=nanavarb, &
11889 nanavard=nanavard, nanavarc=nanavarc,&
11890 nanaattrr=nanaattrr, nanaattri=nanaattri, nanaattrb=nanaattrb,&
11891 nanaattrd=nanaattrd, nanaattrc=nanaattrc,&
11892 nanavarattrr=nanavarattrr, nanavarattri=nanavarattri, nanavarattrb=nanavarattrb, &
11893 nanavarattrd=nanavarattrd, nanavarattrc=nanavarattrc)
11894
11895
11898if (associated(this%level)) read(unit=lunit)this%level
11899if (associated(this%timerange)) read(unit=lunit)this%timerange
11900if (associated(this%network)) read(unit=lunit)this%network
11901
11902if (associated(this%anavar%r)) read(unit=lunit)this%anavar%r
11903if (associated(this%anavar%i)) read(unit=lunit)this%anavar%i
11904if (associated(this%anavar%b)) read(unit=lunit)this%anavar%b
11905if (associated(this%anavar%d)) read(unit=lunit)this%anavar%d
11906if (associated(this%anavar%c)) read(unit=lunit)this%anavar%c
11907
11908if (associated(this%anaattr%r)) read(unit=lunit)this%anaattr%r
11909if (associated(this%anaattr%i)) read(unit=lunit)this%anaattr%i
11910if (associated(this%anaattr%b)) read(unit=lunit)this%anaattr%b
11911if (associated(this%anaattr%d)) read(unit=lunit)this%anaattr%d
11912if (associated(this%anaattr%c)) read(unit=lunit)this%anaattr%c
11913
11914if (associated(this%anavarattr%r)) read(unit=lunit)this%anavarattr%r
11915if (associated(this%anavarattr%i)) read(unit=lunit)this%anavarattr%i
11916if (associated(this%anavarattr%b)) read(unit=lunit)this%anavarattr%b
11917if (associated(this%anavarattr%d)) read(unit=lunit)this%anavarattr%d
11918if (associated(this%anavarattr%c)) read(unit=lunit)this%anavarattr%c
11919
11920if (associated(this%dativar%r)) read(unit=lunit)this%dativar%r
11921if (associated(this%dativar%i)) read(unit=lunit)this%dativar%i
11922if (associated(this%dativar%b)) read(unit=lunit)this%dativar%b
11923if (associated(this%dativar%d)) read(unit=lunit)this%dativar%d
11924if (associated(this%dativar%c)) read(unit=lunit)this%dativar%c
11925
11926if (associated(this%datiattr%r)) read(unit=lunit)this%datiattr%r
11927if (associated(this%datiattr%i)) read(unit=lunit)this%datiattr%i
11928if (associated(this%datiattr%b)) read(unit=lunit)this%datiattr%b
11929if (associated(this%datiattr%d)) read(unit=lunit)this%datiattr%d
11930if (associated(this%datiattr%c)) read(unit=lunit)this%datiattr%c
11931
11932if (associated(this%dativarattr%r)) read(unit=lunit)this%dativarattr%r
11933if (associated(this%dativarattr%i)) read(unit=lunit)this%dativarattr%i
11934if (associated(this%dativarattr%b)) read(unit=lunit)this%dativarattr%b
11935if (associated(this%dativarattr%d)) read(unit=lunit)this%dativarattr%d
11936if (associated(this%dativarattr%c)) read(unit=lunit)this%dativarattr%c
11937
11938call vol7d_alloc_vol (this)
11939
11940!! Volumi di valori e attributi per anagrafica e dati
11941
11942if (associated(this%volanar)) read(unit=lunit)this%volanar
11943if (associated(this%volanaattrr)) read(unit=lunit)this%volanaattrr
11944if (associated(this%voldatir)) read(unit=lunit)this%voldatir
11945if (associated(this%voldatiattrr)) read(unit=lunit)this%voldatiattrr
11946
11947if (associated(this%volanai)) read(unit=lunit)this%volanai
11948if (associated(this%volanaattri)) read(unit=lunit)this%volanaattri
11949if (associated(this%voldatii)) read(unit=lunit)this%voldatii
11950if (associated(this%voldatiattri)) read(unit=lunit)this%voldatiattri
11951
11952if (associated(this%volanab)) read(unit=lunit)this%volanab
11953if (associated(this%volanaattrb)) read(unit=lunit)this%volanaattrb
11954if (associated(this%voldatib)) read(unit=lunit)this%voldatib
11955if (associated(this%voldatiattrb)) read(unit=lunit)this%voldatiattrb
11956
11957if (associated(this%volanad)) read(unit=lunit)this%volanad
11958if (associated(this%volanaattrd)) read(unit=lunit)this%volanaattrd
11959if (associated(this%voldatid)) read(unit=lunit)this%voldatid
11960if (associated(this%voldatiattrd)) read(unit=lunit)this%voldatiattrd
11961
11962if (associated(this%volanac)) read(unit=lunit)this%volanac
11963if (associated(this%volanaattrc)) read(unit=lunit)this%volanaattrc
11964if (associated(this%voldatic)) read(unit=lunit)this%voldatic
11965if (associated(this%voldatiattrc)) read(unit=lunit)this%voldatiattrc
11966
11967if (.not. present(unit)) close(unit=lunit)
11968
11969end subroutine vol7d_read_from_file
11970
11971
11972! to double precision
11973elemental doubleprecision function doubledatd(voldat,var)
11974doubleprecision,intent(in) :: voldat
11975type(vol7d_var),intent(in) :: var
11976
11977doubledatd=voldat
11978
11979end function doubledatd
11980
11981
11982elemental doubleprecision function doubledatr(voldat,var)
11983real,intent(in) :: voldat
11984type(vol7d_var),intent(in) :: var
11985
11987 doubledatr=dble(voldat)
11988else
11989 doubledatr=dmiss
11990end if
11991
11992end function doubledatr
11993
11994
11995elemental doubleprecision function doubledati(voldat,var)
11996integer,intent(in) :: voldat
11997type(vol7d_var),intent(in) :: var
11998
12001 doubledati=dble(voldat)/10.d0**var%scalefactor
12002 else
12003 doubledati=dble(voldat)
12004 endif
12005else
12006 doubledati=dmiss
12007end if
12008
12009end function doubledati
12010
12011
12012elemental doubleprecision function doubledatb(voldat,var)
12013integer(kind=int_b),intent(in) :: voldat
12014type(vol7d_var),intent(in) :: var
12015
12018 doubledatb=dble(voldat)/10.d0**var%scalefactor
12019 else
12020 doubledatb=dble(voldat)
12021 endif
12022else
12023 doubledatb=dmiss
12024end if
12025
12026end function doubledatb
12027
12028
12029elemental doubleprecision function doubledatc(voldat,var)
12030CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
12031type(vol7d_var),intent(in) :: var
12032
12033doubledatc = c2d(voldat)
12035 doubledatc=doubledatc/10.d0**var%scalefactor
12036end if
12037
12038end function doubledatc
12039
12040
12041! to integer
12042elemental integer function integerdatd(voldat,var)
12043doubleprecision,intent(in) :: voldat
12044type(vol7d_var),intent(in) :: var
12045
12048 integerdatd=nint(voldat*10d0**var%scalefactor)
12049 else
12050 integerdatd=nint(voldat)
12051 endif
12052else
12053 integerdatd=imiss
12054end if
12055
12056end function integerdatd
12057
12058
12059elemental integer function integerdatr(voldat,var)
12060real,intent(in) :: voldat
12061type(vol7d_var),intent(in) :: var
12062
12065 integerdatr=nint(voldat*10d0**var%scalefactor)
12066 else
12067 integerdatr=nint(voldat)
12068 endif
12069else
12070 integerdatr=imiss
12071end if
12072
12073end function integerdatr
12074
12075
12076elemental integer function integerdati(voldat,var)
12077integer,intent(in) :: voldat
12078type(vol7d_var),intent(in) :: var
12079
12080integerdati=voldat
12081
12082end function integerdati
12083
12084
12085elemental integer function integerdatb(voldat,var)
12086integer(kind=int_b),intent(in) :: voldat
12087type(vol7d_var),intent(in) :: var
12088
12090 integerdatb=voldat
12091else
12092 integerdatb=imiss
12093end if
12094
12095end function integerdatb
12096
12097
12098elemental integer function integerdatc(voldat,var)
12099CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
12100type(vol7d_var),intent(in) :: var
12101
12102integerdatc=c2i(voldat)
12103
12104end function integerdatc
12105
12106
12107! to real
12108elemental real function realdatd(voldat,var)
12109doubleprecision,intent(in) :: voldat
12110type(vol7d_var),intent(in) :: var
12111
12113 realdatd=real(voldat)
12114else
12115 realdatd=rmiss
12116end if
12117
12118end function realdatd
12119
12120
12121elemental real function realdatr(voldat,var)
12122real,intent(in) :: voldat
12123type(vol7d_var),intent(in) :: var
12124
12125realdatr=voldat
12126
12127end function realdatr
12128
12129
12130elemental real function realdati(voldat,var)
12131integer,intent(in) :: voldat
12132type(vol7d_var),intent(in) :: var
12133
12136 realdati=float(voldat)/10.**var%scalefactor
12137 else
12138 realdati=float(voldat)
12139 endif
12140else
12141 realdati=rmiss
12142end if
12143
12144end function realdati
12145
12146
12147elemental real function realdatb(voldat,var)
12148integer(kind=int_b),intent(in) :: voldat
12149type(vol7d_var),intent(in) :: var
12150
12153 realdatb=float(voldat)/10**var%scalefactor
12154 else
12155 realdatb=float(voldat)
12156 endif
12157else
12158 realdatb=rmiss
12159end if
12160
12161end function realdatb
12162
12163
12164elemental real function realdatc(voldat,var)
12165CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
12166type(vol7d_var),intent(in) :: var
12167
12168realdatc=c2r(voldat)
12170 realdatc=realdatc/10.**var%scalefactor
12171end if
12172
12173end function realdatc
12174
12175
12181FUNCTION realanavol(this, var) RESULT(vol)
12182TYPE(vol7d),INTENT(in) :: this
12183TYPE(vol7d_var),INTENT(in) :: var
12184REAL :: vol(SIZE(this%ana),size(this%network))
12185
12186CHARACTER(len=1) :: dtype
12187INTEGER :: indvar
12188
12189dtype = cmiss
12190indvar = index(this%anavar, var, type=dtype)
12191
12192IF (indvar > 0) THEN
12193 SELECT CASE (dtype)
12194 CASE("d")
12195 vol = realdat(this%volanad(:,indvar,:), var)
12196 CASE("r")
12197 vol = this%volanar(:,indvar,:)
12198 CASE("i")
12199 vol = realdat(this%volanai(:,indvar,:), var)
12200 CASE("b")
12201 vol = realdat(this%volanab(:,indvar,:), var)
12202 CASE("c")
12203 vol = realdat(this%volanac(:,indvar,:), var)
12204 CASE default
12205 vol = rmiss
12206 END SELECT
12207ELSE
12208 vol = rmiss
12209ENDIF
12210
12211END FUNCTION realanavol
12212
12213
12219FUNCTION integeranavol(this, var) RESULT(vol)
12220TYPE(vol7d),INTENT(in) :: this
12221TYPE(vol7d_var),INTENT(in) :: var
12222INTEGER :: vol(SIZE(this%ana),size(this%network))
12223
12224CHARACTER(len=1) :: dtype
12225INTEGER :: indvar
12226
12227dtype = cmiss
12228indvar = index(this%anavar, var, type=dtype)
12229
12230IF (indvar > 0) THEN
12231 SELECT CASE (dtype)
12232 CASE("d")
12233 vol = integerdat(this%volanad(:,indvar,:), var)
12234 CASE("r")
12235 vol = integerdat(this%volanar(:,indvar,:), var)
12236 CASE("i")
12237 vol = this%volanai(:,indvar,:)
12238 CASE("b")
12239 vol = integerdat(this%volanab(:,indvar,:), var)
12240 CASE("c")
12241 vol = integerdat(this%volanac(:,indvar,:), var)
12242 CASE default
12243 vol = imiss
12244 END SELECT
12245ELSE
12246 vol = imiss
12247ENDIF
12248
12249END FUNCTION integeranavol
12250
12251
12257subroutine move_datac (v7d,&
12258 indana,indtime,indlevel,indtimerange,indnetwork,&
12259 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
12260
12261TYPE(vol7d),intent(inout) :: v7d
12262
12263integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
12264integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
12265integer :: inddativar,inddativarattr
12266
12267
12268do inddativar=1,size(v7d%dativar%c)
12269
12271 .not. c_e(v7d%voldatic(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
12272 ) then
12273
12274 ! dati
12275 v7d%voldatic &
12276 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
12277 v7d%voldatic &
12278 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
12279
12280
12281 ! attributi
12282 if (associated (v7d%dativarattr%i)) then
12283 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%c(inddativar))
12284 if (inddativarattr > 0 ) then
12285 v7d%voldatiattri &
12286 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12287 v7d%voldatiattri &
12288 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12289 end if
12290 end if
12291
12292 if (associated (v7d%dativarattr%r)) then
12293 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%c(inddativar))
12294 if (inddativarattr > 0 ) then
12295 v7d%voldatiattrr &
12296 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12297 v7d%voldatiattrr &
12298 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12299 end if
12300 end if
12301
12302 if (associated (v7d%dativarattr%d)) then
12303 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%c(inddativar))
12304 if (inddativarattr > 0 ) then
12305 v7d%voldatiattrd &
12306 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12307 v7d%voldatiattrd &
12308 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12309 end if
12310 end if
12311
12312 if (associated (v7d%dativarattr%b)) then
12313 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%c(inddativar))
12314 if (inddativarattr > 0 ) then
12315 v7d%voldatiattrb &
12316 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12317 v7d%voldatiattrb &
12318 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12319 end if
12320 end if
12321
12322 if (associated (v7d%dativarattr%c)) then
12323 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%c(inddativar))
12324 if (inddativarattr > 0 ) then
12325 v7d%voldatiattrc &
12326 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12327 v7d%voldatiattrc &
12328 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12329 end if
12330 end if
12331
12332 end if
12333
12334end do
12335
12336end subroutine move_datac
12337
12343subroutine move_datar (v7d,&
12344 indana,indtime,indlevel,indtimerange,indnetwork,&
12345 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
12346
12347TYPE(vol7d),intent(inout) :: v7d
12348
12349integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
12350integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
12351integer :: inddativar,inddativarattr
12352
12353
12354do inddativar=1,size(v7d%dativar%r)
12355
12357 .not. c_e(v7d%voldatir(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
12358 ) then
12359
12360 ! dati
12361 v7d%voldatir &
12362 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
12363 v7d%voldatir &
12364 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
12365
12366
12367 ! attributi
12368 if (associated (v7d%dativarattr%i)) then
12369 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%r(inddativar))
12370 if (inddativarattr > 0 ) then
12371 v7d%voldatiattri &
12372 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12373 v7d%voldatiattri &
12374 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12375 end if
12376 end if
12377
12378 if (associated (v7d%dativarattr%r)) then
12379 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%r(inddativar))
12380 if (inddativarattr > 0 ) then
12381 v7d%voldatiattrr &
12382 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12383 v7d%voldatiattrr &
12384 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12385 end if
12386 end if
12387
12388 if (associated (v7d%dativarattr%d)) then
12389 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%r(inddativar))
12390 if (inddativarattr > 0 ) then
12391 v7d%voldatiattrd &
12392 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12393 v7d%voldatiattrd &
12394 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12395 end if
12396 end if
12397
12398 if (associated (v7d%dativarattr%b)) then
12399 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%r(inddativar))
12400 if (inddativarattr > 0 ) then
12401 v7d%voldatiattrb &
12402 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12403 v7d%voldatiattrb &
12404 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12405 end if
12406 end if
12407
12408 if (associated (v7d%dativarattr%c)) then
12409 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%r(inddativar))
12410 if (inddativarattr > 0 ) then
12411 v7d%voldatiattrc &
12412 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12413 v7d%voldatiattrc &
12414 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12415 end if
12416 end if
12417
12418 end if
12419
12420end do
12421
12422end subroutine move_datar
12423
12424
12438subroutine v7d_rounding(v7din,v7dout,level,timerange,nostatproc)
12439type(vol7d),intent(inout) :: v7din
12440type(vol7d),intent(out) :: v7dout
12441type(vol7d_level),intent(in),optional :: level(:)
12442type(vol7d_timerange),intent(in),optional :: timerange(:)
12443!logical,intent(in),optional :: merge !< if there are data on more then one almost equal levels or timeranges
12444!! will be merged POINT BY POINT with priority for the fird data found ordered by icreasing var index
12445logical,intent(in),optional :: nostatproc
12446
12447integer :: nana,nlevel,ntime,ntimerange,nnetwork,nbin
12448integer :: iana,ilevel,itimerange,indl,indt,itime,inetwork
12449type(vol7d_level) :: roundlevel(size(v7din%level))
12450type(vol7d_timerange) :: roundtimerange(size(v7din%timerange))
12451type(vol7d) :: v7d_tmp
12452
12453
12454nbin=0
12455
12456if (associated(v7din%dativar%r)) nbin = nbin + size(v7din%dativar%r)
12457if (associated(v7din%dativar%i)) nbin = nbin + size(v7din%dativar%i)
12458if (associated(v7din%dativar%d)) nbin = nbin + size(v7din%dativar%d)
12459if (associated(v7din%dativar%b)) nbin = nbin + size(v7din%dativar%b)
12460
12462
12463roundlevel=v7din%level
12464
12465if (present(level))then
12466 do ilevel = 1, size(v7din%level)
12467 if ((any(v7din%level(ilevel) .almosteq. level))) then
12468 roundlevel(ilevel)=level(1)
12469 end if
12470 end do
12471end if
12472
12473roundtimerange=v7din%timerange
12474
12475if (present(timerange))then
12476 do itimerange = 1, size(v7din%timerange)
12477 if ((any(v7din%timerange(itimerange) .almosteq. timerange))) then
12478 roundtimerange(itimerange)=timerange(1)
12479 end if
12480 end do
12481end if
12482
12483!set istantaneous values everywere
12484!preserve p1 for forecast time
12485if (optio_log(nostatproc)) then
12486 roundtimerange(:)%timerange=254
12487 roundtimerange(:)%p2=0
12488end if
12489
12490
12491nana=size(v7din%ana)
12492nlevel=count_distinct(roundlevel,back=.true.)
12493ntime=size(v7din%time)
12494ntimerange=count_distinct(roundtimerange,back=.true.)
12495nnetwork=size(v7din%network)
12496
12498
12499if (nbin == 0) then
12501else
12502 call vol7d_convr(v7din,v7d_tmp)
12503end if
12504
12505v7d_tmp%level=roundlevel
12506v7d_tmp%timerange=roundtimerange
12507
12508do ilevel=1, size(v7d_tmp%level)
12509 indl=index(v7d_tmp%level,roundlevel(ilevel))
12510 do itimerange=1,size(v7d_tmp%timerange)
12511 indt=index(v7d_tmp%timerange,roundtimerange(itimerange))
12512
12513 if (indl /= ilevel .or. indt /= itimerange) then
12514
12515 do iana=1, nana
12516 do itime=1,ntime
12517 do inetwork=1,nnetwork
12518
12519 if (nbin > 0) then
12520 call move_datar (v7d_tmp,&
12521 iana,itime,ilevel,itimerange,inetwork,&
12522 iana,itime,indl,indt,inetwork)
12523 else
12524 call move_datac (v7d_tmp,&
12525 iana,itime,ilevel,itimerange,inetwork,&
12526 iana,itime,indl,indt,inetwork)
12527 end if
12528
12529 end do
12530 end do
12531 end do
12532
12533 end if
12534
12535 end do
12536end do
12537
12538! set to missing level and time > nlevel
12539do ilevel=nlevel+1,size(v7d_tmp%level)
12541end do
12542
12543do itimerange=ntimerange+1,size(v7d_tmp%timerange)
12545end do
12546
12547!copy with remove
12550
12551!call display(v7dout)
12552
12553end subroutine v7d_rounding
12554
12555
12557
12563
12564
Set of functions that return a trimmed CHARACTER representation of the input variable. Definition: char_utilities.F90:284 Legge un oggetto datetime/timedelta o un vettore di oggetti datetime/timedelta da un file FORMATTED o... Definition: datetime_class.F90:484 Scrive un oggetto datetime/timedelta o un vettore di oggetti datetime/timedelta su un file FORMATTED ... Definition: datetime_class.F90:491 Generic subroutine for checking OPTIONAL parameters. Definition: optional_values.f90:36 Check for problems return 0 if all check passed print diagnostics with log4f. Definition: vol7d_class.F90:451 Reduce some dimensions (level and timerage) for semplification (rounding). Definition: vol7d_class.F90:468 This module defines usefull general purpose function and subroutine. Definition: array_utilities.F90:218 Definition of constants to be used for declaring variables of a desired type. Definition: kinds.F90:251 Module for quickly interpreting the OPTIONAL parameters passed to a subprogram. Definition: optional_values.f90:28 Classe per la gestione dell'anagrafica di stazioni meteo e affini. Definition: vol7d_ana_class.F90:218 Classe per la gestione di un volume completo di dati osservati. Definition: vol7d_class.F90:279 Classe per la gestione dei livelli verticali in osservazioni meteo e affini. Definition: vol7d_level_class.F90:219 Classe per la gestione delle reti di stazioni per osservazioni meteo e affini. Definition: vol7d_network_class.F90:220 Classe per la gestione degli intervalli temporali di osservazioni meteo e affini. Definition: vol7d_timerange_class.F90:221 Classe per gestire un vettore di oggetti di tipo vol7d_var_class::vol7d_var. Definition: vol7d_varvect_class.f90:22 Definisce un oggetto contenente i volumi anagrafica e dati e tutti i descrittori delle loro dimension... Definition: vol7d_class.F90:318 |