libsim Versione 7.2.1

◆ vol7d_get_voldatic()

subroutine vol7d_get_voldatic ( type(vol7d), intent(in) this,
integer, dimension(:), intent(in) dimlist,
character(len=vol7d_cdatalen), dimension(:), optional, pointer vol1dp,
character(len=vol7d_cdatalen), dimension(:,:), optional, pointer vol2dp,
character(len=vol7d_cdatalen), dimension(:,:,:), optional, pointer vol3dp,
character(len=vol7d_cdatalen), dimension(:,:,:,:), optional, pointer vol4dp,
character(len=vol7d_cdatalen), dimension(:,:,:,:,:), optional, pointer vol5dp,
character(len=vol7d_cdatalen), dimension(:,:,:,:,:,:), optional, pointer vol6dp )

Crea una vista a dimensione ridotta di un volume di dati di tipo CHARACTER(len=vol7d_cdatalen).

È necessario fornire uno solo dei parametri opzionali vol*dp corrispondente al numero di dimensioni richieste. L'ordine delle dimensioni nella vista è quello prefissato in ::vol7d indipendentemente dall'ordine delle dimensioni fornito in dimlist. In caso di fallimento, in particolare se dimlist non contiene tutte le dimensioni non degeneri del volume richiesto oppure se una delle dimensioni è =0, il puntatore vol*dp è restituito in uno stato disassociato, per cui è opportuno controllare sempre in uscita, lo stato del puntatore per evitare che il programma abortisca con un errore di sistema, ad esempio:

CHARACTER(len=vol7d_cdatalen), POINTER :: vol2d(:,:)
...
CALL vol7d_get_voldatic(v7d1, (/vol7d_ana_d, vol7d_time_d/), vol2d)
IF (ASSOCIATED(vol2d)) THEN
print*,vol2d
...
ENDIF
return
Parametri
[in]thisoggetto di cui creare la vista
[in]dimlistlista delle dimensioni da includere nella vista, attenzione tutte le dimensioni non degeneri (cioè con estensione >1) devono essere incluse nella lista; utilizzare le costanti vol7d_ana_d ... vol7d_attr_d, ecc.
vol1dparray che in uscita conterrà la vista 1d
vol2dparray che in uscita conterrà la vista 2d
vol3dparray che in uscita conterrà la vista 3d
vol4dparray che in uscita conterrà la vista 4d
vol5dparray che in uscita conterrà la vista 5d
vol6dparray che in uscita conterrà la vista 6d

Definizione alla linea 6349 del file vol7d_class.F90.

