libsim Versione 7.2.1

◆ move_datac()

subroutine move_datac ( type(vol7d), intent(inout)  v7d,
integer, intent(in)  indana,
integer, intent(in)  indtime,
integer, intent(in)  indlevel,
integer, intent(in)  indtimerange,
integer, intent(in)  indnetwork,
integer, intent(in)  indananew,
integer, intent(in)  indtimenew,
integer, intent(in)  indlevelnew,
integer, intent(in)  indtimerangenew,
integer, intent(in)  indnetworknew 
)

Move data for all variables from one coordinate in the character volume to other.

Only not missing data will be copyed and all attributes will be moved together. Usefull to colapse data spread in more indices (level or time or ....). After the move is possible to set to missing some descriptor and make a copy with miss=.true. to obtain a new vol7d with less data shape.

Parametri
[in,out]v7ddata in form of character in this object will be moved
[in]indnetworksource coordinate of the data
[in]indnetworknewdestination coordinate of data

Definizione alla linea 9186 del file vol7d_class.F90.

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

Generated with Doxygen.