libsim Versione 7.2.1

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

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

Generated with Doxygen.