libsim Versione 7.2.1

◆ vol7d_write_on_file()

subroutine, private vol7d_write_on_file ( type(vol7d), intent(in)  this,
integer, intent(inout), optional  unit,
character(len=*), intent(in), optional  description,
character(len=*), intent(in), optional  filename,
character(len=*), intent(out), optional  filename_auto 
)
private

Scrittura su file di un volume Vol7d.

Scrittura su file unformatted di un intero volume Vol7d. Il volume viene serializzato e scritto su file. Il file puo' essere aperto opzionalmente dall'utente. Si possono passare opzionalmente unità e nome del file altrimenti assegnati internamente a dei default; se assegnati internamente tali parametri saranno in output. Se non viene fornito il nome file viene utilizzato un file di default con nome pari al nome del programma in esecuzione con postfisso ".v7d". Come parametro opzionale c'è la description che insieme alla data corrente viene inserita nell'header del file.

Parametri
[in]thisvolume vol7d da scrivere
[in,out]unitunità su cui scrivere; se passata =0 ritorna il valore rielaborato (default =rielaborato internamente con getlun )
[in]filenamenome del file su cui scrivere
[out]filename_autonome del file generato se "filename" è omesso
[in]descriptiondescrizione del volume

Definizione alla linea 8468 del file vol7d_class.F90.