6351! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
6352! authors:
6353! Davide Cesari <dcesari@arpa.emr.it>
6354! Paolo Patruno <ppatruno@arpa.emr.it>
6355
6356! This program is free software; you can redistribute it and/or
6357! modify it under the terms of the GNU General Public License as
6358! published by the Free Software Foundation; either version 2 of
6359! the License, or (at your option) any later version.
6360
6361! This program is distributed in the hope that it will be useful,
6362! but WITHOUT ANY WARRANTY; without even the implied warranty of
6363! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
6364! GNU General Public License for more details.
6365
6366! You should have received a copy of the GNU General Public License
6367! along with this program. If not, see <http://www.gnu.org/licenses/>.
6368#include "config.h"
6369
6381
6435MODULE vol7d_class
6436USE kinds
6440USE log4fortran
6441USE err_handling
6442USE io_units
6449IMPLICIT NONE
6450
6451
6452INTEGER, PARAMETER :: vol7d_maxdim_a = 3, vol7d_maxdim_aa = 4, &
6453 vol7d_maxdim_d = 6, vol7d_maxdim_ad = 7
6454
6455INTEGER, PARAMETER :: vol7d_ana_a=1
6456INTEGER, PARAMETER :: vol7d_var_a=2
6457INTEGER, PARAMETER :: vol7d_network_a=3
6458INTEGER, PARAMETER :: vol7d_attr_a=4
6459INTEGER, PARAMETER :: vol7d_ana_d=1
6460INTEGER, PARAMETER :: vol7d_time_d=2
6461INTEGER, PARAMETER :: vol7d_level_d=3
6462INTEGER, PARAMETER :: vol7d_timerange_d=4
6463INTEGER, PARAMETER :: vol7d_var_d=5
6464INTEGER, PARAMETER :: vol7d_network_d=6
6465INTEGER, PARAMETER :: vol7d_attr_d=7
6466INTEGER, PARAMETER :: vol7d_cdatalen=32
6467
6468TYPE vol7d_varmap
6469 INTEGER :: r, d, i, b, c
6470END TYPE vol7d_varmap
6471
6474TYPE vol7d
6476 TYPE(vol7d_ana),POINTER :: ana(:)
6478 TYPE(datetime),POINTER :: time(:)
6480 TYPE(vol7d_level),POINTER :: level(:)
6482 TYPE(vol7d_timerange),POINTER :: timerange(:)
6484 TYPE(vol7d_network),POINTER :: network(:)
6486 TYPE(vol7d_varvect) :: anavar
6488 TYPE(vol7d_varvect) :: anaattr
6490 TYPE(vol7d_varvect) :: anavarattr
6492 TYPE(vol7d_varvect) :: dativar
6494 TYPE(vol7d_varvect) :: datiattr
6496 TYPE(vol7d_varvect) :: dativarattr
6497
6499 REAL,POINTER :: volanar(:,:,:)
6501 DOUBLE PRECISION,POINTER :: volanad(:,:,:)
6503 INTEGER,POINTER :: volanai(:,:,:)
6505 INTEGER(kind=int_b),POINTER :: volanab(:,:,:)
6507 CHARACTER(len=vol7d_cdatalen),POINTER :: volanac(:,:,:)
6508
6510 REAL,POINTER :: volanaattrr(:,:,:,:)
6512 DOUBLE PRECISION,POINTER :: volanaattrd(:,:,:,:)
6514 INTEGER,POINTER :: volanaattri(:,:,:,:)
6516 INTEGER(kind=int_b),POINTER :: volanaattrb(:,:,:,:)
6518 CHARACTER(len=vol7d_cdatalen),POINTER :: volanaattrc(:,:,:,:)
6519
6521 REAL,POINTER :: voldatir(:,:,:,:,:,:) ! sono i dati
6523 DOUBLE PRECISION,POINTER :: voldatid(:,:,:,:,:,:)
6525 INTEGER,POINTER :: voldatii(:,:,:,:,:,:)
6527 INTEGER(kind=int_b),POINTER :: voldatib(:,:,:,:,:,:)
6529 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatic(:,:,:,:,:,:)
6530
6532 REAL,POINTER :: voldatiattrr(:,:,:,:,:,:,:)
6534 DOUBLE PRECISION,POINTER :: voldatiattrd(:,:,:,:,:,:,:)
6536 INTEGER,POINTER :: voldatiattri(:,:,:,:,:,:,:)
6538 INTEGER(kind=int_b),POINTER :: voldatiattrb(:,:,:,:,:,:,:)
6540 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatiattrc(:,:,:,:,:,:,:)
6541
6543 integer :: time_definition
6544
6545END TYPE vol7d
6546
6550INTERFACE init
6551 MODULE PROCEDURE vol7d_init
6552END INTERFACE
6553
6555INTERFACE delete
6556 MODULE PROCEDURE vol7d_delete
6557END INTERFACE
6558
6560INTERFACE export
6561 MODULE PROCEDURE vol7d_write_on_file
6562END INTERFACE
6563
6565INTERFACE import
6566 MODULE PROCEDURE vol7d_read_from_file
6567END INTERFACE
6568
6570INTERFACE display
6571 MODULE PROCEDURE vol7d_display, dat_display, dat_vect_display
6572END INTERFACE
6573
6575INTERFACE to_char
6576 MODULE PROCEDURE to_char_dat
6577END INTERFACE
6578
6580INTERFACE doubledat
6581 MODULE PROCEDURE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
6582END INTERFACE
6583
6585INTERFACE realdat
6586 MODULE PROCEDURE realdatd,realdatr,realdati,realdatb,realdatc
6587END INTERFACE
6588
6590INTERFACE integerdat
6591 MODULE PROCEDURE integerdatd,integerdatr,integerdati,integerdatb,integerdatc
6592END INTERFACE
6593
6595INTERFACE copy
6596 MODULE PROCEDURE vol7d_copy
6597END INTERFACE
6598
6600INTERFACE c_e
6601 MODULE PROCEDURE vol7d_c_e
6602END INTERFACE
6603
6607INTERFACE check
6608 MODULE PROCEDURE vol7d_check
6609END INTERFACE
6610
6624INTERFACE rounding
6625 MODULE PROCEDURE v7d_rounding
6626END INTERFACE
6627
6628!!$INTERFACE get_volana
6629!!$ MODULE PROCEDURE vol7d_get_volanar, vol7d_get_volanad, vol7d_get_volanai, &
6630!!$ vol7d_get_volanab, vol7d_get_volanac
6631!!$END INTERFACE
6632!!$
6633!!$INTERFACE get_voldati
6634!!$ MODULE PROCEDURE vol7d_get_voldatir, vol7d_get_voldatid, vol7d_get_voldatii, &
6635!!$ vol7d_get_voldatib, vol7d_get_voldatic
6636!!$END INTERFACE
6637!!$
6638!!$INTERFACE get_volanaattr
6639!!$ MODULE PROCEDURE vol7d_get_volanaattrr, vol7d_get_volanaattrd, &
6640!!$ vol7d_get_volanaattri, vol7d_get_volanaattrb, vol7d_get_volanaattrc
6641!!$END INTERFACE
6642!!$
6643!!$INTERFACE get_voldatiattr
6644!!$ MODULE PROCEDURE vol7d_get_voldatiattrr, vol7d_get_voldatiattrd, &
6645!!$ vol7d_get_voldatiattri, vol7d_get_voldatiattrb, vol7d_get_voldatiattrc
6646!!$END INTERFACE
6647
6648PRIVATE vol7d_get_volr, vol7d_get_vold, vol7d_get_voli, vol7d_get_volb, &
6649 vol7d_get_volc, &
6650 volptr1dr, volptr2dr, volptr3dr, volptr4dr, volptr5dr, volptr6dr, volptr7dr, &
6651 volptr1dd, volptr2dd, volptr3dd, volptr4dd, volptr5dd, volptr6dd, volptr7dd, &
6652 volptr1di, volptr2di, volptr3di, volptr4di, volptr5di, volptr6di, volptr7di, &
6653 volptr1db, volptr2db, volptr3db, volptr4db, volptr5db, volptr6db, volptr7db, &
6654 volptr1dc, volptr2dc, volptr3dc, volptr4dc, volptr5dc, volptr6dc, volptr7dc, &
6655 vol7d_nullifyr, vol7d_nullifyd, vol7d_nullifyi, vol7d_nullifyb, vol7d_nullifyc, &
6656 vol7d_init, vol7d_delete, vol7d_write_on_file, vol7d_read_from_file, &
6657 vol7d_check_alloc_ana, vol7d_force_alloc_ana, &
6658 vol7d_check_alloc_dati, vol7d_force_alloc_dati, vol7d_force_alloc, &
6659 vol7d_display, dat_display, dat_vect_display, &
6660 to_char_dat, vol7d_check
6661
6662PRIVATE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
6663
6664PRIVATE vol7d_c_e
6665
6666CONTAINS
6667
6668
6673SUBROUTINE vol7d_init(this,time_definition)
6674TYPE(vol7d),intent(out) :: this
6675integer,INTENT(IN),OPTIONAL :: time_definition
6676
6677CALL init(this%anavar)
6678CALL init(this%anaattr)
6679CALL init(this%anavarattr)
6680CALL init(this%dativar)
6681CALL init(this%datiattr)
6682CALL init(this%dativarattr)
6683CALL vol7d_var_features_init() ! initialise var features table once
6684
6685NULLIFY(this%ana, this%time, this%level, this%timerange, this%network)
6686
6687NULLIFY(this%volanar, this%volanaattrr, this%voldatir, this%voldatiattrr)
6688NULLIFY(this%volanad, this%volanaattrd, this%voldatid, this%voldatiattrd)
6689NULLIFY(this%volanai, this%volanaattri, this%voldatii, this%voldatiattri)
6690NULLIFY(this%volanab, this%volanaattrb, this%voldatib, this%voldatiattrb)
6691NULLIFY(this%volanac, this%volanaattrc, this%voldatic, this%voldatiattrc)
6692
6693if(present(time_definition)) then
6694 this%time_definition=time_definition
6695else
6696 this%time_definition=1 !default to validity time
6697end if
6698
6699END SUBROUTINE vol7d_init
6700
6701
6705ELEMENTAL SUBROUTINE vol7d_delete(this, dataonly)
6706TYPE(vol7d),intent(inout) :: this
6707LOGICAL, INTENT(in), OPTIONAL :: dataonly
6708
6709
6710IF (.NOT. optio_log(dataonly)) THEN
6711 IF (ASSOCIATED(this%volanar)) DEALLOCATE(this%volanar)
6712 IF (ASSOCIATED(this%volanad)) DEALLOCATE(this%volanad)
6713 IF (ASSOCIATED(this%volanai)) DEALLOCATE(this%volanai)
6714 IF (ASSOCIATED(this%volanab)) DEALLOCATE(this%volanab)
6715 IF (ASSOCIATED(this%volanac)) DEALLOCATE(this%volanac)
6716 IF (ASSOCIATED(this%volanaattrr)) DEALLOCATE(this%volanaattrr)
6717 IF (ASSOCIATED(this%volanaattrd)) DEALLOCATE(this%volanaattrd)
6718 IF (ASSOCIATED(this%volanaattri)) DEALLOCATE(this%volanaattri)
6719 IF (ASSOCIATED(this%volanaattrb)) DEALLOCATE(this%volanaattrb)
6720 IF (ASSOCIATED(this%volanaattrc)) DEALLOCATE(this%volanaattrc)
6721ENDIF
6722IF (ASSOCIATED(this%voldatir)) DEALLOCATE(this%voldatir)
6723IF (ASSOCIATED(this%voldatid)) DEALLOCATE(this%voldatid)
6724IF (ASSOCIATED(this%voldatii)) DEALLOCATE(this%voldatii)
6725IF (ASSOCIATED(this%voldatib)) DEALLOCATE(this%voldatib)
6726IF (ASSOCIATED(this%voldatic)) DEALLOCATE(this%voldatic)
6727IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
6728IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
6729IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
6730IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
6731IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
6732
6733IF (.NOT. optio_log(dataonly)) THEN
6734 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
6735 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
6736ENDIF
6737IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
6738IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
6739IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
6740
6741IF (.NOT. optio_log(dataonly)) THEN
6742 CALL delete(this%anavar)
6743 CALL delete(this%anaattr)
6744 CALL delete(this%anavarattr)
6745ENDIF
6746CALL delete(this%dativar)
6747CALL delete(this%datiattr)
6748CALL delete(this%dativarattr)
6749
6750END SUBROUTINE vol7d_delete
6751
6752
6753
6754integer function vol7d_check(this)
6755TYPE(vol7d),intent(in) :: this
6756integer :: i,j,k,l,m,n
6757
6758vol7d_check=0
6759
6760if (associated(this%voldatii)) then
6761do i = 1,size(this%voldatii,1)
6762 do j = 1,size(this%voldatii,2)
6763 do k = 1,size(this%voldatii,3)
6764 do l = 1,size(this%voldatii,4)
6765 do m = 1,size(this%voldatii,5)
6766 do n = 1,size(this%voldatii,6)
6767 if (this%voldatii(i,j,k,l,m,n) /= this%voldatii(i,j,k,l,m,n) ) then
6768 CALL l4f_log(l4f_warn,"check: abnormal value at voldatii("&
6769 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
6770 vol7d_check=1
6771 end if
6772 end do
6773 end do
6774 end do
6775 end do
6776 end do
6777end do
6778end if
6779
6780
6781if (associated(this%voldatir)) then
6782do i = 1,size(this%voldatir,1)
6783 do j = 1,size(this%voldatir,2)
6784 do k = 1,size(this%voldatir,3)
6785 do l = 1,size(this%voldatir,4)
6786 do m = 1,size(this%voldatir,5)
6787 do n = 1,size(this%voldatir,6)
6788 if (this%voldatir(i,j,k,l,m,n) /= this%voldatir(i,j,k,l,m,n) ) then
6789 CALL l4f_log(l4f_warn,"check: abnormal value at voldatir("&
6790 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
6791 vol7d_check=2
6792 end if
6793 end do
6794 end do
6795 end do
6796 end do
6797 end do
6798end do
6799end if
6800
6801if (associated(this%voldatid)) then
6802do i = 1,size(this%voldatid,1)
6803 do j = 1,size(this%voldatid,2)
6804 do k = 1,size(this%voldatid,3)
6805 do l = 1,size(this%voldatid,4)
6806 do m = 1,size(this%voldatid,5)
6807 do n = 1,size(this%voldatid,6)
6808 if (this%voldatid(i,j,k,l,m,n) /= this%voldatid(i,j,k,l,m,n) ) then
6809 CALL l4f_log(l4f_warn,"check: abnormal value at voldatid("&
6810 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
6811 vol7d_check=3
6812 end if
6813 end do
6814 end do
6815 end do
6816 end do
6817 end do
6818end do
6819end if
6820
6821if (associated(this%voldatib)) then
6822do i = 1,size(this%voldatib,1)
6823 do j = 1,size(this%voldatib,2)
6824 do k = 1,size(this%voldatib,3)
6825 do l = 1,size(this%voldatib,4)
6826 do m = 1,size(this%voldatib,5)
6827 do n = 1,size(this%voldatib,6)
6828 if (this%voldatib(i,j,k,l,m,n) /= this%voldatib(i,j,k,l,m,n) ) then
6829 CALL l4f_log(l4f_warn,"check: abnormal value at voldatib("&
6830 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
6831 vol7d_check=4
6832 end if
6833 end do
6834 end do
6835 end do
6836 end do
6837 end do
6838end do
6839end if
6840
6841end function vol7d_check
6842
6843
6844
6845!TODO da completare ! aborta se i volumi sono allocati a dimensione 0
6847SUBROUTINE vol7d_display(this)
6848TYPE(vol7d),intent(in) :: this
6849integer :: i
6850
6851REAL :: rdat
6852DOUBLE PRECISION :: ddat
6853INTEGER :: idat
6854INTEGER(kind=int_b) :: bdat
6855CHARACTER(len=vol7d_cdatalen) :: cdat
6856
6857
6858print*,"<<<<<<<<<<<<<<<<<<< vol7d object >>>>>>>>>>>>>>>>>>>>"
6859if (this%time_definition == 0) then
6860 print*,"TIME DEFINITION: time is reference time"
6861else if (this%time_definition == 1) then
6862 print*,"TIME DEFINITION: time is validity time"
6863else
6864 print*,"Time definition have a wrong walue:", this%time_definition
6865end if
6866
6867IF (ASSOCIATED(this%network))then
6868 print*,"---- network vector ----"
6869 print*,"elements=",size(this%network)
6870 do i=1, size(this%network)
6871 call display(this%network(i))
6872 end do
6873end IF
6874
6875IF (ASSOCIATED(this%ana))then
6876 print*,"---- ana vector ----"
6877 print*,"elements=",size(this%ana)
6878 do i=1, size(this%ana)
6879 call display(this%ana(i))
6880 end do
6881end IF
6882
6883IF (ASSOCIATED(this%time))then
6884 print*,"---- time vector ----"
6885 print*,"elements=",size(this%time)
6886 do i=1, size(this%time)
6887 call display(this%time(i))
6888 end do
6889end if
6890
6891IF (ASSOCIATED(this%level)) then
6892 print*,"---- level vector ----"
6893 print*,"elements=",size(this%level)
6894 do i =1,size(this%level)
6895 call display(this%level(i))
6896 end do
6897end if
6898
6899IF (ASSOCIATED(this%timerange))then
6900 print*,"---- timerange vector ----"
6901 print*,"elements=",size(this%timerange)
6902 do i =1,size(this%timerange)
6903 call display(this%timerange(i))
6904 end do
6905end if
6906
6907
6908print*,"---- ana vector ----"
6909print*,""
6910print*,"->>>>>>>>> anavar -"
6911call display(this%anavar)
6912print*,""
6913print*,"->>>>>>>>> anaattr -"
6914call display(this%anaattr)
6915print*,""
6916print*,"->>>>>>>>> anavarattr -"
6917call display(this%anavarattr)
6918
6919print*,"-- ana data section (first point) --"
6920
6921idat=imiss
6922rdat=rmiss
6923ddat=dmiss
6924bdat=ibmiss
6925cdat=cmiss
6926
6927!ntime = MIN(SIZE(this%time),nprint)
6928!ntimerange = MIN(SIZE(this%timerange),nprint)
6929!nlevel = MIN(SIZE(this%level),nprint)
6930!nnetwork = MIN(SIZE(this%network),nprint)
6931!nana = MIN(SIZE(this%ana),nprint)
6932
6933IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0) THEN
6934if (associated(this%volanai)) then
6935 do i=1,size(this%anavar%i)
6936 idat=this%volanai(1,i,1)
6937 if (associated(this%anavar%i)) call display(this%anavar%i(i),idat,rdat,ddat,bdat,cdat)
6938 end do
6939end if
6940idat=imiss
6941
6942if (associated(this%volanar)) then
6943 do i=1,size(this%anavar%r)
6944 rdat=this%volanar(1,i,1)
6945 if (associated(this%anavar%r)) call display(this%anavar%r(i),idat,rdat,ddat,bdat,cdat)
6946 end do
6947end if
6948rdat=rmiss
6949
6950if (associated(this%volanad)) then
6951 do i=1,size(this%anavar%d)
6952 ddat=this%volanad(1,i,1)
6953 if (associated(this%anavar%d)) call display(this%anavar%d(i),idat,rdat,ddat,bdat,cdat)
6954 end do
6955end if
6956ddat=dmiss
6957
6958if (associated(this%volanab)) then
6959 do i=1,size(this%anavar%b)
6960 bdat=this%volanab(1,i,1)
6961 if (associated(this%anavar%b)) call display(this%anavar%b(i),idat,rdat,ddat,bdat,cdat)
6962 end do
6963end if
6964bdat=ibmiss
6965
6966if (associated(this%volanac)) then
6967 do i=1,size(this%anavar%c)
6968 cdat=this%volanac(1,i,1)
6969 if (associated(this%anavar%c)) call display(this%anavar%c(i),idat,rdat,ddat,bdat,cdat)
6970 end do
6971end if
6972cdat=cmiss
6973ENDIF
6974
6975print*,"---- data vector ----"
6976print*,""
6977print*,"->>>>>>>>> dativar -"
6978call display(this%dativar)
6979print*,""
6980print*,"->>>>>>>>> datiattr -"
6981call display(this%datiattr)
6982print*,""
6983print*,"->>>>>>>>> dativarattr -"
6984call display(this%dativarattr)
6985
6986print*,"-- data data section (first point) --"
6987
6988idat=imiss
6989rdat=rmiss
6990ddat=dmiss
6991bdat=ibmiss
6992cdat=cmiss
6993
6994IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0 .AND. size(this%time) > 0 &
6995 .AND. size(this%level) > 0 .AND. size(this%timerange) > 0) THEN
6996if (associated(this%voldatii)) then
6997 do i=1,size(this%dativar%i)
6998 idat=this%voldatii(1,1,1,1,i,1)
6999 if (associated(this%dativar%i)) call display(this%dativar%i(i),idat,rdat,ddat,bdat,cdat)
7000 end do
7001end if
7002idat=imiss
7003
7004if (associated(this%voldatir)) then
7005 do i=1,size(this%dativar%r)
7006 rdat=this%voldatir(1,1,1,1,i,1)
7007 if (associated(this%dativar%r)) call display(this%dativar%r(i),idat,rdat,ddat,bdat,cdat)
7008 end do
7009end if
7010rdat=rmiss
7011
7012if (associated(this%voldatid)) then
7013 do i=1,size(this%dativar%d)
7014 ddat=this%voldatid(1,1,1,1,i,1)
7015 if (associated(this%dativar%d)) call display(this%dativar%d(i),idat,rdat,ddat,bdat,cdat)
7016 end do
7017end if
7018ddat=dmiss
7019
7020if (associated(this%voldatib)) then
7021 do i=1,size(this%dativar%b)
7022 bdat=this%voldatib(1,1,1,1,i,1)
7023 if (associated(this%dativar%b)) call display(this%dativar%b(i),idat,rdat,ddat,bdat,cdat)
7024 end do
7025end if
7026bdat=ibmiss
7027
7028if (associated(this%voldatic)) then
7029 do i=1,size(this%dativar%c)
7030 cdat=this%voldatic(1,1,1,1,i,1)
7031 if (associated(this%dativar%c)) call display(this%dativar%c(i),idat,rdat,ddat,bdat,cdat)
7032 end do
7033end if
7034cdat=cmiss
7035ENDIF
7036
7037print*,"<<<<<<<<<<<<<<<<<<< END vol7d object >>>>>>>>>>>>>>>>>>>>"
7038
7039END SUBROUTINE vol7d_display
7040
7041
7043SUBROUTINE dat_display(this,idat,rdat,ddat,bdat,cdat)
7044TYPE(vol7d_var),intent(in) :: this
7046REAL :: rdat
7048DOUBLE PRECISION :: ddat
7050INTEGER :: idat
7052INTEGER(kind=int_b) :: bdat
7054CHARACTER(len=*) :: cdat
7055
7056print *, to_char_dat(this,idat,rdat,ddat,bdat,cdat)
7057
7058end SUBROUTINE dat_display
7059
7061SUBROUTINE dat_vect_display(this,idat,rdat,ddat,bdat,cdat)
7062
7063TYPE(vol7d_var),intent(in) :: this(:)
7065REAL :: rdat(:)
7067DOUBLE PRECISION :: ddat(:)
7069INTEGER :: idat(:)
7071INTEGER(kind=int_b) :: bdat(:)
7073CHARACTER(len=*):: cdat(:)
7074
7075integer :: i
7076
7077do i =1,size(this)
7078 call display(this(i),idat(i),rdat(i),ddat(i),bdat(i),cdat(i))
7079end do
7080
7081end SUBROUTINE dat_vect_display
7082
7083
7084FUNCTION to_char_dat(this,idat,rdat,ddat,bdat,cdat)
7085#ifdef HAVE_DBALLE
7086USE dballef
7087#endif
7088TYPE(vol7d_var),INTENT(in) :: this
7090REAL :: rdat
7092DOUBLE PRECISION :: ddat
7094INTEGER :: idat
7096INTEGER(kind=int_b) :: bdat
7098CHARACTER(len=*) :: cdat
7099CHARACTER(len=80) :: to_char_dat
7100
7101CHARACTER(len=LEN(to_char_dat)) :: to_char_tmp
7102
7103
7104#ifdef HAVE_DBALLE
7105INTEGER :: handle, ier
7106
7107handle = 0
7108to_char_dat="VALUE: "
7109
7110if (c_e(idat)) to_char_dat=trim(to_char_dat)//" ;int> "//trim(to_char(idat))
7111if (c_e(rdat)) to_char_dat=trim(to_char_dat)//" ;real> "//trim(to_char(rdat))
7112if (c_e(ddat)) to_char_dat=trim(to_char_dat)//" ;double> "//trim(to_char(ddat))
7113if (c_e(bdat)) to_char_dat=trim(to_char_dat)//" ;byte> "//trim(to_char(bdat))
7114
7115if ( c_e(cdat))then
7116 ier = idba_messaggi(handle,"/dev/null", "w", "BUFR")
7117 ier = idba_spiegab(handle,this%btable,cdat,to_char_tmp)
7118 ier = idba_fatto(handle)
7119 to_char_dat=trim(to_char_dat)//" ;char> "//trim(to_char_tmp)
7120endif
7121
7122#else
7123
7124to_char_dat="VALUE: "
7125if (c_e(idat)) to_char_dat=trim(to_char_dat)//" ;int> "//trim(to_char(idat))
7126if (c_e(rdat)) to_char_dat=trim(to_char_dat)//" ;real> "//trim(to_char(rdat))
7127if (c_e(ddat)) to_char_dat=trim(to_char_dat)//" ;double> "//trim(to_char(ddat))
7128if (c_e(bdat)) to_char_dat=trim(to_char_dat)//" ;byte> "//trim(to_char(bdat))
7129if (c_e(cdat)) to_char_dat=trim(to_char_dat)//" ;char> "//trim(cdat)
7130
7131#endif
7132
7133END FUNCTION to_char_dat
7134
7135
7138FUNCTION vol7d_c_e(this) RESULT(c_e)
7139TYPE(vol7d), INTENT(in) :: this
7140
7141LOGICAL :: c_e
7142
7143c_e = ASSOCIATED(this%ana) .OR. ASSOCIATED(this%time) .OR. &
7144 ASSOCIATED(this%level) .OR. ASSOCIATED(this%timerange) .OR. &
7145 ASSOCIATED(this%network) .OR. &
7146 ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
7147 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
7148 ASSOCIATED(this%anavar%c) .OR. &
7149 ASSOCIATED(this%anaattr%r) .OR. ASSOCIATED(this%anaattr%d) .OR. &
7150 ASSOCIATED(this%anaattr%i) .OR. ASSOCIATED(this%anaattr%b) .OR. &
7151 ASSOCIATED(this%anaattr%c) .OR. &
7152 ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
7153 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
7154 ASSOCIATED(this%dativar%c) .OR. &
7155 ASSOCIATED(this%datiattr%r) .OR. ASSOCIATED(this%datiattr%d) .OR. &
7156 ASSOCIATED(this%datiattr%i) .OR. ASSOCIATED(this%datiattr%b) .OR. &
7157 ASSOCIATED(this%datiattr%c)
7158
7159END FUNCTION vol7d_c_e
7160
7161
7200SUBROUTINE vol7d_alloc(this, nana, ntime, nlevel, ntimerange, nnetwork, &
7201 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
7202 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
7203 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
7204 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
7205 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
7206 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc, &
7207 ini)
7208TYPE(vol7d),INTENT(inout) :: this
7209INTEGER,INTENT(in),OPTIONAL :: nana
7210INTEGER,INTENT(in),OPTIONAL :: ntime
7211INTEGER,INTENT(in),OPTIONAL :: nlevel
7212INTEGER,INTENT(in),OPTIONAL :: ntimerange
7213INTEGER,INTENT(in),OPTIONAL :: nnetwork
7215INTEGER,INTENT(in),OPTIONAL :: &
7216 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
7217 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
7218 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
7219 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
7220 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
7221 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc
7222LOGICAL,INTENT(in),OPTIONAL :: ini
7223
7224INTEGER :: i
7225LOGICAL :: linit
7226
7227IF (PRESENT(ini)) THEN
7228 linit = ini
7229ELSE
7230 linit = .false.
7231ENDIF
7232
7233! Dimensioni principali
7234IF (PRESENT(nana)) THEN
7235 IF (nana >= 0) THEN
7236 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
7237 ALLOCATE(this%ana(nana))
7238 IF (linit) THEN
7239 DO i = 1, nana
7240 CALL init(this%ana(i))
7241 ENDDO
7242 ENDIF
7243 ENDIF
7244ENDIF
7245IF (PRESENT(ntime)) THEN
7246 IF (ntime >= 0) THEN
7247 IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
7248 ALLOCATE(this%time(ntime))
7249 IF (linit) THEN
7250 DO i = 1, ntime
7251 CALL init(this%time(i))
7252 ENDDO
7253 ENDIF
7254 ENDIF
7255ENDIF
7256IF (PRESENT(nlevel)) THEN
7257 IF (nlevel >= 0) THEN
7258 IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
7259 ALLOCATE(this%level(nlevel))
7260 IF (linit) THEN
7261 DO i = 1, nlevel
7262 CALL init(this%level(i))
7263 ENDDO
7264 ENDIF
7265 ENDIF
7266ENDIF
7267IF (PRESENT(ntimerange)) THEN
7268 IF (ntimerange >= 0) THEN
7269 IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
7270 ALLOCATE(this%timerange(ntimerange))
7271 IF (linit) THEN
7272 DO i = 1, ntimerange
7273 CALL init(this%timerange(i))
7274 ENDDO
7275 ENDIF
7276 ENDIF
7277ENDIF
7278IF (PRESENT(nnetwork)) THEN
7279 IF (nnetwork >= 0) THEN
7280 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
7281 ALLOCATE(this%network(nnetwork))
7282 IF (linit) THEN
7283 DO i = 1, nnetwork
7284 CALL init(this%network(i))
7285 ENDDO
7286 ENDIF
7287 ENDIF
7288ENDIF
7289! Dimensioni dei tipi delle variabili
7290CALL vol7d_varvect_alloc(this%anavar, nanavarr, nanavard, &
7291 nanavari, nanavarb, nanavarc, ini)
7292CALL vol7d_varvect_alloc(this%anaattr, nanaattrr, nanaattrd, &
7293 nanaattri, nanaattrb, nanaattrc, ini)
7294CALL vol7d_varvect_alloc(this%anavarattr, nanavarattrr, nanavarattrd, &
7295 nanavarattri, nanavarattrb, nanavarattrc, ini)
7296CALL vol7d_varvect_alloc(this%dativar, ndativarr, ndativard, &
7297 ndativari, ndativarb, ndativarc, ini)
7298CALL vol7d_varvect_alloc(this%datiattr, ndatiattrr, ndatiattrd, &
7299 ndatiattri, ndatiattrb, ndatiattrc, ini)
7300CALL vol7d_varvect_alloc(this%dativarattr, ndativarattrr, ndativarattrd, &
7301 ndativarattri, ndativarattrb, ndativarattrc, ini)
7302
7303END SUBROUTINE vol7d_alloc
7304
7305
7306FUNCTION vol7d_check_alloc_ana(this)
7307TYPE(vol7d),INTENT(in) :: this
7308LOGICAL :: vol7d_check_alloc_ana
7309
7310vol7d_check_alloc_ana = ASSOCIATED(this%ana) .AND. ASSOCIATED(this%network)
7311
7312END FUNCTION vol7d_check_alloc_ana
7313
7314SUBROUTINE vol7d_force_alloc_ana(this, ini)
7315TYPE(vol7d),INTENT(inout) :: this
7316LOGICAL,INTENT(in),OPTIONAL :: ini
7317
7318! Alloco i descrittori minimi per avere un volume di anagrafica
7319IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=1, ini=ini)
7320IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=1, ini=ini)
7321
7322END SUBROUTINE vol7d_force_alloc_ana
7323
7324
7325FUNCTION vol7d_check_alloc_dati(this)
7326TYPE(vol7d),INTENT(in) :: this
7327LOGICAL :: vol7d_check_alloc_dati
7328
7329vol7d_check_alloc_dati = vol7d_check_alloc_ana(this) .AND. &
7330 ASSOCIATED(this%time) .AND. ASSOCIATED(this%level) .AND. &
7331 ASSOCIATED(this%timerange)
7332
7333END FUNCTION vol7d_check_alloc_dati
7334
7335SUBROUTINE vol7d_force_alloc_dati(this, ini)
7336TYPE(vol7d),INTENT(inout) :: this
7337LOGICAL,INTENT(in),OPTIONAL :: ini
7338
7339! Alloco i descrittori minimi per avere un volume di dati
7340CALL vol7d_force_alloc_ana(this, ini)
7341IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=1, ini=ini)
7342IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=1, ini=ini)
7343IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=1, ini=ini)
7344
7345END SUBROUTINE vol7d_force_alloc_dati
7346
7347
7348SUBROUTINE vol7d_force_alloc(this)
7349TYPE(vol7d),INTENT(inout) :: this
7350
7351! If anything really not allocated yet, allocate with size 0
7352IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=0)
7353IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=0)
7354IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=0)
7355IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=0)
7356IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=0)
7357
7358END SUBROUTINE vol7d_force_alloc
7359
7360
7361FUNCTION vol7d_check_vol(this)
7362TYPE(vol7d),INTENT(in) :: this
7363LOGICAL :: vol7d_check_vol
7364
7365vol7d_check_vol = c_e(this)
7366
7367! Anagrafica
7368IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
7369 vol7d_check_vol = .false.
7370ENDIF
7371
7372IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
7373 vol7d_check_vol = .false.
7374ENDIF
7375
7376IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
7377 vol7d_check_vol = .false.
7378ENDIF
7379
7380IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
7381 vol7d_check_vol = .false.
7382ENDIF
7383
7384IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
7385 vol7d_check_vol = .false.
7386ENDIF
7387IF (ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
7388 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
7389 ASSOCIATED(this%anavar%c)) THEN
7390 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_ana(this)
7391ENDIF
7392
7393! Attributi dell'anagrafica
7394IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
7395 .NOT.ASSOCIATED(this%volanaattrr)) THEN
7396 vol7d_check_vol = .false.
7397ENDIF
7398
7399IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
7400 .NOT.ASSOCIATED(this%volanaattrd)) THEN
7401 vol7d_check_vol = .false.
7402ENDIF
7403
7404IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
7405 .NOT.ASSOCIATED(this%volanaattri)) THEN
7406 vol7d_check_vol = .false.
7407ENDIF
7408
7409IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
7410 .NOT.ASSOCIATED(this%volanaattrb)) THEN
7411 vol7d_check_vol = .false.
7412ENDIF
7413
7414IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
7415 .NOT.ASSOCIATED(this%volanaattrc)) THEN
7416 vol7d_check_vol = .false.
7417ENDIF
7418
7419! Dati
7420IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
7421 vol7d_check_vol = .false.
7422ENDIF
7423
7424IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
7425 vol7d_check_vol = .false.
7426ENDIF
7427
7428IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
7429 vol7d_check_vol = .false.
7430ENDIF
7431
7432IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
7433 vol7d_check_vol = .false.
7434ENDIF
7435
7436IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
7437 vol7d_check_vol = .false.
7438ENDIF
7439
7440! Attributi dei dati
7441IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
7442 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
7443 vol7d_check_vol = .false.
7444ENDIF
7445
7446IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
7447 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
7448 vol7d_check_vol = .false.
7449ENDIF
7450
7451IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
7452 .NOT.ASSOCIATED(this%voldatiattri)) THEN
7453 vol7d_check_vol = .false.
7454ENDIF
7455
7456IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
7457 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
7458 vol7d_check_vol = .false.
7459ENDIF
7460
7461IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
7462 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
7463 vol7d_check_vol = .false.
7464ENDIF
7465IF (ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
7466 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
7467 ASSOCIATED(this%dativar%c)) THEN
7468 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_dati(this)
7469ENDIF
7470
7471END FUNCTION vol7d_check_vol
7472
7473
7488SUBROUTINE vol7d_alloc_vol(this, ini, inivol)
7489TYPE(vol7d),INTENT(inout) :: this
7490LOGICAL,INTENT(in),OPTIONAL :: ini
7491LOGICAL,INTENT(in),OPTIONAL :: inivol
7492
7493LOGICAL :: linivol
7494
7495IF (PRESENT(inivol)) THEN
7496 linivol = inivol
7497ELSE
7498 linivol = .true.
7499ENDIF
7500
7501! Anagrafica
7502IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
7503 CALL vol7d_force_alloc_ana(this, ini)
7504 ALLOCATE(this%volanar(SIZE(this%ana), SIZE(this%anavar%r), SIZE(this%network)))
7505 IF (linivol) this%volanar(:,:,:) = rmiss
7506ENDIF
7507
7508IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
7509 CALL vol7d_force_alloc_ana(this, ini)
7510 ALLOCATE(this%volanad(SIZE(this%ana), SIZE(this%anavar%d), SIZE(this%network)))
7511 IF (linivol) this%volanad(:,:,:) = rdmiss
7512ENDIF
7513
7514IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
7515 CALL vol7d_force_alloc_ana(this, ini)
7516 ALLOCATE(this%volanai(SIZE(this%ana), SIZE(this%anavar%i), SIZE(this%network)))
7517 IF (linivol) this%volanai(:,:,:) = imiss
7518ENDIF
7519
7520IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
7521 CALL vol7d_force_alloc_ana(this, ini)
7522 ALLOCATE(this%volanab(SIZE(this%ana), SIZE(this%anavar%b), SIZE(this%network)))
7523 IF (linivol) this%volanab(:,:,:) = ibmiss
7524ENDIF
7525
7526IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
7527 CALL vol7d_force_alloc_ana(this, ini)
7528 ALLOCATE(this%volanac(SIZE(this%ana), SIZE(this%anavar%c), SIZE(this%network)))
7529 IF (linivol) this%volanac(:,:,:) = cmiss
7530ENDIF
7531
7532! Attributi dell'anagrafica
7533IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
7534 .NOT.ASSOCIATED(this%volanaattrr)) THEN
7535 CALL vol7d_force_alloc_ana(this, ini)
7536 ALLOCATE(this%volanaattrr(SIZE(this%ana), SIZE(this%anavarattr%r), &
7537 SIZE(this%network), SIZE(this%anaattr%r)))
7538 IF (linivol) this%volanaattrr(:,:,:,:) = rmiss
7539ENDIF
7540
7541IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
7542 .NOT.ASSOCIATED(this%volanaattrd)) THEN
7543 CALL vol7d_force_alloc_ana(this, ini)
7544 ALLOCATE(this%volanaattrd(SIZE(this%ana), SIZE(this%anavarattr%d), &
7545 SIZE(this%network), SIZE(this%anaattr%d)))
7546 IF (linivol) this%volanaattrd(:,:,:,:) = rdmiss
7547ENDIF
7548
7549IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
7550 .NOT.ASSOCIATED(this%volanaattri)) THEN
7551 CALL vol7d_force_alloc_ana(this, ini)
7552 ALLOCATE(this%volanaattri(SIZE(this%ana), SIZE(this%anavarattr%i), &
7553 SIZE(this%network), SIZE(this%anaattr%i)))
7554 IF (linivol) this%volanaattri(:,:,:,:) = imiss
7555ENDIF
7556
7557IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
7558 .NOT.ASSOCIATED(this%volanaattrb)) THEN
7559 CALL vol7d_force_alloc_ana(this, ini)
7560 ALLOCATE(this%volanaattrb(SIZE(this%ana), SIZE(this%anavarattr%b), &
7561 SIZE(this%network), SIZE(this%anaattr%b)))
7562 IF (linivol) this%volanaattrb(:,:,:,:) = ibmiss
7563ENDIF
7564
7565IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
7566 .NOT.ASSOCIATED(this%volanaattrc)) THEN
7567 CALL vol7d_force_alloc_ana(this, ini)
7568 ALLOCATE(this%volanaattrc(SIZE(this%ana), SIZE(this%anavarattr%c), &
7569 SIZE(this%network), SIZE(this%anaattr%c)))
7570 IF (linivol) this%volanaattrc(:,:,:,:) = cmiss
7571ENDIF
7572
7573! Dati
7574IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
7575 CALL vol7d_force_alloc_dati(this, ini)
7576 ALLOCATE(this%voldatir(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
7577 SIZE(this%timerange), SIZE(this%dativar%r), SIZE(this%network)))
7578 IF (linivol) this%voldatir(:,:,:,:,:,:) = rmiss
7579ENDIF
7580
7581IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
7582 CALL vol7d_force_alloc_dati(this, ini)
7583 ALLOCATE(this%voldatid(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
7584 SIZE(this%timerange), SIZE(this%dativar%d), SIZE(this%network)))
7585 IF (linivol) this%voldatid(:,:,:,:,:,:) = rdmiss
7586ENDIF
7587
7588IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
7589 CALL vol7d_force_alloc_dati(this, ini)
7590 ALLOCATE(this%voldatii(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
7591 SIZE(this%timerange), SIZE(this%dativar%i), SIZE(this%network)))
7592 IF (linivol) this%voldatii(:,:,:,:,:,:) = imiss
7593ENDIF
7594
7595IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
7596 CALL vol7d_force_alloc_dati(this, ini)
7597 ALLOCATE(this%voldatib(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
7598 SIZE(this%timerange), SIZE(this%dativar%b), SIZE(this%network)))
7599 IF (linivol) this%voldatib(:,:,:,:,:,:) = ibmiss
7600ENDIF
7601
7602IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
7603 CALL vol7d_force_alloc_dati(this, ini)
7604 ALLOCATE(this%voldatic(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
7605 SIZE(this%timerange), SIZE(this%dativar%c), SIZE(this%network)))
7606 IF (linivol) this%voldatic(:,:,:,:,:,:) = cmiss
7607ENDIF
7608
7609! Attributi dei dati
7610IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
7611 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
7612 CALL vol7d_force_alloc_dati(this, ini)
7613 ALLOCATE(this%voldatiattrr(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
7614 SIZE(this%timerange), SIZE(this%dativarattr%r), SIZE(this%network), &
7615 SIZE(this%datiattr%r)))
7616 IF (linivol) this%voldatiattrr(:,:,:,:,:,:,:) = rmiss
7617ENDIF
7618
7619IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
7620 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
7621 CALL vol7d_force_alloc_dati(this, ini)
7622 ALLOCATE(this%voldatiattrd(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
7623 SIZE(this%timerange), SIZE(this%dativarattr%d), SIZE(this%network), &
7624 SIZE(this%datiattr%d)))
7625 IF (linivol) this%voldatiattrd(:,:,:,:,:,:,:) = rdmiss
7626ENDIF
7627
7628IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
7629 .NOT.ASSOCIATED(this%voldatiattri)) THEN
7630 CALL vol7d_force_alloc_dati(this, ini)
7631 ALLOCATE(this%voldatiattri(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
7632 SIZE(this%timerange), SIZE(this%dativarattr%i), SIZE(this%network), &
7633 SIZE(this%datiattr%i)))
7634 IF (linivol) this%voldatiattri(:,:,:,:,:,:,:) = imiss
7635ENDIF
7636
7637IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
7638 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
7639 CALL vol7d_force_alloc_dati(this, ini)
7640 ALLOCATE(this%voldatiattrb(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
7641 SIZE(this%timerange), SIZE(this%dativarattr%b), SIZE(this%network), &
7642 SIZE(this%datiattr%b)))
7643 IF (linivol) this%voldatiattrb(:,:,:,:,:,:,:) = ibmiss
7644ENDIF
7645
7646IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
7647 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
7648 CALL vol7d_force_alloc_dati(this, ini)
7649 ALLOCATE(this%voldatiattrc(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
7650 SIZE(this%timerange), SIZE(this%dativarattr%c), SIZE(this%network), &
7651 SIZE(this%datiattr%c)))
7652 IF (linivol) this%voldatiattrc(:,:,:,:,:,:,:) = cmiss
7653ENDIF
7654
7655! Catch-all method
7656CALL vol7d_force_alloc(this)
7657
7658! Creo gli indici var-attr
7659
7660#ifdef DEBUG
7661CALL l4f_log(l4f_debug,"calling: vol7d_set_attr_ind")
7662#endif
7663
7664CALL vol7d_set_attr_ind(this)
7665
7666
7667
7668END SUBROUTINE vol7d_alloc_vol
7669
7670
7677SUBROUTINE vol7d_set_attr_ind(this)
7678TYPE(vol7d),INTENT(inout) :: this
7679
7680INTEGER :: i
7681
7682! real
7683IF (ASSOCIATED(this%dativar%r)) THEN
7684 IF (ASSOCIATED(this%dativarattr%r)) THEN
7685 DO i = 1, SIZE(this%dativar%r)
7686 this%dativar%r(i)%r = &
7687 firsttrue(this%dativar%r(i)%btable == this%dativarattr%r(:)%btable)
7688 ENDDO
7689 ENDIF
7690
7691 IF (ASSOCIATED(this%dativarattr%d)) THEN
7692 DO i = 1, SIZE(this%dativar%r)
7693 this%dativar%r(i)%d = &
7694 firsttrue(this%dativar%r(i)%btable == this%dativarattr%d(:)%btable)
7695 ENDDO
7696 ENDIF
7697
7698 IF (ASSOCIATED(this%dativarattr%i)) THEN
7699 DO i = 1, SIZE(this%dativar%r)
7700 this%dativar%r(i)%i = &
7701 firsttrue(this%dativar%r(i)%btable == this%dativarattr%i(:)%btable)
7702 ENDDO
7703 ENDIF
7704
7705 IF (ASSOCIATED(this%dativarattr%b)) THEN
7706 DO i = 1, SIZE(this%dativar%r)
7707 this%dativar%r(i)%b = &
7708 firsttrue(this%dativar%r(i)%btable == this%dativarattr%b(:)%btable)
7709 ENDDO
7710 ENDIF
7711
7712 IF (ASSOCIATED(this%dativarattr%c)) THEN
7713 DO i = 1, SIZE(this%dativar%r)
7714 this%dativar%r(i)%c = &
7715 firsttrue(this%dativar%r(i)%btable == this%dativarattr%c(:)%btable)
7716 ENDDO
7717 ENDIF
7718ENDIF
7719! double
7720IF (ASSOCIATED(this%dativar%d)) THEN
7721 IF (ASSOCIATED(this%dativarattr%r)) THEN
7722 DO i = 1, SIZE(this%dativar%d)
7723 this%dativar%d(i)%r = &
7724 firsttrue(this%dativar%d(i)%btable == this%dativarattr%r(:)%btable)
7725 ENDDO
7726 ENDIF
7727
7728 IF (ASSOCIATED(this%dativarattr%d)) THEN
7729 DO i = 1, SIZE(this%dativar%d)
7730 this%dativar%d(i)%d = &
7731 firsttrue(this%dativar%d(i)%btable == this%dativarattr%d(:)%btable)
7732 ENDDO
7733 ENDIF
7734
7735 IF (ASSOCIATED(this%dativarattr%i)) THEN
7736 DO i = 1, SIZE(this%dativar%d)
7737 this%dativar%d(i)%i = &
7738 firsttrue(this%dativar%d(i)%btable == this%dativarattr%i(:)%btable)
7739 ENDDO
7740 ENDIF
7741
7742 IF (ASSOCIATED(this%dativarattr%b)) THEN
7743 DO i = 1, SIZE(this%dativar%d)
7744 this%dativar%d(i)%b = &
7745 firsttrue(this%dativar%d(i)%btable == this%dativarattr%b(:)%btable)
7746 ENDDO
7747 ENDIF
7748
7749 IF (ASSOCIATED(this%dativarattr%c)) THEN
7750 DO i = 1, SIZE(this%dativar%d)
7751 this%dativar%d(i)%c = &
7752 firsttrue(this%dativar%d(i)%btable == this%dativarattr%c(:)%btable)
7753 ENDDO
7754 ENDIF
7755ENDIF
7756! integer
7757IF (ASSOCIATED(this%dativar%i)) THEN
7758 IF (ASSOCIATED(this%dativarattr%r)) THEN
7759 DO i = 1, SIZE(this%dativar%i)
7760 this%dativar%i(i)%r = &
7761 firsttrue(this%dativar%i(i)%btable == this%dativarattr%r(:)%btable)
7762 ENDDO
7763 ENDIF
7764
7765 IF (ASSOCIATED(this%dativarattr%d)) THEN
7766 DO i = 1, SIZE(this%dativar%i)
7767 this%dativar%i(i)%d = &
7768 firsttrue(this%dativar%i(i)%btable == this%dativarattr%d(:)%btable)
7769 ENDDO
7770 ENDIF
7771
7772 IF (ASSOCIATED(this%dativarattr%i)) THEN
7773 DO i = 1, SIZE(this%dativar%i)
7774 this%dativar%i(i)%i = &
7775 firsttrue(this%dativar%i(i)%btable == this%dativarattr%i(:)%btable)
7776 ENDDO
7777 ENDIF
7778
7779 IF (ASSOCIATED(this%dativarattr%b)) THEN
7780 DO i = 1, SIZE(this%dativar%i)
7781 this%dativar%i(i)%b = &
7782 firsttrue(this%dativar%i(i)%btable == this%dativarattr%b(:)%btable)
7783 ENDDO
7784 ENDIF
7785
7786 IF (ASSOCIATED(this%dativarattr%c)) THEN
7787 DO i = 1, SIZE(this%dativar%i)
7788 this%dativar%i(i)%c = &
7789 firsttrue(this%dativar%i(i)%btable == this%dativarattr%c(:)%btable)
7790 ENDDO
7791 ENDIF
7792ENDIF
7793! byte
7794IF (ASSOCIATED(this%dativar%b)) THEN
7795 IF (ASSOCIATED(this%dativarattr%r)) THEN
7796 DO i = 1, SIZE(this%dativar%b)
7797 this%dativar%b(i)%r = &
7798 firsttrue(this%dativar%b(i)%btable == this%dativarattr%r(:)%btable)
7799 ENDDO
7800 ENDIF
7801
7802 IF (ASSOCIATED(this%dativarattr%d)) THEN
7803 DO i = 1, SIZE(this%dativar%b)
7804 this%dativar%b(i)%d = &
7805 firsttrue(this%dativar%b(i)%btable == this%dativarattr%d(:)%btable)
7806 ENDDO
7807 ENDIF
7808
7809 IF (ASSOCIATED(this%dativarattr%i)) THEN
7810 DO i = 1, SIZE(this%dativar%b)
7811 this%dativar%b(i)%i = &
7812 firsttrue(this%dativar%b(i)%btable == this%dativarattr%i(:)%btable)
7813 ENDDO
7814 ENDIF
7815
7816 IF (ASSOCIATED(this%dativarattr%b)) THEN
7817 DO i = 1, SIZE(this%dativar%b)
7818 this%dativar%b(i)%b = &
7819 firsttrue(this%dativar%b(i)%btable == this%dativarattr%b(:)%btable)
7820 ENDDO
7821 ENDIF
7822
7823 IF (ASSOCIATED(this%dativarattr%c)) THEN
7824 DO i = 1, SIZE(this%dativar%b)
7825 this%dativar%b(i)%c = &
7826 firsttrue(this%dativar%b(i)%btable == this%dativarattr%c(:)%btable)
7827 ENDDO
7828 ENDIF
7829ENDIF
7830! character
7831IF (ASSOCIATED(this%dativar%c)) THEN
7832 IF (ASSOCIATED(this%dativarattr%r)) THEN
7833 DO i = 1, SIZE(this%dativar%c)
7834 this%dativar%c(i)%r = &
7835 firsttrue(this%dativar%c(i)%btable == this%dativarattr%r(:)%btable)
7836 ENDDO
7837 ENDIF
7838
7839 IF (ASSOCIATED(this%dativarattr%d)) THEN
7840 DO i = 1, SIZE(this%dativar%c)
7841 this%dativar%c(i)%d = &
7842 firsttrue(this%dativar%c(i)%btable == this%dativarattr%d(:)%btable)
7843 ENDDO
7844 ENDIF
7845
7846 IF (ASSOCIATED(this%dativarattr%i)) THEN
7847 DO i = 1, SIZE(this%dativar%c)
7848 this%dativar%c(i)%i = &
7849 firsttrue(this%dativar%c(i)%btable == this%dativarattr%i(:)%btable)
7850 ENDDO
7851 ENDIF
7852
7853 IF (ASSOCIATED(this%dativarattr%b)) THEN
7854 DO i = 1, SIZE(this%dativar%c)
7855 this%dativar%c(i)%b = &
7856 firsttrue(this%dativar%c(i)%btable == this%dativarattr%b(:)%btable)
7857 ENDDO
7858 ENDIF
7859
7860 IF (ASSOCIATED(this%dativarattr%c)) THEN
7861 DO i = 1, SIZE(this%dativar%c)
7862 this%dativar%c(i)%c = &
7863 firsttrue(this%dativar%c(i)%btable == this%dativarattr%c(:)%btable)
7864 ENDDO
7865 ENDIF
7866ENDIF
7867
7868END SUBROUTINE vol7d_set_attr_ind
7869
7870
7875SUBROUTINE vol7d_merge(this, that, sort, bestdata, &
7876 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
7877TYPE(vol7d),INTENT(INOUT) :: this
7878TYPE(vol7d),INTENT(INOUT) :: that
7879LOGICAL,INTENT(IN),OPTIONAL :: sort
7880LOGICAL,INTENT(in),OPTIONAL :: bestdata
7881LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple ! experimental, please do not use outside the library now
7882
7883TYPE(vol7d) :: v7d_clean
7884
7885
7886IF (.NOT.c_e(this)) THEN ! speedup
7887 this = that
7888 CALL init(v7d_clean)
7889 that = v7d_clean ! destroy that without deallocating
7890ELSE ! Append that to this and destroy that
7891 CALL vol7d_append(this, that, sort, bestdata, &
7892 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
7893 CALL delete(that)
7894ENDIF
7895
7896END SUBROUTINE vol7d_merge
7897
7898
7927SUBROUTINE vol7d_append(this, that, sort, bestdata, &
7928 ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple)
7929TYPE(vol7d),INTENT(INOUT) :: this
7930TYPE(vol7d),INTENT(IN) :: that
7931LOGICAL,INTENT(IN),OPTIONAL :: sort
7932! experimental, please do not use outside the library now, they force the use
7933! of a simplified mapping algorithm which is valid only whene the dimension
7934! content is the same in both volumes , or when one of them is empty
7935LOGICAL,INTENT(in),OPTIONAL :: bestdata
7936LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple
7937
7938
7939TYPE(vol7d) :: v7dtmp
7940LOGICAL :: lsort, lbestdata
7941INTEGER,POINTER :: remapt1(:), remapt2(:), remaptr1(:), remaptr2(:), &
7942 remapl1(:), remapl2(:), remapa1(:), remapa2(:), remapn1(:), remapn2(:)
7943
7944IF (.NOT.c_e(that)) RETURN ! speedup, nothing to do
7945IF (.NOT.vol7d_check_vol(that)) RETURN ! be safe
7946IF (.NOT.c_e(this)) THEN ! this case is like a vol7d_copy, more efficient to copy?
7947 CALL vol7d_copy(that, this, sort=sort)
7948 RETURN
7949ENDIF
7950
7951IF (this%time_definition /= that%time_definition) THEN
7952 CALL l4f_log(l4f_fatal, &
7953 'in vol7d_append, cannot append volumes with different &
7954 &time definition')
7955 CALL raise_fatal_error()
7956ENDIF
7957
7958! Completo l'allocazione per avere volumi a norma
7959CALL vol7d_alloc_vol(this)
7960
7961CALL init(v7dtmp, time_definition=this%time_definition)
7962CALL optio(sort, lsort)
7963CALL optio(bestdata, lbestdata)
7964
7965! Calcolo le mappature tra volumi vecchi e volume nuovo
7966! I puntatori remap* vengono tutti o allocati o nullificati
7967IF (optio_log(ltimesimple)) THEN
7968 CALL vol7d_remap2simple_datetime(this%time, that%time, v7dtmp%time, &
7969 lsort, remapt1, remapt2)
7970ELSE
7971 CALL vol7d_remap2_datetime(this%time, that%time, v7dtmp%time, &
7972 lsort, remapt1, remapt2)
7973ENDIF
7974IF (optio_log(ltimerangesimple)) THEN
7975 CALL vol7d_remap2simple_vol7d_timerange(this%timerange, that%timerange, &
7976 v7dtmp%timerange, lsort, remaptr1, remaptr2)
7977ELSE
7978 CALL vol7d_remap2_vol7d_timerange(this%timerange, that%timerange, &
7979 v7dtmp%timerange, lsort, remaptr1, remaptr2)
7980ENDIF
7981IF (optio_log(llevelsimple)) THEN
7982 CALL vol7d_remap2simple_vol7d_level(this%level, that%level, v7dtmp%level, &
7983 lsort, remapl1, remapl2)
7984ELSE
7985 CALL vol7d_remap2_vol7d_level(this%level, that%level, v7dtmp%level, &
7986 lsort, remapl1, remapl2)
7987ENDIF
7988IF (optio_log(lanasimple)) THEN
7989 CALL vol7d_remap2simple_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
7990 .false., remapa1, remapa2)
7991ELSE
7992 CALL vol7d_remap2_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
7993 .false., remapa1, remapa2)
7994ENDIF
7995IF (optio_log(lnetworksimple)) THEN
7996 CALL vol7d_remap2simple_vol7d_network(this%network, that%network, v7dtmp%network, &
7997 .false., remapn1, remapn2)
7998ELSE
7999 CALL vol7d_remap2_vol7d_network(this%network, that%network, v7dtmp%network, &
8000 .false., remapn1, remapn2)
8001ENDIF
8002
8003! Faccio la fusione fisica dei volumi
8004CALL vol7d_merge_finalr(this, that, v7dtmp, &
8005 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
8006 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
8007CALL vol7d_merge_finald(this, that, v7dtmp, &
8008 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
8009 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
8010CALL vol7d_merge_finali(this, that, v7dtmp, &
8011 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
8012 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
8013CALL vol7d_merge_finalb(this, that, v7dtmp, &
8014 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
8015 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
8016CALL vol7d_merge_finalc(this, that, v7dtmp, &
8017 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
8018 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
8019
8020! Dealloco i vettori di rimappatura
8021IF (ASSOCIATED(remapt1)) DEALLOCATE(remapt1)
8022IF (ASSOCIATED(remapt2)) DEALLOCATE(remapt2)
8023IF (ASSOCIATED(remaptr1)) DEALLOCATE(remaptr1)
8024IF (ASSOCIATED(remaptr2)) DEALLOCATE(remaptr2)
8025IF (ASSOCIATED(remapl1)) DEALLOCATE(remapl1)
8026IF (ASSOCIATED(remapl2)) DEALLOCATE(remapl2)
8027IF (ASSOCIATED(remapa1)) DEALLOCATE(remapa1)
8028IF (ASSOCIATED(remapa2)) DEALLOCATE(remapa2)
8029IF (ASSOCIATED(remapn1)) DEALLOCATE(remapn1)
8030IF (ASSOCIATED(remapn2)) DEALLOCATE(remapn2)
8031
8032! Distruggo il vecchio volume e assegno il nuovo a this
8033CALL delete(this)
8034this = v7dtmp
8035! Ricreo gli indici var-attr
8036CALL vol7d_set_attr_ind(this)
8037
8038END SUBROUTINE vol7d_append
8039
8040
8073SUBROUTINE vol7d_copy(this, that, sort, unique, miss, &
8074 lsort_time, lsort_timerange, lsort_level, &
8075 ltime, ltimerange, llevel, lana, lnetwork, &
8076 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
8077 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
8078 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
8079 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
8080 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
8081 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
8082TYPE(vol7d),INTENT(IN) :: this
8083TYPE(vol7d),INTENT(INOUT) :: that
8084LOGICAL,INTENT(IN),OPTIONAL :: sort
8085LOGICAL,INTENT(IN),OPTIONAL :: unique
8086LOGICAL,INTENT(IN),OPTIONAL :: miss
8087LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
8088LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
8089LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
8097LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
8099LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
8101LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
8103LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
8105LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
8107LOGICAL,INTENT(in),OPTIONAL :: &
8108 lanavarr(:), lanavard(:), lanavari(:), lanavarb(:), lanavarc(:), &
8109 lanaattrr(:), lanaattrd(:), lanaattri(:), lanaattrb(:), lanaattrc(:), &
8110 lanavarattrr(:), lanavarattrd(:), lanavarattri(:), lanavarattrb(:), lanavarattrc(:), &
8111 ldativarr(:), ldativard(:), ldativari(:), ldativarb(:), ldativarc(:), &
8112 ldatiattrr(:), ldatiattrd(:), ldatiattri(:), ldatiattrb(:), ldatiattrc(:), &
8113 ldativarattrr(:), ldativarattrd(:), ldativarattri(:), ldativarattrb(:), ldativarattrc(:)
8114
8115LOGICAL :: lsort, lunique, lmiss
8116INTEGER,POINTER :: remapt(:), remaptr(:), remapl(:), remapa(:), remapn(:)
8117
8118CALL init(that)
8119IF (.NOT.c_e(this)) RETURN ! speedup, nothing to do
8120IF (.NOT.vol7d_check_vol(this)) RETURN ! be safe
8121
8122CALL optio(sort, lsort)
8123CALL optio(unique, lunique)
8124CALL optio(miss, lmiss)
8125
8126! Calcolo le mappature tra volume vecchio e volume nuovo
8127! I puntatori remap* vengono tutti o allocati o nullificati
8128CALL vol7d_remap1_datetime(this%time, that%time, &
8129 lsort.OR.optio_log(lsort_time), lunique, lmiss, remapt, ltime)
8130CALL vol7d_remap1_vol7d_timerange(this%timerange, that%timerange, &
8131 lsort.OR.optio_log(lsort_timerange), lunique, lmiss, remaptr, ltimerange)
8132CALL vol7d_remap1_vol7d_level(this%level, that%level, &
8133 lsort.OR.optio_log(lsort_level), lunique, lmiss, remapl, llevel)
8134CALL vol7d_remap1_vol7d_ana(this%ana, that%ana, &
8135 lsort, lunique, lmiss, remapa, lana)
8136CALL vol7d_remap1_vol7d_network(this%network, that%network, &
8137 lsort, lunique, lmiss, remapn, lnetwork)
8138
8139! lanavari, lanavarb, lanavarc, &
8140! lanaattri, lanaattrb, lanaattrc, &
8141! lanavarattri, lanavarattrb, lanavarattrc, &
8142! ldativari, ldativarb, ldativarc, &
8143! ldatiattri, ldatiattrb, ldatiattrc, &
8144! ldativarattri, ldativarattrb, ldativarattrc
8145! Faccio la riforma fisica dei volumi
8146CALL vol7d_reform_finalr(this, that, &
8147 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
8148 lanavarr, lanaattrr, lanavarattrr, ldativarr, ldatiattrr, ldativarattrr)
8149CALL vol7d_reform_finald(this, that, &
8150 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
8151 lanavard, lanaattrd, lanavarattrd, ldativard, ldatiattrd, ldativarattrd)
8152CALL vol7d_reform_finali(this, that, &
8153 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
8154 lanavari, lanaattri, lanavarattri, ldativari, ldatiattri, ldativarattri)
8155CALL vol7d_reform_finalb(this, that, &
8156 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
8157 lanavarb, lanaattrb, lanavarattrb, ldativarb, ldatiattrb, ldativarattrb)
8158CALL vol7d_reform_finalc(this, that, &
8159 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
8160 lanavarc, lanaattrc, lanavarattrc, ldativarc, ldatiattrc, ldativarattrc)
8161
8162! Dealloco i vettori di rimappatura
8163IF (ASSOCIATED(remapt)) DEALLOCATE(remapt)
8164IF (ASSOCIATED(remaptr)) DEALLOCATE(remaptr)
8165IF (ASSOCIATED(remapl)) DEALLOCATE(remapl)
8166IF (ASSOCIATED(remapa)) DEALLOCATE(remapa)
8167IF (ASSOCIATED(remapn)) DEALLOCATE(remapn)
8168
8169! Ricreo gli indici var-attr
8170CALL vol7d_set_attr_ind(that)
8171that%time_definition = this%time_definition
8172
8173END SUBROUTINE vol7d_copy
8174
8175
8186SUBROUTINE vol7d_reform(this, sort, unique, miss, &
8187 lsort_time, lsort_timerange, lsort_level, &
8188 ltime, ltimerange, llevel, lana, lnetwork, &
8189 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
8190 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
8191 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
8192 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
8193 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
8194 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc&
8195 ,purgeana)
8196TYPE(vol7d),INTENT(INOUT) :: this
8197LOGICAL,INTENT(IN),OPTIONAL :: sort
8198LOGICAL,INTENT(IN),OPTIONAL :: unique
8199LOGICAL,INTENT(IN),OPTIONAL :: miss
8200LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
8201LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
8202LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
8210LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
8211LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
8212LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
8213LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
8214LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
8216LOGICAL,INTENT(in),OPTIONAL :: &
8217 lanavarr(:), lanavard(:), lanavari(:), lanavarb(:), lanavarc(:), &
8218 lanaattrr(:), lanaattrd(:), lanaattri(:), lanaattrb(:), lanaattrc(:), &
8219 lanavarattrr(:), lanavarattrd(:), lanavarattri(:), lanavarattrb(:), lanavarattrc(:), &
8220 ldativarr(:), ldativard(:), ldativari(:), ldativarb(:), ldativarc(:), &
8221 ldatiattrr(:), ldatiattrd(:), ldatiattri(:), ldatiattrb(:), ldatiattrc(:), &
8222 ldativarattrr(:), ldativarattrd(:), ldativarattri(:), ldativarattrb(:), ldativarattrc(:)
8223LOGICAL,INTENT(IN),OPTIONAL :: purgeana
8224
8225TYPE(vol7d) :: v7dtmp
8226logical,allocatable :: llana(:)
8227integer :: i
8228
8229CALL vol7d_copy(this, v7dtmp, sort, unique, miss, &
8230 lsort_time, lsort_timerange, lsort_level, &
8231 ltime, ltimerange, llevel, lana, lnetwork, &
8232 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
8233 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
8234 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
8235 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
8236 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
8237 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
8238
8239! destroy old volume
8240CALL delete(this)
8241
8242if (optio_log(purgeana)) then
8243 allocate(llana(size(v7dtmp%ana)))
8244 llana =.false.
8245 do i =1,size(v7dtmp%ana)
8246 if (associated(v7dtmp%voldatii)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatii(i,:,:,:,:,:)))
8247 if (associated(v7dtmp%voldatir)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatir(i,:,:,:,:,:)))
8248 if (associated(v7dtmp%voldatid)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatid(i,:,:,:,:,:)))
8249 if (associated(v7dtmp%voldatib)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatib(i,:,:,:,:,:)))
8250 if (associated(v7dtmp%voldatic)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatic(i,:,:,:,:,:)))
8251 end do
8252 CALL vol7d_copy(v7dtmp, this,lana=llana)
8253 CALL delete(v7dtmp)
8254 deallocate(llana)
8255else
8256 this=v7dtmp
8257end if
8258
8259END SUBROUTINE vol7d_reform
8260
8261
8269SUBROUTINE vol7d_smart_sort(this, lsort_time, lsort_timerange, lsort_level)
8270TYPE(vol7d),INTENT(INOUT) :: this
8271LOGICAL,OPTIONAL,INTENT(in) :: lsort_time
8272LOGICAL,OPTIONAL,INTENT(in) :: lsort_timerange
8273LOGICAL,OPTIONAL,INTENT(in) :: lsort_level
8274
8275INTEGER :: i
8276LOGICAL :: to_be_sorted
8277
8278to_be_sorted = .false.
8279CALL vol7d_alloc_vol(this) ! usual safety check
8280
8281IF (optio_log(lsort_time)) THEN
8282 DO i = 2, SIZE(this%time)
8283 IF (this%time(i) < this%time(i-1)) THEN
8284 to_be_sorted = .true.
8285 EXIT
8286 ENDIF
8287 ENDDO
8288ENDIF
8289IF (optio_log(lsort_timerange)) THEN
8290 DO i = 2, SIZE(this%timerange)
8291 IF (this%timerange(i) < this%timerange(i-1)) THEN
8292 to_be_sorted = .true.
8293 EXIT
8294 ENDIF
8295 ENDDO
8296ENDIF
8297IF (optio_log(lsort_level)) THEN
8298 DO i = 2, SIZE(this%level)
8299 IF (this%level(i) < this%level(i-1)) THEN
8300 to_be_sorted = .true.
8301 EXIT
8302 ENDIF
8303 ENDDO
8304ENDIF
8305
8306IF (to_be_sorted) CALL vol7d_reform(this, &
8307 lsort_time=lsort_time, lsort_timerange=lsort_timerange, lsort_level=lsort_level )
8308
8309END SUBROUTINE vol7d_smart_sort
8310
8318SUBROUTINE vol7d_filter(this, avl, vl, nl, s_d, e_d)
8319TYPE(vol7d),INTENT(inout) :: this
8320CHARACTER(len=*),INTENT(in),OPTIONAL :: avl(:)
8321CHARACTER(len=*),INTENT(in),OPTIONAL :: vl(:)
8322TYPE(vol7d_network),OPTIONAL :: nl(:)
8323TYPE(datetime),INTENT(in),OPTIONAL :: s_d
8324TYPE(datetime),INTENT(in),OPTIONAL :: e_d
8325
8326INTEGER :: i
8327
8328IF (PRESENT(avl)) THEN
8329 IF (SIZE(avl) > 0) THEN
8330
8331 IF (ASSOCIATED(this%anavar%r)) THEN
8332 DO i = 1, SIZE(this%anavar%r)
8333 IF (all(this%anavar%r(i)%btable /= avl)) this%anavar%r(i) = vol7d_var_miss
8334 ENDDO
8335 ENDIF
8336
8337 IF (ASSOCIATED(this%anavar%i)) THEN
8338 DO i = 1, SIZE(this%anavar%i)
8339 IF (all(this%anavar%i(i)%btable /= avl)) this%anavar%i(i) = vol7d_var_miss
8340 ENDDO
8341 ENDIF
8342
8343 IF (ASSOCIATED(this%anavar%b)) THEN
8344 DO i = 1, SIZE(this%anavar%b)
8345 IF (all(this%anavar%b(i)%btable /= avl)) this%anavar%b(i) = vol7d_var_miss
8346 ENDDO
8347 ENDIF
8348
8349 IF (ASSOCIATED(this%anavar%d)) THEN
8350 DO i = 1, SIZE(this%anavar%d)
8351 IF (all(this%anavar%d(i)%btable /= avl)) this%anavar%d(i) = vol7d_var_miss
8352 ENDDO
8353 ENDIF
8354
8355 IF (ASSOCIATED(this%anavar%c)) THEN
8356 DO i = 1, SIZE(this%anavar%c)
8357 IF (all(this%anavar%c(i)%btable /= avl)) this%anavar%c(i) = vol7d_var_miss
8358 ENDDO
8359 ENDIF
8360
8361 ENDIF
8362ENDIF
8363
8364
8365IF (PRESENT(vl)) THEN
8366 IF (size(vl) > 0) THEN
8367 IF (ASSOCIATED(this%dativar%r)) THEN
8368 DO i = 1, SIZE(this%dativar%r)
8369 IF (all(this%dativar%r(i)%btable /= vl)) this%dativar%r(i) = vol7d_var_miss
8370 ENDDO
8371 ENDIF
8372
8373 IF (ASSOCIATED(this%dativar%i)) THEN
8374 DO i = 1, SIZE(this%dativar%i)
8375 IF (all(this%dativar%i(i)%btable /= vl)) this%dativar%i(i) = vol7d_var_miss
8376 ENDDO
8377 ENDIF
8378
8379 IF (ASSOCIATED(this%dativar%b)) THEN
8380 DO i = 1, SIZE(this%dativar%b)
8381 IF (all(this%dativar%b(i)%btable /= vl)) this%dativar%b(i) = vol7d_var_miss
8382 ENDDO
8383 ENDIF
8384
8385 IF (ASSOCIATED(this%dativar%d)) THEN
8386 DO i = 1, SIZE(this%dativar%d)
8387 IF (all(this%dativar%d(i)%btable /= vl)) this%dativar%d(i) = vol7d_var_miss
8388 ENDDO
8389 ENDIF
8390
8391 IF (ASSOCIATED(this%dativar%c)) THEN
8392 DO i = 1, SIZE(this%dativar%c)
8393 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
8394 ENDDO
8395 ENDIF
8396
8397 IF (ASSOCIATED(this%dativar%c)) THEN
8398 DO i = 1, SIZE(this%dativar%c)
8399 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
8400 ENDDO
8401 ENDIF
8402
8403 ENDIF
8404ENDIF
8405
8406IF (PRESENT(nl)) THEN
8407 IF (SIZE(nl) > 0) THEN
8408 DO i = 1, SIZE(this%network)
8409 IF (all(this%network(i) /= nl)) this%network(i) = vol7d_network_miss
8410 ENDDO
8411 ENDIF
8412ENDIF
8413
8414IF (PRESENT(s_d)) THEN
8415 IF (c_e(s_d)) THEN
8416 WHERE (this%time < s_d)
8417 this%time = datetime_miss
8418 END WHERE
8419 ENDIF
8420ENDIF
8421
8422IF (PRESENT(e_d)) THEN
8423 IF (c_e(e_d)) THEN
8424 WHERE (this%time > e_d)
8425 this%time = datetime_miss
8426 END WHERE
8427 ENDIF
8428ENDIF
8429
8430CALL vol7d_reform(this, miss=.true.)
8431
8432END SUBROUTINE vol7d_filter
8433
8434
8441SUBROUTINE vol7d_convr(this, that, anaconv)
8442TYPE(vol7d),INTENT(IN) :: this
8443TYPE(vol7d),INTENT(INOUT) :: that
8444LOGICAL,OPTIONAL,INTENT(in) :: anaconv
8445INTEGER :: i
8446LOGICAL :: fv(1)=(/.false./), tv(1)=(/.true./), acp(1), acn(1)
8447TYPE(vol7d) :: v7d_tmp
8448
8449IF (optio_log(anaconv)) THEN
8450 acp=fv
8451 acn=tv
8452ELSE
8453 acp=tv
8454 acn=fv
8455ENDIF
8456
8457! Volume con solo i dati reali e tutti gli attributi
8458! l'anagrafica e` copiata interamente se necessario
8459CALL vol7d_copy(this, that, &
8460 lanavarr=tv, lanavard=acp, lanavari=acp, lanavarb=acp, lanavarc=acp, &
8461 ldativarr=tv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=fv)
8462
8463! Volume solo di dati double
8464CALL vol7d_copy(this, v7d_tmp, &
8465 lanavarr=fv, lanavard=acn, lanavari=fv, lanavarb=fv, lanavarc=fv, &
8466 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
8467 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
8468 ldativarr=fv, ldativard=tv, ldativari=fv, ldativarb=fv, ldativarc=fv, &
8469 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
8470 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
8471
8472! converto a dati reali
8473IF (ASSOCIATED(v7d_tmp%anavar%d) .OR. ASSOCIATED(v7d_tmp%dativar%d)) THEN
8474
8475 IF (ASSOCIATED(v7d_tmp%anavar%d)) THEN
8476! alloco i dati reali e vi trasferisco i double
8477 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanad, 1), SIZE(v7d_tmp%volanad, 2), &
8478 SIZE(v7d_tmp%volanad, 3)))
8479 DO i = 1, SIZE(v7d_tmp%anavar%d)
8480 v7d_tmp%volanar(:,i,:) = &
8481 realdat(v7d_tmp%volanad(:,i,:), v7d_tmp%anavar%d(i))
8482 ENDDO
8483 DEALLOCATE(v7d_tmp%volanad)
8484! trasferisco le variabili
8485 v7d_tmp%anavar%r => v7d_tmp%anavar%d
8486 NULLIFY(v7d_tmp%anavar%d)
8487 ENDIF
8488
8489 IF (ASSOCIATED(v7d_tmp%dativar%d)) THEN
8490! alloco i dati reali e vi trasferisco i double
8491 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatid, 1), SIZE(v7d_tmp%voldatid, 2), &
8492 SIZE(v7d_tmp%voldatid, 3), SIZE(v7d_tmp%voldatid, 4), SIZE(v7d_tmp%voldatid, 5), &
8493 SIZE(v7d_tmp%voldatid, 6)))
8494 DO i = 1, SIZE(v7d_tmp%dativar%d)
8495 v7d_tmp%voldatir(:,:,:,:,i,:) = &
8496 realdat(v7d_tmp%voldatid(:,:,:,:,i,:), v7d_tmp%dativar%d(i))
8497 ENDDO
8498 DEALLOCATE(v7d_tmp%voldatid)
8499! trasferisco le variabili
8500 v7d_tmp%dativar%r => v7d_tmp%dativar%d
8501 NULLIFY(v7d_tmp%dativar%d)
8502 ENDIF
8503
8504! fondo con il volume definitivo
8505 CALL vol7d_merge(that, v7d_tmp)
8506ELSE
8507 CALL delete(v7d_tmp)
8508ENDIF
8509
8510
8511! Volume solo di dati interi
8512CALL vol7d_copy(this, v7d_tmp, &
8513 lanavarr=fv, lanavard=fv, lanavari=acn, lanavarb=fv, lanavarc=fv, &
8514 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
8515 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
8516 ldativarr=fv, ldativard=fv, ldativari=tv, ldativarb=fv, ldativarc=fv, &
8517 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
8518 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
8519
8520! converto a dati reali
8521IF (ASSOCIATED(v7d_tmp%anavar%i) .OR. ASSOCIATED(v7d_tmp%dativar%i)) THEN
8522
8523 IF (ASSOCIATED(v7d_tmp%anavar%i)) THEN
8524! alloco i dati reali e vi trasferisco gli interi
8525 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanai, 1), SIZE(v7d_tmp%volanai, 2), &
8526 SIZE(v7d_tmp%volanai, 3)))
8527 DO i = 1, SIZE(v7d_tmp%anavar%i)
8528 v7d_tmp%volanar(:,i,:) = &
8529 realdat(v7d_tmp%volanai(:,i,:), v7d_tmp%anavar%i(i))
8530 ENDDO
8531 DEALLOCATE(v7d_tmp%volanai)
8532! trasferisco le variabili
8533 v7d_tmp%anavar%r => v7d_tmp%anavar%i
8534 NULLIFY(v7d_tmp%anavar%i)
8535 ENDIF
8536
8537 IF (ASSOCIATED(v7d_tmp%dativar%i)) THEN
8538! alloco i dati reali e vi trasferisco gli interi
8539 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatii, 1), SIZE(v7d_tmp%voldatii, 2), &
8540 SIZE(v7d_tmp%voldatii, 3), SIZE(v7d_tmp%voldatii, 4), SIZE(v7d_tmp%voldatii, 5), &
8541 SIZE(v7d_tmp%voldatii, 6)))
8542 DO i = 1, SIZE(v7d_tmp%dativar%i)
8543 v7d_tmp%voldatir(:,:,:,:,i,:) = &
8544 realdat(v7d_tmp%voldatii(:,:,:,:,i,:), v7d_tmp%dativar%i(i))
8545 ENDDO
8546 DEALLOCATE(v7d_tmp%voldatii)
8547! trasferisco le variabili
8548 v7d_tmp%dativar%r => v7d_tmp%dativar%i
8549 NULLIFY(v7d_tmp%dativar%i)
8550 ENDIF
8551
8552! fondo con il volume definitivo
8553 CALL vol7d_merge(that, v7d_tmp)
8554ELSE
8555 CALL delete(v7d_tmp)
8556ENDIF
8557
8558
8559! Volume solo di dati byte
8560CALL vol7d_copy(this, v7d_tmp, &
8561 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=acn, lanavarc=fv, &
8562 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
8563 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
8564 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=tv, ldativarc=fv, &
8565 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
8566 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
8567
8568! converto a dati reali
8569IF (ASSOCIATED(v7d_tmp%anavar%b) .OR. ASSOCIATED(v7d_tmp%dativar%b)) THEN
8570
8571 IF (ASSOCIATED(v7d_tmp%anavar%b)) THEN
8572! alloco i dati reali e vi trasferisco i byte
8573 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanab, 1), SIZE(v7d_tmp%volanab, 2), &
8574 SIZE(v7d_tmp%volanab, 3)))
8575 DO i = 1, SIZE(v7d_tmp%anavar%b)
8576 v7d_tmp%volanar(:,i,:) = &
8577 realdat(v7d_tmp%volanab(:,i,:), v7d_tmp%anavar%b(i))
8578 ENDDO
8579 DEALLOCATE(v7d_tmp%volanab)
8580! trasferisco le variabili
8581 v7d_tmp%anavar%r => v7d_tmp%anavar%b
8582 NULLIFY(v7d_tmp%anavar%b)
8583 ENDIF
8584
8585 IF (ASSOCIATED(v7d_tmp%dativar%b)) THEN
8586! alloco i dati reali e vi trasferisco i byte
8587 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatib, 1), SIZE(v7d_tmp%voldatib, 2), &
8588 SIZE(v7d_tmp%voldatib, 3), SIZE(v7d_tmp%voldatib, 4), SIZE(v7d_tmp%voldatib, 5), &
8589 SIZE(v7d_tmp%voldatib, 6)))
8590 DO i = 1, SIZE(v7d_tmp%dativar%b)
8591 v7d_tmp%voldatir(:,:,:,:,i,:) = &
8592 realdat(v7d_tmp%voldatib(:,:,:,:,i,:), v7d_tmp%dativar%b(i))
8593 ENDDO
8594 DEALLOCATE(v7d_tmp%voldatib)
8595! trasferisco le variabili
8596 v7d_tmp%dativar%r => v7d_tmp%dativar%b
8597 NULLIFY(v7d_tmp%dativar%b)
8598 ENDIF
8599
8600! fondo con il volume definitivo
8601 CALL vol7d_merge(that, v7d_tmp)
8602ELSE
8603 CALL delete(v7d_tmp)
8604ENDIF
8605
8606
8607! Volume solo di dati character
8608CALL vol7d_copy(this, v7d_tmp, &
8609 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=fv, lanavarc=acn, &
8610 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
8611 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
8612 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=tv, &
8613 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
8614 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
8615
8616! converto a dati reali
8617IF (ASSOCIATED(v7d_tmp%anavar%c) .OR. ASSOCIATED(v7d_tmp%dativar%c)) THEN
8618
8619 IF (ASSOCIATED(v7d_tmp%anavar%c)) THEN
8620! alloco i dati reali e vi trasferisco i character
8621 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanac, 1), SIZE(v7d_tmp%volanac, 2), &
8622 SIZE(v7d_tmp%volanac, 3)))
8623 DO i = 1, SIZE(v7d_tmp%anavar%c)
8624 v7d_tmp%volanar(:,i,:) = &
8625 realdat(v7d_tmp%volanac(:,i,:), v7d_tmp%anavar%c(i))
8626 ENDDO
8627 DEALLOCATE(v7d_tmp%volanac)
8628! trasferisco le variabili
8629 v7d_tmp%anavar%r => v7d_tmp%anavar%c
8630 NULLIFY(v7d_tmp%anavar%c)
8631 ENDIF
8632
8633 IF (ASSOCIATED(v7d_tmp%dativar%c)) THEN
8634! alloco i dati reali e vi trasferisco i character
8635 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatic, 1), SIZE(v7d_tmp%voldatic, 2), &
8636 SIZE(v7d_tmp%voldatic, 3), SIZE(v7d_tmp%voldatic, 4), SIZE(v7d_tmp%voldatic, 5), &
8637 SIZE(v7d_tmp%voldatic, 6)))
8638 DO i = 1, SIZE(v7d_tmp%dativar%c)
8639 v7d_tmp%voldatir(:,:,:,:,i,:) = &
8640 realdat(v7d_tmp%voldatic(:,:,:,:,i,:), v7d_tmp%dativar%c(i))
8641 ENDDO
8642 DEALLOCATE(v7d_tmp%voldatic)
8643! trasferisco le variabili
8644 v7d_tmp%dativar%r => v7d_tmp%dativar%c
8645 NULLIFY(v7d_tmp%dativar%c)
8646 ENDIF
8647
8648! fondo con il volume definitivo
8649 CALL vol7d_merge(that, v7d_tmp)
8650ELSE
8651 CALL delete(v7d_tmp)
8652ENDIF
8653
8654END SUBROUTINE vol7d_convr
8655
8656
8660SUBROUTINE vol7d_diff_only (this, that, data_only,ana)
8661TYPE(vol7d),INTENT(IN) :: this
8662TYPE(vol7d),INTENT(OUT) :: that
8663logical , optional, intent(in) :: data_only
8664logical , optional, intent(in) :: ana
8665logical :: ldata_only,lana
8666
8667IF (PRESENT(data_only)) THEN
8668 ldata_only = data_only
8669ELSE
8670 ldata_only = .false.
8671ENDIF
8672
8673IF (PRESENT(ana)) THEN
8674 lana = ana
8675ELSE
8676 lana = .false.
8677ENDIF
8678
8679
8680#undef VOL7D_POLY_ARRAY
8681#define VOL7D_POLY_ARRAY voldati
8682#include "vol7d_class_diff.F90"
8683#undef VOL7D_POLY_ARRAY
8684#define VOL7D_POLY_ARRAY voldatiattr
8685#include "vol7d_class_diff.F90"
8686#undef VOL7D_POLY_ARRAY
8687
8688if ( .not. ldata_only) then
8689
8690#define VOL7D_POLY_ARRAY volana
8691#include "vol7d_class_diff.F90"
8692#undef VOL7D_POLY_ARRAY
8693#define VOL7D_POLY_ARRAY volanaattr
8694#include "vol7d_class_diff.F90"
8695#undef VOL7D_POLY_ARRAY
8696
8697 if(lana)then
8698 where ( this%ana == that%ana )
8699 that%ana = vol7d_ana_miss
8700 end where
8701 end if
8702
8703end if
8704
8705
8706
8707END SUBROUTINE vol7d_diff_only
8708
8709
8710
8711! Creo le routine da ripetere per i vari tipi di dati di v7d
8712! tramite un template e il preprocessore
8713#undef VOL7D_POLY_TYPE
8714#undef VOL7D_POLY_TYPES
8715#define VOL7D_POLY_TYPE REAL
8716#define VOL7D_POLY_TYPES r
8717#include "vol7d_class_type_templ.F90"
8718#undef VOL7D_POLY_TYPE
8719#undef VOL7D_POLY_TYPES
8720#define VOL7D_POLY_TYPE DOUBLE PRECISION
8721#define VOL7D_POLY_TYPES d
8722#include "vol7d_class_type_templ.F90"
8723#undef VOL7D_POLY_TYPE
8724#undef VOL7D_POLY_TYPES
8725#define VOL7D_POLY_TYPE INTEGER
8726#define VOL7D_POLY_TYPES i
8727#include "vol7d_class_type_templ.F90"
8728#undef VOL7D_POLY_TYPE
8729#undef VOL7D_POLY_TYPES
8730#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
8731#define VOL7D_POLY_TYPES b
8732#include "vol7d_class_type_templ.F90"
8733#undef VOL7D_POLY_TYPE
8734#undef VOL7D_POLY_TYPES
8735#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
8736#define VOL7D_POLY_TYPES c
8737#include "vol7d_class_type_templ.F90"
8738
8739! Creo le routine da ripetere per i vari descrittori di dimensioni di v7d
8740! tramite un template e il preprocessore
8741#define VOL7D_SORT
8742#undef VOL7D_NO_ZERO_ALLOC
8743#undef VOL7D_POLY_TYPE
8744#define VOL7D_POLY_TYPE datetime
8745#include "vol7d_class_desc_templ.F90"
8746#undef VOL7D_POLY_TYPE
8747#define VOL7D_POLY_TYPE vol7d_timerange
8748#include "vol7d_class_desc_templ.F90"
8749#undef VOL7D_POLY_TYPE
8750#define VOL7D_POLY_TYPE vol7d_level
8751#include "vol7d_class_desc_templ.F90"
8752#undef VOL7D_SORT
8753#undef VOL7D_POLY_TYPE
8754#define VOL7D_POLY_TYPE vol7d_network
8755#include "vol7d_class_desc_templ.F90"
8756#undef VOL7D_POLY_TYPE
8757#define VOL7D_POLY_TYPE vol7d_ana
8758#include "vol7d_class_desc_templ.F90"
8759#define VOL7D_NO_ZERO_ALLOC
8760#undef VOL7D_POLY_TYPE
8761#define VOL7D_POLY_TYPE vol7d_var
8762#include "vol7d_class_desc_templ.F90"
8763
8773subroutine vol7d_write_on_file (this,unit,description,filename,filename_auto)
8774
8775TYPE(vol7d),INTENT(IN) :: this
8776integer,optional,intent(inout) :: unit
8777character(len=*),intent(in),optional :: filename
8778character(len=*),intent(out),optional :: filename_auto
8779character(len=*),INTENT(IN),optional :: description
8780
8781integer :: lunit
8782character(len=254) :: ldescription,arg,lfilename
8783integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
8784 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
8785 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
8786 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
8787 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
8788 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
8789 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
8790!integer :: im,id,iy
8791integer :: tarray(8)
8792logical :: opened,exist
8793
8794 nana=0
8795 ntime=0
8796 ntimerange=0
8797 nlevel=0
8798 nnetwork=0
8799 ndativarr=0
8800 ndativari=0
8801 ndativarb=0
8802 ndativard=0
8803 ndativarc=0
8804 ndatiattrr=0
8805 ndatiattri=0
8806 ndatiattrb=0
8807 ndatiattrd=0
8808 ndatiattrc=0
8809 ndativarattrr=0
8810 ndativarattri=0
8811 ndativarattrb=0
8812 ndativarattrd=0
8813 ndativarattrc=0
8814 nanavarr=0
8815 nanavari=0
8816 nanavarb=0
8817 nanavard=0
8818 nanavarc=0
8819 nanaattrr=0
8820 nanaattri=0
8821 nanaattrb=0
8822 nanaattrd=0
8823 nanaattrc=0
8824 nanavarattrr=0
8825 nanavarattri=0
8826 nanavarattrb=0
8827 nanavarattrd=0
8828 nanavarattrc=0
8829
8830
8831!call idate(im,id,iy)
8832call date_and_time(values=tarray)
8833call getarg(0,arg)
8834
8835if (present(description))then
8836 ldescription=description
8837else
8838 ldescription="Vol7d generated by: "//trim(arg)
8839end if
8840
8841if (.not. present(unit))then
8842 lunit=getunit()
8843else
8844 if (unit==0)then
8845 lunit=getunit()
8846 unit=lunit
8847 else
8848 lunit=unit
8849 end if
8850end if
8851
8852lfilename=trim(arg)//".v7d"
8853if (index(arg,'/',back=.true.) > 0) lfilename=lfilename(index(arg,'/',back=.true.)+1 : )
8854
8855if (present(filename))then
8856 if (filename /= "")then
8857 lfilename=filename
8858 end if
8859end if
8860
8861if (present(filename_auto))filename_auto=lfilename
8862
8863
8864inquire(unit=lunit,opened=opened)
8865if (.not. opened) then
8866! inquire(file=lfilename, EXIST=exist)
8867! IF (exist) THEN
8868! CALL l4f_log(L4F_FATAL, &
8869! 'in vol7d_write_on_file, file exists, cannot open file '//TRIM(lfilename))
8870! CALL raise_fatal_error()
8871! ENDIF
8872 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM')
8873 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
8874end if
8875
8876if (associated(this%ana)) nana=size(this%ana)
8877if (associated(this%time)) ntime=size(this%time)
8878if (associated(this%timerange)) ntimerange=size(this%timerange)
8879if (associated(this%level)) nlevel=size(this%level)
8880if (associated(this%network)) nnetwork=size(this%network)
8881
8882if (associated(this%dativar%r)) ndativarr=size(this%dativar%r)
8883if (associated(this%dativar%i)) ndativari=size(this%dativar%i)
8884if (associated(this%dativar%b)) ndativarb=size(this%dativar%b)
8885if (associated(this%dativar%d)) ndativard=size(this%dativar%d)
8886if (associated(this%dativar%c)) ndativarc=size(this%dativar%c)
8887
8888if (associated(this%datiattr%r)) ndatiattrr=size(this%datiattr%r)
8889if (associated(this%datiattr%i)) ndatiattri=size(this%datiattr%i)
8890if (associated(this%datiattr%b)) ndatiattrb=size(this%datiattr%b)
8891if (associated(this%datiattr%d)) ndatiattrd=size(this%datiattr%d)
8892if (associated(this%datiattr%c)) ndatiattrc=size(this%datiattr%c)
8893
8894if (associated(this%dativarattr%r)) ndativarattrr=size(this%dativarattr%r)
8895if (associated(this%dativarattr%i)) ndativarattri=size(this%dativarattr%i)
8896if (associated(this%dativarattr%b)) ndativarattrb=size(this%dativarattr%b)
8897if (associated(this%dativarattr%d)) ndativarattrd=size(this%dativarattr%d)
8898if (associated(this%dativarattr%c)) ndativarattrc=size(this%dativarattr%c)
8899
8900if (associated(this%anavar%r)) nanavarr=size(this%anavar%r)
8901if (associated(this%anavar%i)) nanavari=size(this%anavar%i)
8902if (associated(this%anavar%b)) nanavarb=size(this%anavar%b)
8903if (associated(this%anavar%d)) nanavard=size(this%anavar%d)
8904if (associated(this%anavar%c)) nanavarc=size(this%anavar%c)
8905
8906if (associated(this%anaattr%r)) nanaattrr=size(this%anaattr%r)
8907if (associated(this%anaattr%i)) nanaattri=size(this%anaattr%i)
8908if (associated(this%anaattr%b)) nanaattrb=size(this%anaattr%b)
8909if (associated(this%anaattr%d)) nanaattrd=size(this%anaattr%d)
8910if (associated(this%anaattr%c)) nanaattrc=size(this%anaattr%c)
8911
8912if (associated(this%anavarattr%r)) nanavarattrr=size(this%anavarattr%r)
8913if (associated(this%anavarattr%i)) nanavarattri=size(this%anavarattr%i)
8914if (associated(this%anavarattr%b)) nanavarattrb=size(this%anavarattr%b)
8915if (associated(this%anavarattr%d)) nanavarattrd=size(this%anavarattr%d)
8916if (associated(this%anavarattr%c)) nanavarattrc=size(this%anavarattr%c)
8917
8918write(unit=lunit)ldescription
8919write(unit=lunit)tarray
8920
8921write(unit=lunit)&
8922 nana, ntime, ntimerange, nlevel, nnetwork, &
8923 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
8924 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
8925 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
8926 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
8927 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
8928 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
8929 this%time_definition
8930
8931
8932!write(unit=lunit)this
8933
8934
8935!! prime 5 dimensioni
8936if (associated(this%ana)) call write_unit(this%ana, lunit)
8937if (associated(this%time)) call write_unit(this%time, lunit)
8938if (associated(this%level)) write(unit=lunit)this%level
8939if (associated(this%timerange)) write(unit=lunit)this%timerange
8940if (associated(this%network)) write(unit=lunit)this%network
8941
8942 !! 6a dimensione: variabile dell'anagrafica e dei dati
8943 !! con relativi attributi e in 5 tipi diversi
8944
8945if (associated(this%anavar%r)) write(unit=lunit)this%anavar%r
8946if (associated(this%anavar%i)) write(unit=lunit)this%anavar%i
8947if (associated(this%anavar%b)) write(unit=lunit)this%anavar%b
8948if (associated(this%anavar%d)) write(unit=lunit)this%anavar%d
8949if (associated(this%anavar%c)) write(unit=lunit)this%anavar%c
8950
8951if (associated(this%anaattr%r)) write(unit=lunit)this%anaattr%r
8952if (associated(this%anaattr%i)) write(unit=lunit)this%anaattr%i
8953if (associated(this%anaattr%b)) write(unit=lunit)this%anaattr%b
8954if (associated(this%anaattr%d)) write(unit=lunit)this%anaattr%d
8955if (associated(this%anaattr%c)) write(unit=lunit)this%anaattr%c
8956
8957if (associated(this%anavarattr%r)) write(unit=lunit)this%anavarattr%r
8958if (associated(this%anavarattr%i)) write(unit=lunit)this%anavarattr%i
8959if (associated(this%anavarattr%b)) write(unit=lunit)this%anavarattr%b
8960if (associated(this%anavarattr%d)) write(unit=lunit)this%anavarattr%d
8961if (associated(this%anavarattr%c)) write(unit=lunit)this%anavarattr%c
8962
8963if (associated(this%dativar%r)) write(unit=lunit)this%dativar%r
8964if (associated(this%dativar%i)) write(unit=lunit)this%dativar%i
8965if (associated(this%dativar%b)) write(unit=lunit)this%dativar%b
8966if (associated(this%dativar%d)) write(unit=lunit)this%dativar%d
8967if (associated(this%dativar%c)) write(unit=lunit)this%dativar%c
8968
8969if (associated(this%datiattr%r)) write(unit=lunit)this%datiattr%r
8970if (associated(this%datiattr%i)) write(unit=lunit)this%datiattr%i
8971if (associated(this%datiattr%b)) write(unit=lunit)this%datiattr%b
8972if (associated(this%datiattr%d)) write(unit=lunit)this%datiattr%d
8973if (associated(this%datiattr%c)) write(unit=lunit)this%datiattr%c
8974
8975if (associated(this%dativarattr%r)) write(unit=lunit)this%dativarattr%r
8976if (associated(this%dativarattr%i)) write(unit=lunit)this%dativarattr%i
8977if (associated(this%dativarattr%b)) write(unit=lunit)this%dativarattr%b
8978if (associated(this%dativarattr%d)) write(unit=lunit)this%dativarattr%d
8979if (associated(this%dativarattr%c)) write(unit=lunit)this%dativarattr%c
8980
8981!! Volumi di valori e attributi per anagrafica e dati
8982
8983if (associated(this%volanar)) write(unit=lunit)this%volanar
8984if (associated(this%volanaattrr)) write(unit=lunit)this%volanaattrr
8985if (associated(this%voldatir)) write(unit=lunit)this%voldatir
8986if (associated(this%voldatiattrr)) write(unit=lunit)this%voldatiattrr
8987
8988if (associated(this%volanai)) write(unit=lunit)this%volanai
8989if (associated(this%volanaattri)) write(unit=lunit)this%volanaattri
8990if (associated(this%voldatii)) write(unit=lunit)this%voldatii
8991if (associated(this%voldatiattri)) write(unit=lunit)this%voldatiattri
8992
8993if (associated(this%volanab)) write(unit=lunit)this%volanab
8994if (associated(this%volanaattrb)) write(unit=lunit)this%volanaattrb
8995if (associated(this%voldatib)) write(unit=lunit)this%voldatib
8996if (associated(this%voldatiattrb)) write(unit=lunit)this%voldatiattrb
8997
8998if (associated(this%volanad)) write(unit=lunit)this%volanad
8999if (associated(this%volanaattrd)) write(unit=lunit)this%volanaattrd
9000if (associated(this%voldatid)) write(unit=lunit)this%voldatid
9001if (associated(this%voldatiattrd)) write(unit=lunit)this%voldatiattrd
9002
9003if (associated(this%volanac)) write(unit=lunit)this%volanac
9004if (associated(this%volanaattrc)) write(unit=lunit)this%volanaattrc
9005if (associated(this%voldatic)) write(unit=lunit)this%voldatic
9006if (associated(this%voldatiattrc)) write(unit=lunit)this%voldatiattrc
9007
9008if (.not. present(unit)) close(unit=lunit)
9009
9010end subroutine vol7d_write_on_file
9011
9012
9019
9020
9021subroutine vol7d_read_from_file (this,unit,filename,description,tarray,filename_auto)
9022
9023TYPE(vol7d),INTENT(OUT) :: this
9024integer,intent(inout),optional :: unit
9025character(len=*),INTENT(in),optional :: filename
9026character(len=*),intent(out),optional :: filename_auto
9027character(len=*),INTENT(out),optional :: description
9028integer,intent(out),optional :: tarray(8)
9029
9030
9031integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
9032 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
9033 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
9034 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
9035 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
9036 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
9037 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
9038
9039character(len=254) :: ldescription,lfilename,arg
9040integer :: ltarray(8),lunit,ios
9041logical :: opened,exist
9042
9043
9044call getarg(0,arg)
9045
9046if (.not. present(unit))then
9047 lunit=getunit()
9048else
9049 if (unit==0)then
9050 lunit=getunit()
9051 unit=lunit
9052 else
9053 lunit=unit
9054 end if
9055end if
9056
9057lfilename=trim(arg)//".v7d"
9058if (index(arg,'/',back=.true.) > 0) lfilename=lfilename(index(arg,'/',back=.true.)+1 : )
9059
9060if (present(filename))then
9061 if (filename /= "")then
9062 lfilename=filename
9063 end if
9064end if
9065
9066if (present(filename_auto))filename_auto=lfilename
9067
9068
9069inquire(unit=lunit,opened=opened)
9070IF (.NOT. opened) THEN
9071 inquire(file=lfilename,exist=exist)
9072 IF (.NOT.exist) THEN
9073 CALL l4f_log(l4f_fatal, &
9074 'in vol7d_read_from_file, file does not exists, cannot open')
9075 CALL raise_fatal_error()
9076 ENDIF
9077 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM', &
9078 status='OLD', action='READ')
9079 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
9080end if
9081
9082
9083call init(this)
9084read(unit=lunit,iostat=ios)ldescription
9085
9086if (ios < 0) then ! A negative value indicates that the End of File or End of Record
9087 call vol7d_alloc (this)
9088 call vol7d_alloc_vol (this)
9089 if (present(description))description=ldescription
9090 if (present(tarray))tarray=ltarray
9091 if (.not. present(unit)) close(unit=lunit)
9092end if
9093
9094read(unit=lunit)ltarray
9095
9096CALL l4f_log(l4f_info, 'Reading vol7d from file')
9097CALL l4f_log(l4f_info, 'description: '//trim(ldescription))
9098CALL l4f_log(l4f_info, 'written on '//trim(to_char(ltarray(1)))//' '// &
9099 trim(to_char(ltarray(2)))//' '//trim(to_char(ltarray(3))))
9100
9101if (present(description))description=ldescription
9102if (present(tarray))tarray=ltarray
9103
9104read(unit=lunit)&
9105 nana, ntime, ntimerange, nlevel, nnetwork, &
9106 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
9107 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
9108 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
9109 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
9110 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
9111 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
9112 this%time_definition
9113
9114call vol7d_alloc (this, &
9115 nana=nana, ntime=ntime, ntimerange=ntimerange, nlevel=nlevel, nnetwork=nnetwork,&
9116 ndativarr=ndativarr, ndativari=ndativari, ndativarb=ndativarb,&
9117 ndativard=ndativard, ndativarc=ndativarc,&
9118 ndatiattrr=ndatiattrr, ndatiattri=ndatiattri, ndatiattrb=ndatiattrb,&
9119 ndatiattrd=ndatiattrd, ndatiattrc=ndatiattrc,&
9120 ndativarattrr=ndativarattrr, ndativarattri=ndativarattri, ndativarattrb=ndativarattrb, &
9121 ndativarattrd=ndativarattrd, ndativarattrc=ndativarattrc,&
9122 nanavarr=nanavarr, nanavari=nanavari, nanavarb=nanavarb, &
9123 nanavard=nanavard, nanavarc=nanavarc,&
9124 nanaattrr=nanaattrr, nanaattri=nanaattri, nanaattrb=nanaattrb,&
9125 nanaattrd=nanaattrd, nanaattrc=nanaattrc,&
9126 nanavarattrr=nanavarattrr, nanavarattri=nanavarattri, nanavarattrb=nanavarattrb, &
9127 nanavarattrd=nanavarattrd, nanavarattrc=nanavarattrc)
9128
9129
9130if (associated(this%ana)) call read_unit(this%ana, lunit)
9131if (associated(this%time)) call read_unit(this%time, lunit)
9132if (associated(this%level)) read(unit=lunit)this%level
9133if (associated(this%timerange)) read(unit=lunit)this%timerange
9134if (associated(this%network)) read(unit=lunit)this%network
9135
9136if (associated(this%anavar%r)) read(unit=lunit)this%anavar%r
9137if (associated(this%anavar%i)) read(unit=lunit)this%anavar%i
9138if (associated(this%anavar%b)) read(unit=lunit)this%anavar%b
9139if (associated(this%anavar%d)) read(unit=lunit)this%anavar%d
9140if (associated(this%anavar%c)) read(unit=lunit)this%anavar%c
9141
9142if (associated(this%anaattr%r)) read(unit=lunit)this%anaattr%r
9143if (associated(this%anaattr%i)) read(unit=lunit)this%anaattr%i
9144if (associated(this%anaattr%b)) read(unit=lunit)this%anaattr%b
9145if (associated(this%anaattr%d)) read(unit=lunit)this%anaattr%d
9146if (associated(this%anaattr%c)) read(unit=lunit)this%anaattr%c
9147
9148if (associated(this%anavarattr%r)) read(unit=lunit)this%anavarattr%r
9149if (associated(this%anavarattr%i)) read(unit=lunit)this%anavarattr%i
9150if (associated(this%anavarattr%b)) read(unit=lunit)this%anavarattr%b
9151if (associated(this%anavarattr%d)) read(unit=lunit)this%anavarattr%d
9152if (associated(this%anavarattr%c)) read(unit=lunit)this%anavarattr%c
9153
9154if (associated(this%dativar%r)) read(unit=lunit)this%dativar%r
9155if (associated(this%dativar%i)) read(unit=lunit)this%dativar%i
9156if (associated(this%dativar%b)) read(unit=lunit)this%dativar%b
9157if (associated(this%dativar%d)) read(unit=lunit)this%dativar%d
9158if (associated(this%dativar%c)) read(unit=lunit)this%dativar%c
9159
9160if (associated(this%datiattr%r)) read(unit=lunit)this%datiattr%r
9161if (associated(this%datiattr%i)) read(unit=lunit)this%datiattr%i
9162if (associated(this%datiattr%b)) read(unit=lunit)this%datiattr%b
9163if (associated(this%datiattr%d)) read(unit=lunit)this%datiattr%d
9164if (associated(this%datiattr%c)) read(unit=lunit)this%datiattr%c
9165
9166if (associated(this%dativarattr%r)) read(unit=lunit)this%dativarattr%r
9167if (associated(this%dativarattr%i)) read(unit=lunit)this%dativarattr%i
9168if (associated(this%dativarattr%b)) read(unit=lunit)this%dativarattr%b
9169if (associated(this%dativarattr%d)) read(unit=lunit)this%dativarattr%d
9170if (associated(this%dativarattr%c)) read(unit=lunit)this%dativarattr%c
9171
9172call vol7d_alloc_vol (this)
9173
9174!! Volumi di valori e attributi per anagrafica e dati
9175
9176if (associated(this%volanar)) read(unit=lunit)this%volanar
9177if (associated(this%volanaattrr)) read(unit=lunit)this%volanaattrr
9178if (associated(this%voldatir)) read(unit=lunit)this%voldatir
9179if (associated(this%voldatiattrr)) read(unit=lunit)this%voldatiattrr
9180
9181if (associated(this%volanai)) read(unit=lunit)this%volanai
9182if (associated(this%volanaattri)) read(unit=lunit)this%volanaattri
9183if (associated(this%voldatii)) read(unit=lunit)this%voldatii
9184if (associated(this%voldatiattri)) read(unit=lunit)this%voldatiattri
9185
9186if (associated(this%volanab)) read(unit=lunit)this%volanab
9187if (associated(this%volanaattrb)) read(unit=lunit)this%volanaattrb
9188if (associated(this%voldatib)) read(unit=lunit)this%voldatib
9189if (associated(this%voldatiattrb)) read(unit=lunit)this%voldatiattrb
9190
9191if (associated(this%volanad)) read(unit=lunit)this%volanad
9192if (associated(this%volanaattrd)) read(unit=lunit)this%volanaattrd
9193if (associated(this%voldatid)) read(unit=lunit)this%voldatid
9194if (associated(this%voldatiattrd)) read(unit=lunit)this%voldatiattrd
9195
9196if (associated(this%volanac)) read(unit=lunit)this%volanac
9197if (associated(this%volanaattrc)) read(unit=lunit)this%volanaattrc
9198if (associated(this%voldatic)) read(unit=lunit)this%voldatic
9199if (associated(this%voldatiattrc)) read(unit=lunit)this%voldatiattrc
9200
9201if (.not. present(unit)) close(unit=lunit)
9202
9203end subroutine vol7d_read_from_file
9204
9205
9206! to double precision
9207elemental doubleprecision function doubledatd(voldat,var)
9208doubleprecision,intent(in) :: voldat
9209type(vol7d_var),intent(in) :: var
9210
9211doubledatd=voldat
9212
9213end function doubledatd
9214
9215
9216elemental doubleprecision function doubledatr(voldat,var)
9217real,intent(in) :: voldat
9218type(vol7d_var),intent(in) :: var
9219
9220if (c_e(voldat))then
9221 doubledatr=dble(voldat)
9222else
9223 doubledatr=dmiss
9224end if
9225
9226end function doubledatr
9227
9228
9229elemental doubleprecision function doubledati(voldat,var)
9230integer,intent(in) :: voldat
9231type(vol7d_var),intent(in) :: var
9232
9233if (c_e(voldat)) then
9234 if (c_e(var%scalefactor))then
9235 doubledati=dble(voldat)/10.d0**var%scalefactor
9236 else
9237 doubledati=dble(voldat)
9238 endif
9239else
9240 doubledati=dmiss
9241end if
9242
9243end function doubledati
9244
9245
9246elemental doubleprecision function doubledatb(voldat,var)
9247integer(kind=int_b),intent(in) :: voldat
9248type(vol7d_var),intent(in) :: var
9249
9250if (c_e(voldat)) then
9251 if (c_e(var%scalefactor))then
9252 doubledatb=dble(voldat)/10.d0**var%scalefactor
9253 else
9254 doubledatb=dble(voldat)
9255 endif
9256else
9257 doubledatb=dmiss
9258end if
9259
9260end function doubledatb
9261
9262
9263elemental doubleprecision function doubledatc(voldat,var)
9264CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
9265type(vol7d_var),intent(in) :: var
9266
9267doubledatc = c2d(voldat)
9268if (c_e(doubledatc) .and. c_e(var%scalefactor))then
9269 doubledatc=doubledatc/10.d0**var%scalefactor
9270end if
9271
9272end function doubledatc
9273
9274
9275! to integer
9276elemental integer function integerdatd(voldat,var)
9277doubleprecision,intent(in) :: voldat
9278type(vol7d_var),intent(in) :: var
9279
9280if (c_e(voldat))then
9281 if (c_e(var%scalefactor)) then
9282 integerdatd=nint(voldat*10d0**var%scalefactor)
9283 else
9284 integerdatd=nint(voldat)
9285 endif
9286else
9287 integerdatd=imiss
9288end if
9289
9290end function integerdatd
9291
9292
9293elemental integer function integerdatr(voldat,var)
9294real,intent(in) :: voldat
9295type(vol7d_var),intent(in) :: var
9296
9297if (c_e(voldat))then
9298 if (c_e(var%scalefactor)) then
9299 integerdatr=nint(voldat*10d0**var%scalefactor)
9300 else
9301 integerdatr=nint(voldat)
9302 endif
9303else
9304 integerdatr=imiss
9305end if
9306
9307end function integerdatr
9308
9309
9310elemental integer function integerdati(voldat,var)
9311integer,intent(in) :: voldat
9312type(vol7d_var),intent(in) :: var
9313
9314integerdati=voldat
9315
9316end function integerdati
9317
9318
9319elemental integer function integerdatb(voldat,var)
9320integer(kind=int_b),intent(in) :: voldat
9321type(vol7d_var),intent(in) :: var
9322
9323if (c_e(voldat))then
9324 integerdatb=voldat
9325else
9326 integerdatb=imiss
9327end if
9328
9329end function integerdatb
9330
9331
9332elemental integer function integerdatc(voldat,var)
9333CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
9334type(vol7d_var),intent(in) :: var
9335
9336integerdatc=c2i(voldat)
9337
9338end function integerdatc
9339
9340
9341! to real
9342elemental real function realdatd(voldat,var)
9343doubleprecision,intent(in) :: voldat
9344type(vol7d_var),intent(in) :: var
9345
9346if (c_e(voldat))then
9347 realdatd=real(voldat)
9348else
9349 realdatd=rmiss
9350end if
9351
9352end function realdatd
9353
9354
9355elemental real function realdatr(voldat,var)
9356real,intent(in) :: voldat
9357type(vol7d_var),intent(in) :: var
9358
9359realdatr=voldat
9360
9361end function realdatr
9362
9363
9364elemental real function realdati(voldat,var)
9365integer,intent(in) :: voldat
9366type(vol7d_var),intent(in) :: var
9367
9368if (c_e(voldat)) then
9369 if (c_e(var%scalefactor))then
9370 realdati=float(voldat)/10.**var%scalefactor
9371 else
9372 realdati=float(voldat)
9373 endif
9374else
9375 realdati=rmiss
9376end if
9377
9378end function realdati
9379
9380
9381elemental real function realdatb(voldat,var)
9382integer(kind=int_b),intent(in) :: voldat
9383type(vol7d_var),intent(in) :: var
9384
9385if (c_e(voldat)) then
9386 if (c_e(var%scalefactor))then
9387 realdatb=float(voldat)/10**var%scalefactor
9388 else
9389 realdatb=float(voldat)
9390 endif
9391else
9392 realdatb=rmiss
9393end if
9394
9395end function realdatb
9396
9397
9398elemental real function realdatc(voldat,var)
9399CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
9400type(vol7d_var),intent(in) :: var
9401
9402realdatc=c2r(voldat)
9403if (c_e(realdatc) .and. c_e(var%scalefactor))then
9404 realdatc=realdatc/10.**var%scalefactor
9405end if
9406
9407end function realdatc
9408
9409
9415FUNCTION realanavol(this, var) RESULT(vol)
9416TYPE(vol7d),INTENT(in) :: this
9417TYPE(vol7d_var),INTENT(in) :: var
9418REAL :: vol(SIZE(this%ana),size(this%network))
9419
9420CHARACTER(len=1) :: dtype
9421INTEGER :: indvar
9422
9423dtype = cmiss
9424indvar = index(this%anavar, var, type=dtype)
9425
9426IF (indvar > 0) THEN
9427 SELECT CASE (dtype)
9428 CASE("d")
9429 vol = realdat(this%volanad(:,indvar,:), var)
9430 CASE("r")
9431 vol = this%volanar(:,indvar,:)
9432 CASE("i")
9433 vol = realdat(this%volanai(:,indvar,:), var)
9434 CASE("b")
9435 vol = realdat(this%volanab(:,indvar,:), var)
9436 CASE("c")
9437 vol = realdat(this%volanac(:,indvar,:), var)
9438 CASE default
9439 vol = rmiss
9440 END SELECT
9441ELSE
9442 vol = rmiss
9443ENDIF
9444
9445END FUNCTION realanavol
9446
9447
9453FUNCTION integeranavol(this, var) RESULT(vol)
9454TYPE(vol7d),INTENT(in) :: this
9455TYPE(vol7d_var),INTENT(in) :: var
9456INTEGER :: vol(SIZE(this%ana),size(this%network))
9457
9458CHARACTER(len=1) :: dtype
9459INTEGER :: indvar
9460
9461dtype = cmiss
9462indvar = index(this%anavar, var, type=dtype)
9463
9464IF (indvar > 0) THEN
9465 SELECT CASE (dtype)
9466 CASE("d")
9467 vol = integerdat(this%volanad(:,indvar,:), var)
9468 CASE("r")
9469 vol = integerdat(this%volanar(:,indvar,:), var)
9470 CASE("i")
9471 vol = this%volanai(:,indvar,:)
9472 CASE("b")
9473 vol = integerdat(this%volanab(:,indvar,:), var)
9474 CASE("c")
9475 vol = integerdat(this%volanac(:,indvar,:), var)
9476 CASE default
9477 vol = imiss
9478 END SELECT
9479ELSE
9480 vol = imiss
9481ENDIF
9482
9483END FUNCTION integeranavol
9484
9485
9491subroutine move_datac (v7d,&
9492 indana,indtime,indlevel,indtimerange,indnetwork,&
9493 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
9494
9495TYPE(vol7d),intent(inout) :: v7d
9496
9497integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
9498integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
9499integer :: inddativar,inddativarattr
9500
9501
9502do inddativar=1,size(v7d%dativar%c)
9503
9504 if (c_e(v7d%voldatic(indana,indtime,indlevel,indtimerange,inddativar,indnetwork)) .and. &
9505 .not. c_e(v7d%voldatic(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
9506 ) then
9507
9508 ! dati
9509 v7d%voldatic &
9510 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
9511 v7d%voldatic &
9512 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
9513
9514
9515 ! attributi
9516 if (associated (v7d%dativarattr%i)) then
9517 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%c(inddativar))
9518 if (inddativarattr > 0 ) then
9519 v7d%voldatiattri &
9520 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
9521 v7d%voldatiattri &
9522 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
9523 end if
9524 end if
9525
9526 if (associated (v7d%dativarattr%r)) then
9527 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%c(inddativar))
9528 if (inddativarattr > 0 ) then
9529 v7d%voldatiattrr &
9530 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
9531 v7d%voldatiattrr &
9532 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
9533 end if
9534 end if
9535
9536 if (associated (v7d%dativarattr%d)) then
9537 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%c(inddativar))
9538 if (inddativarattr > 0 ) then
9539 v7d%voldatiattrd &
9540 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
9541 v7d%voldatiattrd &
9542 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
9543 end if
9544 end if
9545
9546 if (associated (v7d%dativarattr%b)) then
9547 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%c(inddativar))
9548 if (inddativarattr > 0 ) then
9549 v7d%voldatiattrb &
9550 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
9551 v7d%voldatiattrb &
9552 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
9553 end if
9554 end if
9555
9556 if (associated (v7d%dativarattr%c)) then
9557 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%c(inddativar))
9558 if (inddativarattr > 0 ) then
9559 v7d%voldatiattrc &
9560 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
9561 v7d%voldatiattrc &
9562 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
9563 end if
9564 end if
9565
9566 end if
9567
9568end do
9569
9570end subroutine move_datac
9571
9577subroutine move_datar (v7d,&
9578 indana,indtime,indlevel,indtimerange,indnetwork,&
9579 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
9580
9581TYPE(vol7d),intent(inout) :: v7d
9582
9583integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
9584integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
9585integer :: inddativar,inddativarattr
9586
9587
9588do inddativar=1,size(v7d%dativar%r)
9589
9590 if (c_e(v7d%voldatir(indana,indtime,indlevel,indtimerange,inddativar,indnetwork)) .and. &
9591 .not. c_e(v7d%voldatir(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
9592 ) then
9593
9594 ! dati
9595 v7d%voldatir &
9596 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
9597 v7d%voldatir &
9598 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
9599
9600
9601 ! attributi
9602 if (associated (v7d%dativarattr%i)) then
9603 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%r(inddativar))
9604 if (inddativarattr > 0 ) then
9605 v7d%voldatiattri &
9606 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
9607 v7d%voldatiattri &
9608 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
9609 end if
9610 end if
9611
9612 if (associated (v7d%dativarattr%r)) then
9613 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%r(inddativar))
9614 if (inddativarattr > 0 ) then
9615 v7d%voldatiattrr &
9616 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
9617 v7d%voldatiattrr &
9618 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
9619 end if
9620 end if
9621
9622 if (associated (v7d%dativarattr%d)) then
9623 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%r(inddativar))
9624 if (inddativarattr > 0 ) then
9625 v7d%voldatiattrd &
9626 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
9627 v7d%voldatiattrd &
9628 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
9629 end if
9630 end if
9631
9632 if (associated (v7d%dativarattr%b)) then
9633 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%r(inddativar))
9634 if (inddativarattr > 0 ) then
9635 v7d%voldatiattrb &
9636 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
9637 v7d%voldatiattrb &
9638 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
9639 end if
9640 end if
9641
9642 if (associated (v7d%dativarattr%c)) then
9643 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%r(inddativar))
9644 if (inddativarattr > 0 ) then
9645 v7d%voldatiattrc &
9646 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
9647 v7d%voldatiattrc &
9648 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
9649 end if
9650 end if
9651
9652 end if
9653
9654end do
9655
9656end subroutine move_datar
9657
9658
9672subroutine v7d_rounding(v7din,v7dout,level,timerange,nostatproc)
9673type(vol7d),intent(inout) :: v7din
9674type(vol7d),intent(out) :: v7dout
9675type(vol7d_level),intent(in),optional :: level(:)
9676type(vol7d_timerange),intent(in),optional :: timerange(:)
9677!logical,intent(in),optional :: merge !< if there are data on more then one almost equal levels or timeranges
9678!! will be merged POINT BY POINT with priority for the fird data found ordered by icreasing var index
9679logical,intent(in),optional :: nostatproc
9680
9681integer :: nana,nlevel,ntime,ntimerange,nnetwork,nbin
9682integer :: iana,ilevel,itimerange,indl,indt,itime,inetwork
9683type(vol7d_level) :: roundlevel(size(v7din%level))
9684type(vol7d_timerange) :: roundtimerange(size(v7din%timerange))
9685type(vol7d) :: v7d_tmp
9686
9687
9688nbin=0
9689
9690if (associated(v7din%dativar%r)) nbin = nbin + size(v7din%dativar%r)
9691if (associated(v7din%dativar%i)) nbin = nbin + size(v7din%dativar%i)
9692if (associated(v7din%dativar%d)) nbin = nbin + size(v7din%dativar%d)
9693if (associated(v7din%dativar%b)) nbin = nbin + size(v7din%dativar%b)
9694
9695call init(v7d_tmp)
9696
9697roundlevel=v7din%level
9698
9699if (present(level))then
9700 do ilevel = 1, size(v7din%level)
9701 if ((any(v7din%level(ilevel) .almosteq. level))) then
9702 roundlevel(ilevel)=level(1)
9703 end if
9704 end do
9705end if
9706
9707roundtimerange=v7din%timerange
9708
9709if (present(timerange))then
9710 do itimerange = 1, size(v7din%timerange)
9711 if ((any(v7din%timerange(itimerange) .almosteq. timerange))) then
9712 roundtimerange(itimerange)=timerange(1)
9713 end if
9714 end do
9715end if
9716
9717!set istantaneous values everywere
9718!preserve p1 for forecast time
9719if (optio_log(nostatproc)) then
9720 roundtimerange(:)%timerange=254
9721 roundtimerange(:)%p2=0
9722end if
9723
9724
9725nana=size(v7din%ana)
9726nlevel=count_distinct(roundlevel,back=.true.)
9727ntime=size(v7din%time)
9728ntimerange=count_distinct(roundtimerange,back=.true.)
9729nnetwork=size(v7din%network)
9730
9731call init(v7d_tmp)
9732
9733if (nbin == 0) then
9734 call copy(v7din,v7d_tmp)
9735else
9736 call vol7d_convr(v7din,v7d_tmp)
9737end if
9738
9739v7d_tmp%level=roundlevel
9740v7d_tmp%timerange=roundtimerange
9741
9742do ilevel=1, size(v7d_tmp%level)
9743 indl=index(v7d_tmp%level,roundlevel(ilevel))
9744 do itimerange=1,size(v7d_tmp%timerange)
9745 indt=index(v7d_tmp%timerange,roundtimerange(itimerange))
9746
9747 if (indl /= ilevel .or. indt /= itimerange) then
9748
9749 do iana=1, nana
9750 do itime=1,ntime
9751 do inetwork=1,nnetwork
9752
9753 if (nbin > 0) then
9754 call move_datar (v7d_tmp,&
9755 iana,itime,ilevel,itimerange,inetwork,&
9756 iana,itime,indl,indt,inetwork)
9757 else
9758 call move_datac (v7d_tmp,&
9759 iana,itime,ilevel,itimerange,inetwork,&
9760 iana,itime,indl,indt,inetwork)
9761 end if
9762
9763 end do
9764 end do
9765 end do
9766
9767 end if
9768
9769 end do
9770end do
9771
9772! set to missing level and time > nlevel
9773do ilevel=nlevel+1,size(v7d_tmp%level)
9774 call init (v7d_tmp%level(ilevel))
9775end do
9776
9777do itimerange=ntimerange+1,size(v7d_tmp%timerange)
9778 call init (v7d_tmp%timerange(itimerange))
9779end do
9780
9781!copy with remove
9782CALL copy(v7d_tmp,v7dout,miss=.true.,lsort_timerange=.true.,lsort_level=.true.)
9783CALL delete(v7d_tmp)
9784
9785!call display(v7dout)
9786
9787end subroutine v7d_rounding
9788
9789
9790END MODULE vol7d_class
9791
9797
9798
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.