libsim Versione 7.2.4

◆ vol7d_get_volanac()

subroutine vol7d_get_volanac ( 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 )

Crea una vista a dimensione ridotta di un volume di anagrafica 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 :: vol1d(:)
...
CALL vol7d_get_volanac(v7d1, (/vol7d_ana_d/), vol1d)
IF (ASSOCIATED(vol1d)) THEN
print*,vol1d
...
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_a ... vol7d_attr_a, 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

Definizione alla linea 6244 del file vol7d_class.F90.

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