8469! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
8470! authors:
8471! Davide Cesari <dcesari@arpa.emr.it>
8472! Paolo Patruno <ppatruno@arpa.emr.it>
8473
8474! This program is free software; you can redistribute it and/or
8475! modify it under the terms of the GNU General Public License as
8476! published by the Free Software Foundation; either version 2 of
8477! the License, or (at your option) any later version.
8478
8479! This program is distributed in the hope that it will be useful,
8480! but WITHOUT ANY WARRANTY; without even the implied warranty of
8481! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
8482! GNU General Public License for more details.
8483
8484! You should have received a copy of the GNU General Public License
8485! along with this program. If not, see <http://www.gnu.org/licenses/>.
8486#include "config.h"
8487
8499
8553MODULE vol7d_class
8554USE kinds
8558USE log4fortran
8559USE err_handling
8560USE io_units
8567IMPLICIT NONE
8568
8569
8570INTEGER, PARAMETER :: vol7d_maxdim_a = 3, vol7d_maxdim_aa = 4, &
8571 vol7d_maxdim_d = 6, vol7d_maxdim_ad = 7
8572
8573INTEGER, PARAMETER :: vol7d_ana_a=1
8574INTEGER, PARAMETER :: vol7d_var_a=2
8575INTEGER, PARAMETER :: vol7d_network_a=3
8576INTEGER, PARAMETER :: vol7d_attr_a=4
8577INTEGER, PARAMETER :: vol7d_ana_d=1
8578INTEGER, PARAMETER :: vol7d_time_d=2
8579INTEGER, PARAMETER :: vol7d_level_d=3
8580INTEGER, PARAMETER :: vol7d_timerange_d=4
8581INTEGER, PARAMETER :: vol7d_var_d=5
8582INTEGER, PARAMETER :: vol7d_network_d=6
8583INTEGER, PARAMETER :: vol7d_attr_d=7
8584INTEGER, PARAMETER :: vol7d_cdatalen=32
8585
8586TYPE vol7d_varmap
8587 INTEGER :: r, d, i, b, c
8588END TYPE vol7d_varmap
8589
8592TYPE vol7d
8594 TYPE(vol7d_ana),POINTER :: ana(:)
8596 TYPE(datetime),POINTER :: time(:)
8598 TYPE(vol7d_level),POINTER :: level(:)
8600 TYPE(vol7d_timerange),POINTER :: timerange(:)
8602 TYPE(vol7d_network),POINTER :: network(:)
8604 TYPE(vol7d_varvect) :: anavar
8606 TYPE(vol7d_varvect) :: anaattr
8608 TYPE(vol7d_varvect) :: anavarattr
8610 TYPE(vol7d_varvect) :: dativar
8612 TYPE(vol7d_varvect) :: datiattr
8614 TYPE(vol7d_varvect) :: dativarattr
8615
8617 REAL,POINTER :: volanar(:,:,:)
8619 DOUBLE PRECISION,POINTER :: volanad(:,:,:)
8621 INTEGER,POINTER :: volanai(:,:,:)
8623 INTEGER(kind=int_b),POINTER :: volanab(:,:,:)
8625 CHARACTER(len=vol7d_cdatalen),POINTER :: volanac(:,:,:)
8626
8628 REAL,POINTER :: volanaattrr(:,:,:,:)
8630 DOUBLE PRECISION,POINTER :: volanaattrd(:,:,:,:)
8632 INTEGER,POINTER :: volanaattri(:,:,:,:)
8634 INTEGER(kind=int_b),POINTER :: volanaattrb(:,:,:,:)
8636 CHARACTER(len=vol7d_cdatalen),POINTER :: volanaattrc(:,:,:,:)
8637
8639 REAL,POINTER :: voldatir(:,:,:,:,:,:) ! sono i dati
8641 DOUBLE PRECISION,POINTER :: voldatid(:,:,:,:,:,:)
8643 INTEGER,POINTER :: voldatii(:,:,:,:,:,:)
8645 INTEGER(kind=int_b),POINTER :: voldatib(:,:,:,:,:,:)
8647 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatic(:,:,:,:,:,:)
8648
8650 REAL,POINTER :: voldatiattrr(:,:,:,:,:,:,:)
8652 DOUBLE PRECISION,POINTER :: voldatiattrd(:,:,:,:,:,:,:)
8654 INTEGER,POINTER :: voldatiattri(:,:,:,:,:,:,:)
8656 INTEGER(kind=int_b),POINTER :: voldatiattrb(:,:,:,:,:,:,:)
8658 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatiattrc(:,:,:,:,:,:,:)
8659
8661 integer :: time_definition
8662
8663END TYPE vol7d
8664
8668INTERFACE init
8669 MODULE PROCEDURE vol7d_init
8670END INTERFACE
8671
8673INTERFACE delete
8674 MODULE PROCEDURE vol7d_delete
8675END INTERFACE
8676
8678INTERFACE export
8679 MODULE PROCEDURE vol7d_write_on_file
8680END INTERFACE
8681
8683INTERFACE import
8684 MODULE PROCEDURE vol7d_read_from_file
8685END INTERFACE
8686
8688INTERFACE display
8689 MODULE PROCEDURE vol7d_display, dat_display, dat_vect_display
8690END INTERFACE
8691
8693INTERFACE to_char
8694 MODULE PROCEDURE to_char_dat
8695END INTERFACE
8696
8698INTERFACE doubledat
8699 MODULE PROCEDURE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
8700END INTERFACE
8701
8703INTERFACE realdat
8704 MODULE PROCEDURE realdatd,realdatr,realdati,realdatb,realdatc
8705END INTERFACE
8706
8708INTERFACE integerdat
8709 MODULE PROCEDURE integerdatd,integerdatr,integerdati,integerdatb,integerdatc
8710END INTERFACE
8711
8713INTERFACE copy
8714 MODULE PROCEDURE vol7d_copy
8715END INTERFACE
8716
8718INTERFACE c_e
8719 MODULE PROCEDURE vol7d_c_e
8720END INTERFACE
8721
8725INTERFACE check
8726 MODULE PROCEDURE vol7d_check
8727END INTERFACE
8728
8742INTERFACE rounding
8743 MODULE PROCEDURE v7d_rounding
8744END INTERFACE
8745
8746!!$INTERFACE get_volana
8747!!$ MODULE PROCEDURE vol7d_get_volanar, vol7d_get_volanad, vol7d_get_volanai, &
8748!!$ vol7d_get_volanab, vol7d_get_volanac
8749!!$END INTERFACE
8750!!$
8751!!$INTERFACE get_voldati
8752!!$ MODULE PROCEDURE vol7d_get_voldatir, vol7d_get_voldatid, vol7d_get_voldatii, &
8753!!$ vol7d_get_voldatib, vol7d_get_voldatic
8754!!$END INTERFACE
8755!!$
8756!!$INTERFACE get_volanaattr
8757!!$ MODULE PROCEDURE vol7d_get_volanaattrr, vol7d_get_volanaattrd, &
8758!!$ vol7d_get_volanaattri, vol7d_get_volanaattrb, vol7d_get_volanaattrc
8759!!$END INTERFACE
8760!!$
8761!!$INTERFACE get_voldatiattr
8762!!$ MODULE PROCEDURE vol7d_get_voldatiattrr, vol7d_get_voldatiattrd, &
8763!!$ vol7d_get_voldatiattri, vol7d_get_voldatiattrb, vol7d_get_voldatiattrc
8764!!$END INTERFACE
8765
8766PRIVATE vol7d_get_volr, vol7d_get_vold, vol7d_get_voli, vol7d_get_volb, &
8767 vol7d_get_volc, &
8768 volptr1dr, volptr2dr, volptr3dr, volptr4dr, volptr5dr, volptr6dr, volptr7dr, &
8769 volptr1dd, volptr2dd, volptr3dd, volptr4dd, volptr5dd, volptr6dd, volptr7dd, &
8770 volptr1di, volptr2di, volptr3di, volptr4di, volptr5di, volptr6di, volptr7di, &
8771 volptr1db, volptr2db, volptr3db, volptr4db, volptr5db, volptr6db, volptr7db, &
8772 volptr1dc, volptr2dc, volptr3dc, volptr4dc, volptr5dc, volptr6dc, volptr7dc, &
8773 vol7d_nullifyr, vol7d_nullifyd, vol7d_nullifyi, vol7d_nullifyb, vol7d_nullifyc, &
8774 vol7d_init, vol7d_delete, vol7d_write_on_file, vol7d_read_from_file, &
8775 vol7d_check_alloc_ana, vol7d_force_alloc_ana, &
8776 vol7d_check_alloc_dati, vol7d_force_alloc_dati, vol7d_force_alloc, &
8777 vol7d_display, dat_display, dat_vect_display, &
8778 to_char_dat, vol7d_check
8779
8780PRIVATE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
8781
8782PRIVATE vol7d_c_e
8783
8784CONTAINS
8785
8786
8791SUBROUTINE vol7d_init(this,time_definition)
8792TYPE(vol7d),intent(out) :: this
8793integer,INTENT(IN),OPTIONAL :: time_definition
8794
8795CALL init(this%anavar)
8796CALL init(this%anaattr)
8797CALL init(this%anavarattr)
8798CALL init(this%dativar)
8799CALL init(this%datiattr)
8800CALL init(this%dativarattr)
8801CALL vol7d_var_features_init() ! initialise var features table once
8802
8803NULLIFY(this%ana, this%time, this%level, this%timerange, this%network)
8804
8805NULLIFY(this%volanar, this%volanaattrr, this%voldatir, this%voldatiattrr)
8806NULLIFY(this%volanad, this%volanaattrd, this%voldatid, this%voldatiattrd)
8807NULLIFY(this%volanai, this%volanaattri, this%voldatii, this%voldatiattri)
8808NULLIFY(this%volanab, this%volanaattrb, this%voldatib, this%voldatiattrb)
8809NULLIFY(this%volanac, this%volanaattrc, this%voldatic, this%voldatiattrc)
8810
8811if(present(time_definition)) then
8812 this%time_definition=time_definition
8813else
8814 this%time_definition=1 !default to validity time
8815end if
8816
8817END SUBROUTINE vol7d_init
8818
8819
8823ELEMENTAL SUBROUTINE vol7d_delete(this, dataonly)
8824TYPE(vol7d),intent(inout) :: this
8825LOGICAL, INTENT(in), OPTIONAL :: dataonly
8826
8827
8828IF (.NOT. optio_log(dataonly)) THEN
8829 IF (ASSOCIATED(this%volanar)) DEALLOCATE(this%volanar)
8830 IF (ASSOCIATED(this%volanad)) DEALLOCATE(this%volanad)
8831 IF (ASSOCIATED(this%volanai)) DEALLOCATE(this%volanai)
8832 IF (ASSOCIATED(this%volanab)) DEALLOCATE(this%volanab)
8833 IF (ASSOCIATED(this%volanac)) DEALLOCATE(this%volanac)
8834 IF (ASSOCIATED(this%volanaattrr)) DEALLOCATE(this%volanaattrr)
8835 IF (ASSOCIATED(this%volanaattrd)) DEALLOCATE(this%volanaattrd)
8836 IF (ASSOCIATED(this%volanaattri)) DEALLOCATE(this%volanaattri)
8837 IF (ASSOCIATED(this%volanaattrb)) DEALLOCATE(this%volanaattrb)
8838 IF (ASSOCIATED(this%volanaattrc)) DEALLOCATE(this%volanaattrc)
8839ENDIF
8840IF (ASSOCIATED(this%voldatir)) DEALLOCATE(this%voldatir)
8841IF (ASSOCIATED(this%voldatid)) DEALLOCATE(this%voldatid)
8842IF (ASSOCIATED(this%voldatii)) DEALLOCATE(this%voldatii)
8843IF (ASSOCIATED(this%voldatib)) DEALLOCATE(this%voldatib)
8844IF (ASSOCIATED(this%voldatic)) DEALLOCATE(this%voldatic)
8845IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
8846IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
8847IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
8848IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
8849IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
8850
8851IF (.NOT. optio_log(dataonly)) THEN
8852 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
8853 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
8854ENDIF
8855IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
8856IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
8857IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
8858
8859IF (.NOT. optio_log(dataonly)) THEN
8860 CALL delete(this%anavar)
8861 CALL delete(this%anaattr)
8862 CALL delete(this%anavarattr)
8863ENDIF
8864CALL delete(this%dativar)
8865CALL delete(this%datiattr)
8866CALL delete(this%dativarattr)
8867
8868END SUBROUTINE vol7d_delete
8869
8870
8871
8872integer function vol7d_check(this)
8873TYPE(vol7d),intent(in) :: this
8874integer :: i,j,k,l,m,n
8875
8876vol7d_check=0
8877
8878if (associated(this%voldatii)) then
8879do i = 1,size(this%voldatii,1)
8880 do j = 1,size(this%voldatii,2)
8881 do k = 1,size(this%voldatii,3)
8882 do l = 1,size(this%voldatii,4)
8883 do m = 1,size(this%voldatii,5)
8884 do n = 1,size(this%voldatii,6)
8885 if (this%voldatii(i,j,k,l,m,n) /= this%voldatii(i,j,k,l,m,n) ) then
8886 CALL l4f_log(l4f_warn,"check: abnormal value at voldatii("&
8887 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
8888 vol7d_check=1
8889 end if
8890 end do
8891 end do
8892 end do
8893 end do
8894 end do
8895end do
8896end if
8897
8898
8899if (associated(this%voldatir)) then
8900do i = 1,size(this%voldatir,1)
8901 do j = 1,size(this%voldatir,2)
8902 do k = 1,size(this%voldatir,3)
8903 do l = 1,size(this%voldatir,4)
8904 do m = 1,size(this%voldatir,5)
8905 do n = 1,size(this%voldatir,6)
8906 if (this%voldatir(i,j,k,l,m,n) /= this%voldatir(i,j,k,l,m,n) ) then
8907 CALL l4f_log(l4f_warn,"check: abnormal value at voldatir("&
8908 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
8909 vol7d_check=2
8910 end if
8911 end do
8912 end do
8913 end do
8914 end do
8915 end do
8916end do
8917end if
8918
8919if (associated(this%voldatid)) then
8920do i = 1,size(this%voldatid,1)
8921 do j = 1,size(this%voldatid,2)
8922 do k = 1,size(this%voldatid,3)
8923 do l = 1,size(this%voldatid,4)
8924 do m = 1,size(this%voldatid,5)
8925 do n = 1,size(this%voldatid,6)
8926 if (this%voldatid(i,j,k,l,m,n) /= this%voldatid(i,j,k,l,m,n) ) then
8927 CALL l4f_log(l4f_warn,"check: abnormal value at voldatid("&
8928 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
8929 vol7d_check=3
8930 end if
8931 end do
8932 end do
8933 end do
8934 end do
8935 end do
8936end do
8937end if
8938
8939if (associated(this%voldatib)) then
8940do i = 1,size(this%voldatib,1)
8941 do j = 1,size(this%voldatib,2)
8942 do k = 1,size(this%voldatib,3)
8943 do l = 1,size(this%voldatib,4)
8944 do m = 1,size(this%voldatib,5)
8945 do n = 1,size(this%voldatib,6)
8946 if (this%voldatib(i,j,k,l,m,n) /= this%voldatib(i,j,k,l,m,n) ) then
8947 CALL l4f_log(l4f_warn,"check: abnormal value at voldatib("&
8948 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
8949 vol7d_check=4
8950 end if
8951 end do
8952 end do
8953 end do
8954 end do
8955 end do
8956end do
8957end if
8958
8959end function vol7d_check
8960
8961
8962
8963!TODO da completare ! aborta se i volumi sono allocati a dimensione 0
8965SUBROUTINE vol7d_display(this)
8966TYPE(vol7d),intent(in) :: this
8967integer :: i
8968
8969REAL :: rdat
8970DOUBLE PRECISION :: ddat
8971INTEGER :: idat
8972INTEGER(kind=int_b) :: bdat
8973CHARACTER(len=vol7d_cdatalen) :: cdat
8974
8975
8976print*,"<<<<<<<<<<<<<<<<<<< vol7d object >>>>>>>>>>>>>>>>>>>>"
8977if (this%time_definition == 0) then
8978 print*,"TIME DEFINITION: time is reference time"
8979else if (this%time_definition == 1) then
8980 print*,"TIME DEFINITION: time is validity time"
8981else
8982 print*,"Time definition have a wrong walue:", this%time_definition
8983end if
8984
8985IF (ASSOCIATED(this%network))then
8986 print*,"---- network vector ----"
8987 print*,"elements=",size(this%network)
8988 do i=1, size(this%network)
8989 call display(this%network(i))
8990 end do
8991end IF
8992
8993IF (ASSOCIATED(this%ana))then
8994 print*,"---- ana vector ----"
8995 print*,"elements=",size(this%ana)
8996 do i=1, size(this%ana)
8997 call display(this%ana(i))
8998 end do
8999end IF
9000
9001IF (ASSOCIATED(this%time))then
9002 print*,"---- time vector ----"
9003 print*,"elements=",size(this%time)
9004 do i=1, size(this%time)
9005 call display(this%time(i))
9006 end do
9007end if
9008
9009IF (ASSOCIATED(this%level)) then
9010 print*,"---- level vector ----"
9011 print*,"elements=",size(this%level)
9012 do i =1,size(this%level)
9013 call display(this%level(i))
9014 end do
9015end if
9016
9017IF (ASSOCIATED(this%timerange))then
9018 print*,"---- timerange vector ----"
9019 print*,"elements=",size(this%timerange)
9020 do i =1,size(this%timerange)
9021 call display(this%timerange(i))
9022 end do
9023end if
9024
9025
9026print*,"---- ana vector ----"
9027print*,""
9028print*,"->>>>>>>>> anavar -"
9029call display(this%anavar)
9030print*,""
9031print*,"->>>>>>>>> anaattr -"
9032call display(this%anaattr)
9033print*,""
9034print*,"->>>>>>>>> anavarattr -"
9035call display(this%anavarattr)
9036
9037print*,"-- ana data section (first point) --"
9038
9039idat=imiss
9040rdat=rmiss
9041ddat=dmiss
9042bdat=ibmiss
9043cdat=cmiss
9044
9045!ntime = MIN(SIZE(this%time),nprint)
9046!ntimerange = MIN(SIZE(this%timerange),nprint)
9047!nlevel = MIN(SIZE(this%level),nprint)
9048!nnetwork = MIN(SIZE(this%network),nprint)
9049!nana = MIN(SIZE(this%ana),nprint)
9050
9051IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0) THEN
9052if (associated(this%volanai)) then
9053 do i=1,size(this%anavar%i)
9054 idat=this%volanai(1,i,1)
9055 if (associated(this%anavar%i)) call display(this%anavar%i(i),idat,rdat,ddat,bdat,cdat)
9056 end do
9057end if
9058idat=imiss
9059
9060if (associated(this%volanar)) then
9061 do i=1,size(this%anavar%r)
9062 rdat=this%volanar(1,i,1)
9063 if (associated(this%anavar%r)) call display(this%anavar%r(i),idat,rdat,ddat,bdat,cdat)
9064 end do
9065end if
9066rdat=rmiss
9067
9068if (associated(this%volanad)) then
9069 do i=1,size(this%anavar%d)
9070 ddat=this%volanad(1,i,1)
9071 if (associated(this%anavar%d)) call display(this%anavar%d(i),idat,rdat,ddat,bdat,cdat)
9072 end do
9073end if
9074ddat=dmiss
9075
9076if (associated(this%volanab)) then
9077 do i=1,size(this%anavar%b)
9078 bdat=this%volanab(1,i,1)
9079 if (associated(this%anavar%b)) call display(this%anavar%b(i),idat,rdat,ddat,bdat,cdat)
9080 end do
9081end if
9082bdat=ibmiss
9083
9084if (associated(this%volanac)) then
9085 do i=1,size(this%anavar%c)
9086 cdat=this%volanac(1,i,1)
9087 if (associated(this%anavar%c)) call display(this%anavar%c(i),idat,rdat,ddat,bdat,cdat)
9088 end do
9089end if
9090cdat=cmiss
9091ENDIF
9092
9093print*,"---- data vector ----"
9094print*,""
9095print*,"->>>>>>>>> dativar -"
9096call display(this%dativar)
9097print*,""
9098print*,"->>>>>>>>> datiattr -"
9099call display(this%datiattr)
9100print*,""
9101print*,"->>>>>>>>> dativarattr -"
9102call display(this%dativarattr)
9103
9104print*,"-- data data section (first point) --"
9105
9106idat=imiss
9107rdat=rmiss
9108ddat=dmiss
9109bdat=ibmiss
9110cdat=cmiss
9111
9112IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0 .AND. size(this%time) > 0 &
9113 .AND. size(this%level) > 0 .AND. size(this%timerange) > 0) THEN
9114if (associated(this%voldatii)) then
9115 do i=1,size(this%dativar%i)
9116 idat=this%voldatii(1,1,1,1,i,1)
9117 if (associated(this%dativar%i)) call display(this%dativar%i(i),idat,rdat,ddat,bdat,cdat)
9118 end do
9119end if
9120idat=imiss
9121
9122if (associated(this%voldatir)) then
9123 do i=1,size(this%dativar%r)
9124 rdat=this%voldatir(1,1,1,1,i,1)
9125 if (associated(this%dativar%r)) call display(this%dativar%r(i),idat,rdat,ddat,bdat,cdat)
9126 end do
9127end if
9128rdat=rmiss
9129
9130if (associated(this%voldatid)) then
9131 do i=1,size(this%dativar%d)
9132 ddat=this%voldatid(1,1,1,1,i,1)
9133 if (associated(this%dativar%d)) call display(this%dativar%d(i),idat,rdat,ddat,bdat,cdat)
9134 end do
9135end if
9136ddat=dmiss
9137
9138if (associated(this%voldatib)) then
9139 do i=1,size(this%dativar%b)
9140 bdat=this%voldatib(1,1,1,1,i,1)
9141 if (associated(this%dativar%b)) call display(this%dativar%b(i),idat,rdat,ddat,bdat,cdat)
9142 end do
9143end if
9144bdat=ibmiss
9145
9146if (associated(this%voldatic)) then
9147 do i=1,size(this%dativar%c)
9148 cdat=this%voldatic(1,1,1,1,i,1)
9149 if (associated(this%dativar%c)) call display(this%dativar%c(i),idat,rdat,ddat,bdat,cdat)
9150 end do
9151end if
9152cdat=cmiss
9153ENDIF
9154
9155print*,"<<<<<<<<<<<<<<<<<<< END vol7d object >>>>>>>>>>>>>>>>>>>>"
9156
9157END SUBROUTINE vol7d_display
9158
9159
9161SUBROUTINE dat_display(this,idat,rdat,ddat,bdat,cdat)
9162TYPE(vol7d_var),intent(in) :: this
9164REAL :: rdat
9166DOUBLE PRECISION :: ddat
9168INTEGER :: idat
9170INTEGER(kind=int_b) :: bdat
9172CHARACTER(len=*) :: cdat
9173
9174print *, to_char_dat(this,idat,rdat,ddat,bdat,cdat)
9175
9176end SUBROUTINE dat_display
9177
9179SUBROUTINE dat_vect_display(this,idat,rdat,ddat,bdat,cdat)
9180
9181TYPE(vol7d_var),intent(in) :: this(:)
9183REAL :: rdat(:)
9185DOUBLE PRECISION :: ddat(:)
9187INTEGER :: idat(:)
9189INTEGER(kind=int_b) :: bdat(:)
9191CHARACTER(len=*):: cdat(:)
9192
9193integer :: i
9194
9195do i =1,size(this)
9196 call display(this(i),idat(i),rdat(i),ddat(i),bdat(i),cdat(i))
9197end do
9198
9199end SUBROUTINE dat_vect_display
9200
9201
9202FUNCTION to_char_dat(this,idat,rdat,ddat,bdat,cdat)
9203#ifdef HAVE_DBALLE
9204USE dballef
9205#endif
9206TYPE(vol7d_var),INTENT(in) :: this
9208REAL :: rdat
9210DOUBLE PRECISION :: ddat
9212INTEGER :: idat
9214INTEGER(kind=int_b) :: bdat
9216CHARACTER(len=*) :: cdat
9217CHARACTER(len=80) :: to_char_dat
9218
9219CHARACTER(len=LEN(to_char_dat)) :: to_char_tmp
9220
9221
9222#ifdef HAVE_DBALLE
9223INTEGER :: handle, ier
9224
9225handle = 0
9226to_char_dat="VALUE: "
9227
9228if (c_e(idat)) to_char_dat=trim(to_char_dat)//" ;int> "//trim(to_char(idat))
9229if (c_e(rdat)) to_char_dat=trim(to_char_dat)//" ;real> "//trim(to_char(rdat))
9230if (c_e(ddat)) to_char_dat=trim(to_char_dat)//" ;double> "//trim(to_char(ddat))
9231if (c_e(bdat)) to_char_dat=trim(to_char_dat)//" ;byte> "//trim(to_char(bdat))
9232
9233if ( c_e(cdat))then
9234 ier = idba_messaggi(handle,"/dev/null", "w", "BUFR")
9235 ier = idba_spiegab(handle,this%btable,cdat,to_char_tmp)
9236 ier = idba_fatto(handle)
9237 to_char_dat=trim(to_char_dat)//" ;char> "//trim(to_char_tmp)
9238endif
9239
9240#else
9241
9242to_char_dat="VALUE: "
9243if (c_e(idat)) to_char_dat=trim(to_char_dat)//" ;int> "//trim(to_char(idat))
9244if (c_e(rdat)) to_char_dat=trim(to_char_dat)//" ;real> "//trim(to_char(rdat))
9245if (c_e(ddat)) to_char_dat=trim(to_char_dat)//" ;double> "//trim(to_char(ddat))
9246if (c_e(bdat)) to_char_dat=trim(to_char_dat)//" ;byte> "//trim(to_char(bdat))
9247if (c_e(cdat)) to_char_dat=trim(to_char_dat)//" ;char> "//trim(cdat)
9248
9249#endif
9250
9251END FUNCTION to_char_dat
9252
9253
9256FUNCTION vol7d_c_e(this) RESULT(c_e)
9257TYPE(vol7d), INTENT(in) :: this
9258
9259LOGICAL :: c_e
9260
9261c_e = ASSOCIATED(this%ana) .OR. ASSOCIATED(this%time) .OR. &
9262 ASSOCIATED(this%level) .OR. ASSOCIATED(this%timerange) .OR. &
9263 ASSOCIATED(this%network) .OR. &
9264 ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
9265 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
9266 ASSOCIATED(this%anavar%c) .OR. &
9267 ASSOCIATED(this%anaattr%r) .OR. ASSOCIATED(this%anaattr%d) .OR. &
9268 ASSOCIATED(this%anaattr%i) .OR. ASSOCIATED(this%anaattr%b) .OR. &
9269 ASSOCIATED(this%anaattr%c) .OR. &
9270 ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
9271 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
9272 ASSOCIATED(this%dativar%c) .OR. &
9273 ASSOCIATED(this%datiattr%r) .OR. ASSOCIATED(this%datiattr%d) .OR. &
9274 ASSOCIATED(this%datiattr%i) .OR. ASSOCIATED(this%datiattr%b) .OR. &
9275 ASSOCIATED(this%datiattr%c)
9276
9277END FUNCTION vol7d_c_e
9278
9279
9318SUBROUTINE vol7d_alloc(this, nana, ntime, nlevel, ntimerange, nnetwork, &
9319 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
9320 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
9321 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
9322 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
9323 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
9324 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc, &
9325 ini)
9326TYPE(vol7d),INTENT(inout) :: this
9327INTEGER,INTENT(in),OPTIONAL :: nana
9328INTEGER,INTENT(in),OPTIONAL :: ntime
9329INTEGER,INTENT(in),OPTIONAL :: nlevel
9330INTEGER,INTENT(in),OPTIONAL :: ntimerange
9331INTEGER,INTENT(in),OPTIONAL :: nnetwork
9333INTEGER,INTENT(in),OPTIONAL :: &
9334 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
9335 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
9336 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
9337 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
9338 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
9339 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc
9340LOGICAL,INTENT(in),OPTIONAL :: ini
9341
9342INTEGER :: i
9343LOGICAL :: linit
9344
9345IF (PRESENT(ini)) THEN
9346 linit = ini
9347ELSE
9348 linit = .false.
9349ENDIF
9350
9351! Dimensioni principali
9352IF (PRESENT(nana)) THEN
9353 IF (nana >= 0) THEN
9354 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
9355 ALLOCATE(this%ana(nana))
9356 IF (linit) THEN
9357 DO i = 1, nana
9358 CALL init(this%ana(i))
9359 ENDDO
9360 ENDIF
9361 ENDIF
9362ENDIF
9363IF (PRESENT(ntime)) THEN
9364 IF (ntime >= 0) THEN
9365 IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
9366 ALLOCATE(this%time(ntime))
9367 IF (linit) THEN
9368 DO i = 1, ntime
9369 CALL init(this%time(i))
9370 ENDDO
9371 ENDIF
9372 ENDIF
9373ENDIF
9374IF (PRESENT(nlevel)) THEN
9375 IF (nlevel >= 0) THEN
9376 IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
9377 ALLOCATE(this%level(nlevel))
9378 IF (linit) THEN
9379 DO i = 1, nlevel
9380 CALL init(this%level(i))
9381 ENDDO
9382 ENDIF
9383 ENDIF
9384ENDIF
9385IF (PRESENT(ntimerange)) THEN
9386 IF (ntimerange >= 0) THEN
9387 IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
9388 ALLOCATE(this%timerange(ntimerange))
9389 IF (linit) THEN
9390 DO i = 1, ntimerange
9391 CALL init(this%timerange(i))
9392 ENDDO
9393 ENDIF
9394 ENDIF
9395ENDIF
9396IF (PRESENT(nnetwork)) THEN
9397 IF (nnetwork >= 0) THEN
9398 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
9399 ALLOCATE(this%network(nnetwork))
9400 IF (linit) THEN
9401 DO i = 1, nnetwork
9402 CALL init(this%network(i))
9403 ENDDO
9404 ENDIF
9405 ENDIF
9406ENDIF
9407! Dimensioni dei tipi delle variabili
9408CALL vol7d_varvect_alloc(this%anavar, nanavarr, nanavard, &
9409 nanavari, nanavarb, nanavarc, ini)
9410CALL vol7d_varvect_alloc(this%anaattr, nanaattrr, nanaattrd, &
9411 nanaattri, nanaattrb, nanaattrc, ini)
9412CALL vol7d_varvect_alloc(this%anavarattr, nanavarattrr, nanavarattrd, &
9413 nanavarattri, nanavarattrb, nanavarattrc, ini)
9414CALL vol7d_varvect_alloc(this%dativar, ndativarr, ndativard, &
9415 ndativari, ndativarb, ndativarc, ini)
9416CALL vol7d_varvect_alloc(this%datiattr, ndatiattrr, ndatiattrd, &
9417 ndatiattri, ndatiattrb, ndatiattrc, ini)
9418CALL vol7d_varvect_alloc(this%dativarattr, ndativarattrr, ndativarattrd, &
9419 ndativarattri, ndativarattrb, ndativarattrc, ini)
9420
9421END SUBROUTINE vol7d_alloc
9422
9423
9424FUNCTION vol7d_check_alloc_ana(this)
9425TYPE(vol7d),INTENT(in) :: this
9426LOGICAL :: vol7d_check_alloc_ana
9427
9428vol7d_check_alloc_ana = ASSOCIATED(this%ana) .AND. ASSOCIATED(this%network)
9429
9430END FUNCTION vol7d_check_alloc_ana
9431
9432SUBROUTINE vol7d_force_alloc_ana(this, ini)
9433TYPE(vol7d),INTENT(inout) :: this
9434LOGICAL,INTENT(in),OPTIONAL :: ini
9435
9436! Alloco i descrittori minimi per avere un volume di anagrafica
9437IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=1, ini=ini)
9438IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=1, ini=ini)
9439
9440END SUBROUTINE vol7d_force_alloc_ana
9441
9442
9443FUNCTION vol7d_check_alloc_dati(this)
9444TYPE(vol7d),INTENT(in) :: this
9445LOGICAL :: vol7d_check_alloc_dati
9446
9447vol7d_check_alloc_dati = vol7d_check_alloc_ana(this) .AND. &
9448 ASSOCIATED(this%time) .AND. ASSOCIATED(this%level) .AND. &
9449 ASSOCIATED(this%timerange)
9450
9451END FUNCTION vol7d_check_alloc_dati
9452
9453SUBROUTINE vol7d_force_alloc_dati(this, ini)
9454TYPE(vol7d),INTENT(inout) :: this
9455LOGICAL,INTENT(in),OPTIONAL :: ini
9456
9457! Alloco i descrittori minimi per avere un volume di dati
9458CALL vol7d_force_alloc_ana(this, ini)
9459IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=1, ini=ini)
9460IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=1, ini=ini)
9461IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=1, ini=ini)
9462
9463END SUBROUTINE vol7d_force_alloc_dati
9464
9465
9466SUBROUTINE vol7d_force_alloc(this)
9467TYPE(vol7d),INTENT(inout) :: this
9468
9469! If anything really not allocated yet, allocate with size 0
9470IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=0)
9471IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=0)
9472IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=0)
9473IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=0)
9474IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=0)
9475
9476END SUBROUTINE vol7d_force_alloc
9477
9478
9479FUNCTION vol7d_check_vol(this)
9480TYPE(vol7d),INTENT(in) :: this
9481LOGICAL :: vol7d_check_vol
9482
9483vol7d_check_vol = c_e(this)
9484
9485! Anagrafica
9486IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
9487 vol7d_check_vol = .false.
9488ENDIF
9489
9490IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
9491 vol7d_check_vol = .false.
9492ENDIF
9493
9494IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
9495 vol7d_check_vol = .false.
9496ENDIF
9497
9498IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
9499 vol7d_check_vol = .false.
9500ENDIF
9501
9502IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
9503 vol7d_check_vol = .false.
9504ENDIF
9505IF (ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
9506 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
9507 ASSOCIATED(this%anavar%c)) THEN
9508 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_ana(this)
9509ENDIF
9510
9511! Attributi dell'anagrafica
9512IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
9513 .NOT.ASSOCIATED(this%volanaattrr)) THEN
9514 vol7d_check_vol = .false.
9515ENDIF
9516
9517IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
9518 .NOT.ASSOCIATED(this%volanaattrd)) THEN
9519 vol7d_check_vol = .false.
9520ENDIF
9521
9522IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
9523 .NOT.ASSOCIATED(this%volanaattri)) THEN
9524 vol7d_check_vol = .false.
9525ENDIF
9526
9527IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
9528 .NOT.ASSOCIATED(this%volanaattrb)) THEN
9529 vol7d_check_vol = .false.
9530ENDIF
9531
9532IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
9533 .NOT.ASSOCIATED(this%volanaattrc)) THEN
9534 vol7d_check_vol = .false.
9535ENDIF
9536
9537! Dati
9538IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
9539 vol7d_check_vol = .false.
9540ENDIF
9541
9542IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
9543 vol7d_check_vol = .false.
9544ENDIF
9545
9546IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
9547 vol7d_check_vol = .false.
9548ENDIF
9549
9550IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
9551 vol7d_check_vol = .false.
9552ENDIF
9553
9554IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
9555 vol7d_check_vol = .false.
9556ENDIF
9557
9558! Attributi dei dati
9559IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
9560 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
9561 vol7d_check_vol = .false.
9562ENDIF
9563
9564IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
9565 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
9566 vol7d_check_vol = .false.
9567ENDIF
9568
9569IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
9570 .NOT.ASSOCIATED(this%voldatiattri)) THEN
9571 vol7d_check_vol = .false.
9572ENDIF
9573
9574IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
9575 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
9576 vol7d_check_vol = .false.
9577ENDIF
9578
9579IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
9580 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
9581 vol7d_check_vol = .false.
9582ENDIF
9583IF (ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
9584 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
9585 ASSOCIATED(this%dativar%c)) THEN
9586 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_dati(this)
9587ENDIF
9588
9589END FUNCTION vol7d_check_vol
9590
9591
9606SUBROUTINE vol7d_alloc_vol(this, ini, inivol)
9607TYPE(vol7d),INTENT(inout) :: this
9608LOGICAL,INTENT(in),OPTIONAL :: ini
9609LOGICAL,INTENT(in),OPTIONAL :: inivol
9610
9611LOGICAL :: linivol
9612
9613IF (PRESENT(inivol)) THEN
9614 linivol = inivol
9615ELSE
9616 linivol = .true.
9617ENDIF
9618
9619! Anagrafica
9620IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
9621 CALL vol7d_force_alloc_ana(this, ini)
9622 ALLOCATE(this%volanar(SIZE(this%ana), SIZE(this%anavar%r), SIZE(this%network)))
9623 IF (linivol) this%volanar(:,:,:) = rmiss
9624ENDIF
9625
9626IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
9627 CALL vol7d_force_alloc_ana(this, ini)
9628 ALLOCATE(this%volanad(SIZE(this%ana), SIZE(this%anavar%d), SIZE(this%network)))
9629 IF (linivol) this%volanad(:,:,:) = rdmiss
9630ENDIF
9631
9632IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
9633 CALL vol7d_force_alloc_ana(this, ini)
9634 ALLOCATE(this%volanai(SIZE(this%ana), SIZE(this%anavar%i), SIZE(this%network)))
9635 IF (linivol) this%volanai(:,:,:) = imiss
9636ENDIF
9637
9638IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
9639 CALL vol7d_force_alloc_ana(this, ini)
9640 ALLOCATE(this%volanab(SIZE(this%ana), SIZE(this%anavar%b), SIZE(this%network)))
9641 IF (linivol) this%volanab(:,:,:) = ibmiss
9642ENDIF
9643
9644IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
9645 CALL vol7d_force_alloc_ana(this, ini)
9646 ALLOCATE(this%volanac(SIZE(this%ana), SIZE(this%anavar%c), SIZE(this%network)))
9647 IF (linivol) this%volanac(:,:,:) = cmiss
9648ENDIF
9649
9650! Attributi dell'anagrafica
9651IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
9652 .NOT.ASSOCIATED(this%volanaattrr)) THEN
9653 CALL vol7d_force_alloc_ana(this, ini)
9654 ALLOCATE(this%volanaattrr(SIZE(this%ana), SIZE(this%anavarattr%r), &
9655 SIZE(this%network), SIZE(this%anaattr%r)))
9656 IF (linivol) this%volanaattrr(:,:,:,:) = rmiss
9657ENDIF
9658
9659IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
9660 .NOT.ASSOCIATED(this%volanaattrd)) THEN
9661 CALL vol7d_force_alloc_ana(this, ini)
9662 ALLOCATE(this%volanaattrd(SIZE(this%ana), SIZE(this%anavarattr%d), &
9663 SIZE(this%network), SIZE(this%anaattr%d)))
9664 IF (linivol) this%volanaattrd(:,:,:,:) = rdmiss
9665ENDIF
9666
9667IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
9668 .NOT.ASSOCIATED(this%volanaattri)) THEN
9669 CALL vol7d_force_alloc_ana(this, ini)
9670 ALLOCATE(this%volanaattri(SIZE(this%ana), SIZE(this%anavarattr%i), &
9671 SIZE(this%network), SIZE(this%anaattr%i)))
9672 IF (linivol) this%volanaattri(:,:,:,:) = imiss
9673ENDIF
9674
9675IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
9676 .NOT.ASSOCIATED(this%volanaattrb)) THEN
9677 CALL vol7d_force_alloc_ana(this, ini)
9678 ALLOCATE(this%volanaattrb(SIZE(this%ana), SIZE(this%anavarattr%b), &
9679 SIZE(this%network), SIZE(this%anaattr%b)))
9680 IF (linivol) this%volanaattrb(:,:,:,:) = ibmiss
9681ENDIF
9682
9683IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
9684 .NOT.ASSOCIATED(this%volanaattrc)) THEN
9685 CALL vol7d_force_alloc_ana(this, ini)
9686 ALLOCATE(this%volanaattrc(SIZE(this%ana), SIZE(this%anavarattr%c), &
9687 SIZE(this%network), SIZE(this%anaattr%c)))
9688 IF (linivol) this%volanaattrc(:,:,:,:) = cmiss
9689ENDIF
9690
9691! Dati
9692IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
9693 CALL vol7d_force_alloc_dati(this, ini)
9694 ALLOCATE(this%voldatir(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
9695 SIZE(this%timerange), SIZE(this%dativar%r), SIZE(this%network)))
9696 IF (linivol) this%voldatir(:,:,:,:,:,:) = rmiss
9697ENDIF
9698
9699IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
9700 CALL vol7d_force_alloc_dati(this, ini)
9701 ALLOCATE(this%voldatid(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
9702 SIZE(this%timerange), SIZE(this%dativar%d), SIZE(this%network)))
9703 IF (linivol) this%voldatid(:,:,:,:,:,:) = rdmiss
9704ENDIF
9705
9706IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
9707 CALL vol7d_force_alloc_dati(this, ini)
9708 ALLOCATE(this%voldatii(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
9709 SIZE(this%timerange), SIZE(this%dativar%i), SIZE(this%network)))
9710 IF (linivol) this%voldatii(:,:,:,:,:,:) = imiss
9711ENDIF
9712
9713IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
9714 CALL vol7d_force_alloc_dati(this, ini)
9715 ALLOCATE(this%voldatib(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
9716 SIZE(this%timerange), SIZE(this%dativar%b), SIZE(this%network)))
9717 IF (linivol) this%voldatib(:,:,:,:,:,:) = ibmiss
9718ENDIF
9719
9720IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
9721 CALL vol7d_force_alloc_dati(this, ini)
9722 ALLOCATE(this%voldatic(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
9723 SIZE(this%timerange), SIZE(this%dativar%c), SIZE(this%network)))
9724 IF (linivol) this%voldatic(:,:,:,:,:,:) = cmiss
9725ENDIF
9726
9727! Attributi dei dati
9728IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
9729 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
9730 CALL vol7d_force_alloc_dati(this, ini)
9731 ALLOCATE(this%voldatiattrr(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
9732 SIZE(this%timerange), SIZE(this%dativarattr%r), SIZE(this%network), &
9733 SIZE(this%datiattr%r)))
9734 IF (linivol) this%voldatiattrr(:,:,:,:,:,:,:) = rmiss
9735ENDIF
9736
9737IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
9738 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
9739 CALL vol7d_force_alloc_dati(this, ini)
9740 ALLOCATE(this%voldatiattrd(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
9741 SIZE(this%timerange), SIZE(this%dativarattr%d), SIZE(this%network), &
9742 SIZE(this%datiattr%d)))
9743 IF (linivol) this%voldatiattrd(:,:,:,:,:,:,:) = rdmiss
9744ENDIF
9745
9746IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
9747 .NOT.ASSOCIATED(this%voldatiattri)) THEN
9748 CALL vol7d_force_alloc_dati(this, ini)
9749 ALLOCATE(this%voldatiattri(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
9750 SIZE(this%timerange), SIZE(this%dativarattr%i), SIZE(this%network), &
9751 SIZE(this%datiattr%i)))
9752 IF (linivol) this%voldatiattri(:,:,:,:,:,:,:) = imiss
9753ENDIF
9754
9755IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
9756 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
9757 CALL vol7d_force_alloc_dati(this, ini)
9758 ALLOCATE(this%voldatiattrb(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
9759 SIZE(this%timerange), SIZE(this%dativarattr%b), SIZE(this%network), &
9760 SIZE(this%datiattr%b)))
9761 IF (linivol) this%voldatiattrb(:,:,:,:,:,:,:) = ibmiss
9762ENDIF
9763
9764IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
9765 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
9766 CALL vol7d_force_alloc_dati(this, ini)
9767 ALLOCATE(this%voldatiattrc(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
9768 SIZE(this%timerange), SIZE(this%dativarattr%c), SIZE(this%network), &
9769 SIZE(this%datiattr%c)))
9770 IF (linivol) this%voldatiattrc(:,:,:,:,:,:,:) = cmiss
9771ENDIF
9772
9773! Catch-all method
9774CALL vol7d_force_alloc(this)
9775
9776! Creo gli indici var-attr
9777
9778#ifdef DEBUG
9779CALL l4f_log(l4f_debug,"calling: vol7d_set_attr_ind")
9780#endif
9781
9782CALL vol7d_set_attr_ind(this)
9783
9784
9785
9786END SUBROUTINE vol7d_alloc_vol
9787
9788
9795SUBROUTINE vol7d_set_attr_ind(this)
9796TYPE(vol7d),INTENT(inout) :: this
9797
9798INTEGER :: i
9799
9800! real
9801IF (ASSOCIATED(this%dativar%r)) THEN
9802 IF (ASSOCIATED(this%dativarattr%r)) THEN
9803 DO i = 1, SIZE(this%dativar%r)
9804 this%dativar%r(i)%r = &
9805 firsttrue(this%dativar%r(i)%btable == this%dativarattr%r(:)%btable)
9806 ENDDO
9807 ENDIF
9808
9809 IF (ASSOCIATED(this%dativarattr%d)) THEN
9810 DO i = 1, SIZE(this%dativar%r)
9811 this%dativar%r(i)%d = &
9812 firsttrue(this%dativar%r(i)%btable == this%dativarattr%d(:)%btable)
9813 ENDDO
9814 ENDIF
9815
9816 IF (ASSOCIATED(this%dativarattr%i)) THEN
9817 DO i = 1, SIZE(this%dativar%r)
9818 this%dativar%r(i)%i = &
9819 firsttrue(this%dativar%r(i)%btable == this%dativarattr%i(:)%btable)
9820 ENDDO
9821 ENDIF
9822
9823 IF (ASSOCIATED(this%dativarattr%b)) THEN
9824 DO i = 1, SIZE(this%dativar%r)
9825 this%dativar%r(i)%b = &
9826 firsttrue(this%dativar%r(i)%btable == this%dativarattr%b(:)%btable)
9827 ENDDO
9828 ENDIF
9829
9830 IF (ASSOCIATED(this%dativarattr%c)) THEN
9831 DO i = 1, SIZE(this%dativar%r)
9832 this%dativar%r(i)%c = &
9833 firsttrue(this%dativar%r(i)%btable == this%dativarattr%c(:)%btable)
9834 ENDDO
9835 ENDIF
9836ENDIF
9837! double
9838IF (ASSOCIATED(this%dativar%d)) THEN
9839 IF (ASSOCIATED(this%dativarattr%r)) THEN
9840 DO i = 1, SIZE(this%dativar%d)
9841 this%dativar%d(i)%r = &
9842 firsttrue(this%dativar%d(i)%btable == this%dativarattr%r(:)%btable)
9843 ENDDO
9844 ENDIF
9845
9846 IF (ASSOCIATED(this%dativarattr%d)) THEN
9847 DO i = 1, SIZE(this%dativar%d)
9848 this%dativar%d(i)%d = &
9849 firsttrue(this%dativar%d(i)%btable == this%dativarattr%d(:)%btable)
9850 ENDDO
9851 ENDIF
9852
9853 IF (ASSOCIATED(this%dativarattr%i)) THEN
9854 DO i = 1, SIZE(this%dativar%d)
9855 this%dativar%d(i)%i = &
9856 firsttrue(this%dativar%d(i)%btable == this%dativarattr%i(:)%btable)
9857 ENDDO
9858 ENDIF
9859
9860 IF (ASSOCIATED(this%dativarattr%b)) THEN
9861 DO i = 1, SIZE(this%dativar%d)
9862 this%dativar%d(i)%b = &
9863 firsttrue(this%dativar%d(i)%btable == this%dativarattr%b(:)%btable)
9864 ENDDO
9865 ENDIF
9866
9867 IF (ASSOCIATED(this%dativarattr%c)) THEN
9868 DO i = 1, SIZE(this%dativar%d)
9869 this%dativar%d(i)%c = &
9870 firsttrue(this%dativar%d(i)%btable == this%dativarattr%c(:)%btable)
9871 ENDDO
9872 ENDIF
9873ENDIF
9874! integer
9875IF (ASSOCIATED(this%dativar%i)) THEN
9876 IF (ASSOCIATED(this%dativarattr%r)) THEN
9877 DO i = 1, SIZE(this%dativar%i)
9878 this%dativar%i(i)%r = &
9879 firsttrue(this%dativar%i(i)%btable == this%dativarattr%r(:)%btable)
9880 ENDDO
9881 ENDIF
9882
9883 IF (ASSOCIATED(this%dativarattr%d)) THEN
9884 DO i = 1, SIZE(this%dativar%i)
9885 this%dativar%i(i)%d = &
9886 firsttrue(this%dativar%i(i)%btable == this%dativarattr%d(:)%btable)
9887 ENDDO
9888 ENDIF
9889
9890 IF (ASSOCIATED(this%dativarattr%i)) THEN
9891 DO i = 1, SIZE(this%dativar%i)
9892 this%dativar%i(i)%i = &
9893 firsttrue(this%dativar%i(i)%btable == this%dativarattr%i(:)%btable)
9894 ENDDO
9895 ENDIF
9896
9897 IF (ASSOCIATED(this%dativarattr%b)) THEN
9898 DO i = 1, SIZE(this%dativar%i)
9899 this%dativar%i(i)%b = &
9900 firsttrue(this%dativar%i(i)%btable == this%dativarattr%b(:)%btable)
9901 ENDDO
9902 ENDIF
9903
9904 IF (ASSOCIATED(this%dativarattr%c)) THEN
9905 DO i = 1, SIZE(this%dativar%i)
9906 this%dativar%i(i)%c = &
9907 firsttrue(this%dativar%i(i)%btable == this%dativarattr%c(:)%btable)
9908 ENDDO
9909 ENDIF
9910ENDIF
9911! byte
9912IF (ASSOCIATED(this%dativar%b)) THEN
9913 IF (ASSOCIATED(this%dativarattr%r)) THEN
9914 DO i = 1, SIZE(this%dativar%b)
9915 this%dativar%b(i)%r = &
9916 firsttrue(this%dativar%b(i)%btable == this%dativarattr%r(:)%btable)
9917 ENDDO
9918 ENDIF
9919
9920 IF (ASSOCIATED(this%dativarattr%d)) THEN
9921 DO i = 1, SIZE(this%dativar%b)
9922 this%dativar%b(i)%d = &
9923 firsttrue(this%dativar%b(i)%btable == this%dativarattr%d(:)%btable)
9924 ENDDO
9925 ENDIF
9926
9927 IF (ASSOCIATED(this%dativarattr%i)) THEN
9928 DO i = 1, SIZE(this%dativar%b)
9929 this%dativar%b(i)%i = &
9930 firsttrue(this%dativar%b(i)%btable == this%dativarattr%i(:)%btable)
9931 ENDDO
9932 ENDIF
9933
9934 IF (ASSOCIATED(this%dativarattr%b)) THEN
9935 DO i = 1, SIZE(this%dativar%b)
9936 this%dativar%b(i)%b = &
9937 firsttrue(this%dativar%b(i)%btable == this%dativarattr%b(:)%btable)
9938 ENDDO
9939 ENDIF
9940
9941 IF (ASSOCIATED(this%dativarattr%c)) THEN
9942 DO i = 1, SIZE(this%dativar%b)
9943 this%dativar%b(i)%c = &
9944 firsttrue(this%dativar%b(i)%btable == this%dativarattr%c(:)%btable)
9945 ENDDO
9946 ENDIF
9947ENDIF
9948! character
9949IF (ASSOCIATED(this%dativar%c)) THEN
9950 IF (ASSOCIATED(this%dativarattr%r)) THEN
9951 DO i = 1, SIZE(this%dativar%c)
9952 this%dativar%c(i)%r = &
9953 firsttrue(this%dativar%c(i)%btable == this%dativarattr%r(:)%btable)
9954 ENDDO
9955 ENDIF
9956
9957 IF (ASSOCIATED(this%dativarattr%d)) THEN
9958 DO i = 1, SIZE(this%dativar%c)
9959 this%dativar%c(i)%d = &
9960 firsttrue(this%dativar%c(i)%btable == this%dativarattr%d(:)%btable)
9961 ENDDO
9962 ENDIF
9963
9964 IF (ASSOCIATED(this%dativarattr%i)) THEN
9965 DO i = 1, SIZE(this%dativar%c)
9966 this%dativar%c(i)%i = &
9967 firsttrue(this%dativar%c(i)%btable == this%dativarattr%i(:)%btable)
9968 ENDDO
9969 ENDIF
9970
9971 IF (ASSOCIATED(this%dativarattr%b)) THEN
9972 DO i = 1, SIZE(this%dativar%c)
9973 this%dativar%c(i)%b = &
9974 firsttrue(this%dativar%c(i)%btable == this%dativarattr%b(:)%btable)
9975 ENDDO
9976 ENDIF
9977
9978 IF (ASSOCIATED(this%dativarattr%c)) THEN
9979 DO i = 1, SIZE(this%dativar%c)
9980 this%dativar%c(i)%c = &
9981 firsttrue(this%dativar%c(i)%btable == this%dativarattr%c(:)%btable)
9982 ENDDO
9983 ENDIF
9984ENDIF
9985
9986END SUBROUTINE vol7d_set_attr_ind
9987
9988
9993SUBROUTINE vol7d_merge(this, that, sort, bestdata, &
9994 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
9995TYPE(vol7d),INTENT(INOUT) :: this
9996TYPE(vol7d),INTENT(INOUT) :: that
9997LOGICAL,INTENT(IN),OPTIONAL :: sort
9998LOGICAL,INTENT(in),OPTIONAL :: bestdata
9999LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple ! experimental, please do not use outside the library now
10000
10001TYPE(vol7d) :: v7d_clean
10002
10003
10004IF (.NOT.c_e(this)) THEN ! speedup
10005 this = that
10006 CALL init(v7d_clean)
10007 that = v7d_clean ! destroy that without deallocating
10008ELSE ! Append that to this and destroy that
10009 CALL vol7d_append(this, that, sort, bestdata, &
10010 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
10011 CALL delete(that)
10012ENDIF
10013
10014END SUBROUTINE vol7d_merge
10015
10016
10045SUBROUTINE vol7d_append(this, that, sort, bestdata, &
10046 ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple)
10047TYPE(vol7d),INTENT(INOUT) :: this
10048TYPE(vol7d),INTENT(IN) :: that
10049LOGICAL,INTENT(IN),OPTIONAL :: sort
10050! experimental, please do not use outside the library now, they force the use
10051! of a simplified mapping algorithm which is valid only whene the dimension
10052! content is the same in both volumes , or when one of them is empty
10053LOGICAL,INTENT(in),OPTIONAL :: bestdata
10054LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple
10055
10056
10057TYPE(vol7d) :: v7dtmp
10058LOGICAL :: lsort, lbestdata
10059INTEGER,POINTER :: remapt1(:), remapt2(:), remaptr1(:), remaptr2(:), &
10060 remapl1(:), remapl2(:), remapa1(:), remapa2(:), remapn1(:), remapn2(:)
10061
10062IF (.NOT.c_e(that)) RETURN ! speedup, nothing to do
10063IF (.NOT.vol7d_check_vol(that)) RETURN ! be safe
10064IF (.NOT.c_e(this)) THEN ! this case is like a vol7d_copy, more efficient to copy?
10065 CALL vol7d_copy(that, this, sort=sort)
10066 RETURN
10067ENDIF
10068
10069IF (this%time_definition /= that%time_definition) THEN
10070 CALL l4f_log(l4f_fatal, &
10071 'in vol7d_append, cannot append volumes with different &
10072 &time definition')
10073 CALL raise_fatal_error()
10074ENDIF
10075
10076! Completo l'allocazione per avere volumi a norma
10077CALL vol7d_alloc_vol(this)
10078
10079CALL init(v7dtmp, time_definition=this%time_definition)
10080CALL optio(sort, lsort)
10081CALL optio(bestdata, lbestdata)
10082
10083! Calcolo le mappature tra volumi vecchi e volume nuovo
10084! I puntatori remap* vengono tutti o allocati o nullificati
10085IF (optio_log(ltimesimple)) THEN
10086 CALL vol7d_remap2simple_datetime(this%time, that%time, v7dtmp%time, &
10087 lsort, remapt1, remapt2)
10088ELSE
10089 CALL vol7d_remap2_datetime(this%time, that%time, v7dtmp%time, &
10090 lsort, remapt1, remapt2)
10091ENDIF
10092IF (optio_log(ltimerangesimple)) THEN
10093 CALL vol7d_remap2simple_vol7d_timerange(this%timerange, that%timerange, &
10094 v7dtmp%timerange, lsort, remaptr1, remaptr2)
10095ELSE
10096 CALL vol7d_remap2_vol7d_timerange(this%timerange, that%timerange, &
10097 v7dtmp%timerange, lsort, remaptr1, remaptr2)
10098ENDIF
10099IF (optio_log(llevelsimple)) THEN
10100 CALL vol7d_remap2simple_vol7d_level(this%level, that%level, v7dtmp%level, &
10101 lsort, remapl1, remapl2)
10102ELSE
10103 CALL vol7d_remap2_vol7d_level(this%level, that%level, v7dtmp%level, &
10104 lsort, remapl1, remapl2)
10105ENDIF
10106IF (optio_log(lanasimple)) THEN
10107 CALL vol7d_remap2simple_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
10108 .false., remapa1, remapa2)
10109ELSE
10110 CALL vol7d_remap2_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
10111 .false., remapa1, remapa2)
10112ENDIF
10113IF (optio_log(lnetworksimple)) THEN
10114 CALL vol7d_remap2simple_vol7d_network(this%network, that%network, v7dtmp%network, &
10115 .false., remapn1, remapn2)
10116ELSE
10117 CALL vol7d_remap2_vol7d_network(this%network, that%network, v7dtmp%network, &
10118 .false., remapn1, remapn2)
10119ENDIF
10120
10121! Faccio la fusione fisica dei volumi
10122CALL vol7d_merge_finalr(this, that, v7dtmp, &
10123 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10124 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10125CALL vol7d_merge_finald(this, that, v7dtmp, &
10126 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10127 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10128CALL vol7d_merge_finali(this, that, v7dtmp, &
10129 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10130 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10131CALL vol7d_merge_finalb(this, that, v7dtmp, &
10132 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10133 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10134CALL vol7d_merge_finalc(this, that, v7dtmp, &
10135 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10136 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10137
10138! Dealloco i vettori di rimappatura
10139IF (ASSOCIATED(remapt1)) DEALLOCATE(remapt1)
10140IF (ASSOCIATED(remapt2)) DEALLOCATE(remapt2)
10141IF (ASSOCIATED(remaptr1)) DEALLOCATE(remaptr1)
10142IF (ASSOCIATED(remaptr2)) DEALLOCATE(remaptr2)
10143IF (ASSOCIATED(remapl1)) DEALLOCATE(remapl1)
10144IF (ASSOCIATED(remapl2)) DEALLOCATE(remapl2)
10145IF (ASSOCIATED(remapa1)) DEALLOCATE(remapa1)
10146IF (ASSOCIATED(remapa2)) DEALLOCATE(remapa2)
10147IF (ASSOCIATED(remapn1)) DEALLOCATE(remapn1)
10148IF (ASSOCIATED(remapn2)) DEALLOCATE(remapn2)
10149
10150! Distruggo il vecchio volume e assegno il nuovo a this
10151CALL delete(this)
10152this = v7dtmp
10153! Ricreo gli indici var-attr
10154CALL vol7d_set_attr_ind(this)
10155
10156END SUBROUTINE vol7d_append
10157
10158
10191SUBROUTINE vol7d_copy(this, that, sort, unique, miss, &
10192 lsort_time, lsort_timerange, lsort_level, &
10193 ltime, ltimerange, llevel, lana, lnetwork, &
10194 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
10195 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
10196 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
10197 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
10198 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
10199 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
10200TYPE(vol7d),INTENT(IN) :: this
10201TYPE(vol7d),INTENT(INOUT) :: that
10202LOGICAL,INTENT(IN),OPTIONAL :: sort
10203LOGICAL,INTENT(IN),OPTIONAL :: unique
10204LOGICAL,INTENT(IN),OPTIONAL :: miss
10205LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
10206LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
10207LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
10215LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
10217LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
10219LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
10221LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
10223LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
10225LOGICAL,INTENT(in),OPTIONAL :: &
10226 lanavarr(:), lanavard(:), lanavari(:), lanavarb(:), lanavarc(:), &
10227 lanaattrr(:), lanaattrd(:), lanaattri(:), lanaattrb(:), lanaattrc(:), &
10228 lanavarattrr(:), lanavarattrd(:), lanavarattri(:), lanavarattrb(:), lanavarattrc(:), &
10229 ldativarr(:), ldativard(:), ldativari(:), ldativarb(:), ldativarc(:), &
10230 ldatiattrr(:), ldatiattrd(:), ldatiattri(:), ldatiattrb(:), ldatiattrc(:), &
10231 ldativarattrr(:), ldativarattrd(:), ldativarattri(:), ldativarattrb(:), ldativarattrc(:)
10232
10233LOGICAL :: lsort, lunique, lmiss
10234INTEGER,POINTER :: remapt(:), remaptr(:), remapl(:), remapa(:), remapn(:)
10235
10236CALL init(that)
10237IF (.NOT.c_e(this)) RETURN ! speedup, nothing to do
10238IF (.NOT.vol7d_check_vol(this)) RETURN ! be safe
10239
10240CALL optio(sort, lsort)
10241CALL optio(unique, lunique)
10242CALL optio(miss, lmiss)
10243
10244! Calcolo le mappature tra volume vecchio e volume nuovo
10245! I puntatori remap* vengono tutti o allocati o nullificati
10246CALL vol7d_remap1_datetime(this%time, that%time, &
10247 lsort.OR.optio_log(lsort_time), lunique, lmiss, remapt, ltime)
10248CALL vol7d_remap1_vol7d_timerange(this%timerange, that%timerange, &
10249 lsort.OR.optio_log(lsort_timerange), lunique, lmiss, remaptr, ltimerange)
10250CALL vol7d_remap1_vol7d_level(this%level, that%level, &
10251 lsort.OR.optio_log(lsort_level), lunique, lmiss, remapl, llevel)
10252CALL vol7d_remap1_vol7d_ana(this%ana, that%ana, &
10253 lsort, lunique, lmiss, remapa, lana)
10254CALL vol7d_remap1_vol7d_network(this%network, that%network, &
10255 lsort, lunique, lmiss, remapn, lnetwork)
10256
10257! lanavari, lanavarb, lanavarc, &
10258! lanaattri, lanaattrb, lanaattrc, &
10259! lanavarattri, lanavarattrb, lanavarattrc, &
10260! ldativari, ldativarb, ldativarc, &
10261! ldatiattri, ldatiattrb, ldatiattrc, &
10262! ldativarattri, ldativarattrb, ldativarattrc
10263! Faccio la riforma fisica dei volumi
10264CALL vol7d_reform_finalr(this, that, &
10265 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10266 lanavarr, lanaattrr, lanavarattrr, ldativarr, ldatiattrr, ldativarattrr)
10267CALL vol7d_reform_finald(this, that, &
10268 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10269 lanavard, lanaattrd, lanavarattrd, ldativard, ldatiattrd, ldativarattrd)
10270CALL vol7d_reform_finali(this, that, &
10271 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10272 lanavari, lanaattri, lanavarattri, ldativari, ldatiattri, ldativarattri)
10273CALL vol7d_reform_finalb(this, that, &
10274 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10275 lanavarb, lanaattrb, lanavarattrb, ldativarb, ldatiattrb, ldativarattrb)
10276CALL vol7d_reform_finalc(this, that, &
10277 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10278 lanavarc, lanaattrc, lanavarattrc, ldativarc, ldatiattrc, ldativarattrc)
10279
10280! Dealloco i vettori di rimappatura
10281IF (ASSOCIATED(remapt)) DEALLOCATE(remapt)
10282IF (ASSOCIATED(remaptr)) DEALLOCATE(remaptr)
10283IF (ASSOCIATED(remapl)) DEALLOCATE(remapl)
10284IF (ASSOCIATED(remapa)) DEALLOCATE(remapa)
10285IF (ASSOCIATED(remapn)) DEALLOCATE(remapn)
10286
10287! Ricreo gli indici var-attr
10288CALL vol7d_set_attr_ind(that)
10289that%time_definition = this%time_definition
10290
10291END SUBROUTINE vol7d_copy
10292
10293
10304SUBROUTINE vol7d_reform(this, sort, unique, miss, &
10305 lsort_time, lsort_timerange, lsort_level, &
10306 ltime, ltimerange, llevel, lana, lnetwork, &
10307 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
10308 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
10309 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
10310 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
10311 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
10312 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc&
10313 ,purgeana)
10314TYPE(vol7d),INTENT(INOUT) :: this
10315LOGICAL,INTENT(IN),OPTIONAL :: sort
10316LOGICAL,INTENT(IN),OPTIONAL :: unique
10317LOGICAL,INTENT(IN),OPTIONAL :: miss
10318LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
10319LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
10320LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
10328LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
10329LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
10330LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
10331LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
10332LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
10334LOGICAL,INTENT(in),OPTIONAL :: &
10335 lanavarr(:), lanavard(:), lanavari(:), lanavarb(:), lanavarc(:), &
10336 lanaattrr(:), lanaattrd(:), lanaattri(:), lanaattrb(:), lanaattrc(:), &
10337 lanavarattrr(:), lanavarattrd(:), lanavarattri(:), lanavarattrb(:), lanavarattrc(:), &
10338 ldativarr(:), ldativard(:), ldativari(:), ldativarb(:), ldativarc(:), &
10339 ldatiattrr(:), ldatiattrd(:), ldatiattri(:), ldatiattrb(:), ldatiattrc(:), &
10340 ldativarattrr(:), ldativarattrd(:), ldativarattri(:), ldativarattrb(:), ldativarattrc(:)
10341LOGICAL,INTENT(IN),OPTIONAL :: purgeana
10342
10343TYPE(vol7d) :: v7dtmp
10344logical,allocatable :: llana(:)
10345integer :: i
10346
10347CALL vol7d_copy(this, v7dtmp, sort, unique, miss, &
10348 lsort_time, lsort_timerange, lsort_level, &
10349 ltime, ltimerange, llevel, lana, lnetwork, &
10350 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
10351 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
10352 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
10353 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
10354 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
10355 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
10356
10357! destroy old volume
10358CALL delete(this)
10359
10360if (optio_log(purgeana)) then
10361 allocate(llana(size(v7dtmp%ana)))
10362 llana =.false.
10363 do i =1,size(v7dtmp%ana)
10364 if (associated(v7dtmp%voldatii)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatii(i,:,:,:,:,:)))
10365 if (associated(v7dtmp%voldatir)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatir(i,:,:,:,:,:)))
10366 if (associated(v7dtmp%voldatid)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatid(i,:,:,:,:,:)))
10367 if (associated(v7dtmp%voldatib)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatib(i,:,:,:,:,:)))
10368 if (associated(v7dtmp%voldatic)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatic(i,:,:,:,:,:)))
10369 end do
10370 CALL vol7d_copy(v7dtmp, this,lana=llana)
10371 CALL delete(v7dtmp)
10372 deallocate(llana)
10373else
10374 this=v7dtmp
10375end if
10376
10377END SUBROUTINE vol7d_reform
10378
10379
10387SUBROUTINE vol7d_smart_sort(this, lsort_time, lsort_timerange, lsort_level)
10388TYPE(vol7d),INTENT(INOUT) :: this
10389LOGICAL,OPTIONAL,INTENT(in) :: lsort_time
10390LOGICAL,OPTIONAL,INTENT(in) :: lsort_timerange
10391LOGICAL,OPTIONAL,INTENT(in) :: lsort_level
10392
10393INTEGER :: i
10394LOGICAL :: to_be_sorted
10395
10396to_be_sorted = .false.
10397CALL vol7d_alloc_vol(this) ! usual safety check
10398
10399IF (optio_log(lsort_time)) THEN
10400 DO i = 2, SIZE(this%time)
10401 IF (this%time(i) < this%time(i-1)) THEN
10402 to_be_sorted = .true.
10403 EXIT
10404 ENDIF
10405 ENDDO
10406ENDIF
10407IF (optio_log(lsort_timerange)) THEN
10408 DO i = 2, SIZE(this%timerange)
10409 IF (this%timerange(i) < this%timerange(i-1)) THEN
10410 to_be_sorted = .true.
10411 EXIT
10412 ENDIF
10413 ENDDO
10414ENDIF
10415IF (optio_log(lsort_level)) THEN
10416 DO i = 2, SIZE(this%level)
10417 IF (this%level(i) < this%level(i-1)) THEN
10418 to_be_sorted = .true.
10419 EXIT
10420 ENDIF
10421 ENDDO
10422ENDIF
10423
10424IF (to_be_sorted) CALL vol7d_reform(this, &
10425 lsort_time=lsort_time, lsort_timerange=lsort_timerange, lsort_level=lsort_level )
10426
10427END SUBROUTINE vol7d_smart_sort
10428
10436SUBROUTINE vol7d_filter(this, avl, vl, nl, s_d, e_d)
10437TYPE(vol7d),INTENT(inout) :: this
10438CHARACTER(len=*),INTENT(in),OPTIONAL :: avl(:)
10439CHARACTER(len=*),INTENT(in),OPTIONAL :: vl(:)
10440TYPE(vol7d_network),OPTIONAL :: nl(:)
10441TYPE(datetime),INTENT(in),OPTIONAL :: s_d
10442TYPE(datetime),INTENT(in),OPTIONAL :: e_d
10443
10444INTEGER :: i
10445
10446IF (PRESENT(avl)) THEN
10447 IF (SIZE(avl) > 0) THEN
10448
10449 IF (ASSOCIATED(this%anavar%r)) THEN
10450 DO i = 1, SIZE(this%anavar%r)
10451 IF (all(this%anavar%r(i)%btable /= avl)) this%anavar%r(i) = vol7d_var_miss
10452 ENDDO
10453 ENDIF
10454
10455 IF (ASSOCIATED(this%anavar%i)) THEN
10456 DO i = 1, SIZE(this%anavar%i)
10457 IF (all(this%anavar%i(i)%btable /= avl)) this%anavar%i(i) = vol7d_var_miss
10458 ENDDO
10459 ENDIF
10460
10461 IF (ASSOCIATED(this%anavar%b)) THEN
10462 DO i = 1, SIZE(this%anavar%b)
10463 IF (all(this%anavar%b(i)%btable /= avl)) this%anavar%b(i) = vol7d_var_miss
10464 ENDDO
10465 ENDIF
10466
10467 IF (ASSOCIATED(this%anavar%d)) THEN
10468 DO i = 1, SIZE(this%anavar%d)
10469 IF (all(this%anavar%d(i)%btable /= avl)) this%anavar%d(i) = vol7d_var_miss
10470 ENDDO
10471 ENDIF
10472
10473 IF (ASSOCIATED(this%anavar%c)) THEN
10474 DO i = 1, SIZE(this%anavar%c)
10475 IF (all(this%anavar%c(i)%btable /= avl)) this%anavar%c(i) = vol7d_var_miss
10476 ENDDO
10477 ENDIF
10478
10479 ENDIF
10480ENDIF
10481
10482
10483IF (PRESENT(vl)) THEN
10484 IF (size(vl) > 0) THEN
10485 IF (ASSOCIATED(this%dativar%r)) THEN
10486 DO i = 1, SIZE(this%dativar%r)
10487 IF (all(this%dativar%r(i)%btable /= vl)) this%dativar%r(i) = vol7d_var_miss
10488 ENDDO
10489 ENDIF
10490
10491 IF (ASSOCIATED(this%dativar%i)) THEN
10492 DO i = 1, SIZE(this%dativar%i)
10493 IF (all(this%dativar%i(i)%btable /= vl)) this%dativar%i(i) = vol7d_var_miss
10494 ENDDO
10495 ENDIF
10496
10497 IF (ASSOCIATED(this%dativar%b)) THEN
10498 DO i = 1, SIZE(this%dativar%b)
10499 IF (all(this%dativar%b(i)%btable /= vl)) this%dativar%b(i) = vol7d_var_miss
10500 ENDDO
10501 ENDIF
10502
10503 IF (ASSOCIATED(this%dativar%d)) THEN
10504 DO i = 1, SIZE(this%dativar%d)
10505 IF (all(this%dativar%d(i)%btable /= vl)) this%dativar%d(i) = vol7d_var_miss
10506 ENDDO
10507 ENDIF
10508
10509 IF (ASSOCIATED(this%dativar%c)) THEN
10510 DO i = 1, SIZE(this%dativar%c)
10511 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
10512 ENDDO
10513 ENDIF
10514
10515 IF (ASSOCIATED(this%dativar%c)) THEN
10516 DO i = 1, SIZE(this%dativar%c)
10517 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
10518 ENDDO
10519 ENDIF
10520
10521 ENDIF
10522ENDIF
10523
10524IF (PRESENT(nl)) THEN
10525 IF (SIZE(nl) > 0) THEN
10526 DO i = 1, SIZE(this%network)
10527 IF (all(this%network(i) /= nl)) this%network(i) = vol7d_network_miss
10528 ENDDO
10529 ENDIF
10530ENDIF
10531
10532IF (PRESENT(s_d)) THEN
10533 IF (c_e(s_d)) THEN
10534 WHERE (this%time < s_d)
10535 this%time = datetime_miss
10536 END WHERE
10537 ENDIF
10538ENDIF
10539
10540IF (PRESENT(e_d)) THEN
10541 IF (c_e(e_d)) THEN
10542 WHERE (this%time > e_d)
10543 this%time = datetime_miss
10544 END WHERE
10545 ENDIF
10546ENDIF
10547
10548CALL vol7d_reform(this, miss=.true.)
10549
10550END SUBROUTINE vol7d_filter
10551
10552
10559SUBROUTINE vol7d_convr(this, that, anaconv)
10560TYPE(vol7d),INTENT(IN) :: this
10561TYPE(vol7d),INTENT(INOUT) :: that
10562LOGICAL,OPTIONAL,INTENT(in) :: anaconv
10563INTEGER :: i
10564LOGICAL :: fv(1)=(/.false./), tv(1)=(/.true./), acp(1), acn(1)
10565TYPE(vol7d) :: v7d_tmp
10566
10567IF (optio_log(anaconv)) THEN
10568 acp=fv
10569 acn=tv
10570ELSE
10571 acp=tv
10572 acn=fv
10573ENDIF
10574
10575! Volume con solo i dati reali e tutti gli attributi
10576! l'anagrafica e` copiata interamente se necessario
10577CALL vol7d_copy(this, that, &
10578 lanavarr=tv, lanavard=acp, lanavari=acp, lanavarb=acp, lanavarc=acp, &
10579 ldativarr=tv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=fv)
10580
10581! Volume solo di dati double
10582CALL vol7d_copy(this, v7d_tmp, &
10583 lanavarr=fv, lanavard=acn, lanavari=fv, lanavarb=fv, lanavarc=fv, &
10584 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
10585 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
10586 ldativarr=fv, ldativard=tv, ldativari=fv, ldativarb=fv, ldativarc=fv, &
10587 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
10588 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
10589
10590! converto a dati reali
10591IF (ASSOCIATED(v7d_tmp%anavar%d) .OR. ASSOCIATED(v7d_tmp%dativar%d)) THEN
10592
10593 IF (ASSOCIATED(v7d_tmp%anavar%d)) THEN
10594! alloco i dati reali e vi trasferisco i double
10595 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanad, 1), SIZE(v7d_tmp%volanad, 2), &
10596 SIZE(v7d_tmp%volanad, 3)))
10597 DO i = 1, SIZE(v7d_tmp%anavar%d)
10598 v7d_tmp%volanar(:,i,:) = &
10599 realdat(v7d_tmp%volanad(:,i,:), v7d_tmp%anavar%d(i))
10600 ENDDO
10601 DEALLOCATE(v7d_tmp%volanad)
10602! trasferisco le variabili
10603 v7d_tmp%anavar%r => v7d_tmp%anavar%d
10604 NULLIFY(v7d_tmp%anavar%d)
10605 ENDIF
10606
10607 IF (ASSOCIATED(v7d_tmp%dativar%d)) THEN
10608! alloco i dati reali e vi trasferisco i double
10609 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatid, 1), SIZE(v7d_tmp%voldatid, 2), &
10610 SIZE(v7d_tmp%voldatid, 3), SIZE(v7d_tmp%voldatid, 4), SIZE(v7d_tmp%voldatid, 5), &
10611 SIZE(v7d_tmp%voldatid, 6)))
10612 DO i = 1, SIZE(v7d_tmp%dativar%d)
10613 v7d_tmp%voldatir(:,:,:,:,i,:) = &
10614 realdat(v7d_tmp%voldatid(:,:,:,:,i,:), v7d_tmp%dativar%d(i))
10615 ENDDO
10616 DEALLOCATE(v7d_tmp%voldatid)
10617! trasferisco le variabili
10618 v7d_tmp%dativar%r => v7d_tmp%dativar%d
10619 NULLIFY(v7d_tmp%dativar%d)
10620 ENDIF
10621
10622! fondo con il volume definitivo
10623 CALL vol7d_merge(that, v7d_tmp)
10624ELSE
10625 CALL delete(v7d_tmp)
10626ENDIF
10627
10628
10629! Volume solo di dati interi
10630CALL vol7d_copy(this, v7d_tmp, &
10631 lanavarr=fv, lanavard=fv, lanavari=acn, lanavarb=fv, lanavarc=fv, &
10632 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
10633 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
10634 ldativarr=fv, ldativard=fv, ldativari=tv, ldativarb=fv, ldativarc=fv, &
10635 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
10636 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
10637
10638! converto a dati reali
10639IF (ASSOCIATED(v7d_tmp%anavar%i) .OR. ASSOCIATED(v7d_tmp%dativar%i)) THEN
10640
10641 IF (ASSOCIATED(v7d_tmp%anavar%i)) THEN
10642! alloco i dati reali e vi trasferisco gli interi
10643 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanai, 1), SIZE(v7d_tmp%volanai, 2), &
10644 SIZE(v7d_tmp%volanai, 3)))
10645 DO i = 1, SIZE(v7d_tmp%anavar%i)
10646 v7d_tmp%volanar(:,i,:) = &
10647 realdat(v7d_tmp%volanai(:,i,:), v7d_tmp%anavar%i(i))
10648 ENDDO
10649 DEALLOCATE(v7d_tmp%volanai)
10650! trasferisco le variabili
10651 v7d_tmp%anavar%r => v7d_tmp%anavar%i
10652 NULLIFY(v7d_tmp%anavar%i)
10653 ENDIF
10654
10655 IF (ASSOCIATED(v7d_tmp%dativar%i)) THEN
10656! alloco i dati reali e vi trasferisco gli interi
10657 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatii, 1), SIZE(v7d_tmp%voldatii, 2), &
10658 SIZE(v7d_tmp%voldatii, 3), SIZE(v7d_tmp%voldatii, 4), SIZE(v7d_tmp%voldatii, 5), &
10659 SIZE(v7d_tmp%voldatii, 6)))
10660 DO i = 1, SIZE(v7d_tmp%dativar%i)
10661 v7d_tmp%voldatir(:,:,:,:,i,:) = &
10662 realdat(v7d_tmp%voldatii(:,:,:,:,i,:), v7d_tmp%dativar%i(i))
10663 ENDDO
10664 DEALLOCATE(v7d_tmp%voldatii)
10665! trasferisco le variabili
10666 v7d_tmp%dativar%r => v7d_tmp%dativar%i
10667 NULLIFY(v7d_tmp%dativar%i)
10668 ENDIF
10669
10670! fondo con il volume definitivo
10671 CALL vol7d_merge(that, v7d_tmp)
10672ELSE
10673 CALL delete(v7d_tmp)
10674ENDIF
10675
10676
10677! Volume solo di dati byte
10678CALL vol7d_copy(this, v7d_tmp, &
10679 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=acn, lanavarc=fv, &
10680 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
10681 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
10682 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=tv, ldativarc=fv, &
10683 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
10684 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
10685
10686! converto a dati reali
10687IF (ASSOCIATED(v7d_tmp%anavar%b) .OR. ASSOCIATED(v7d_tmp%dativar%b)) THEN
10688
10689 IF (ASSOCIATED(v7d_tmp%anavar%b)) THEN
10690! alloco i dati reali e vi trasferisco i byte
10691 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanab, 1), SIZE(v7d_tmp%volanab, 2), &
10692 SIZE(v7d_tmp%volanab, 3)))
10693 DO i = 1, SIZE(v7d_tmp%anavar%b)
10694 v7d_tmp%volanar(:,i,:) = &
10695 realdat(v7d_tmp%volanab(:,i,:), v7d_tmp%anavar%b(i))
10696 ENDDO
10697 DEALLOCATE(v7d_tmp%volanab)
10698! trasferisco le variabili
10699 v7d_tmp%anavar%r => v7d_tmp%anavar%b
10700 NULLIFY(v7d_tmp%anavar%b)
10701 ENDIF
10702
10703 IF (ASSOCIATED(v7d_tmp%dativar%b)) THEN
10704! alloco i dati reali e vi trasferisco i byte
10705 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatib, 1), SIZE(v7d_tmp%voldatib, 2), &
10706 SIZE(v7d_tmp%voldatib, 3), SIZE(v7d_tmp%voldatib, 4), SIZE(v7d_tmp%voldatib, 5), &
10707 SIZE(v7d_tmp%voldatib, 6)))
10708 DO i = 1, SIZE(v7d_tmp%dativar%b)
10709 v7d_tmp%voldatir(:,:,:,:,i,:) = &
10710 realdat(v7d_tmp%voldatib(:,:,:,:,i,:), v7d_tmp%dativar%b(i))
10711 ENDDO
10712 DEALLOCATE(v7d_tmp%voldatib)
10713! trasferisco le variabili
10714 v7d_tmp%dativar%r => v7d_tmp%dativar%b
10715 NULLIFY(v7d_tmp%dativar%b)
10716 ENDIF
10717
10718! fondo con il volume definitivo
10719 CALL vol7d_merge(that, v7d_tmp)
10720ELSE
10721 CALL delete(v7d_tmp)
10722ENDIF
10723
10724
10725! Volume solo di dati character
10726CALL vol7d_copy(this, v7d_tmp, &
10727 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=fv, lanavarc=acn, &
10728 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
10729 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
10730 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=tv, &
10731 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
10732 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
10733
10734! converto a dati reali
10735IF (ASSOCIATED(v7d_tmp%anavar%c) .OR. ASSOCIATED(v7d_tmp%dativar%c)) THEN
10736
10737 IF (ASSOCIATED(v7d_tmp%anavar%c)) THEN
10738! alloco i dati reali e vi trasferisco i character
10739 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanac, 1), SIZE(v7d_tmp%volanac, 2), &
10740 SIZE(v7d_tmp%volanac, 3)))
10741 DO i = 1, SIZE(v7d_tmp%anavar%c)
10742 v7d_tmp%volanar(:,i,:) = &
10743 realdat(v7d_tmp%volanac(:,i,:), v7d_tmp%anavar%c(i))
10744 ENDDO
10745 DEALLOCATE(v7d_tmp%volanac)
10746! trasferisco le variabili
10747 v7d_tmp%anavar%r => v7d_tmp%anavar%c
10748 NULLIFY(v7d_tmp%anavar%c)
10749 ENDIF
10750
10751 IF (ASSOCIATED(v7d_tmp%dativar%c)) THEN
10752! alloco i dati reali e vi trasferisco i character
10753 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatic, 1), SIZE(v7d_tmp%voldatic, 2), &
10754 SIZE(v7d_tmp%voldatic, 3), SIZE(v7d_tmp%voldatic, 4), SIZE(v7d_tmp%voldatic, 5), &
10755 SIZE(v7d_tmp%voldatic, 6)))
10756 DO i = 1, SIZE(v7d_tmp%dativar%c)
10757 v7d_tmp%voldatir(:,:,:,:,i,:) = &
10758 realdat(v7d_tmp%voldatic(:,:,:,:,i,:), v7d_tmp%dativar%c(i))
10759 ENDDO
10760 DEALLOCATE(v7d_tmp%voldatic)
10761! trasferisco le variabili
10762 v7d_tmp%dativar%r => v7d_tmp%dativar%c
10763 NULLIFY(v7d_tmp%dativar%c)
10764 ENDIF
10765
10766! fondo con il volume definitivo
10767 CALL vol7d_merge(that, v7d_tmp)
10768ELSE
10769 CALL delete(v7d_tmp)
10770ENDIF
10771
10772END SUBROUTINE vol7d_convr
10773
10774
10778SUBROUTINE vol7d_diff_only (this, that, data_only,ana)
10779TYPE(vol7d),INTENT(IN) :: this
10780TYPE(vol7d),INTENT(OUT) :: that
10781logical , optional, intent(in) :: data_only
10782logical , optional, intent(in) :: ana
10783logical :: ldata_only,lana
10784
10785IF (PRESENT(data_only)) THEN
10786 ldata_only = data_only
10787ELSE
10788 ldata_only = .false.
10789ENDIF
10790
10791IF (PRESENT(ana)) THEN
10792 lana = ana
10793ELSE
10794 lana = .false.
10795ENDIF
10796
10797
10798#undef VOL7D_POLY_ARRAY
10799#define VOL7D_POLY_ARRAY voldati
10800#include "vol7d_class_diff.F90"
10801#undef VOL7D_POLY_ARRAY
10802#define VOL7D_POLY_ARRAY voldatiattr
10803#include "vol7d_class_diff.F90"
10804#undef VOL7D_POLY_ARRAY
10805
10806if ( .not. ldata_only) then
10807
10808#define VOL7D_POLY_ARRAY volana
10809#include "vol7d_class_diff.F90"
10810#undef VOL7D_POLY_ARRAY
10811#define VOL7D_POLY_ARRAY volanaattr
10812#include "vol7d_class_diff.F90"
10813#undef VOL7D_POLY_ARRAY
10814
10815 if(lana)then
10816 where ( this%ana == that%ana )
10817 that%ana = vol7d_ana_miss
10818 end where
10819 end if
10820
10821end if
10822
10823
10824
10825END SUBROUTINE vol7d_diff_only
10826
10827
10828
10829! Creo le routine da ripetere per i vari tipi di dati di v7d
10830! tramite un template e il preprocessore
10831#undef VOL7D_POLY_TYPE
10832#undef VOL7D_POLY_TYPES
10833#define VOL7D_POLY_TYPE REAL
10834#define VOL7D_POLY_TYPES r
10835#include "vol7d_class_type_templ.F90"
10836#undef VOL7D_POLY_TYPE
10837#undef VOL7D_POLY_TYPES
10838#define VOL7D_POLY_TYPE DOUBLE PRECISION
10839#define VOL7D_POLY_TYPES d
10840#include "vol7d_class_type_templ.F90"
10841#undef VOL7D_POLY_TYPE
10842#undef VOL7D_POLY_TYPES
10843#define VOL7D_POLY_TYPE INTEGER
10844#define VOL7D_POLY_TYPES i
10845#include "vol7d_class_type_templ.F90"
10846#undef VOL7D_POLY_TYPE
10847#undef VOL7D_POLY_TYPES
10848#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
10849#define VOL7D_POLY_TYPES b
10850#include "vol7d_class_type_templ.F90"
10851#undef VOL7D_POLY_TYPE
10852#undef VOL7D_POLY_TYPES
10853#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
10854#define VOL7D_POLY_TYPES c
10855#include "vol7d_class_type_templ.F90"
10856
10857! Creo le routine da ripetere per i vari descrittori di dimensioni di v7d
10858! tramite un template e il preprocessore
10859#define VOL7D_SORT
10860#undef VOL7D_NO_ZERO_ALLOC
10861#undef VOL7D_POLY_TYPE
10862#define VOL7D_POLY_TYPE datetime
10863#include "vol7d_class_desc_templ.F90"
10864#undef VOL7D_POLY_TYPE
10865#define VOL7D_POLY_TYPE vol7d_timerange
10866#include "vol7d_class_desc_templ.F90"
10867#undef VOL7D_POLY_TYPE
10868#define VOL7D_POLY_TYPE vol7d_level
10869#include "vol7d_class_desc_templ.F90"
10870#undef VOL7D_SORT
10871#undef VOL7D_POLY_TYPE
10872#define VOL7D_POLY_TYPE vol7d_network
10873#include "vol7d_class_desc_templ.F90"
10874#undef VOL7D_POLY_TYPE
10875#define VOL7D_POLY_TYPE vol7d_ana
10876#include "vol7d_class_desc_templ.F90"
10877#define VOL7D_NO_ZERO_ALLOC
10878#undef VOL7D_POLY_TYPE
10879#define VOL7D_POLY_TYPE vol7d_var
10880#include "vol7d_class_desc_templ.F90"
10881
10891subroutine vol7d_write_on_file (this,unit,description,filename,filename_auto)
10892
10893TYPE(vol7d),INTENT(IN) :: this
10894integer,optional,intent(inout) :: unit
10895character(len=*),intent(in),optional :: filename
10896character(len=*),intent(out),optional :: filename_auto
10897character(len=*),INTENT(IN),optional :: description
10898
10899integer :: lunit
10900character(len=254) :: ldescription,arg,lfilename
10901integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
10902 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
10903 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
10904 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
10905 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
10906 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
10907 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
10908!integer :: im,id,iy
10909integer :: tarray(8)
10910logical :: opened,exist
10911
10912 nana=0
10913 ntime=0
10914 ntimerange=0
10915 nlevel=0
10916 nnetwork=0
10917 ndativarr=0
10918 ndativari=0
10919 ndativarb=0
10920 ndativard=0
10921 ndativarc=0
10922 ndatiattrr=0
10923 ndatiattri=0
10924 ndatiattrb=0
10925 ndatiattrd=0
10926 ndatiattrc=0
10927 ndativarattrr=0
10928 ndativarattri=0
10929 ndativarattrb=0
10930 ndativarattrd=0
10931 ndativarattrc=0
10932 nanavarr=0
10933 nanavari=0
10934 nanavarb=0
10935 nanavard=0
10936 nanavarc=0
10937 nanaattrr=0
10938 nanaattri=0
10939 nanaattrb=0
10940 nanaattrd=0
10941 nanaattrc=0
10942 nanavarattrr=0
10943 nanavarattri=0
10944 nanavarattrb=0
10945 nanavarattrd=0
10946 nanavarattrc=0
10947
10948
10949!call idate(im,id,iy)
10950call date_and_time(values=tarray)
10951call getarg(0,arg)
10952
10953if (present(description))then
10954 ldescription=description
10955else
10956 ldescription="Vol7d generated by: "//trim(arg)
10957end if
10958
10959if (.not. present(unit))then
10960 lunit=getunit()
10961else
10962 if (unit==0)then
10963 lunit=getunit()
10964 unit=lunit
10965 else
10966 lunit=unit
10967 end if
10968end if
10969
10970lfilename=trim(arg)//".v7d"
10971if (index(arg,'/',back=.true.) > 0) lfilename=lfilename(index(arg,'/',back=.true.)+1 : )
10972
10973if (present(filename))then
10974 if (filename /= "")then
10975 lfilename=filename
10976 end if
10977end if
10978
10979if (present(filename_auto))filename_auto=lfilename
10980
10981
10982inquire(unit=lunit,opened=opened)
10983if (.not. opened) then
10984! inquire(file=lfilename, EXIST=exist)
10985! IF (exist) THEN
10986! CALL l4f_log(L4F_FATAL, &
10987! 'in vol7d_write_on_file, file exists, cannot open file '//TRIM(lfilename))
10988! CALL raise_fatal_error()
10989! ENDIF
10990 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM')
10991 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
10992end if
10993
10994if (associated(this%ana)) nana=size(this%ana)
10995if (associated(this%time)) ntime=size(this%time)
10996if (associated(this%timerange)) ntimerange=size(this%timerange)
10997if (associated(this%level)) nlevel=size(this%level)
10998if (associated(this%network)) nnetwork=size(this%network)
10999
11000if (associated(this%dativar%r)) ndativarr=size(this%dativar%r)
11001if (associated(this%dativar%i)) ndativari=size(this%dativar%i)
11002if (associated(this%dativar%b)) ndativarb=size(this%dativar%b)
11003if (associated(this%dativar%d)) ndativard=size(this%dativar%d)
11004if (associated(this%dativar%c)) ndativarc=size(this%dativar%c)
11005
11006if (associated(this%datiattr%r)) ndatiattrr=size(this%datiattr%r)
11007if (associated(this%datiattr%i)) ndatiattri=size(this%datiattr%i)
11008if (associated(this%datiattr%b)) ndatiattrb=size(this%datiattr%b)
11009if (associated(this%datiattr%d)) ndatiattrd=size(this%datiattr%d)
11010if (associated(this%datiattr%c)) ndatiattrc=size(this%datiattr%c)
11011
11012if (associated(this%dativarattr%r)) ndativarattrr=size(this%dativarattr%r)
11013if (associated(this%dativarattr%i)) ndativarattri=size(this%dativarattr%i)
11014if (associated(this%dativarattr%b)) ndativarattrb=size(this%dativarattr%b)
11015if (associated(this%dativarattr%d)) ndativarattrd=size(this%dativarattr%d)
11016if (associated(this%dativarattr%c)) ndativarattrc=size(this%dativarattr%c)
11017
11018if (associated(this%anavar%r)) nanavarr=size(this%anavar%r)
11019if (associated(this%anavar%i)) nanavari=size(this%anavar%i)
11020if (associated(this%anavar%b)) nanavarb=size(this%anavar%b)
11021if (associated(this%anavar%d)) nanavard=size(this%anavar%d)
11022if (associated(this%anavar%c)) nanavarc=size(this%anavar%c)
11023
11024if (associated(this%anaattr%r)) nanaattrr=size(this%anaattr%r)
11025if (associated(this%anaattr%i)) nanaattri=size(this%anaattr%i)
11026if (associated(this%anaattr%b)) nanaattrb=size(this%anaattr%b)
11027if (associated(this%anaattr%d)) nanaattrd=size(this%anaattr%d)
11028if (associated(this%anaattr%c)) nanaattrc=size(this%anaattr%c)
11029
11030if (associated(this%anavarattr%r)) nanavarattrr=size(this%anavarattr%r)
11031if (associated(this%anavarattr%i)) nanavarattri=size(this%anavarattr%i)
11032if (associated(this%anavarattr%b)) nanavarattrb=size(this%anavarattr%b)
11033if (associated(this%anavarattr%d)) nanavarattrd=size(this%anavarattr%d)
11034if (associated(this%anavarattr%c)) nanavarattrc=size(this%anavarattr%c)
11035
11036write(unit=lunit)ldescription
11037write(unit=lunit)tarray
11038
11039write(unit=lunit)&
11040 nana, ntime, ntimerange, nlevel, nnetwork, &
11041 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11042 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11043 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11044 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11045 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11046 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
11047 this%time_definition
11048
11049
11050!write(unit=lunit)this
11051
11052
11053!! prime 5 dimensioni
11054if (associated(this%ana)) call write_unit(this%ana, lunit)
11055if (associated(this%time)) call write_unit(this%time, lunit)
11056if (associated(this%level)) write(unit=lunit)this%level
11057if (associated(this%timerange)) write(unit=lunit)this%timerange
11058if (associated(this%network)) write(unit=lunit)this%network
11059
11060 !! 6a dimensione: variabile dell'anagrafica e dei dati
11061 !! con relativi attributi e in 5 tipi diversi
11062
11063if (associated(this%anavar%r)) write(unit=lunit)this%anavar%r
11064if (associated(this%anavar%i)) write(unit=lunit)this%anavar%i
11065if (associated(this%anavar%b)) write(unit=lunit)this%anavar%b
11066if (associated(this%anavar%d)) write(unit=lunit)this%anavar%d
11067if (associated(this%anavar%c)) write(unit=lunit)this%anavar%c
11068
11069if (associated(this%anaattr%r)) write(unit=lunit)this%anaattr%r
11070if (associated(this%anaattr%i)) write(unit=lunit)this%anaattr%i
11071if (associated(this%anaattr%b)) write(unit=lunit)this%anaattr%b
11072if (associated(this%anaattr%d)) write(unit=lunit)this%anaattr%d
11073if (associated(this%anaattr%c)) write(unit=lunit)this%anaattr%c
11074
11075if (associated(this%anavarattr%r)) write(unit=lunit)this%anavarattr%r
11076if (associated(this%anavarattr%i)) write(unit=lunit)this%anavarattr%i
11077if (associated(this%anavarattr%b)) write(unit=lunit)this%anavarattr%b
11078if (associated(this%anavarattr%d)) write(unit=lunit)this%anavarattr%d
11079if (associated(this%anavarattr%c)) write(unit=lunit)this%anavarattr%c
11080
11081if (associated(this%dativar%r)) write(unit=lunit)this%dativar%r
11082if (associated(this%dativar%i)) write(unit=lunit)this%dativar%i
11083if (associated(this%dativar%b)) write(unit=lunit)this%dativar%b
11084if (associated(this%dativar%d)) write(unit=lunit)this%dativar%d
11085if (associated(this%dativar%c)) write(unit=lunit)this%dativar%c
11086
11087if (associated(this%datiattr%r)) write(unit=lunit)this%datiattr%r
11088if (associated(this%datiattr%i)) write(unit=lunit)this%datiattr%i
11089if (associated(this%datiattr%b)) write(unit=lunit)this%datiattr%b
11090if (associated(this%datiattr%d)) write(unit=lunit)this%datiattr%d
11091if (associated(this%datiattr%c)) write(unit=lunit)this%datiattr%c
11092
11093if (associated(this%dativarattr%r)) write(unit=lunit)this%dativarattr%r
11094if (associated(this%dativarattr%i)) write(unit=lunit)this%dativarattr%i
11095if (associated(this%dativarattr%b)) write(unit=lunit)this%dativarattr%b
11096if (associated(this%dativarattr%d)) write(unit=lunit)this%dativarattr%d
11097if (associated(this%dativarattr%c)) write(unit=lunit)this%dativarattr%c
11098
11099!! Volumi di valori e attributi per anagrafica e dati
11100
11101if (associated(this%volanar)) write(unit=lunit)this%volanar
11102if (associated(this%volanaattrr)) write(unit=lunit)this%volanaattrr
11103if (associated(this%voldatir)) write(unit=lunit)this%voldatir
11104if (associated(this%voldatiattrr)) write(unit=lunit)this%voldatiattrr
11105
11106if (associated(this%volanai)) write(unit=lunit)this%volanai
11107if (associated(this%volanaattri)) write(unit=lunit)this%volanaattri
11108if (associated(this%voldatii)) write(unit=lunit)this%voldatii
11109if (associated(this%voldatiattri)) write(unit=lunit)this%voldatiattri
11110
11111if (associated(this%volanab)) write(unit=lunit)this%volanab
11112if (associated(this%volanaattrb)) write(unit=lunit)this%volanaattrb
11113if (associated(this%voldatib)) write(unit=lunit)this%voldatib
11114if (associated(this%voldatiattrb)) write(unit=lunit)this%voldatiattrb
11115
11116if (associated(this%volanad)) write(unit=lunit)this%volanad
11117if (associated(this%volanaattrd)) write(unit=lunit)this%volanaattrd
11118if (associated(this%voldatid)) write(unit=lunit)this%voldatid
11119if (associated(this%voldatiattrd)) write(unit=lunit)this%voldatiattrd
11120
11121if (associated(this%volanac)) write(unit=lunit)this%volanac
11122if (associated(this%volanaattrc)) write(unit=lunit)this%volanaattrc
11123if (associated(this%voldatic)) write(unit=lunit)this%voldatic
11124if (associated(this%voldatiattrc)) write(unit=lunit)this%voldatiattrc
11125
11126if (.not. present(unit)) close(unit=lunit)
11127
11128end subroutine vol7d_write_on_file
11129
11130
11137
11138
11139subroutine vol7d_read_from_file (this,unit,filename,description,tarray,filename_auto)
11140
11141TYPE(vol7d),INTENT(OUT) :: this
11142integer,intent(inout),optional :: unit
11143character(len=*),INTENT(in),optional :: filename
11144character(len=*),intent(out),optional :: filename_auto
11145character(len=*),INTENT(out),optional :: description
11146integer,intent(out),optional :: tarray(8)
11147
11148
11149integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
11150 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11151 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11152 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11153 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11154 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11155 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
11156
11157character(len=254) :: ldescription,lfilename,arg
11158integer :: ltarray(8),lunit,ios
11159logical :: opened,exist
11160
11161
11162call getarg(0,arg)
11163
11164if (.not. present(unit))then
11165 lunit=getunit()
11166else
11167 if (unit==0)then
11168 lunit=getunit()
11169 unit=lunit
11170 else
11171 lunit=unit
11172 end if
11173end if
11174
11175lfilename=trim(arg)//".v7d"
11176if (index(arg,'/',back=.true.) > 0) lfilename=lfilename(index(arg,'/',back=.true.)+1 : )
11177
11178if (present(filename))then
11179 if (filename /= "")then
11180 lfilename=filename
11181 end if
11182end if
11183
11184if (present(filename_auto))filename_auto=lfilename
11185
11186
11187inquire(unit=lunit,opened=opened)
11188IF (.NOT. opened) THEN
11189 inquire(file=lfilename,exist=exist)
11190 IF (.NOT.exist) THEN
11191 CALL l4f_log(l4f_fatal, &
11192 'in vol7d_read_from_file, file does not exists, cannot open')
11193 CALL raise_fatal_error()
11194 ENDIF
11195 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM', &
11196 status='OLD', action='READ')
11197 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
11198end if
11199
11200
11201call init(this)
11202read(unit=lunit,iostat=ios)ldescription
11203
11204if (ios < 0) then ! A negative value indicates that the End of File or End of Record
11205 call vol7d_alloc (this)
11206 call vol7d_alloc_vol (this)
11207 if (present(description))description=ldescription
11208 if (present(tarray))tarray=ltarray
11209 if (.not. present(unit)) close(unit=lunit)
11210end if
11211
11212read(unit=lunit)ltarray
11213
11214CALL l4f_log(l4f_info, 'Reading vol7d from file')
11215CALL l4f_log(l4f_info, 'description: '//trim(ldescription))
11216CALL l4f_log(l4f_info, 'written on '//trim(to_char(ltarray(1)))//' '// &
11217 trim(to_char(ltarray(2)))//' '//trim(to_char(ltarray(3))))
11218
11219if (present(description))description=ldescription
11220if (present(tarray))tarray=ltarray
11221
11222read(unit=lunit)&
11223 nana, ntime, ntimerange, nlevel, nnetwork, &
11224 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11225 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11226 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11227 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11228 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11229 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
11230 this%time_definition
11231
11232call vol7d_alloc (this, &
11233 nana=nana, ntime=ntime, ntimerange=ntimerange, nlevel=nlevel, nnetwork=nnetwork,&
11234 ndativarr=ndativarr, ndativari=ndativari, ndativarb=ndativarb,&
11235 ndativard=ndativard, ndativarc=ndativarc,&
11236 ndatiattrr=ndatiattrr, ndatiattri=ndatiattri, ndatiattrb=ndatiattrb,&
11237 ndatiattrd=ndatiattrd, ndatiattrc=ndatiattrc,&
11238 ndativarattrr=ndativarattrr, ndativarattri=ndativarattri, ndativarattrb=ndativarattrb, &
11239 ndativarattrd=ndativarattrd, ndativarattrc=ndativarattrc,&
11240 nanavarr=nanavarr, nanavari=nanavari, nanavarb=nanavarb, &
11241 nanavard=nanavard, nanavarc=nanavarc,&
11242 nanaattrr=nanaattrr, nanaattri=nanaattri, nanaattrb=nanaattrb,&
11243 nanaattrd=nanaattrd, nanaattrc=nanaattrc,&
11244 nanavarattrr=nanavarattrr, nanavarattri=nanavarattri, nanavarattrb=nanavarattrb, &
11245 nanavarattrd=nanavarattrd, nanavarattrc=nanavarattrc)
11246
11247
11248if (associated(this%ana)) call read_unit(this%ana, lunit)
11249if (associated(this%time)) call read_unit(this%time, lunit)
11250if (associated(this%level)) read(unit=lunit)this%level
11251if (associated(this%timerange)) read(unit=lunit)this%timerange
11252if (associated(this%network)) read(unit=lunit)this%network
11253
11254if (associated(this%anavar%r)) read(unit=lunit)this%anavar%r
11255if (associated(this%anavar%i)) read(unit=lunit)this%anavar%i
11256if (associated(this%anavar%b)) read(unit=lunit)this%anavar%b
11257if (associated(this%anavar%d)) read(unit=lunit)this%anavar%d
11258if (associated(this%anavar%c)) read(unit=lunit)this%anavar%c
11259
11260if (associated(this%anaattr%r)) read(unit=lunit)this%anaattr%r
11261if (associated(this%anaattr%i)) read(unit=lunit)this%anaattr%i
11262if (associated(this%anaattr%b)) read(unit=lunit)this%anaattr%b
11263if (associated(this%anaattr%d)) read(unit=lunit)this%anaattr%d
11264if (associated(this%anaattr%c)) read(unit=lunit)this%anaattr%c
11265
11266if (associated(this%anavarattr%r)) read(unit=lunit)this%anavarattr%r
11267if (associated(this%anavarattr%i)) read(unit=lunit)this%anavarattr%i
11268if (associated(this%anavarattr%b)) read(unit=lunit)this%anavarattr%b
11269if (associated(this%anavarattr%d)) read(unit=lunit)this%anavarattr%d
11270if (associated(this%anavarattr%c)) read(unit=lunit)this%anavarattr%c
11271
11272if (associated(this%dativar%r)) read(unit=lunit)this%dativar%r
11273if (associated(this%dativar%i)) read(unit=lunit)this%dativar%i
11274if (associated(this%dativar%b)) read(unit=lunit)this%dativar%b
11275if (associated(this%dativar%d)) read(unit=lunit)this%dativar%d
11276if (associated(this%dativar%c)) read(unit=lunit)this%dativar%c
11277
11278if (associated(this%datiattr%r)) read(unit=lunit)this%datiattr%r
11279if (associated(this%datiattr%i)) read(unit=lunit)this%datiattr%i
11280if (associated(this%datiattr%b)) read(unit=lunit)this%datiattr%b
11281if (associated(this%datiattr%d)) read(unit=lunit)this%datiattr%d
11282if (associated(this%datiattr%c)) read(unit=lunit)this%datiattr%c
11283
11284if (associated(this%dativarattr%r)) read(unit=lunit)this%dativarattr%r
11285if (associated(this%dativarattr%i)) read(unit=lunit)this%dativarattr%i
11286if (associated(this%dativarattr%b)) read(unit=lunit)this%dativarattr%b
11287if (associated(this%dativarattr%d)) read(unit=lunit)this%dativarattr%d
11288if (associated(this%dativarattr%c)) read(unit=lunit)this%dativarattr%c
11289
11290call vol7d_alloc_vol (this)
11291
11292!! Volumi di valori e attributi per anagrafica e dati
11293
11294if (associated(this%volanar)) read(unit=lunit)this%volanar
11295if (associated(this%volanaattrr)) read(unit=lunit)this%volanaattrr
11296if (associated(this%voldatir)) read(unit=lunit)this%voldatir
11297if (associated(this%voldatiattrr)) read(unit=lunit)this%voldatiattrr
11298
11299if (associated(this%volanai)) read(unit=lunit)this%volanai
11300if (associated(this%volanaattri)) read(unit=lunit)this%volanaattri
11301if (associated(this%voldatii)) read(unit=lunit)this%voldatii
11302if (associated(this%voldatiattri)) read(unit=lunit)this%voldatiattri
11303
11304if (associated(this%volanab)) read(unit=lunit)this%volanab
11305if (associated(this%volanaattrb)) read(unit=lunit)this%volanaattrb
11306if (associated(this%voldatib)) read(unit=lunit)this%voldatib
11307if (associated(this%voldatiattrb)) read(unit=lunit)this%voldatiattrb
11308
11309if (associated(this%volanad)) read(unit=lunit)this%volanad
11310if (associated(this%volanaattrd)) read(unit=lunit)this%volanaattrd
11311if (associated(this%voldatid)) read(unit=lunit)this%voldatid
11312if (associated(this%voldatiattrd)) read(unit=lunit)this%voldatiattrd
11313
11314if (associated(this%volanac)) read(unit=lunit)this%volanac
11315if (associated(this%volanaattrc)) read(unit=lunit)this%volanaattrc
11316if (associated(this%voldatic)) read(unit=lunit)this%voldatic
11317if (associated(this%voldatiattrc)) read(unit=lunit)this%voldatiattrc
11318
11319if (.not. present(unit)) close(unit=lunit)
11320
11321end subroutine vol7d_read_from_file
11322
11323
11324! to double precision
11325elemental doubleprecision function doubledatd(voldat,var)
11326doubleprecision,intent(in) :: voldat
11327type(vol7d_var),intent(in) :: var
11328
11329doubledatd=voldat
11330
11331end function doubledatd
11332
11333
11334elemental doubleprecision function doubledatr(voldat,var)
11335real,intent(in) :: voldat
11336type(vol7d_var),intent(in) :: var
11337
11338if (c_e(voldat))then
11339 doubledatr=dble(voldat)
11340else
11341 doubledatr=dmiss
11342end if
11343
11344end function doubledatr
11345
11346
11347elemental doubleprecision function doubledati(voldat,var)
11348integer,intent(in) :: voldat
11349type(vol7d_var),intent(in) :: var
11350
11351if (c_e(voldat)) then
11352 if (c_e(var%scalefactor))then
11353 doubledati=dble(voldat)/10.d0**var%scalefactor
11354 else
11355 doubledati=dble(voldat)
11356 endif
11357else
11358 doubledati=dmiss
11359end if
11360
11361end function doubledati
11362
11363
11364elemental doubleprecision function doubledatb(voldat,var)
11365integer(kind=int_b),intent(in) :: voldat
11366type(vol7d_var),intent(in) :: var
11367
11368if (c_e(voldat)) then
11369 if (c_e(var%scalefactor))then
11370 doubledatb=dble(voldat)/10.d0**var%scalefactor
11371 else
11372 doubledatb=dble(voldat)
11373 endif
11374else
11375 doubledatb=dmiss
11376end if
11377
11378end function doubledatb
11379
11380
11381elemental doubleprecision function doubledatc(voldat,var)
11382CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
11383type(vol7d_var),intent(in) :: var
11384
11385doubledatc = c2d(voldat)
11386if (c_e(doubledatc) .and. c_e(var%scalefactor))then
11387 doubledatc=doubledatc/10.d0**var%scalefactor
11388end if
11389
11390end function doubledatc
11391
11392
11393! to integer
11394elemental integer function integerdatd(voldat,var)
11395doubleprecision,intent(in) :: voldat
11396type(vol7d_var),intent(in) :: var
11397
11398if (c_e(voldat))then
11399 if (c_e(var%scalefactor)) then
11400 integerdatd=nint(voldat*10d0**var%scalefactor)
11401 else
11402 integerdatd=nint(voldat)
11403 endif
11404else
11405 integerdatd=imiss
11406end if
11407
11408end function integerdatd
11409
11410
11411elemental integer function integerdatr(voldat,var)
11412real,intent(in) :: voldat
11413type(vol7d_var),intent(in) :: var
11414
11415if (c_e(voldat))then
11416 if (c_e(var%scalefactor)) then
11417 integerdatr=nint(voldat*10d0**var%scalefactor)
11418 else
11419 integerdatr=nint(voldat)
11420 endif
11421else
11422 integerdatr=imiss
11423end if
11424
11425end function integerdatr
11426
11427
11428elemental integer function integerdati(voldat,var)
11429integer,intent(in) :: voldat
11430type(vol7d_var),intent(in) :: var
11431
11432integerdati=voldat
11433
11434end function integerdati
11435
11436
11437elemental integer function integerdatb(voldat,var)
11438integer(kind=int_b),intent(in) :: voldat
11439type(vol7d_var),intent(in) :: var
11440
11441if (c_e(voldat))then
11442 integerdatb=voldat
11443else
11444 integerdatb=imiss
11445end if
11446
11447end function integerdatb
11448
11449
11450elemental integer function integerdatc(voldat,var)
11451CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
11452type(vol7d_var),intent(in) :: var
11453
11454integerdatc=c2i(voldat)
11455
11456end function integerdatc
11457
11458
11459! to real
11460elemental real function realdatd(voldat,var)
11461doubleprecision,intent(in) :: voldat
11462type(vol7d_var),intent(in) :: var
11463
11464if (c_e(voldat))then
11465 realdatd=real(voldat)
11466else
11467 realdatd=rmiss
11468end if
11469
11470end function realdatd
11471
11472
11473elemental real function realdatr(voldat,var)
11474real,intent(in) :: voldat
11475type(vol7d_var),intent(in) :: var
11476
11477realdatr=voldat
11478
11479end function realdatr
11480
11481
11482elemental real function realdati(voldat,var)
11483integer,intent(in) :: voldat
11484type(vol7d_var),intent(in) :: var
11485
11486if (c_e(voldat)) then
11487 if (c_e(var%scalefactor))then
11488 realdati=float(voldat)/10.**var%scalefactor
11489 else
11490 realdati=float(voldat)
11491 endif
11492else
11493 realdati=rmiss
11494end if
11495
11496end function realdati
11497
11498
11499elemental real function realdatb(voldat,var)
11500integer(kind=int_b),intent(in) :: voldat
11501type(vol7d_var),intent(in) :: var
11502
11503if (c_e(voldat)) then
11504 if (c_e(var%scalefactor))then
11505 realdatb=float(voldat)/10**var%scalefactor
11506 else
11507 realdatb=float(voldat)
11508 endif
11509else
11510 realdatb=rmiss
11511end if
11512
11513end function realdatb
11514
11515
11516elemental real function realdatc(voldat,var)
11517CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
11518type(vol7d_var),intent(in) :: var
11519
11520realdatc=c2r(voldat)
11521if (c_e(realdatc) .and. c_e(var%scalefactor))then
11522 realdatc=realdatc/10.**var%scalefactor
11523end if
11524
11525end function realdatc
11526
11527
11533FUNCTION realanavol(this, var) RESULT(vol)
11534TYPE(vol7d),INTENT(in) :: this
11535TYPE(vol7d_var),INTENT(in) :: var
11536REAL :: vol(SIZE(this%ana),size(this%network))
11537
11538CHARACTER(len=1) :: dtype
11539INTEGER :: indvar
11540
11541dtype = cmiss
11542indvar = index(this%anavar, var, type=dtype)
11543
11544IF (indvar > 0) THEN
11545 SELECT CASE (dtype)
11546 CASE("d")
11547 vol = realdat(this%volanad(:,indvar,:), var)
11548 CASE("r")
11549 vol = this%volanar(:,indvar,:)
11550 CASE("i")
11551 vol = realdat(this%volanai(:,indvar,:), var)
11552 CASE("b")
11553 vol = realdat(this%volanab(:,indvar,:), var)
11554 CASE("c")
11555 vol = realdat(this%volanac(:,indvar,:), var)
11556 CASE default
11557 vol = rmiss
11558 END SELECT
11559ELSE
11560 vol = rmiss
11561ENDIF
11562
11563END FUNCTION realanavol
11564
11565
11571FUNCTION integeranavol(this, var) RESULT(vol)
11572TYPE(vol7d),INTENT(in) :: this
11573TYPE(vol7d_var),INTENT(in) :: var
11574INTEGER :: vol(SIZE(this%ana),size(this%network))
11575
11576CHARACTER(len=1) :: dtype
11577INTEGER :: indvar
11578
11579dtype = cmiss
11580indvar = index(this%anavar, var, type=dtype)
11581
11582IF (indvar > 0) THEN
11583 SELECT CASE (dtype)
11584 CASE("d")
11585 vol = integerdat(this%volanad(:,indvar,:), var)
11586 CASE("r")
11587 vol = integerdat(this%volanar(:,indvar,:), var)
11588 CASE("i")
11589 vol = this%volanai(:,indvar,:)
11590 CASE("b")
11591 vol = integerdat(this%volanab(:,indvar,:), var)
11592 CASE("c")
11593 vol = integerdat(this%volanac(:,indvar,:), var)
11594 CASE default
11595 vol = imiss
11596 END SELECT
11597ELSE
11598 vol = imiss
11599ENDIF
11600
11601END FUNCTION integeranavol
11602
11603
11609subroutine move_datac (v7d,&
11610 indana,indtime,indlevel,indtimerange,indnetwork,&
11611 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
11612
11613TYPE(vol7d),intent(inout) :: v7d
11614
11615integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
11616integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
11617integer :: inddativar,inddativarattr
11618
11619
11620do inddativar=1,size(v7d%dativar%c)
11621
11622 if (c_e(v7d%voldatic(indana,indtime,indlevel,indtimerange,inddativar,indnetwork)) .and. &
11623 .not. c_e(v7d%voldatic(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
11624 ) then
11625
11626 ! dati
11627 v7d%voldatic &
11628 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
11629 v7d%voldatic &
11630 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
11631
11632
11633 ! attributi
11634 if (associated (v7d%dativarattr%i)) then
11635 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%c(inddativar))
11636 if (inddativarattr > 0 ) then
11637 v7d%voldatiattri &
11638 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
11639 v7d%voldatiattri &
11640 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
11641 end if
11642 end if
11643
11644 if (associated (v7d%dativarattr%r)) then
11645 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%c(inddativar))
11646 if (inddativarattr > 0 ) then
11647 v7d%voldatiattrr &
11648 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
11649 v7d%voldatiattrr &
11650 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
11651 end if
11652 end if
11653
11654 if (associated (v7d%dativarattr%d)) then
11655 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%c(inddativar))
11656 if (inddativarattr > 0 ) then
11657 v7d%voldatiattrd &
11658 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
11659 v7d%voldatiattrd &
11660 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
11661 end if
11662 end if
11663
11664 if (associated (v7d%dativarattr%b)) then
11665 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%c(inddativar))
11666 if (inddativarattr > 0 ) then
11667 v7d%voldatiattrb &
11668 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
11669 v7d%voldatiattrb &
11670 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
11671 end if
11672 end if
11673
11674 if (associated (v7d%dativarattr%c)) then
11675 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%c(inddativar))
11676 if (inddativarattr > 0 ) then
11677 v7d%voldatiattrc &
11678 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
11679 v7d%voldatiattrc &
11680 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
11681 end if
11682 end if
11683
11684 end if
11685
11686end do
11687
11688end subroutine move_datac
11689
11695subroutine move_datar (v7d,&
11696 indana,indtime,indlevel,indtimerange,indnetwork,&
11697 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
11698
11699TYPE(vol7d),intent(inout) :: v7d
11700
11701integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
11702integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
11703integer :: inddativar,inddativarattr
11704
11705
11706do inddativar=1,size(v7d%dativar%r)
11707
11708 if (c_e(v7d%voldatir(indana,indtime,indlevel,indtimerange,inddativar,indnetwork)) .and. &
11709 .not. c_e(v7d%voldatir(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
11710 ) then
11711
11712 ! dati
11713 v7d%voldatir &
11714 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
11715 v7d%voldatir &
11716 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
11717
11718
11719 ! attributi
11720 if (associated (v7d%dativarattr%i)) then
11721 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%r(inddativar))
11722 if (inddativarattr > 0 ) then
11723 v7d%voldatiattri &
11724 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
11725 v7d%voldatiattri &
11726 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
11727 end if
11728 end if
11729
11730 if (associated (v7d%dativarattr%r)) then
11731 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%r(inddativar))
11732 if (inddativarattr > 0 ) then
11733 v7d%voldatiattrr &
11734 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
11735 v7d%voldatiattrr &
11736 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
11737 end if
11738 end if
11739
11740 if (associated (v7d%dativarattr%d)) then
11741 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%r(inddativar))
11742 if (inddativarattr > 0 ) then
11743 v7d%voldatiattrd &
11744 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
11745 v7d%voldatiattrd &
11746 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
11747 end if
11748 end if
11749
11750 if (associated (v7d%dativarattr%b)) then
11751 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%r(inddativar))
11752 if (inddativarattr > 0 ) then
11753 v7d%voldatiattrb &
11754 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
11755 v7d%voldatiattrb &
11756 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
11757 end if
11758 end if
11759
11760 if (associated (v7d%dativarattr%c)) then
11761 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%r(inddativar))
11762 if (inddativarattr > 0 ) then
11763 v7d%voldatiattrc &
11764 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
11765 v7d%voldatiattrc &
11766 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
11767 end if
11768 end if
11769
11770 end if
11771
11772end do
11773
11774end subroutine move_datar
11775
11776
11790subroutine v7d_rounding(v7din,v7dout,level,timerange,nostatproc)
11791type(vol7d),intent(inout) :: v7din
11792type(vol7d),intent(out) :: v7dout
11793type(vol7d_level),intent(in),optional :: level(:)
11794type(vol7d_timerange),intent(in),optional :: timerange(:)
11795!logical,intent(in),optional :: merge !< if there are data on more then one almost equal levels or timeranges
11796!! will be merged POINT BY POINT with priority for the fird data found ordered by icreasing var index
11797logical,intent(in),optional :: nostatproc
11798
11799integer :: nana,nlevel,ntime,ntimerange,nnetwork,nbin
11800integer :: iana,ilevel,itimerange,indl,indt,itime,inetwork
11801type(vol7d_level) :: roundlevel(size(v7din%level))
11802type(vol7d_timerange) :: roundtimerange(size(v7din%timerange))
11803type(vol7d) :: v7d_tmp
11804
11805
11806nbin=0
11807
11808if (associated(v7din%dativar%r)) nbin = nbin + size(v7din%dativar%r)
11809if (associated(v7din%dativar%i)) nbin = nbin + size(v7din%dativar%i)
11810if (associated(v7din%dativar%d)) nbin = nbin + size(v7din%dativar%d)
11811if (associated(v7din%dativar%b)) nbin = nbin + size(v7din%dativar%b)
11812
11813call init(v7d_tmp)
11814
11815roundlevel=v7din%level
11816
11817if (present(level))then
11818 do ilevel = 1, size(v7din%level)
11819 if ((any(v7din%level(ilevel) .almosteq. level))) then
11820 roundlevel(ilevel)=level(1)
11821 end if
11822 end do
11823end if
11824
11825roundtimerange=v7din%timerange
11826
11827if (present(timerange))then
11828 do itimerange = 1, size(v7din%timerange)
11829 if ((any(v7din%timerange(itimerange) .almosteq. timerange))) then
11830 roundtimerange(itimerange)=timerange(1)
11831 end if
11832 end do
11833end if
11834
11835!set istantaneous values everywere
11836!preserve p1 for forecast time
11837if (optio_log(nostatproc)) then
11838 roundtimerange(:)%timerange=254
11839 roundtimerange(:)%p2=0
11840end if
11841
11842
11843nana=size(v7din%ana)
11844nlevel=count_distinct(roundlevel,back=.true.)
11845ntime=size(v7din%time)
11846ntimerange=count_distinct(roundtimerange,back=.true.)
11847nnetwork=size(v7din%network)
11848
11849call init(v7d_tmp)
11850
11851if (nbin == 0) then
11852 call copy(v7din,v7d_tmp)
11853else
11854 call vol7d_convr(v7din,v7d_tmp)
11855end if
11856
11857v7d_tmp%level=roundlevel
11858v7d_tmp%timerange=roundtimerange
11859
11860do ilevel=1, size(v7d_tmp%level)
11861 indl=index(v7d_tmp%level,roundlevel(ilevel))
11862 do itimerange=1,size(v7d_tmp%timerange)
11863 indt=index(v7d_tmp%timerange,roundtimerange(itimerange))
11864
11865 if (indl /= ilevel .or. indt /= itimerange) then
11866
11867 do iana=1, nana
11868 do itime=1,ntime
11869 do inetwork=1,nnetwork
11870
11871 if (nbin > 0) then
11872 call move_datar (v7d_tmp,&
11873 iana,itime,ilevel,itimerange,inetwork,&
11874 iana,itime,indl,indt,inetwork)
11875 else
11876 call move_datac (v7d_tmp,&
11877 iana,itime,ilevel,itimerange,inetwork,&
11878 iana,itime,indl,indt,inetwork)
11879 end if
11880
11881 end do
11882 end do
11883 end do
11884
11885 end if
11886
11887 end do
11888end do
11889
11890! set to missing level and time > nlevel
11891do ilevel=nlevel+1,size(v7d_tmp%level)
11892 call init (v7d_tmp%level(ilevel))
11893end do
11894
11895do itimerange=ntimerange+1,size(v7d_tmp%timerange)
11896 call init (v7d_tmp%timerange(itimerange))
11897end do
11898
11899!copy with remove
11900CALL copy(v7d_tmp,v7dout,miss=.true.,lsort_timerange=.true.,lsort_level=.true.)
11901CALL delete(v7d_tmp)
11902
11903!call display(v7dout)
11904
11905end subroutine v7d_rounding
11906
11907
11908END MODULE vol7d_class
11909
11915
11916
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.