libsim Versione 7.1.11

◆ 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 9192 del file vol7d_class.F90.

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