libsim Versione 7.1.11

◆ v7d_rounding()

subroutine v7d_rounding ( type(vol7d), intent(inout)  v7din,
type(vol7d), intent(out)  v7dout,
type(vol7d_level), dimension(:), intent(in), optional  level,
type(vol7d_timerange), dimension(:), intent(in), optional  timerange,
logical, intent(in), optional  nostatproc 
)

Reduce some dimensions (level and timerage) for semplification (rounding).

You can use this for simplify and use variables in computation like alchimia where fields have to be on the same coordinate It return real or character data only: if input is charcter data only it return character otherwise il return all the data converted to real. examples: means in time for short periods and istantaneous values 2 meter and surface levels If there are data on more then one almost equal levels or timeranges, the first var present (at least one point) will be taken (order is by icreasing var index). You can use predefined values for classic semplification almost_equal_levels and almost_equal_timeranges The level or timerange in output will be defined by the first element of level and timerange list

Parametri
[in,out]v7dininput volume
[in]levelalmost equal level list
[in]timerangealmost equal timerange list
[in]nostatprocdo not take in account statistical processing code in timerange and P2

Definizione alla linea 9373 del file vol7d_class.F90.

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