libsim Versione 7.2.1

◆ realanavol()

real function, dimension(size(this%ana),size(this%network)) realanavol ( type(vol7d), intent(in) this,
type(vol7d_var), intent(in) var )

Return an ana volume of a requested variable as real data.

It returns a 2-d array of the proper shape (ana x network) for the ana variable requested, converted to real type. If the conversion fails or if the variable is not contained in the ana volume, missing data are returned.

Parametri
[in]thisthe vol7d object to query, the method vol7d_alloc_vol must have been called for it otherwise progam may abort
[in]varthe ana variable to be returned

Definizione alla linea 9110 del file vol7d_class.F90.

9111! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
9112! authors:
9113! Davide Cesari <dcesari@arpa.emr.it>
9114! Paolo Patruno <ppatruno@arpa.emr.it>
9115
9116! This program is free software; you can redistribute it and/or
9117! modify it under the terms of the GNU General Public License as
9118! published by the Free Software Foundation; either version 2 of
9119! the License, or (at your option) any later version.
9120
9121! This program is distributed in the hope that it will be useful,
9122! but WITHOUT ANY WARRANTY; without even the implied warranty of
9123! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
9124! GNU General Public License for more details.
9125
9126! You should have received a copy of the GNU General Public License
9127! along with this program. If not, see <http://www.gnu.org/licenses/>.
9128#include "config.h"
9129
9141
9195MODULE vol7d_class
9196USE kinds
9200USE log4fortran
9201USE err_handling
9202USE io_units
9209IMPLICIT NONE
9210
9211
9212INTEGER, PARAMETER :: vol7d_maxdim_a = 3, vol7d_maxdim_aa = 4, &
9213 vol7d_maxdim_d = 6, vol7d_maxdim_ad = 7
9214
9215INTEGER, PARAMETER :: vol7d_ana_a=1
9216INTEGER, PARAMETER :: vol7d_var_a=2
9217INTEGER, PARAMETER :: vol7d_network_a=3
9218INTEGER, PARAMETER :: vol7d_attr_a=4
9219INTEGER, PARAMETER :: vol7d_ana_d=1
9220INTEGER, PARAMETER :: vol7d_time_d=2
9221INTEGER, PARAMETER :: vol7d_level_d=3
9222INTEGER, PARAMETER :: vol7d_timerange_d=4
9223INTEGER, PARAMETER :: vol7d_var_d=5
9224INTEGER, PARAMETER :: vol7d_network_d=6
9225INTEGER, PARAMETER :: vol7d_attr_d=7
9226INTEGER, PARAMETER :: vol7d_cdatalen=32
9227
9228TYPE vol7d_varmap
9229 INTEGER :: r, d, i, b, c
9230END TYPE vol7d_varmap
9231
9234TYPE vol7d
9236 TYPE(vol7d_ana),POINTER :: ana(:)
9238 TYPE(datetime),POINTER :: time(:)
9240 TYPE(vol7d_level),POINTER :: level(:)
9242 TYPE(vol7d_timerange),POINTER :: timerange(:)
9244 TYPE(vol7d_network),POINTER :: network(:)
9246 TYPE(vol7d_varvect) :: anavar
9248 TYPE(vol7d_varvect) :: anaattr
9250 TYPE(vol7d_varvect) :: anavarattr
9252 TYPE(vol7d_varvect) :: dativar
9254 TYPE(vol7d_varvect) :: datiattr
9256 TYPE(vol7d_varvect) :: dativarattr
9257
9259 REAL,POINTER :: volanar(:,:,:)
9261 DOUBLE PRECISION,POINTER :: volanad(:,:,:)
9263 INTEGER,POINTER :: volanai(:,:,:)
9265 INTEGER(kind=int_b),POINTER :: volanab(:,:,:)
9267 CHARACTER(len=vol7d_cdatalen),POINTER :: volanac(:,:,:)
9268
9270 REAL,POINTER :: volanaattrr(:,:,:,:)
9272 DOUBLE PRECISION,POINTER :: volanaattrd(:,:,:,:)
9274 INTEGER,POINTER :: volanaattri(:,:,:,:)
9276 INTEGER(kind=int_b),POINTER :: volanaattrb(:,:,:,:)
9278 CHARACTER(len=vol7d_cdatalen),POINTER :: volanaattrc(:,:,:,:)
9279
9281 REAL,POINTER :: voldatir(:,:,:,:,:,:) ! sono i dati
9283 DOUBLE PRECISION,POINTER :: voldatid(:,:,:,:,:,:)
9285 INTEGER,POINTER :: voldatii(:,:,:,:,:,:)
9287 INTEGER(kind=int_b),POINTER :: voldatib(:,:,:,:,:,:)
9289 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatic(:,:,:,:,:,:)
9290
9292 REAL,POINTER :: voldatiattrr(:,:,:,:,:,:,:)
9294 DOUBLE PRECISION,POINTER :: voldatiattrd(:,:,:,:,:,:,:)
9296 INTEGER,POINTER :: voldatiattri(:,:,:,:,:,:,:)
9298 INTEGER(kind=int_b),POINTER :: voldatiattrb(:,:,:,:,:,:,:)
9300 CHARACTER(len=vol7d_cdatalen),POINTER :: voldatiattrc(:,:,:,:,:,:,:)
9301
9303 integer :: time_definition
9304
9305END TYPE vol7d
9306
9310INTERFACE init
9311 MODULE PROCEDURE vol7d_init
9312END INTERFACE
9313
9315INTERFACE delete
9316 MODULE PROCEDURE vol7d_delete
9317END INTERFACE
9318
9320INTERFACE export
9321 MODULE PROCEDURE vol7d_write_on_file
9322END INTERFACE
9323
9325INTERFACE import
9326 MODULE PROCEDURE vol7d_read_from_file
9327END INTERFACE
9328
9330INTERFACE display
9331 MODULE PROCEDURE vol7d_display, dat_display, dat_vect_display
9332END INTERFACE
9333
9335INTERFACE to_char
9336 MODULE PROCEDURE to_char_dat
9337END INTERFACE
9338
9340INTERFACE doubledat
9341 MODULE PROCEDURE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
9342END INTERFACE
9343
9345INTERFACE realdat
9346 MODULE PROCEDURE realdatd,realdatr,realdati,realdatb,realdatc
9347END INTERFACE
9348
9350INTERFACE integerdat
9351 MODULE PROCEDURE integerdatd,integerdatr,integerdati,integerdatb,integerdatc
9352END INTERFACE
9353
9355INTERFACE copy
9356 MODULE PROCEDURE vol7d_copy
9357END INTERFACE
9358
9360INTERFACE c_e
9361 MODULE PROCEDURE vol7d_c_e
9362END INTERFACE
9363
9367INTERFACE check
9368 MODULE PROCEDURE vol7d_check
9369END INTERFACE
9370
9384INTERFACE rounding
9385 MODULE PROCEDURE v7d_rounding
9386END INTERFACE
9387
9388!!$INTERFACE get_volana
9389!!$ MODULE PROCEDURE vol7d_get_volanar, vol7d_get_volanad, vol7d_get_volanai, &
9390!!$ vol7d_get_volanab, vol7d_get_volanac
9391!!$END INTERFACE
9392!!$
9393!!$INTERFACE get_voldati
9394!!$ MODULE PROCEDURE vol7d_get_voldatir, vol7d_get_voldatid, vol7d_get_voldatii, &
9395!!$ vol7d_get_voldatib, vol7d_get_voldatic
9396!!$END INTERFACE
9397!!$
9398!!$INTERFACE get_volanaattr
9399!!$ MODULE PROCEDURE vol7d_get_volanaattrr, vol7d_get_volanaattrd, &
9400!!$ vol7d_get_volanaattri, vol7d_get_volanaattrb, vol7d_get_volanaattrc
9401!!$END INTERFACE
9402!!$
9403!!$INTERFACE get_voldatiattr
9404!!$ MODULE PROCEDURE vol7d_get_voldatiattrr, vol7d_get_voldatiattrd, &
9405!!$ vol7d_get_voldatiattri, vol7d_get_voldatiattrb, vol7d_get_voldatiattrc
9406!!$END INTERFACE
9407
9408PRIVATE vol7d_get_volr, vol7d_get_vold, vol7d_get_voli, vol7d_get_volb, &
9409 vol7d_get_volc, &
9410 volptr1dr, volptr2dr, volptr3dr, volptr4dr, volptr5dr, volptr6dr, volptr7dr, &
9411 volptr1dd, volptr2dd, volptr3dd, volptr4dd, volptr5dd, volptr6dd, volptr7dd, &
9412 volptr1di, volptr2di, volptr3di, volptr4di, volptr5di, volptr6di, volptr7di, &
9413 volptr1db, volptr2db, volptr3db, volptr4db, volptr5db, volptr6db, volptr7db, &
9414 volptr1dc, volptr2dc, volptr3dc, volptr4dc, volptr5dc, volptr6dc, volptr7dc, &
9415 vol7d_nullifyr, vol7d_nullifyd, vol7d_nullifyi, vol7d_nullifyb, vol7d_nullifyc, &
9416 vol7d_init, vol7d_delete, vol7d_write_on_file, vol7d_read_from_file, &
9417 vol7d_check_alloc_ana, vol7d_force_alloc_ana, &
9418 vol7d_check_alloc_dati, vol7d_force_alloc_dati, vol7d_force_alloc, &
9419 vol7d_display, dat_display, dat_vect_display, &
9420 to_char_dat, vol7d_check
9421
9422PRIVATE doubledatd,doubledatr,doubledati,doubledatb,doubledatc
9423
9424PRIVATE vol7d_c_e
9425
9426CONTAINS
9427
9428
9433SUBROUTINE vol7d_init(this,time_definition)
9434TYPE(vol7d),intent(out) :: this
9435integer,INTENT(IN),OPTIONAL :: time_definition
9436
9437CALL init(this%anavar)
9438CALL init(this%anaattr)
9439CALL init(this%anavarattr)
9440CALL init(this%dativar)
9441CALL init(this%datiattr)
9442CALL init(this%dativarattr)
9443CALL vol7d_var_features_init() ! initialise var features table once
9444
9445NULLIFY(this%ana, this%time, this%level, this%timerange, this%network)
9446
9447NULLIFY(this%volanar, this%volanaattrr, this%voldatir, this%voldatiattrr)
9448NULLIFY(this%volanad, this%volanaattrd, this%voldatid, this%voldatiattrd)
9449NULLIFY(this%volanai, this%volanaattri, this%voldatii, this%voldatiattri)
9450NULLIFY(this%volanab, this%volanaattrb, this%voldatib, this%voldatiattrb)
9451NULLIFY(this%volanac, this%volanaattrc, this%voldatic, this%voldatiattrc)
9452
9453if(present(time_definition)) then
9454 this%time_definition=time_definition
9455else
9456 this%time_definition=1 !default to validity time
9457end if
9458
9459END SUBROUTINE vol7d_init
9460
9461
9465ELEMENTAL SUBROUTINE vol7d_delete(this, dataonly)
9466TYPE(vol7d),intent(inout) :: this
9467LOGICAL, INTENT(in), OPTIONAL :: dataonly
9468
9469
9470IF (.NOT. optio_log(dataonly)) THEN
9471 IF (ASSOCIATED(this%volanar)) DEALLOCATE(this%volanar)
9472 IF (ASSOCIATED(this%volanad)) DEALLOCATE(this%volanad)
9473 IF (ASSOCIATED(this%volanai)) DEALLOCATE(this%volanai)
9474 IF (ASSOCIATED(this%volanab)) DEALLOCATE(this%volanab)
9475 IF (ASSOCIATED(this%volanac)) DEALLOCATE(this%volanac)
9476 IF (ASSOCIATED(this%volanaattrr)) DEALLOCATE(this%volanaattrr)
9477 IF (ASSOCIATED(this%volanaattrd)) DEALLOCATE(this%volanaattrd)
9478 IF (ASSOCIATED(this%volanaattri)) DEALLOCATE(this%volanaattri)
9479 IF (ASSOCIATED(this%volanaattrb)) DEALLOCATE(this%volanaattrb)
9480 IF (ASSOCIATED(this%volanaattrc)) DEALLOCATE(this%volanaattrc)
9481ENDIF
9482IF (ASSOCIATED(this%voldatir)) DEALLOCATE(this%voldatir)
9483IF (ASSOCIATED(this%voldatid)) DEALLOCATE(this%voldatid)
9484IF (ASSOCIATED(this%voldatii)) DEALLOCATE(this%voldatii)
9485IF (ASSOCIATED(this%voldatib)) DEALLOCATE(this%voldatib)
9486IF (ASSOCIATED(this%voldatic)) DEALLOCATE(this%voldatic)
9487IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
9488IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
9489IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
9490IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
9491IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
9492
9493IF (.NOT. optio_log(dataonly)) THEN
9494 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
9495 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
9496ENDIF
9497IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
9498IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
9499IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
9500
9501IF (.NOT. optio_log(dataonly)) THEN
9502 CALL delete(this%anavar)
9503 CALL delete(this%anaattr)
9504 CALL delete(this%anavarattr)
9505ENDIF
9506CALL delete(this%dativar)
9507CALL delete(this%datiattr)
9508CALL delete(this%dativarattr)
9509
9510END SUBROUTINE vol7d_delete
9511
9512
9513
9514integer function vol7d_check(this)
9515TYPE(vol7d),intent(in) :: this
9516integer :: i,j,k,l,m,n
9517
9518vol7d_check=0
9519
9520if (associated(this%voldatii)) then
9521do i = 1,size(this%voldatii,1)
9522 do j = 1,size(this%voldatii,2)
9523 do k = 1,size(this%voldatii,3)
9524 do l = 1,size(this%voldatii,4)
9525 do m = 1,size(this%voldatii,5)
9526 do n = 1,size(this%voldatii,6)
9527 if (this%voldatii(i,j,k,l,m,n) /= this%voldatii(i,j,k,l,m,n) ) then
9528 CALL l4f_log(l4f_warn,"check: abnormal value at voldatii("&
9529 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
9530 vol7d_check=1
9531 end if
9532 end do
9533 end do
9534 end do
9535 end do
9536 end do
9537end do
9538end if
9539
9540
9541if (associated(this%voldatir)) then
9542do i = 1,size(this%voldatir,1)
9543 do j = 1,size(this%voldatir,2)
9544 do k = 1,size(this%voldatir,3)
9545 do l = 1,size(this%voldatir,4)
9546 do m = 1,size(this%voldatir,5)
9547 do n = 1,size(this%voldatir,6)
9548 if (this%voldatir(i,j,k,l,m,n) /= this%voldatir(i,j,k,l,m,n) ) then
9549 CALL l4f_log(l4f_warn,"check: abnormal value at voldatir("&
9550 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
9551 vol7d_check=2
9552 end if
9553 end do
9554 end do
9555 end do
9556 end do
9557 end do
9558end do
9559end if
9560
9561if (associated(this%voldatid)) then
9562do i = 1,size(this%voldatid,1)
9563 do j = 1,size(this%voldatid,2)
9564 do k = 1,size(this%voldatid,3)
9565 do l = 1,size(this%voldatid,4)
9566 do m = 1,size(this%voldatid,5)
9567 do n = 1,size(this%voldatid,6)
9568 if (this%voldatid(i,j,k,l,m,n) /= this%voldatid(i,j,k,l,m,n) ) then
9569 CALL l4f_log(l4f_warn,"check: abnormal value at voldatid("&
9570 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
9571 vol7d_check=3
9572 end if
9573 end do
9574 end do
9575 end do
9576 end do
9577 end do
9578end do
9579end if
9580
9581if (associated(this%voldatib)) then
9582do i = 1,size(this%voldatib,1)
9583 do j = 1,size(this%voldatib,2)
9584 do k = 1,size(this%voldatib,3)
9585 do l = 1,size(this%voldatib,4)
9586 do m = 1,size(this%voldatib,5)
9587 do n = 1,size(this%voldatib,6)
9588 if (this%voldatib(i,j,k,l,m,n) /= this%voldatib(i,j,k,l,m,n) ) then
9589 CALL l4f_log(l4f_warn,"check: abnormal value at voldatib("&
9590 //t2c(i)//","//t2c(j)//","//t2c(k)//","//t2c(l)//","//t2c(m)//","//t2c(n)//",)")
9591 vol7d_check=4
9592 end if
9593 end do
9594 end do
9595 end do
9596 end do
9597 end do
9598end do
9599end if
9600
9601end function vol7d_check
9602
9603
9604
9605!TODO da completare ! aborta se i volumi sono allocati a dimensione 0
9607SUBROUTINE vol7d_display(this)
9608TYPE(vol7d),intent(in) :: this
9609integer :: i
9610
9611REAL :: rdat
9612DOUBLE PRECISION :: ddat
9613INTEGER :: idat
9614INTEGER(kind=int_b) :: bdat
9615CHARACTER(len=vol7d_cdatalen) :: cdat
9616
9617
9618print*,"<<<<<<<<<<<<<<<<<<< vol7d object >>>>>>>>>>>>>>>>>>>>"
9619if (this%time_definition == 0) then
9620 print*,"TIME DEFINITION: time is reference time"
9621else if (this%time_definition == 1) then
9622 print*,"TIME DEFINITION: time is validity time"
9623else
9624 print*,"Time definition have a wrong walue:", this%time_definition
9625end if
9626
9627IF (ASSOCIATED(this%network))then
9628 print*,"---- network vector ----"
9629 print*,"elements=",size(this%network)
9630 do i=1, size(this%network)
9631 call display(this%network(i))
9632 end do
9633end IF
9634
9635IF (ASSOCIATED(this%ana))then
9636 print*,"---- ana vector ----"
9637 print*,"elements=",size(this%ana)
9638 do i=1, size(this%ana)
9639 call display(this%ana(i))
9640 end do
9641end IF
9642
9643IF (ASSOCIATED(this%time))then
9644 print*,"---- time vector ----"
9645 print*,"elements=",size(this%time)
9646 do i=1, size(this%time)
9647 call display(this%time(i))
9648 end do
9649end if
9650
9651IF (ASSOCIATED(this%level)) then
9652 print*,"---- level vector ----"
9653 print*,"elements=",size(this%level)
9654 do i =1,size(this%level)
9655 call display(this%level(i))
9656 end do
9657end if
9658
9659IF (ASSOCIATED(this%timerange))then
9660 print*,"---- timerange vector ----"
9661 print*,"elements=",size(this%timerange)
9662 do i =1,size(this%timerange)
9663 call display(this%timerange(i))
9664 end do
9665end if
9666
9667
9668print*,"---- ana vector ----"
9669print*,""
9670print*,"->>>>>>>>> anavar -"
9671call display(this%anavar)
9672print*,""
9673print*,"->>>>>>>>> anaattr -"
9674call display(this%anaattr)
9675print*,""
9676print*,"->>>>>>>>> anavarattr -"
9677call display(this%anavarattr)
9678
9679print*,"-- ana data section (first point) --"
9680
9681idat=imiss
9682rdat=rmiss
9683ddat=dmiss
9684bdat=ibmiss
9685cdat=cmiss
9686
9687!ntime = MIN(SIZE(this%time),nprint)
9688!ntimerange = MIN(SIZE(this%timerange),nprint)
9689!nlevel = MIN(SIZE(this%level),nprint)
9690!nnetwork = MIN(SIZE(this%network),nprint)
9691!nana = MIN(SIZE(this%ana),nprint)
9692
9693IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0) THEN
9694if (associated(this%volanai)) then
9695 do i=1,size(this%anavar%i)
9696 idat=this%volanai(1,i,1)
9697 if (associated(this%anavar%i)) call display(this%anavar%i(i),idat,rdat,ddat,bdat,cdat)
9698 end do
9699end if
9700idat=imiss
9701
9702if (associated(this%volanar)) then
9703 do i=1,size(this%anavar%r)
9704 rdat=this%volanar(1,i,1)
9705 if (associated(this%anavar%r)) call display(this%anavar%r(i),idat,rdat,ddat,bdat,cdat)
9706 end do
9707end if
9708rdat=rmiss
9709
9710if (associated(this%volanad)) then
9711 do i=1,size(this%anavar%d)
9712 ddat=this%volanad(1,i,1)
9713 if (associated(this%anavar%d)) call display(this%anavar%d(i),idat,rdat,ddat,bdat,cdat)
9714 end do
9715end if
9716ddat=dmiss
9717
9718if (associated(this%volanab)) then
9719 do i=1,size(this%anavar%b)
9720 bdat=this%volanab(1,i,1)
9721 if (associated(this%anavar%b)) call display(this%anavar%b(i),idat,rdat,ddat,bdat,cdat)
9722 end do
9723end if
9724bdat=ibmiss
9725
9726if (associated(this%volanac)) then
9727 do i=1,size(this%anavar%c)
9728 cdat=this%volanac(1,i,1)
9729 if (associated(this%anavar%c)) call display(this%anavar%c(i),idat,rdat,ddat,bdat,cdat)
9730 end do
9731end if
9732cdat=cmiss
9733ENDIF
9734
9735print*,"---- data vector ----"
9736print*,""
9737print*,"->>>>>>>>> dativar -"
9738call display(this%dativar)
9739print*,""
9740print*,"->>>>>>>>> datiattr -"
9741call display(this%datiattr)
9742print*,""
9743print*,"->>>>>>>>> dativarattr -"
9744call display(this%dativarattr)
9745
9746print*,"-- data data section (first point) --"
9747
9748idat=imiss
9749rdat=rmiss
9750ddat=dmiss
9751bdat=ibmiss
9752cdat=cmiss
9753
9754IF (SIZE(this%ana) > 0 .AND. SIZE(this%network) > 0 .AND. size(this%time) > 0 &
9755 .AND. size(this%level) > 0 .AND. size(this%timerange) > 0) THEN
9756if (associated(this%voldatii)) then
9757 do i=1,size(this%dativar%i)
9758 idat=this%voldatii(1,1,1,1,i,1)
9759 if (associated(this%dativar%i)) call display(this%dativar%i(i),idat,rdat,ddat,bdat,cdat)
9760 end do
9761end if
9762idat=imiss
9763
9764if (associated(this%voldatir)) then
9765 do i=1,size(this%dativar%r)
9766 rdat=this%voldatir(1,1,1,1,i,1)
9767 if (associated(this%dativar%r)) call display(this%dativar%r(i),idat,rdat,ddat,bdat,cdat)
9768 end do
9769end if
9770rdat=rmiss
9771
9772if (associated(this%voldatid)) then
9773 do i=1,size(this%dativar%d)
9774 ddat=this%voldatid(1,1,1,1,i,1)
9775 if (associated(this%dativar%d)) call display(this%dativar%d(i),idat,rdat,ddat,bdat,cdat)
9776 end do
9777end if
9778ddat=dmiss
9779
9780if (associated(this%voldatib)) then
9781 do i=1,size(this%dativar%b)
9782 bdat=this%voldatib(1,1,1,1,i,1)
9783 if (associated(this%dativar%b)) call display(this%dativar%b(i),idat,rdat,ddat,bdat,cdat)
9784 end do
9785end if
9786bdat=ibmiss
9787
9788if (associated(this%voldatic)) then
9789 do i=1,size(this%dativar%c)
9790 cdat=this%voldatic(1,1,1,1,i,1)
9791 if (associated(this%dativar%c)) call display(this%dativar%c(i),idat,rdat,ddat,bdat,cdat)
9792 end do
9793end if
9794cdat=cmiss
9795ENDIF
9796
9797print*,"<<<<<<<<<<<<<<<<<<< END vol7d object >>>>>>>>>>>>>>>>>>>>"
9798
9799END SUBROUTINE vol7d_display
9800
9801
9803SUBROUTINE dat_display(this,idat,rdat,ddat,bdat,cdat)
9804TYPE(vol7d_var),intent(in) :: this
9806REAL :: rdat
9808DOUBLE PRECISION :: ddat
9810INTEGER :: idat
9812INTEGER(kind=int_b) :: bdat
9814CHARACTER(len=*) :: cdat
9815
9816print *, to_char_dat(this,idat,rdat,ddat,bdat,cdat)
9817
9818end SUBROUTINE dat_display
9819
9821SUBROUTINE dat_vect_display(this,idat,rdat,ddat,bdat,cdat)
9822
9823TYPE(vol7d_var),intent(in) :: this(:)
9825REAL :: rdat(:)
9827DOUBLE PRECISION :: ddat(:)
9829INTEGER :: idat(:)
9831INTEGER(kind=int_b) :: bdat(:)
9833CHARACTER(len=*):: cdat(:)
9834
9835integer :: i
9836
9837do i =1,size(this)
9838 call display(this(i),idat(i),rdat(i),ddat(i),bdat(i),cdat(i))
9839end do
9840
9841end SUBROUTINE dat_vect_display
9842
9843
9844FUNCTION to_char_dat(this,idat,rdat,ddat,bdat,cdat)
9845#ifdef HAVE_DBALLE
9846USE dballef
9847#endif
9848TYPE(vol7d_var),INTENT(in) :: this
9850REAL :: rdat
9852DOUBLE PRECISION :: ddat
9854INTEGER :: idat
9856INTEGER(kind=int_b) :: bdat
9858CHARACTER(len=*) :: cdat
9859CHARACTER(len=80) :: to_char_dat
9860
9861CHARACTER(len=LEN(to_char_dat)) :: to_char_tmp
9862
9863
9864#ifdef HAVE_DBALLE
9865INTEGER :: handle, ier
9866
9867handle = 0
9868to_char_dat="VALUE: "
9869
9870if (c_e(idat)) to_char_dat=trim(to_char_dat)//" ;int> "//trim(to_char(idat))
9871if (c_e(rdat)) to_char_dat=trim(to_char_dat)//" ;real> "//trim(to_char(rdat))
9872if (c_e(ddat)) to_char_dat=trim(to_char_dat)//" ;double> "//trim(to_char(ddat))
9873if (c_e(bdat)) to_char_dat=trim(to_char_dat)//" ;byte> "//trim(to_char(bdat))
9874
9875if ( c_e(cdat))then
9876 ier = idba_messaggi(handle,"/dev/null", "w", "BUFR")
9877 ier = idba_spiegab(handle,this%btable,cdat,to_char_tmp)
9878 ier = idba_fatto(handle)
9879 to_char_dat=trim(to_char_dat)//" ;char> "//trim(to_char_tmp)
9880endif
9881
9882#else
9883
9884to_char_dat="VALUE: "
9885if (c_e(idat)) to_char_dat=trim(to_char_dat)//" ;int> "//trim(to_char(idat))
9886if (c_e(rdat)) to_char_dat=trim(to_char_dat)//" ;real> "//trim(to_char(rdat))
9887if (c_e(ddat)) to_char_dat=trim(to_char_dat)//" ;double> "//trim(to_char(ddat))
9888if (c_e(bdat)) to_char_dat=trim(to_char_dat)//" ;byte> "//trim(to_char(bdat))
9889if (c_e(cdat)) to_char_dat=trim(to_char_dat)//" ;char> "//trim(cdat)
9890
9891#endif
9892
9893END FUNCTION to_char_dat
9894
9895
9898FUNCTION vol7d_c_e(this) RESULT(c_e)
9899TYPE(vol7d), INTENT(in) :: this
9900
9901LOGICAL :: c_e
9902
9903c_e = ASSOCIATED(this%ana) .OR. ASSOCIATED(this%time) .OR. &
9904 ASSOCIATED(this%level) .OR. ASSOCIATED(this%timerange) .OR. &
9905 ASSOCIATED(this%network) .OR. &
9906 ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
9907 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
9908 ASSOCIATED(this%anavar%c) .OR. &
9909 ASSOCIATED(this%anaattr%r) .OR. ASSOCIATED(this%anaattr%d) .OR. &
9910 ASSOCIATED(this%anaattr%i) .OR. ASSOCIATED(this%anaattr%b) .OR. &
9911 ASSOCIATED(this%anaattr%c) .OR. &
9912 ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
9913 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
9914 ASSOCIATED(this%dativar%c) .OR. &
9915 ASSOCIATED(this%datiattr%r) .OR. ASSOCIATED(this%datiattr%d) .OR. &
9916 ASSOCIATED(this%datiattr%i) .OR. ASSOCIATED(this%datiattr%b) .OR. &
9917 ASSOCIATED(this%datiattr%c)
9918
9919END FUNCTION vol7d_c_e
9920
9921
9960SUBROUTINE vol7d_alloc(this, nana, ntime, nlevel, ntimerange, nnetwork, &
9961 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
9962 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
9963 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
9964 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
9965 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
9966 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc, &
9967 ini)
9968TYPE(vol7d),INTENT(inout) :: this
9969INTEGER,INTENT(in),OPTIONAL :: nana
9970INTEGER,INTENT(in),OPTIONAL :: ntime
9971INTEGER,INTENT(in),OPTIONAL :: nlevel
9972INTEGER,INTENT(in),OPTIONAL :: ntimerange
9973INTEGER,INTENT(in),OPTIONAL :: nnetwork
9975INTEGER,INTENT(in),OPTIONAL :: &
9976 nanavarr, nanavard, nanavari, nanavarb, nanavarc, &
9977 nanaattrr, nanaattrd, nanaattri, nanaattrb, nanaattrc, &
9978 nanavarattrr, nanavarattrd, nanavarattri, nanavarattrb, nanavarattrc, &
9979 ndativarr, ndativard, ndativari, ndativarb, ndativarc, &
9980 ndatiattrr, ndatiattrd, ndatiattri, ndatiattrb, ndatiattrc, &
9981 ndativarattrr, ndativarattrd, ndativarattri, ndativarattrb, ndativarattrc
9982LOGICAL,INTENT(in),OPTIONAL :: ini
9983
9984INTEGER :: i
9985LOGICAL :: linit
9986
9987IF (PRESENT(ini)) THEN
9988 linit = ini
9989ELSE
9990 linit = .false.
9991ENDIF
9992
9993! Dimensioni principali
9994IF (PRESENT(nana)) THEN
9995 IF (nana >= 0) THEN
9996 IF (ASSOCIATED(this%ana)) DEALLOCATE(this%ana)
9997 ALLOCATE(this%ana(nana))
9998 IF (linit) THEN
9999 DO i = 1, nana
10000 CALL init(this%ana(i))
10001 ENDDO
10002 ENDIF
10003 ENDIF
10004ENDIF
10005IF (PRESENT(ntime)) THEN
10006 IF (ntime >= 0) THEN
10007 IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
10008 ALLOCATE(this%time(ntime))
10009 IF (linit) THEN
10010 DO i = 1, ntime
10011 CALL init(this%time(i))
10012 ENDDO
10013 ENDIF
10014 ENDIF
10015ENDIF
10016IF (PRESENT(nlevel)) THEN
10017 IF (nlevel >= 0) THEN
10018 IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
10019 ALLOCATE(this%level(nlevel))
10020 IF (linit) THEN
10021 DO i = 1, nlevel
10022 CALL init(this%level(i))
10023 ENDDO
10024 ENDIF
10025 ENDIF
10026ENDIF
10027IF (PRESENT(ntimerange)) THEN
10028 IF (ntimerange >= 0) THEN
10029 IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
10030 ALLOCATE(this%timerange(ntimerange))
10031 IF (linit) THEN
10032 DO i = 1, ntimerange
10033 CALL init(this%timerange(i))
10034 ENDDO
10035 ENDIF
10036 ENDIF
10037ENDIF
10038IF (PRESENT(nnetwork)) THEN
10039 IF (nnetwork >= 0) THEN
10040 IF (ASSOCIATED(this%network)) DEALLOCATE(this%network)
10041 ALLOCATE(this%network(nnetwork))
10042 IF (linit) THEN
10043 DO i = 1, nnetwork
10044 CALL init(this%network(i))
10045 ENDDO
10046 ENDIF
10047 ENDIF
10048ENDIF
10049! Dimensioni dei tipi delle variabili
10050CALL vol7d_varvect_alloc(this%anavar, nanavarr, nanavard, &
10051 nanavari, nanavarb, nanavarc, ini)
10052CALL vol7d_varvect_alloc(this%anaattr, nanaattrr, nanaattrd, &
10053 nanaattri, nanaattrb, nanaattrc, ini)
10054CALL vol7d_varvect_alloc(this%anavarattr, nanavarattrr, nanavarattrd, &
10055 nanavarattri, nanavarattrb, nanavarattrc, ini)
10056CALL vol7d_varvect_alloc(this%dativar, ndativarr, ndativard, &
10057 ndativari, ndativarb, ndativarc, ini)
10058CALL vol7d_varvect_alloc(this%datiattr, ndatiattrr, ndatiattrd, &
10059 ndatiattri, ndatiattrb, ndatiattrc, ini)
10060CALL vol7d_varvect_alloc(this%dativarattr, ndativarattrr, ndativarattrd, &
10061 ndativarattri, ndativarattrb, ndativarattrc, ini)
10062
10063END SUBROUTINE vol7d_alloc
10064
10065
10066FUNCTION vol7d_check_alloc_ana(this)
10067TYPE(vol7d),INTENT(in) :: this
10068LOGICAL :: vol7d_check_alloc_ana
10069
10070vol7d_check_alloc_ana = ASSOCIATED(this%ana) .AND. ASSOCIATED(this%network)
10071
10072END FUNCTION vol7d_check_alloc_ana
10073
10074SUBROUTINE vol7d_force_alloc_ana(this, ini)
10075TYPE(vol7d),INTENT(inout) :: this
10076LOGICAL,INTENT(in),OPTIONAL :: ini
10077
10078! Alloco i descrittori minimi per avere un volume di anagrafica
10079IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=1, ini=ini)
10080IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=1, ini=ini)
10081
10082END SUBROUTINE vol7d_force_alloc_ana
10083
10084
10085FUNCTION vol7d_check_alloc_dati(this)
10086TYPE(vol7d),INTENT(in) :: this
10087LOGICAL :: vol7d_check_alloc_dati
10088
10089vol7d_check_alloc_dati = vol7d_check_alloc_ana(this) .AND. &
10090 ASSOCIATED(this%time) .AND. ASSOCIATED(this%level) .AND. &
10091 ASSOCIATED(this%timerange)
10092
10093END FUNCTION vol7d_check_alloc_dati
10094
10095SUBROUTINE vol7d_force_alloc_dati(this, ini)
10096TYPE(vol7d),INTENT(inout) :: this
10097LOGICAL,INTENT(in),OPTIONAL :: ini
10098
10099! Alloco i descrittori minimi per avere un volume di dati
10100CALL vol7d_force_alloc_ana(this, ini)
10101IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=1, ini=ini)
10102IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=1, ini=ini)
10103IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=1, ini=ini)
10104
10105END SUBROUTINE vol7d_force_alloc_dati
10106
10107
10108SUBROUTINE vol7d_force_alloc(this)
10109TYPE(vol7d),INTENT(inout) :: this
10110
10111! If anything really not allocated yet, allocate with size 0
10112IF (.NOT. ASSOCIATED(this%ana)) CALL vol7d_alloc(this, nana=0)
10113IF (.NOT. ASSOCIATED(this%network)) CALL vol7d_alloc(this, nnetwork=0)
10114IF (.NOT. ASSOCIATED(this%time)) CALL vol7d_alloc(this, ntime=0)
10115IF (.NOT. ASSOCIATED(this%level)) CALL vol7d_alloc(this, nlevel=0)
10116IF (.NOT. ASSOCIATED(this%timerange)) CALL vol7d_alloc(this, ntimerange=0)
10117
10118END SUBROUTINE vol7d_force_alloc
10119
10120
10121FUNCTION vol7d_check_vol(this)
10122TYPE(vol7d),INTENT(in) :: this
10123LOGICAL :: vol7d_check_vol
10124
10125vol7d_check_vol = c_e(this)
10126
10127! Anagrafica
10128IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
10129 vol7d_check_vol = .false.
10130ENDIF
10131
10132IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
10133 vol7d_check_vol = .false.
10134ENDIF
10135
10136IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
10137 vol7d_check_vol = .false.
10138ENDIF
10139
10140IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
10141 vol7d_check_vol = .false.
10142ENDIF
10143
10144IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
10145 vol7d_check_vol = .false.
10146ENDIF
10147IF (ASSOCIATED(this%anavar%r) .OR. ASSOCIATED(this%anavar%d) .OR. &
10148 ASSOCIATED(this%anavar%i) .OR. ASSOCIATED(this%anavar%b) .OR. &
10149 ASSOCIATED(this%anavar%c)) THEN
10150 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_ana(this)
10151ENDIF
10152
10153! Attributi dell'anagrafica
10154IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
10155 .NOT.ASSOCIATED(this%volanaattrr)) THEN
10156 vol7d_check_vol = .false.
10157ENDIF
10158
10159IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
10160 .NOT.ASSOCIATED(this%volanaattrd)) THEN
10161 vol7d_check_vol = .false.
10162ENDIF
10163
10164IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
10165 .NOT.ASSOCIATED(this%volanaattri)) THEN
10166 vol7d_check_vol = .false.
10167ENDIF
10168
10169IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
10170 .NOT.ASSOCIATED(this%volanaattrb)) THEN
10171 vol7d_check_vol = .false.
10172ENDIF
10173
10174IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
10175 .NOT.ASSOCIATED(this%volanaattrc)) THEN
10176 vol7d_check_vol = .false.
10177ENDIF
10178
10179! Dati
10180IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
10181 vol7d_check_vol = .false.
10182ENDIF
10183
10184IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
10185 vol7d_check_vol = .false.
10186ENDIF
10187
10188IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
10189 vol7d_check_vol = .false.
10190ENDIF
10191
10192IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
10193 vol7d_check_vol = .false.
10194ENDIF
10195
10196IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
10197 vol7d_check_vol = .false.
10198ENDIF
10199
10200! Attributi dei dati
10201IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
10202 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
10203 vol7d_check_vol = .false.
10204ENDIF
10205
10206IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
10207 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
10208 vol7d_check_vol = .false.
10209ENDIF
10210
10211IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
10212 .NOT.ASSOCIATED(this%voldatiattri)) THEN
10213 vol7d_check_vol = .false.
10214ENDIF
10215
10216IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
10217 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
10218 vol7d_check_vol = .false.
10219ENDIF
10220
10221IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
10222 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
10223 vol7d_check_vol = .false.
10224ENDIF
10225IF (ASSOCIATED(this%dativar%r) .OR. ASSOCIATED(this%dativar%d) .OR. &
10226 ASSOCIATED(this%dativar%i) .OR. ASSOCIATED(this%dativar%b) .OR. &
10227 ASSOCIATED(this%dativar%c)) THEN
10228 vol7d_check_vol = vol7d_check_vol .AND. vol7d_check_alloc_dati(this)
10229ENDIF
10230
10231END FUNCTION vol7d_check_vol
10232
10233
10248SUBROUTINE vol7d_alloc_vol(this, ini, inivol)
10249TYPE(vol7d),INTENT(inout) :: this
10250LOGICAL,INTENT(in),OPTIONAL :: ini
10251LOGICAL,INTENT(in),OPTIONAL :: inivol
10252
10253LOGICAL :: linivol
10254
10255IF (PRESENT(inivol)) THEN
10256 linivol = inivol
10257ELSE
10258 linivol = .true.
10259ENDIF
10260
10261! Anagrafica
10262IF (ASSOCIATED(this%anavar%r) .AND. .NOT.ASSOCIATED(this%volanar)) THEN
10263 CALL vol7d_force_alloc_ana(this, ini)
10264 ALLOCATE(this%volanar(SIZE(this%ana), SIZE(this%anavar%r), SIZE(this%network)))
10265 IF (linivol) this%volanar(:,:,:) = rmiss
10266ENDIF
10267
10268IF (ASSOCIATED(this%anavar%d) .AND. .NOT.ASSOCIATED(this%volanad)) THEN
10269 CALL vol7d_force_alloc_ana(this, ini)
10270 ALLOCATE(this%volanad(SIZE(this%ana), SIZE(this%anavar%d), SIZE(this%network)))
10271 IF (linivol) this%volanad(:,:,:) = rdmiss
10272ENDIF
10273
10274IF (ASSOCIATED(this%anavar%i) .AND. .NOT.ASSOCIATED(this%volanai)) THEN
10275 CALL vol7d_force_alloc_ana(this, ini)
10276 ALLOCATE(this%volanai(SIZE(this%ana), SIZE(this%anavar%i), SIZE(this%network)))
10277 IF (linivol) this%volanai(:,:,:) = imiss
10278ENDIF
10279
10280IF (ASSOCIATED(this%anavar%b) .AND. .NOT.ASSOCIATED(this%volanab)) THEN
10281 CALL vol7d_force_alloc_ana(this, ini)
10282 ALLOCATE(this%volanab(SIZE(this%ana), SIZE(this%anavar%b), SIZE(this%network)))
10283 IF (linivol) this%volanab(:,:,:) = ibmiss
10284ENDIF
10285
10286IF (ASSOCIATED(this%anavar%c) .AND. .NOT.ASSOCIATED(this%volanac)) THEN
10287 CALL vol7d_force_alloc_ana(this, ini)
10288 ALLOCATE(this%volanac(SIZE(this%ana), SIZE(this%anavar%c), SIZE(this%network)))
10289 IF (linivol) this%volanac(:,:,:) = cmiss
10290ENDIF
10291
10292! Attributi dell'anagrafica
10293IF (ASSOCIATED(this%anaattr%r) .AND. ASSOCIATED(this%anavarattr%r) .AND. &
10294 .NOT.ASSOCIATED(this%volanaattrr)) THEN
10295 CALL vol7d_force_alloc_ana(this, ini)
10296 ALLOCATE(this%volanaattrr(SIZE(this%ana), SIZE(this%anavarattr%r), &
10297 SIZE(this%network), SIZE(this%anaattr%r)))
10298 IF (linivol) this%volanaattrr(:,:,:,:) = rmiss
10299ENDIF
10300
10301IF (ASSOCIATED(this%anaattr%d) .AND. ASSOCIATED(this%anavarattr%d) .AND. &
10302 .NOT.ASSOCIATED(this%volanaattrd)) THEN
10303 CALL vol7d_force_alloc_ana(this, ini)
10304 ALLOCATE(this%volanaattrd(SIZE(this%ana), SIZE(this%anavarattr%d), &
10305 SIZE(this%network), SIZE(this%anaattr%d)))
10306 IF (linivol) this%volanaattrd(:,:,:,:) = rdmiss
10307ENDIF
10308
10309IF (ASSOCIATED(this%anaattr%i) .AND. ASSOCIATED(this%anavarattr%i) .AND. &
10310 .NOT.ASSOCIATED(this%volanaattri)) THEN
10311 CALL vol7d_force_alloc_ana(this, ini)
10312 ALLOCATE(this%volanaattri(SIZE(this%ana), SIZE(this%anavarattr%i), &
10313 SIZE(this%network), SIZE(this%anaattr%i)))
10314 IF (linivol) this%volanaattri(:,:,:,:) = imiss
10315ENDIF
10316
10317IF (ASSOCIATED(this%anaattr%b) .AND. ASSOCIATED(this%anavarattr%b) .AND. &
10318 .NOT.ASSOCIATED(this%volanaattrb)) THEN
10319 CALL vol7d_force_alloc_ana(this, ini)
10320 ALLOCATE(this%volanaattrb(SIZE(this%ana), SIZE(this%anavarattr%b), &
10321 SIZE(this%network), SIZE(this%anaattr%b)))
10322 IF (linivol) this%volanaattrb(:,:,:,:) = ibmiss
10323ENDIF
10324
10325IF (ASSOCIATED(this%anaattr%c) .AND. ASSOCIATED(this%anavarattr%c) .AND. &
10326 .NOT.ASSOCIATED(this%volanaattrc)) THEN
10327 CALL vol7d_force_alloc_ana(this, ini)
10328 ALLOCATE(this%volanaattrc(SIZE(this%ana), SIZE(this%anavarattr%c), &
10329 SIZE(this%network), SIZE(this%anaattr%c)))
10330 IF (linivol) this%volanaattrc(:,:,:,:) = cmiss
10331ENDIF
10332
10333! Dati
10334IF (ASSOCIATED(this%dativar%r) .AND. .NOT.ASSOCIATED(this%voldatir)) THEN
10335 CALL vol7d_force_alloc_dati(this, ini)
10336 ALLOCATE(this%voldatir(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10337 SIZE(this%timerange), SIZE(this%dativar%r), SIZE(this%network)))
10338 IF (linivol) this%voldatir(:,:,:,:,:,:) = rmiss
10339ENDIF
10340
10341IF (ASSOCIATED(this%dativar%d) .AND. .NOT.ASSOCIATED(this%voldatid)) THEN
10342 CALL vol7d_force_alloc_dati(this, ini)
10343 ALLOCATE(this%voldatid(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10344 SIZE(this%timerange), SIZE(this%dativar%d), SIZE(this%network)))
10345 IF (linivol) this%voldatid(:,:,:,:,:,:) = rdmiss
10346ENDIF
10347
10348IF (ASSOCIATED(this%dativar%i) .AND. .NOT.ASSOCIATED(this%voldatii)) THEN
10349 CALL vol7d_force_alloc_dati(this, ini)
10350 ALLOCATE(this%voldatii(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10351 SIZE(this%timerange), SIZE(this%dativar%i), SIZE(this%network)))
10352 IF (linivol) this%voldatii(:,:,:,:,:,:) = imiss
10353ENDIF
10354
10355IF (ASSOCIATED(this%dativar%b) .AND. .NOT.ASSOCIATED(this%voldatib)) THEN
10356 CALL vol7d_force_alloc_dati(this, ini)
10357 ALLOCATE(this%voldatib(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10358 SIZE(this%timerange), SIZE(this%dativar%b), SIZE(this%network)))
10359 IF (linivol) this%voldatib(:,:,:,:,:,:) = ibmiss
10360ENDIF
10361
10362IF (ASSOCIATED(this%dativar%c) .AND. .NOT.ASSOCIATED(this%voldatic)) THEN
10363 CALL vol7d_force_alloc_dati(this, ini)
10364 ALLOCATE(this%voldatic(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10365 SIZE(this%timerange), SIZE(this%dativar%c), SIZE(this%network)))
10366 IF (linivol) this%voldatic(:,:,:,:,:,:) = cmiss
10367ENDIF
10368
10369! Attributi dei dati
10370IF (ASSOCIATED(this%datiattr%r) .AND. ASSOCIATED(this%dativarattr%r) .AND. &
10371 .NOT.ASSOCIATED(this%voldatiattrr)) THEN
10372 CALL vol7d_force_alloc_dati(this, ini)
10373 ALLOCATE(this%voldatiattrr(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10374 SIZE(this%timerange), SIZE(this%dativarattr%r), SIZE(this%network), &
10375 SIZE(this%datiattr%r)))
10376 IF (linivol) this%voldatiattrr(:,:,:,:,:,:,:) = rmiss
10377ENDIF
10378
10379IF (ASSOCIATED(this%datiattr%d) .AND. ASSOCIATED(this%dativarattr%d) .AND. &
10380 .NOT.ASSOCIATED(this%voldatiattrd)) THEN
10381 CALL vol7d_force_alloc_dati(this, ini)
10382 ALLOCATE(this%voldatiattrd(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10383 SIZE(this%timerange), SIZE(this%dativarattr%d), SIZE(this%network), &
10384 SIZE(this%datiattr%d)))
10385 IF (linivol) this%voldatiattrd(:,:,:,:,:,:,:) = rdmiss
10386ENDIF
10387
10388IF (ASSOCIATED(this%datiattr%i) .AND. ASSOCIATED(this%dativarattr%i) .AND. &
10389 .NOT.ASSOCIATED(this%voldatiattri)) THEN
10390 CALL vol7d_force_alloc_dati(this, ini)
10391 ALLOCATE(this%voldatiattri(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10392 SIZE(this%timerange), SIZE(this%dativarattr%i), SIZE(this%network), &
10393 SIZE(this%datiattr%i)))
10394 IF (linivol) this%voldatiattri(:,:,:,:,:,:,:) = imiss
10395ENDIF
10396
10397IF (ASSOCIATED(this%datiattr%b) .AND. ASSOCIATED(this%dativarattr%b) .AND. &
10398 .NOT.ASSOCIATED(this%voldatiattrb)) THEN
10399 CALL vol7d_force_alloc_dati(this, ini)
10400 ALLOCATE(this%voldatiattrb(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10401 SIZE(this%timerange), SIZE(this%dativarattr%b), SIZE(this%network), &
10402 SIZE(this%datiattr%b)))
10403 IF (linivol) this%voldatiattrb(:,:,:,:,:,:,:) = ibmiss
10404ENDIF
10405
10406IF (ASSOCIATED(this%datiattr%c) .AND. ASSOCIATED(this%dativarattr%c) .AND. &
10407 .NOT.ASSOCIATED(this%voldatiattrc)) THEN
10408 CALL vol7d_force_alloc_dati(this, ini)
10409 ALLOCATE(this%voldatiattrc(SIZE(this%ana), SIZE(this%time), SIZE(this%level), &
10410 SIZE(this%timerange), SIZE(this%dativarattr%c), SIZE(this%network), &
10411 SIZE(this%datiattr%c)))
10412 IF (linivol) this%voldatiattrc(:,:,:,:,:,:,:) = cmiss
10413ENDIF
10414
10415! Catch-all method
10416CALL vol7d_force_alloc(this)
10417
10418! Creo gli indici var-attr
10419
10420#ifdef DEBUG
10421CALL l4f_log(l4f_debug,"calling: vol7d_set_attr_ind")
10422#endif
10423
10424CALL vol7d_set_attr_ind(this)
10425
10426
10427
10428END SUBROUTINE vol7d_alloc_vol
10429
10430
10437SUBROUTINE vol7d_set_attr_ind(this)
10438TYPE(vol7d),INTENT(inout) :: this
10439
10440INTEGER :: i
10441
10442! real
10443IF (ASSOCIATED(this%dativar%r)) THEN
10444 IF (ASSOCIATED(this%dativarattr%r)) THEN
10445 DO i = 1, SIZE(this%dativar%r)
10446 this%dativar%r(i)%r = &
10447 firsttrue(this%dativar%r(i)%btable == this%dativarattr%r(:)%btable)
10448 ENDDO
10449 ENDIF
10450
10451 IF (ASSOCIATED(this%dativarattr%d)) THEN
10452 DO i = 1, SIZE(this%dativar%r)
10453 this%dativar%r(i)%d = &
10454 firsttrue(this%dativar%r(i)%btable == this%dativarattr%d(:)%btable)
10455 ENDDO
10456 ENDIF
10457
10458 IF (ASSOCIATED(this%dativarattr%i)) THEN
10459 DO i = 1, SIZE(this%dativar%r)
10460 this%dativar%r(i)%i = &
10461 firsttrue(this%dativar%r(i)%btable == this%dativarattr%i(:)%btable)
10462 ENDDO
10463 ENDIF
10464
10465 IF (ASSOCIATED(this%dativarattr%b)) THEN
10466 DO i = 1, SIZE(this%dativar%r)
10467 this%dativar%r(i)%b = &
10468 firsttrue(this%dativar%r(i)%btable == this%dativarattr%b(:)%btable)
10469 ENDDO
10470 ENDIF
10471
10472 IF (ASSOCIATED(this%dativarattr%c)) THEN
10473 DO i = 1, SIZE(this%dativar%r)
10474 this%dativar%r(i)%c = &
10475 firsttrue(this%dativar%r(i)%btable == this%dativarattr%c(:)%btable)
10476 ENDDO
10477 ENDIF
10478ENDIF
10479! double
10480IF (ASSOCIATED(this%dativar%d)) THEN
10481 IF (ASSOCIATED(this%dativarattr%r)) THEN
10482 DO i = 1, SIZE(this%dativar%d)
10483 this%dativar%d(i)%r = &
10484 firsttrue(this%dativar%d(i)%btable == this%dativarattr%r(:)%btable)
10485 ENDDO
10486 ENDIF
10487
10488 IF (ASSOCIATED(this%dativarattr%d)) THEN
10489 DO i = 1, SIZE(this%dativar%d)
10490 this%dativar%d(i)%d = &
10491 firsttrue(this%dativar%d(i)%btable == this%dativarattr%d(:)%btable)
10492 ENDDO
10493 ENDIF
10494
10495 IF (ASSOCIATED(this%dativarattr%i)) THEN
10496 DO i = 1, SIZE(this%dativar%d)
10497 this%dativar%d(i)%i = &
10498 firsttrue(this%dativar%d(i)%btable == this%dativarattr%i(:)%btable)
10499 ENDDO
10500 ENDIF
10501
10502 IF (ASSOCIATED(this%dativarattr%b)) THEN
10503 DO i = 1, SIZE(this%dativar%d)
10504 this%dativar%d(i)%b = &
10505 firsttrue(this%dativar%d(i)%btable == this%dativarattr%b(:)%btable)
10506 ENDDO
10507 ENDIF
10508
10509 IF (ASSOCIATED(this%dativarattr%c)) THEN
10510 DO i = 1, SIZE(this%dativar%d)
10511 this%dativar%d(i)%c = &
10512 firsttrue(this%dativar%d(i)%btable == this%dativarattr%c(:)%btable)
10513 ENDDO
10514 ENDIF
10515ENDIF
10516! integer
10517IF (ASSOCIATED(this%dativar%i)) THEN
10518 IF (ASSOCIATED(this%dativarattr%r)) THEN
10519 DO i = 1, SIZE(this%dativar%i)
10520 this%dativar%i(i)%r = &
10521 firsttrue(this%dativar%i(i)%btable == this%dativarattr%r(:)%btable)
10522 ENDDO
10523 ENDIF
10524
10525 IF (ASSOCIATED(this%dativarattr%d)) THEN
10526 DO i = 1, SIZE(this%dativar%i)
10527 this%dativar%i(i)%d = &
10528 firsttrue(this%dativar%i(i)%btable == this%dativarattr%d(:)%btable)
10529 ENDDO
10530 ENDIF
10531
10532 IF (ASSOCIATED(this%dativarattr%i)) THEN
10533 DO i = 1, SIZE(this%dativar%i)
10534 this%dativar%i(i)%i = &
10535 firsttrue(this%dativar%i(i)%btable == this%dativarattr%i(:)%btable)
10536 ENDDO
10537 ENDIF
10538
10539 IF (ASSOCIATED(this%dativarattr%b)) THEN
10540 DO i = 1, SIZE(this%dativar%i)
10541 this%dativar%i(i)%b = &
10542 firsttrue(this%dativar%i(i)%btable == this%dativarattr%b(:)%btable)
10543 ENDDO
10544 ENDIF
10545
10546 IF (ASSOCIATED(this%dativarattr%c)) THEN
10547 DO i = 1, SIZE(this%dativar%i)
10548 this%dativar%i(i)%c = &
10549 firsttrue(this%dativar%i(i)%btable == this%dativarattr%c(:)%btable)
10550 ENDDO
10551 ENDIF
10552ENDIF
10553! byte
10554IF (ASSOCIATED(this%dativar%b)) THEN
10555 IF (ASSOCIATED(this%dativarattr%r)) THEN
10556 DO i = 1, SIZE(this%dativar%b)
10557 this%dativar%b(i)%r = &
10558 firsttrue(this%dativar%b(i)%btable == this%dativarattr%r(:)%btable)
10559 ENDDO
10560 ENDIF
10561
10562 IF (ASSOCIATED(this%dativarattr%d)) THEN
10563 DO i = 1, SIZE(this%dativar%b)
10564 this%dativar%b(i)%d = &
10565 firsttrue(this%dativar%b(i)%btable == this%dativarattr%d(:)%btable)
10566 ENDDO
10567 ENDIF
10568
10569 IF (ASSOCIATED(this%dativarattr%i)) THEN
10570 DO i = 1, SIZE(this%dativar%b)
10571 this%dativar%b(i)%i = &
10572 firsttrue(this%dativar%b(i)%btable == this%dativarattr%i(:)%btable)
10573 ENDDO
10574 ENDIF
10575
10576 IF (ASSOCIATED(this%dativarattr%b)) THEN
10577 DO i = 1, SIZE(this%dativar%b)
10578 this%dativar%b(i)%b = &
10579 firsttrue(this%dativar%b(i)%btable == this%dativarattr%b(:)%btable)
10580 ENDDO
10581 ENDIF
10582
10583 IF (ASSOCIATED(this%dativarattr%c)) THEN
10584 DO i = 1, SIZE(this%dativar%b)
10585 this%dativar%b(i)%c = &
10586 firsttrue(this%dativar%b(i)%btable == this%dativarattr%c(:)%btable)
10587 ENDDO
10588 ENDIF
10589ENDIF
10590! character
10591IF (ASSOCIATED(this%dativar%c)) THEN
10592 IF (ASSOCIATED(this%dativarattr%r)) THEN
10593 DO i = 1, SIZE(this%dativar%c)
10594 this%dativar%c(i)%r = &
10595 firsttrue(this%dativar%c(i)%btable == this%dativarattr%r(:)%btable)
10596 ENDDO
10597 ENDIF
10598
10599 IF (ASSOCIATED(this%dativarattr%d)) THEN
10600 DO i = 1, SIZE(this%dativar%c)
10601 this%dativar%c(i)%d = &
10602 firsttrue(this%dativar%c(i)%btable == this%dativarattr%d(:)%btable)
10603 ENDDO
10604 ENDIF
10605
10606 IF (ASSOCIATED(this%dativarattr%i)) THEN
10607 DO i = 1, SIZE(this%dativar%c)
10608 this%dativar%c(i)%i = &
10609 firsttrue(this%dativar%c(i)%btable == this%dativarattr%i(:)%btable)
10610 ENDDO
10611 ENDIF
10612
10613 IF (ASSOCIATED(this%dativarattr%b)) THEN
10614 DO i = 1, SIZE(this%dativar%c)
10615 this%dativar%c(i)%b = &
10616 firsttrue(this%dativar%c(i)%btable == this%dativarattr%b(:)%btable)
10617 ENDDO
10618 ENDIF
10619
10620 IF (ASSOCIATED(this%dativarattr%c)) THEN
10621 DO i = 1, SIZE(this%dativar%c)
10622 this%dativar%c(i)%c = &
10623 firsttrue(this%dativar%c(i)%btable == this%dativarattr%c(:)%btable)
10624 ENDDO
10625 ENDIF
10626ENDIF
10627
10628END SUBROUTINE vol7d_set_attr_ind
10629
10630
10635SUBROUTINE vol7d_merge(this, that, sort, bestdata, &
10636 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
10637TYPE(vol7d),INTENT(INOUT) :: this
10638TYPE(vol7d),INTENT(INOUT) :: that
10639LOGICAL,INTENT(IN),OPTIONAL :: sort
10640LOGICAL,INTENT(in),OPTIONAL :: bestdata
10641LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple ! experimental, please do not use outside the library now
10642
10643TYPE(vol7d) :: v7d_clean
10644
10645
10646IF (.NOT.c_e(this)) THEN ! speedup
10647 this = that
10648 CALL init(v7d_clean)
10649 that = v7d_clean ! destroy that without deallocating
10650ELSE ! Append that to this and destroy that
10651 CALL vol7d_append(this, that, sort, bestdata, &
10652 ltimesimple, ltimerangesimple, llevelsimple, lanasimple)
10653 CALL delete(that)
10654ENDIF
10655
10656END SUBROUTINE vol7d_merge
10657
10658
10687SUBROUTINE vol7d_append(this, that, sort, bestdata, &
10688 ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple)
10689TYPE(vol7d),INTENT(INOUT) :: this
10690TYPE(vol7d),INTENT(IN) :: that
10691LOGICAL,INTENT(IN),OPTIONAL :: sort
10692! experimental, please do not use outside the library now, they force the use
10693! of a simplified mapping algorithm which is valid only whene the dimension
10694! content is the same in both volumes , or when one of them is empty
10695LOGICAL,INTENT(in),OPTIONAL :: bestdata
10696LOGICAL,INTENT(IN),OPTIONAL :: ltimesimple, ltimerangesimple, llevelsimple, lanasimple, lnetworksimple
10697
10698
10699TYPE(vol7d) :: v7dtmp
10700LOGICAL :: lsort, lbestdata
10701INTEGER,POINTER :: remapt1(:), remapt2(:), remaptr1(:), remaptr2(:), &
10702 remapl1(:), remapl2(:), remapa1(:), remapa2(:), remapn1(:), remapn2(:)
10703
10704IF (.NOT.c_e(that)) RETURN ! speedup, nothing to do
10705IF (.NOT.vol7d_check_vol(that)) RETURN ! be safe
10706IF (.NOT.c_e(this)) THEN ! this case is like a vol7d_copy, more efficient to copy?
10707 CALL vol7d_copy(that, this, sort=sort)
10708 RETURN
10709ENDIF
10710
10711IF (this%time_definition /= that%time_definition) THEN
10712 CALL l4f_log(l4f_fatal, &
10713 'in vol7d_append, cannot append volumes with different &
10714 &time definition')
10715 CALL raise_fatal_error()
10716ENDIF
10717
10718! Completo l'allocazione per avere volumi a norma
10719CALL vol7d_alloc_vol(this)
10720
10721CALL init(v7dtmp, time_definition=this%time_definition)
10722CALL optio(sort, lsort)
10723CALL optio(bestdata, lbestdata)
10724
10725! Calcolo le mappature tra volumi vecchi e volume nuovo
10726! I puntatori remap* vengono tutti o allocati o nullificati
10727IF (optio_log(ltimesimple)) THEN
10728 CALL vol7d_remap2simple_datetime(this%time, that%time, v7dtmp%time, &
10729 lsort, remapt1, remapt2)
10730ELSE
10731 CALL vol7d_remap2_datetime(this%time, that%time, v7dtmp%time, &
10732 lsort, remapt1, remapt2)
10733ENDIF
10734IF (optio_log(ltimerangesimple)) THEN
10735 CALL vol7d_remap2simple_vol7d_timerange(this%timerange, that%timerange, &
10736 v7dtmp%timerange, lsort, remaptr1, remaptr2)
10737ELSE
10738 CALL vol7d_remap2_vol7d_timerange(this%timerange, that%timerange, &
10739 v7dtmp%timerange, lsort, remaptr1, remaptr2)
10740ENDIF
10741IF (optio_log(llevelsimple)) THEN
10742 CALL vol7d_remap2simple_vol7d_level(this%level, that%level, v7dtmp%level, &
10743 lsort, remapl1, remapl2)
10744ELSE
10745 CALL vol7d_remap2_vol7d_level(this%level, that%level, v7dtmp%level, &
10746 lsort, remapl1, remapl2)
10747ENDIF
10748IF (optio_log(lanasimple)) THEN
10749 CALL vol7d_remap2simple_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
10750 .false., remapa1, remapa2)
10751ELSE
10752 CALL vol7d_remap2_vol7d_ana(this%ana, that%ana, v7dtmp%ana, &
10753 .false., remapa1, remapa2)
10754ENDIF
10755IF (optio_log(lnetworksimple)) THEN
10756 CALL vol7d_remap2simple_vol7d_network(this%network, that%network, v7dtmp%network, &
10757 .false., remapn1, remapn2)
10758ELSE
10759 CALL vol7d_remap2_vol7d_network(this%network, that%network, v7dtmp%network, &
10760 .false., remapn1, remapn2)
10761ENDIF
10762
10763! Faccio la fusione fisica dei volumi
10764CALL vol7d_merge_finalr(this, that, v7dtmp, &
10765 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10766 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10767CALL vol7d_merge_finald(this, that, v7dtmp, &
10768 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10769 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10770CALL vol7d_merge_finali(this, that, v7dtmp, &
10771 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10772 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10773CALL vol7d_merge_finalb(this, that, v7dtmp, &
10774 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10775 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10776CALL vol7d_merge_finalc(this, that, v7dtmp, &
10777 remapa1, remapa2, remapt1, remapt2, remapl1, remapl2, &
10778 remaptr1, remaptr2, remapn1, remapn2, lbestdata)
10779
10780! Dealloco i vettori di rimappatura
10781IF (ASSOCIATED(remapt1)) DEALLOCATE(remapt1)
10782IF (ASSOCIATED(remapt2)) DEALLOCATE(remapt2)
10783IF (ASSOCIATED(remaptr1)) DEALLOCATE(remaptr1)
10784IF (ASSOCIATED(remaptr2)) DEALLOCATE(remaptr2)
10785IF (ASSOCIATED(remapl1)) DEALLOCATE(remapl1)
10786IF (ASSOCIATED(remapl2)) DEALLOCATE(remapl2)
10787IF (ASSOCIATED(remapa1)) DEALLOCATE(remapa1)
10788IF (ASSOCIATED(remapa2)) DEALLOCATE(remapa2)
10789IF (ASSOCIATED(remapn1)) DEALLOCATE(remapn1)
10790IF (ASSOCIATED(remapn2)) DEALLOCATE(remapn2)
10791
10792! Distruggo il vecchio volume e assegno il nuovo a this
10793CALL delete(this)
10794this = v7dtmp
10795! Ricreo gli indici var-attr
10796CALL vol7d_set_attr_ind(this)
10797
10798END SUBROUTINE vol7d_append
10799
10800
10833SUBROUTINE vol7d_copy(this, that, sort, unique, miss, &
10834 lsort_time, lsort_timerange, lsort_level, &
10835 ltime, ltimerange, llevel, lana, lnetwork, &
10836 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
10837 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
10838 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
10839 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
10840 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
10841 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
10842TYPE(vol7d),INTENT(IN) :: this
10843TYPE(vol7d),INTENT(INOUT) :: that
10844LOGICAL,INTENT(IN),OPTIONAL :: sort
10845LOGICAL,INTENT(IN),OPTIONAL :: unique
10846LOGICAL,INTENT(IN),OPTIONAL :: miss
10847LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
10848LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
10849LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
10857LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
10859LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
10861LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
10863LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
10865LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
10867LOGICAL,INTENT(in),OPTIONAL :: &
10868 lanavarr(:), lanavard(:), lanavari(:), lanavarb(:), lanavarc(:), &
10869 lanaattrr(:), lanaattrd(:), lanaattri(:), lanaattrb(:), lanaattrc(:), &
10870 lanavarattrr(:), lanavarattrd(:), lanavarattri(:), lanavarattrb(:), lanavarattrc(:), &
10871 ldativarr(:), ldativard(:), ldativari(:), ldativarb(:), ldativarc(:), &
10872 ldatiattrr(:), ldatiattrd(:), ldatiattri(:), ldatiattrb(:), ldatiattrc(:), &
10873 ldativarattrr(:), ldativarattrd(:), ldativarattri(:), ldativarattrb(:), ldativarattrc(:)
10874
10875LOGICAL :: lsort, lunique, lmiss
10876INTEGER,POINTER :: remapt(:), remaptr(:), remapl(:), remapa(:), remapn(:)
10877
10878CALL init(that)
10879IF (.NOT.c_e(this)) RETURN ! speedup, nothing to do
10880IF (.NOT.vol7d_check_vol(this)) RETURN ! be safe
10881
10882CALL optio(sort, lsort)
10883CALL optio(unique, lunique)
10884CALL optio(miss, lmiss)
10885
10886! Calcolo le mappature tra volume vecchio e volume nuovo
10887! I puntatori remap* vengono tutti o allocati o nullificati
10888CALL vol7d_remap1_datetime(this%time, that%time, &
10889 lsort.OR.optio_log(lsort_time), lunique, lmiss, remapt, ltime)
10890CALL vol7d_remap1_vol7d_timerange(this%timerange, that%timerange, &
10891 lsort.OR.optio_log(lsort_timerange), lunique, lmiss, remaptr, ltimerange)
10892CALL vol7d_remap1_vol7d_level(this%level, that%level, &
10893 lsort.OR.optio_log(lsort_level), lunique, lmiss, remapl, llevel)
10894CALL vol7d_remap1_vol7d_ana(this%ana, that%ana, &
10895 lsort, lunique, lmiss, remapa, lana)
10896CALL vol7d_remap1_vol7d_network(this%network, that%network, &
10897 lsort, lunique, lmiss, remapn, lnetwork)
10898
10899! lanavari, lanavarb, lanavarc, &
10900! lanaattri, lanaattrb, lanaattrc, &
10901! lanavarattri, lanavarattrb, lanavarattrc, &
10902! ldativari, ldativarb, ldativarc, &
10903! ldatiattri, ldatiattrb, ldatiattrc, &
10904! ldativarattri, ldativarattrb, ldativarattrc
10905! Faccio la riforma fisica dei volumi
10906CALL vol7d_reform_finalr(this, that, &
10907 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10908 lanavarr, lanaattrr, lanavarattrr, ldativarr, ldatiattrr, ldativarattrr)
10909CALL vol7d_reform_finald(this, that, &
10910 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10911 lanavard, lanaattrd, lanavarattrd, ldativard, ldatiattrd, ldativarattrd)
10912CALL vol7d_reform_finali(this, that, &
10913 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10914 lanavari, lanaattri, lanavarattri, ldativari, ldatiattri, ldativarattri)
10915CALL vol7d_reform_finalb(this, that, &
10916 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10917 lanavarb, lanaattrb, lanavarattrb, ldativarb, ldatiattrb, ldativarattrb)
10918CALL vol7d_reform_finalc(this, that, &
10919 remapa, remapt, remapl, remaptr, remapn, lsort, lunique, lmiss, &
10920 lanavarc, lanaattrc, lanavarattrc, ldativarc, ldatiattrc, ldativarattrc)
10921
10922! Dealloco i vettori di rimappatura
10923IF (ASSOCIATED(remapt)) DEALLOCATE(remapt)
10924IF (ASSOCIATED(remaptr)) DEALLOCATE(remaptr)
10925IF (ASSOCIATED(remapl)) DEALLOCATE(remapl)
10926IF (ASSOCIATED(remapa)) DEALLOCATE(remapa)
10927IF (ASSOCIATED(remapn)) DEALLOCATE(remapn)
10928
10929! Ricreo gli indici var-attr
10930CALL vol7d_set_attr_ind(that)
10931that%time_definition = this%time_definition
10932
10933END SUBROUTINE vol7d_copy
10934
10935
10946SUBROUTINE vol7d_reform(this, sort, unique, miss, &
10947 lsort_time, lsort_timerange, lsort_level, &
10948 ltime, ltimerange, llevel, lana, lnetwork, &
10949 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
10950 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
10951 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
10952 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
10953 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
10954 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc&
10955 ,purgeana)
10956TYPE(vol7d),INTENT(INOUT) :: this
10957LOGICAL,INTENT(IN),OPTIONAL :: sort
10958LOGICAL,INTENT(IN),OPTIONAL :: unique
10959LOGICAL,INTENT(IN),OPTIONAL :: miss
10960LOGICAL,INTENT(IN),OPTIONAL :: lsort_time
10961LOGICAL,INTENT(IN),OPTIONAL :: lsort_timerange
10962LOGICAL,INTENT(IN),OPTIONAL :: lsort_level
10970LOGICAL,INTENT(IN),OPTIONAL :: ltime(:)
10971LOGICAL,INTENT(IN),OPTIONAL :: ltimerange(:)
10972LOGICAL,INTENT(IN),OPTIONAL :: llevel(:)
10973LOGICAL,INTENT(IN),OPTIONAL :: lana(:)
10974LOGICAL,INTENT(IN),OPTIONAL :: lnetwork(:)
10976LOGICAL,INTENT(in),OPTIONAL :: &
10977 lanavarr(:), lanavard(:), lanavari(:), lanavarb(:), lanavarc(:), &
10978 lanaattrr(:), lanaattrd(:), lanaattri(:), lanaattrb(:), lanaattrc(:), &
10979 lanavarattrr(:), lanavarattrd(:), lanavarattri(:), lanavarattrb(:), lanavarattrc(:), &
10980 ldativarr(:), ldativard(:), ldativari(:), ldativarb(:), ldativarc(:), &
10981 ldatiattrr(:), ldatiattrd(:), ldatiattri(:), ldatiattrb(:), ldatiattrc(:), &
10982 ldativarattrr(:), ldativarattrd(:), ldativarattri(:), ldativarattrb(:), ldativarattrc(:)
10983LOGICAL,INTENT(IN),OPTIONAL :: purgeana
10984
10985TYPE(vol7d) :: v7dtmp
10986logical,allocatable :: llana(:)
10987integer :: i
10988
10989CALL vol7d_copy(this, v7dtmp, sort, unique, miss, &
10990 lsort_time, lsort_timerange, lsort_level, &
10991 ltime, ltimerange, llevel, lana, lnetwork, &
10992 lanavarr, lanavard, lanavari, lanavarb, lanavarc, &
10993 lanaattrr, lanaattrd, lanaattri, lanaattrb, lanaattrc, &
10994 lanavarattrr, lanavarattrd, lanavarattri, lanavarattrb, lanavarattrc, &
10995 ldativarr, ldativard, ldativari, ldativarb, ldativarc, &
10996 ldatiattrr, ldatiattrd, ldatiattri, ldatiattrb, ldatiattrc, &
10997 ldativarattrr, ldativarattrd, ldativarattri, ldativarattrb, ldativarattrc)
10998
10999! destroy old volume
11000CALL delete(this)
11001
11002if (optio_log(purgeana)) then
11003 allocate(llana(size(v7dtmp%ana)))
11004 llana =.false.
11005 do i =1,size(v7dtmp%ana)
11006 if (associated(v7dtmp%voldatii)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatii(i,:,:,:,:,:)))
11007 if (associated(v7dtmp%voldatir)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatir(i,:,:,:,:,:)))
11008 if (associated(v7dtmp%voldatid)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatid(i,:,:,:,:,:)))
11009 if (associated(v7dtmp%voldatib)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatib(i,:,:,:,:,:)))
11010 if (associated(v7dtmp%voldatic)) llana(i)= llana(i) .or. any(c_e(v7dtmp%voldatic(i,:,:,:,:,:)))
11011 end do
11012 CALL vol7d_copy(v7dtmp, this,lana=llana)
11013 CALL delete(v7dtmp)
11014 deallocate(llana)
11015else
11016 this=v7dtmp
11017end if
11018
11019END SUBROUTINE vol7d_reform
11020
11021
11029SUBROUTINE vol7d_smart_sort(this, lsort_time, lsort_timerange, lsort_level)
11030TYPE(vol7d),INTENT(INOUT) :: this
11031LOGICAL,OPTIONAL,INTENT(in) :: lsort_time
11032LOGICAL,OPTIONAL,INTENT(in) :: lsort_timerange
11033LOGICAL,OPTIONAL,INTENT(in) :: lsort_level
11034
11035INTEGER :: i
11036LOGICAL :: to_be_sorted
11037
11038to_be_sorted = .false.
11039CALL vol7d_alloc_vol(this) ! usual safety check
11040
11041IF (optio_log(lsort_time)) THEN
11042 DO i = 2, SIZE(this%time)
11043 IF (this%time(i) < this%time(i-1)) THEN
11044 to_be_sorted = .true.
11045 EXIT
11046 ENDIF
11047 ENDDO
11048ENDIF
11049IF (optio_log(lsort_timerange)) THEN
11050 DO i = 2, SIZE(this%timerange)
11051 IF (this%timerange(i) < this%timerange(i-1)) THEN
11052 to_be_sorted = .true.
11053 EXIT
11054 ENDIF
11055 ENDDO
11056ENDIF
11057IF (optio_log(lsort_level)) THEN
11058 DO i = 2, SIZE(this%level)
11059 IF (this%level(i) < this%level(i-1)) THEN
11060 to_be_sorted = .true.
11061 EXIT
11062 ENDIF
11063 ENDDO
11064ENDIF
11065
11066IF (to_be_sorted) CALL vol7d_reform(this, &
11067 lsort_time=lsort_time, lsort_timerange=lsort_timerange, lsort_level=lsort_level )
11068
11069END SUBROUTINE vol7d_smart_sort
11070
11078SUBROUTINE vol7d_filter(this, avl, vl, nl, s_d, e_d)
11079TYPE(vol7d),INTENT(inout) :: this
11080CHARACTER(len=*),INTENT(in),OPTIONAL :: avl(:)
11081CHARACTER(len=*),INTENT(in),OPTIONAL :: vl(:)
11082TYPE(vol7d_network),OPTIONAL :: nl(:)
11083TYPE(datetime),INTENT(in),OPTIONAL :: s_d
11084TYPE(datetime),INTENT(in),OPTIONAL :: e_d
11085
11086INTEGER :: i
11087
11088IF (PRESENT(avl)) THEN
11089 IF (SIZE(avl) > 0) THEN
11090
11091 IF (ASSOCIATED(this%anavar%r)) THEN
11092 DO i = 1, SIZE(this%anavar%r)
11093 IF (all(this%anavar%r(i)%btable /= avl)) this%anavar%r(i) = vol7d_var_miss
11094 ENDDO
11095 ENDIF
11096
11097 IF (ASSOCIATED(this%anavar%i)) THEN
11098 DO i = 1, SIZE(this%anavar%i)
11099 IF (all(this%anavar%i(i)%btable /= avl)) this%anavar%i(i) = vol7d_var_miss
11100 ENDDO
11101 ENDIF
11102
11103 IF (ASSOCIATED(this%anavar%b)) THEN
11104 DO i = 1, SIZE(this%anavar%b)
11105 IF (all(this%anavar%b(i)%btable /= avl)) this%anavar%b(i) = vol7d_var_miss
11106 ENDDO
11107 ENDIF
11108
11109 IF (ASSOCIATED(this%anavar%d)) THEN
11110 DO i = 1, SIZE(this%anavar%d)
11111 IF (all(this%anavar%d(i)%btable /= avl)) this%anavar%d(i) = vol7d_var_miss
11112 ENDDO
11113 ENDIF
11114
11115 IF (ASSOCIATED(this%anavar%c)) THEN
11116 DO i = 1, SIZE(this%anavar%c)
11117 IF (all(this%anavar%c(i)%btable /= avl)) this%anavar%c(i) = vol7d_var_miss
11118 ENDDO
11119 ENDIF
11120
11121 ENDIF
11122ENDIF
11123
11124
11125IF (PRESENT(vl)) THEN
11126 IF (size(vl) > 0) THEN
11127 IF (ASSOCIATED(this%dativar%r)) THEN
11128 DO i = 1, SIZE(this%dativar%r)
11129 IF (all(this%dativar%r(i)%btable /= vl)) this%dativar%r(i) = vol7d_var_miss
11130 ENDDO
11131 ENDIF
11132
11133 IF (ASSOCIATED(this%dativar%i)) THEN
11134 DO i = 1, SIZE(this%dativar%i)
11135 IF (all(this%dativar%i(i)%btable /= vl)) this%dativar%i(i) = vol7d_var_miss
11136 ENDDO
11137 ENDIF
11138
11139 IF (ASSOCIATED(this%dativar%b)) THEN
11140 DO i = 1, SIZE(this%dativar%b)
11141 IF (all(this%dativar%b(i)%btable /= vl)) this%dativar%b(i) = vol7d_var_miss
11142 ENDDO
11143 ENDIF
11144
11145 IF (ASSOCIATED(this%dativar%d)) THEN
11146 DO i = 1, SIZE(this%dativar%d)
11147 IF (all(this%dativar%d(i)%btable /= vl)) this%dativar%d(i) = vol7d_var_miss
11148 ENDDO
11149 ENDIF
11150
11151 IF (ASSOCIATED(this%dativar%c)) THEN
11152 DO i = 1, SIZE(this%dativar%c)
11153 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
11154 ENDDO
11155 ENDIF
11156
11157 IF (ASSOCIATED(this%dativar%c)) THEN
11158 DO i = 1, SIZE(this%dativar%c)
11159 IF (all(this%dativar%c(i)%btable /= vl)) this%dativar%c(i) = vol7d_var_miss
11160 ENDDO
11161 ENDIF
11162
11163 ENDIF
11164ENDIF
11165
11166IF (PRESENT(nl)) THEN
11167 IF (SIZE(nl) > 0) THEN
11168 DO i = 1, SIZE(this%network)
11169 IF (all(this%network(i) /= nl)) this%network(i) = vol7d_network_miss
11170 ENDDO
11171 ENDIF
11172ENDIF
11173
11174IF (PRESENT(s_d)) THEN
11175 IF (c_e(s_d)) THEN
11176 WHERE (this%time < s_d)
11177 this%time = datetime_miss
11178 END WHERE
11179 ENDIF
11180ENDIF
11181
11182IF (PRESENT(e_d)) THEN
11183 IF (c_e(e_d)) THEN
11184 WHERE (this%time > e_d)
11185 this%time = datetime_miss
11186 END WHERE
11187 ENDIF
11188ENDIF
11189
11190CALL vol7d_reform(this, miss=.true.)
11191
11192END SUBROUTINE vol7d_filter
11193
11194
11201SUBROUTINE vol7d_convr(this, that, anaconv)
11202TYPE(vol7d),INTENT(IN) :: this
11203TYPE(vol7d),INTENT(INOUT) :: that
11204LOGICAL,OPTIONAL,INTENT(in) :: anaconv
11205INTEGER :: i
11206LOGICAL :: fv(1)=(/.false./), tv(1)=(/.true./), acp(1), acn(1)
11207TYPE(vol7d) :: v7d_tmp
11208
11209IF (optio_log(anaconv)) THEN
11210 acp=fv
11211 acn=tv
11212ELSE
11213 acp=tv
11214 acn=fv
11215ENDIF
11216
11217! Volume con solo i dati reali e tutti gli attributi
11218! l'anagrafica e` copiata interamente se necessario
11219CALL vol7d_copy(this, that, &
11220 lanavarr=tv, lanavard=acp, lanavari=acp, lanavarb=acp, lanavarc=acp, &
11221 ldativarr=tv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=fv)
11222
11223! Volume solo di dati double
11224CALL vol7d_copy(this, v7d_tmp, &
11225 lanavarr=fv, lanavard=acn, lanavari=fv, lanavarb=fv, lanavarc=fv, &
11226 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11227 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11228 ldativarr=fv, ldativard=tv, ldativari=fv, ldativarb=fv, ldativarc=fv, &
11229 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11230 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11231
11232! converto a dati reali
11233IF (ASSOCIATED(v7d_tmp%anavar%d) .OR. ASSOCIATED(v7d_tmp%dativar%d)) THEN
11234
11235 IF (ASSOCIATED(v7d_tmp%anavar%d)) THEN
11236! alloco i dati reali e vi trasferisco i double
11237 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanad, 1), SIZE(v7d_tmp%volanad, 2), &
11238 SIZE(v7d_tmp%volanad, 3)))
11239 DO i = 1, SIZE(v7d_tmp%anavar%d)
11240 v7d_tmp%volanar(:,i,:) = &
11241 realdat(v7d_tmp%volanad(:,i,:), v7d_tmp%anavar%d(i))
11242 ENDDO
11243 DEALLOCATE(v7d_tmp%volanad)
11244! trasferisco le variabili
11245 v7d_tmp%anavar%r => v7d_tmp%anavar%d
11246 NULLIFY(v7d_tmp%anavar%d)
11247 ENDIF
11248
11249 IF (ASSOCIATED(v7d_tmp%dativar%d)) THEN
11250! alloco i dati reali e vi trasferisco i double
11251 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatid, 1), SIZE(v7d_tmp%voldatid, 2), &
11252 SIZE(v7d_tmp%voldatid, 3), SIZE(v7d_tmp%voldatid, 4), SIZE(v7d_tmp%voldatid, 5), &
11253 SIZE(v7d_tmp%voldatid, 6)))
11254 DO i = 1, SIZE(v7d_tmp%dativar%d)
11255 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11256 realdat(v7d_tmp%voldatid(:,:,:,:,i,:), v7d_tmp%dativar%d(i))
11257 ENDDO
11258 DEALLOCATE(v7d_tmp%voldatid)
11259! trasferisco le variabili
11260 v7d_tmp%dativar%r => v7d_tmp%dativar%d
11261 NULLIFY(v7d_tmp%dativar%d)
11262 ENDIF
11263
11264! fondo con il volume definitivo
11265 CALL vol7d_merge(that, v7d_tmp)
11266ELSE
11267 CALL delete(v7d_tmp)
11268ENDIF
11269
11270
11271! Volume solo di dati interi
11272CALL vol7d_copy(this, v7d_tmp, &
11273 lanavarr=fv, lanavard=fv, lanavari=acn, lanavarb=fv, lanavarc=fv, &
11274 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11275 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11276 ldativarr=fv, ldativard=fv, ldativari=tv, ldativarb=fv, ldativarc=fv, &
11277 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11278 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11279
11280! converto a dati reali
11281IF (ASSOCIATED(v7d_tmp%anavar%i) .OR. ASSOCIATED(v7d_tmp%dativar%i)) THEN
11282
11283 IF (ASSOCIATED(v7d_tmp%anavar%i)) THEN
11284! alloco i dati reali e vi trasferisco gli interi
11285 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanai, 1), SIZE(v7d_tmp%volanai, 2), &
11286 SIZE(v7d_tmp%volanai, 3)))
11287 DO i = 1, SIZE(v7d_tmp%anavar%i)
11288 v7d_tmp%volanar(:,i,:) = &
11289 realdat(v7d_tmp%volanai(:,i,:), v7d_tmp%anavar%i(i))
11290 ENDDO
11291 DEALLOCATE(v7d_tmp%volanai)
11292! trasferisco le variabili
11293 v7d_tmp%anavar%r => v7d_tmp%anavar%i
11294 NULLIFY(v7d_tmp%anavar%i)
11295 ENDIF
11296
11297 IF (ASSOCIATED(v7d_tmp%dativar%i)) THEN
11298! alloco i dati reali e vi trasferisco gli interi
11299 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatii, 1), SIZE(v7d_tmp%voldatii, 2), &
11300 SIZE(v7d_tmp%voldatii, 3), SIZE(v7d_tmp%voldatii, 4), SIZE(v7d_tmp%voldatii, 5), &
11301 SIZE(v7d_tmp%voldatii, 6)))
11302 DO i = 1, SIZE(v7d_tmp%dativar%i)
11303 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11304 realdat(v7d_tmp%voldatii(:,:,:,:,i,:), v7d_tmp%dativar%i(i))
11305 ENDDO
11306 DEALLOCATE(v7d_tmp%voldatii)
11307! trasferisco le variabili
11308 v7d_tmp%dativar%r => v7d_tmp%dativar%i
11309 NULLIFY(v7d_tmp%dativar%i)
11310 ENDIF
11311
11312! fondo con il volume definitivo
11313 CALL vol7d_merge(that, v7d_tmp)
11314ELSE
11315 CALL delete(v7d_tmp)
11316ENDIF
11317
11318
11319! Volume solo di dati byte
11320CALL vol7d_copy(this, v7d_tmp, &
11321 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=acn, lanavarc=fv, &
11322 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11323 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11324 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=tv, ldativarc=fv, &
11325 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11326 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11327
11328! converto a dati reali
11329IF (ASSOCIATED(v7d_tmp%anavar%b) .OR. ASSOCIATED(v7d_tmp%dativar%b)) THEN
11330
11331 IF (ASSOCIATED(v7d_tmp%anavar%b)) THEN
11332! alloco i dati reali e vi trasferisco i byte
11333 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanab, 1), SIZE(v7d_tmp%volanab, 2), &
11334 SIZE(v7d_tmp%volanab, 3)))
11335 DO i = 1, SIZE(v7d_tmp%anavar%b)
11336 v7d_tmp%volanar(:,i,:) = &
11337 realdat(v7d_tmp%volanab(:,i,:), v7d_tmp%anavar%b(i))
11338 ENDDO
11339 DEALLOCATE(v7d_tmp%volanab)
11340! trasferisco le variabili
11341 v7d_tmp%anavar%r => v7d_tmp%anavar%b
11342 NULLIFY(v7d_tmp%anavar%b)
11343 ENDIF
11344
11345 IF (ASSOCIATED(v7d_tmp%dativar%b)) THEN
11346! alloco i dati reali e vi trasferisco i byte
11347 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatib, 1), SIZE(v7d_tmp%voldatib, 2), &
11348 SIZE(v7d_tmp%voldatib, 3), SIZE(v7d_tmp%voldatib, 4), SIZE(v7d_tmp%voldatib, 5), &
11349 SIZE(v7d_tmp%voldatib, 6)))
11350 DO i = 1, SIZE(v7d_tmp%dativar%b)
11351 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11352 realdat(v7d_tmp%voldatib(:,:,:,:,i,:), v7d_tmp%dativar%b(i))
11353 ENDDO
11354 DEALLOCATE(v7d_tmp%voldatib)
11355! trasferisco le variabili
11356 v7d_tmp%dativar%r => v7d_tmp%dativar%b
11357 NULLIFY(v7d_tmp%dativar%b)
11358 ENDIF
11359
11360! fondo con il volume definitivo
11361 CALL vol7d_merge(that, v7d_tmp)
11362ELSE
11363 CALL delete(v7d_tmp)
11364ENDIF
11365
11366
11367! Volume solo di dati character
11368CALL vol7d_copy(this, v7d_tmp, &
11369 lanavarr=fv, lanavard=fv, lanavari=fv, lanavarb=fv, lanavarc=acn, &
11370 lanaattrr=fv, lanaattrd=fv, lanaattri=fv, lanaattrb=fv, lanaattrc=fv, &
11371 lanavarattrr=fv, lanavarattrd=fv, lanavarattri=fv, lanavarattrb=fv, lanavarattrc=fv, &
11372 ldativarr=fv, ldativard=fv, ldativari=fv, ldativarb=fv, ldativarc=tv, &
11373 ldatiattrr=fv, ldatiattrd=fv, ldatiattri=fv, ldatiattrb=fv, ldatiattrc=fv, &
11374 ldativarattrr=fv, ldativarattrd=fv, ldativarattri=fv, ldativarattrb=fv, ldativarattrc=fv)
11375
11376! converto a dati reali
11377IF (ASSOCIATED(v7d_tmp%anavar%c) .OR. ASSOCIATED(v7d_tmp%dativar%c)) THEN
11378
11379 IF (ASSOCIATED(v7d_tmp%anavar%c)) THEN
11380! alloco i dati reali e vi trasferisco i character
11381 ALLOCATE(v7d_tmp%volanar(SIZE(v7d_tmp%volanac, 1), SIZE(v7d_tmp%volanac, 2), &
11382 SIZE(v7d_tmp%volanac, 3)))
11383 DO i = 1, SIZE(v7d_tmp%anavar%c)
11384 v7d_tmp%volanar(:,i,:) = &
11385 realdat(v7d_tmp%volanac(:,i,:), v7d_tmp%anavar%c(i))
11386 ENDDO
11387 DEALLOCATE(v7d_tmp%volanac)
11388! trasferisco le variabili
11389 v7d_tmp%anavar%r => v7d_tmp%anavar%c
11390 NULLIFY(v7d_tmp%anavar%c)
11391 ENDIF
11392
11393 IF (ASSOCIATED(v7d_tmp%dativar%c)) THEN
11394! alloco i dati reali e vi trasferisco i character
11395 ALLOCATE(v7d_tmp%voldatir(SIZE(v7d_tmp%voldatic, 1), SIZE(v7d_tmp%voldatic, 2), &
11396 SIZE(v7d_tmp%voldatic, 3), SIZE(v7d_tmp%voldatic, 4), SIZE(v7d_tmp%voldatic, 5), &
11397 SIZE(v7d_tmp%voldatic, 6)))
11398 DO i = 1, SIZE(v7d_tmp%dativar%c)
11399 v7d_tmp%voldatir(:,:,:,:,i,:) = &
11400 realdat(v7d_tmp%voldatic(:,:,:,:,i,:), v7d_tmp%dativar%c(i))
11401 ENDDO
11402 DEALLOCATE(v7d_tmp%voldatic)
11403! trasferisco le variabili
11404 v7d_tmp%dativar%r => v7d_tmp%dativar%c
11405 NULLIFY(v7d_tmp%dativar%c)
11406 ENDIF
11407
11408! fondo con il volume definitivo
11409 CALL vol7d_merge(that, v7d_tmp)
11410ELSE
11411 CALL delete(v7d_tmp)
11412ENDIF
11413
11414END SUBROUTINE vol7d_convr
11415
11416
11420SUBROUTINE vol7d_diff_only (this, that, data_only,ana)
11421TYPE(vol7d),INTENT(IN) :: this
11422TYPE(vol7d),INTENT(OUT) :: that
11423logical , optional, intent(in) :: data_only
11424logical , optional, intent(in) :: ana
11425logical :: ldata_only,lana
11426
11427IF (PRESENT(data_only)) THEN
11428 ldata_only = data_only
11429ELSE
11430 ldata_only = .false.
11431ENDIF
11432
11433IF (PRESENT(ana)) THEN
11434 lana = ana
11435ELSE
11436 lana = .false.
11437ENDIF
11438
11439
11440#undef VOL7D_POLY_ARRAY
11441#define VOL7D_POLY_ARRAY voldati
11442#include "vol7d_class_diff.F90"
11443#undef VOL7D_POLY_ARRAY
11444#define VOL7D_POLY_ARRAY voldatiattr
11445#include "vol7d_class_diff.F90"
11446#undef VOL7D_POLY_ARRAY
11447
11448if ( .not. ldata_only) then
11449
11450#define VOL7D_POLY_ARRAY volana
11451#include "vol7d_class_diff.F90"
11452#undef VOL7D_POLY_ARRAY
11453#define VOL7D_POLY_ARRAY volanaattr
11454#include "vol7d_class_diff.F90"
11455#undef VOL7D_POLY_ARRAY
11456
11457 if(lana)then
11458 where ( this%ana == that%ana )
11459 that%ana = vol7d_ana_miss
11460 end where
11461 end if
11462
11463end if
11464
11465
11466
11467END SUBROUTINE vol7d_diff_only
11468
11469
11470
11471! Creo le routine da ripetere per i vari tipi di dati di v7d
11472! tramite un template e il preprocessore
11473#undef VOL7D_POLY_TYPE
11474#undef VOL7D_POLY_TYPES
11475#define VOL7D_POLY_TYPE REAL
11476#define VOL7D_POLY_TYPES r
11477#include "vol7d_class_type_templ.F90"
11478#undef VOL7D_POLY_TYPE
11479#undef VOL7D_POLY_TYPES
11480#define VOL7D_POLY_TYPE DOUBLE PRECISION
11481#define VOL7D_POLY_TYPES d
11482#include "vol7d_class_type_templ.F90"
11483#undef VOL7D_POLY_TYPE
11484#undef VOL7D_POLY_TYPES
11485#define VOL7D_POLY_TYPE INTEGER
11486#define VOL7D_POLY_TYPES i
11487#include "vol7d_class_type_templ.F90"
11488#undef VOL7D_POLY_TYPE
11489#undef VOL7D_POLY_TYPES
11490#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
11491#define VOL7D_POLY_TYPES b
11492#include "vol7d_class_type_templ.F90"
11493#undef VOL7D_POLY_TYPE
11494#undef VOL7D_POLY_TYPES
11495#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
11496#define VOL7D_POLY_TYPES c
11497#include "vol7d_class_type_templ.F90"
11498
11499! Creo le routine da ripetere per i vari descrittori di dimensioni di v7d
11500! tramite un template e il preprocessore
11501#define VOL7D_SORT
11502#undef VOL7D_NO_ZERO_ALLOC
11503#undef VOL7D_POLY_TYPE
11504#define VOL7D_POLY_TYPE datetime
11505#include "vol7d_class_desc_templ.F90"
11506#undef VOL7D_POLY_TYPE
11507#define VOL7D_POLY_TYPE vol7d_timerange
11508#include "vol7d_class_desc_templ.F90"
11509#undef VOL7D_POLY_TYPE
11510#define VOL7D_POLY_TYPE vol7d_level
11511#include "vol7d_class_desc_templ.F90"
11512#undef VOL7D_SORT
11513#undef VOL7D_POLY_TYPE
11514#define VOL7D_POLY_TYPE vol7d_network
11515#include "vol7d_class_desc_templ.F90"
11516#undef VOL7D_POLY_TYPE
11517#define VOL7D_POLY_TYPE vol7d_ana
11518#include "vol7d_class_desc_templ.F90"
11519#define VOL7D_NO_ZERO_ALLOC
11520#undef VOL7D_POLY_TYPE
11521#define VOL7D_POLY_TYPE vol7d_var
11522#include "vol7d_class_desc_templ.F90"
11523
11533subroutine vol7d_write_on_file (this,unit,description,filename,filename_auto)
11534
11535TYPE(vol7d),INTENT(IN) :: this
11536integer,optional,intent(inout) :: unit
11537character(len=*),intent(in),optional :: filename
11538character(len=*),intent(out),optional :: filename_auto
11539character(len=*),INTENT(IN),optional :: description
11540
11541integer :: lunit
11542character(len=254) :: ldescription,arg,lfilename
11543integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
11544 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11545 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11546 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11547 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11548 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11549 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
11550!integer :: im,id,iy
11551integer :: tarray(8)
11552logical :: opened,exist
11553
11554 nana=0
11555 ntime=0
11556 ntimerange=0
11557 nlevel=0
11558 nnetwork=0
11559 ndativarr=0
11560 ndativari=0
11561 ndativarb=0
11562 ndativard=0
11563 ndativarc=0
11564 ndatiattrr=0
11565 ndatiattri=0
11566 ndatiattrb=0
11567 ndatiattrd=0
11568 ndatiattrc=0
11569 ndativarattrr=0
11570 ndativarattri=0
11571 ndativarattrb=0
11572 ndativarattrd=0
11573 ndativarattrc=0
11574 nanavarr=0
11575 nanavari=0
11576 nanavarb=0
11577 nanavard=0
11578 nanavarc=0
11579 nanaattrr=0
11580 nanaattri=0
11581 nanaattrb=0
11582 nanaattrd=0
11583 nanaattrc=0
11584 nanavarattrr=0
11585 nanavarattri=0
11586 nanavarattrb=0
11587 nanavarattrd=0
11588 nanavarattrc=0
11589
11590
11591!call idate(im,id,iy)
11592call date_and_time(values=tarray)
11593call getarg(0,arg)
11594
11595if (present(description))then
11596 ldescription=description
11597else
11598 ldescription="Vol7d generated by: "//trim(arg)
11599end if
11600
11601if (.not. present(unit))then
11602 lunit=getunit()
11603else
11604 if (unit==0)then
11605 lunit=getunit()
11606 unit=lunit
11607 else
11608 lunit=unit
11609 end if
11610end if
11611
11612lfilename=trim(arg)//".v7d"
11613if (index(arg,'/',back=.true.) > 0) lfilename=lfilename(index(arg,'/',back=.true.)+1 : )
11614
11615if (present(filename))then
11616 if (filename /= "")then
11617 lfilename=filename
11618 end if
11619end if
11620
11621if (present(filename_auto))filename_auto=lfilename
11622
11623
11624inquire(unit=lunit,opened=opened)
11625if (.not. opened) then
11626! inquire(file=lfilename, EXIST=exist)
11627! IF (exist) THEN
11628! CALL l4f_log(L4F_FATAL, &
11629! 'in vol7d_write_on_file, file exists, cannot open file '//TRIM(lfilename))
11630! CALL raise_fatal_error()
11631! ENDIF
11632 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM')
11633 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
11634end if
11635
11636if (associated(this%ana)) nana=size(this%ana)
11637if (associated(this%time)) ntime=size(this%time)
11638if (associated(this%timerange)) ntimerange=size(this%timerange)
11639if (associated(this%level)) nlevel=size(this%level)
11640if (associated(this%network)) nnetwork=size(this%network)
11641
11642if (associated(this%dativar%r)) ndativarr=size(this%dativar%r)
11643if (associated(this%dativar%i)) ndativari=size(this%dativar%i)
11644if (associated(this%dativar%b)) ndativarb=size(this%dativar%b)
11645if (associated(this%dativar%d)) ndativard=size(this%dativar%d)
11646if (associated(this%dativar%c)) ndativarc=size(this%dativar%c)
11647
11648if (associated(this%datiattr%r)) ndatiattrr=size(this%datiattr%r)
11649if (associated(this%datiattr%i)) ndatiattri=size(this%datiattr%i)
11650if (associated(this%datiattr%b)) ndatiattrb=size(this%datiattr%b)
11651if (associated(this%datiattr%d)) ndatiattrd=size(this%datiattr%d)
11652if (associated(this%datiattr%c)) ndatiattrc=size(this%datiattr%c)
11653
11654if (associated(this%dativarattr%r)) ndativarattrr=size(this%dativarattr%r)
11655if (associated(this%dativarattr%i)) ndativarattri=size(this%dativarattr%i)
11656if (associated(this%dativarattr%b)) ndativarattrb=size(this%dativarattr%b)
11657if (associated(this%dativarattr%d)) ndativarattrd=size(this%dativarattr%d)
11658if (associated(this%dativarattr%c)) ndativarattrc=size(this%dativarattr%c)
11659
11660if (associated(this%anavar%r)) nanavarr=size(this%anavar%r)
11661if (associated(this%anavar%i)) nanavari=size(this%anavar%i)
11662if (associated(this%anavar%b)) nanavarb=size(this%anavar%b)
11663if (associated(this%anavar%d)) nanavard=size(this%anavar%d)
11664if (associated(this%anavar%c)) nanavarc=size(this%anavar%c)
11665
11666if (associated(this%anaattr%r)) nanaattrr=size(this%anaattr%r)
11667if (associated(this%anaattr%i)) nanaattri=size(this%anaattr%i)
11668if (associated(this%anaattr%b)) nanaattrb=size(this%anaattr%b)
11669if (associated(this%anaattr%d)) nanaattrd=size(this%anaattr%d)
11670if (associated(this%anaattr%c)) nanaattrc=size(this%anaattr%c)
11671
11672if (associated(this%anavarattr%r)) nanavarattrr=size(this%anavarattr%r)
11673if (associated(this%anavarattr%i)) nanavarattri=size(this%anavarattr%i)
11674if (associated(this%anavarattr%b)) nanavarattrb=size(this%anavarattr%b)
11675if (associated(this%anavarattr%d)) nanavarattrd=size(this%anavarattr%d)
11676if (associated(this%anavarattr%c)) nanavarattrc=size(this%anavarattr%c)
11677
11678write(unit=lunit)ldescription
11679write(unit=lunit)tarray
11680
11681write(unit=lunit)&
11682 nana, ntime, ntimerange, nlevel, nnetwork, &
11683 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11684 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11685 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11686 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11687 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11688 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
11689 this%time_definition
11690
11691
11692!write(unit=lunit)this
11693
11694
11695!! prime 5 dimensioni
11696if (associated(this%ana)) call write_unit(this%ana, lunit)
11697if (associated(this%time)) call write_unit(this%time, lunit)
11698if (associated(this%level)) write(unit=lunit)this%level
11699if (associated(this%timerange)) write(unit=lunit)this%timerange
11700if (associated(this%network)) write(unit=lunit)this%network
11701
11702 !! 6a dimensione: variabile dell'anagrafica e dei dati
11703 !! con relativi attributi e in 5 tipi diversi
11704
11705if (associated(this%anavar%r)) write(unit=lunit)this%anavar%r
11706if (associated(this%anavar%i)) write(unit=lunit)this%anavar%i
11707if (associated(this%anavar%b)) write(unit=lunit)this%anavar%b
11708if (associated(this%anavar%d)) write(unit=lunit)this%anavar%d
11709if (associated(this%anavar%c)) write(unit=lunit)this%anavar%c
11710
11711if (associated(this%anaattr%r)) write(unit=lunit)this%anaattr%r
11712if (associated(this%anaattr%i)) write(unit=lunit)this%anaattr%i
11713if (associated(this%anaattr%b)) write(unit=lunit)this%anaattr%b
11714if (associated(this%anaattr%d)) write(unit=lunit)this%anaattr%d
11715if (associated(this%anaattr%c)) write(unit=lunit)this%anaattr%c
11716
11717if (associated(this%anavarattr%r)) write(unit=lunit)this%anavarattr%r
11718if (associated(this%anavarattr%i)) write(unit=lunit)this%anavarattr%i
11719if (associated(this%anavarattr%b)) write(unit=lunit)this%anavarattr%b
11720if (associated(this%anavarattr%d)) write(unit=lunit)this%anavarattr%d
11721if (associated(this%anavarattr%c)) write(unit=lunit)this%anavarattr%c
11722
11723if (associated(this%dativar%r)) write(unit=lunit)this%dativar%r
11724if (associated(this%dativar%i)) write(unit=lunit)this%dativar%i
11725if (associated(this%dativar%b)) write(unit=lunit)this%dativar%b
11726if (associated(this%dativar%d)) write(unit=lunit)this%dativar%d
11727if (associated(this%dativar%c)) write(unit=lunit)this%dativar%c
11728
11729if (associated(this%datiattr%r)) write(unit=lunit)this%datiattr%r
11730if (associated(this%datiattr%i)) write(unit=lunit)this%datiattr%i
11731if (associated(this%datiattr%b)) write(unit=lunit)this%datiattr%b
11732if (associated(this%datiattr%d)) write(unit=lunit)this%datiattr%d
11733if (associated(this%datiattr%c)) write(unit=lunit)this%datiattr%c
11734
11735if (associated(this%dativarattr%r)) write(unit=lunit)this%dativarattr%r
11736if (associated(this%dativarattr%i)) write(unit=lunit)this%dativarattr%i
11737if (associated(this%dativarattr%b)) write(unit=lunit)this%dativarattr%b
11738if (associated(this%dativarattr%d)) write(unit=lunit)this%dativarattr%d
11739if (associated(this%dativarattr%c)) write(unit=lunit)this%dativarattr%c
11740
11741!! Volumi di valori e attributi per anagrafica e dati
11742
11743if (associated(this%volanar)) write(unit=lunit)this%volanar
11744if (associated(this%volanaattrr)) write(unit=lunit)this%volanaattrr
11745if (associated(this%voldatir)) write(unit=lunit)this%voldatir
11746if (associated(this%voldatiattrr)) write(unit=lunit)this%voldatiattrr
11747
11748if (associated(this%volanai)) write(unit=lunit)this%volanai
11749if (associated(this%volanaattri)) write(unit=lunit)this%volanaattri
11750if (associated(this%voldatii)) write(unit=lunit)this%voldatii
11751if (associated(this%voldatiattri)) write(unit=lunit)this%voldatiattri
11752
11753if (associated(this%volanab)) write(unit=lunit)this%volanab
11754if (associated(this%volanaattrb)) write(unit=lunit)this%volanaattrb
11755if (associated(this%voldatib)) write(unit=lunit)this%voldatib
11756if (associated(this%voldatiattrb)) write(unit=lunit)this%voldatiattrb
11757
11758if (associated(this%volanad)) write(unit=lunit)this%volanad
11759if (associated(this%volanaattrd)) write(unit=lunit)this%volanaattrd
11760if (associated(this%voldatid)) write(unit=lunit)this%voldatid
11761if (associated(this%voldatiattrd)) write(unit=lunit)this%voldatiattrd
11762
11763if (associated(this%volanac)) write(unit=lunit)this%volanac
11764if (associated(this%volanaattrc)) write(unit=lunit)this%volanaattrc
11765if (associated(this%voldatic)) write(unit=lunit)this%voldatic
11766if (associated(this%voldatiattrc)) write(unit=lunit)this%voldatiattrc
11767
11768if (.not. present(unit)) close(unit=lunit)
11769
11770end subroutine vol7d_write_on_file
11771
11772
11779
11780
11781subroutine vol7d_read_from_file (this,unit,filename,description,tarray,filename_auto)
11782
11783TYPE(vol7d),INTENT(OUT) :: this
11784integer,intent(inout),optional :: unit
11785character(len=*),INTENT(in),optional :: filename
11786character(len=*),intent(out),optional :: filename_auto
11787character(len=*),INTENT(out),optional :: description
11788integer,intent(out),optional :: tarray(8)
11789
11790
11791integer :: nana, ntime, ntimerange, nlevel, nnetwork, &
11792 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11793 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11794 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11795 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11796 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11797 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc
11798
11799character(len=254) :: ldescription,lfilename,arg
11800integer :: ltarray(8),lunit,ios
11801logical :: opened,exist
11802
11803
11804call getarg(0,arg)
11805
11806if (.not. present(unit))then
11807 lunit=getunit()
11808else
11809 if (unit==0)then
11810 lunit=getunit()
11811 unit=lunit
11812 else
11813 lunit=unit
11814 end if
11815end if
11816
11817lfilename=trim(arg)//".v7d"
11818if (index(arg,'/',back=.true.) > 0) lfilename=lfilename(index(arg,'/',back=.true.)+1 : )
11819
11820if (present(filename))then
11821 if (filename /= "")then
11822 lfilename=filename
11823 end if
11824end if
11825
11826if (present(filename_auto))filename_auto=lfilename
11827
11828
11829inquire(unit=lunit,opened=opened)
11830IF (.NOT. opened) THEN
11831 inquire(file=lfilename,exist=exist)
11832 IF (.NOT.exist) THEN
11833 CALL l4f_log(l4f_fatal, &
11834 'in vol7d_read_from_file, file does not exists, cannot open')
11835 CALL raise_fatal_error()
11836 ENDIF
11837 OPEN(unit=lunit, file=lfilename, form='UNFORMATTED', access='STREAM', &
11838 status='OLD', action='READ')
11839 CALL l4f_log(l4f_info, 'opened: '//trim(lfilename))
11840end if
11841
11842
11843call init(this)
11844read(unit=lunit,iostat=ios)ldescription
11845
11846if (ios < 0) then ! A negative value indicates that the End of File or End of Record
11847 call vol7d_alloc (this)
11848 call vol7d_alloc_vol (this)
11849 if (present(description))description=ldescription
11850 if (present(tarray))tarray=ltarray
11851 if (.not. present(unit)) close(unit=lunit)
11852end if
11853
11854read(unit=lunit)ltarray
11855
11856CALL l4f_log(l4f_info, 'Reading vol7d from file')
11857CALL l4f_log(l4f_info, 'description: '//trim(ldescription))
11858CALL l4f_log(l4f_info, 'written on '//trim(to_char(ltarray(1)))//' '// &
11859 trim(to_char(ltarray(2)))//' '//trim(to_char(ltarray(3))))
11860
11861if (present(description))description=ldescription
11862if (present(tarray))tarray=ltarray
11863
11864read(unit=lunit)&
11865 nana, ntime, ntimerange, nlevel, nnetwork, &
11866 ndativarr, ndativari, ndativarb, ndativard, ndativarc,&
11867 ndatiattrr, ndatiattri, ndatiattrb, ndatiattrd, ndatiattrc,&
11868 ndativarattrr, ndativarattri, ndativarattrb, ndativarattrd, ndativarattrc,&
11869 nanavarr, nanavari, nanavarb, nanavard, nanavarc,&
11870 nanaattrr, nanaattri, nanaattrb, nanaattrd, nanaattrc,&
11871 nanavarattrr, nanavarattri, nanavarattrb, nanavarattrd, nanavarattrc, &
11872 this%time_definition
11873
11874call vol7d_alloc (this, &
11875 nana=nana, ntime=ntime, ntimerange=ntimerange, nlevel=nlevel, nnetwork=nnetwork,&
11876 ndativarr=ndativarr, ndativari=ndativari, ndativarb=ndativarb,&
11877 ndativard=ndativard, ndativarc=ndativarc,&
11878 ndatiattrr=ndatiattrr, ndatiattri=ndatiattri, ndatiattrb=ndatiattrb,&
11879 ndatiattrd=ndatiattrd, ndatiattrc=ndatiattrc,&
11880 ndativarattrr=ndativarattrr, ndativarattri=ndativarattri, ndativarattrb=ndativarattrb, &
11881 ndativarattrd=ndativarattrd, ndativarattrc=ndativarattrc,&
11882 nanavarr=nanavarr, nanavari=nanavari, nanavarb=nanavarb, &
11883 nanavard=nanavard, nanavarc=nanavarc,&
11884 nanaattrr=nanaattrr, nanaattri=nanaattri, nanaattrb=nanaattrb,&
11885 nanaattrd=nanaattrd, nanaattrc=nanaattrc,&
11886 nanavarattrr=nanavarattrr, nanavarattri=nanavarattri, nanavarattrb=nanavarattrb, &
11887 nanavarattrd=nanavarattrd, nanavarattrc=nanavarattrc)
11888
11889
11890if (associated(this%ana)) call read_unit(this%ana, lunit)
11891if (associated(this%time)) call read_unit(this%time, lunit)
11892if (associated(this%level)) read(unit=lunit)this%level
11893if (associated(this%timerange)) read(unit=lunit)this%timerange
11894if (associated(this%network)) read(unit=lunit)this%network
11895
11896if (associated(this%anavar%r)) read(unit=lunit)this%anavar%r
11897if (associated(this%anavar%i)) read(unit=lunit)this%anavar%i
11898if (associated(this%anavar%b)) read(unit=lunit)this%anavar%b
11899if (associated(this%anavar%d)) read(unit=lunit)this%anavar%d
11900if (associated(this%anavar%c)) read(unit=lunit)this%anavar%c
11901
11902if (associated(this%anaattr%r)) read(unit=lunit)this%anaattr%r
11903if (associated(this%anaattr%i)) read(unit=lunit)this%anaattr%i
11904if (associated(this%anaattr%b)) read(unit=lunit)this%anaattr%b
11905if (associated(this%anaattr%d)) read(unit=lunit)this%anaattr%d
11906if (associated(this%anaattr%c)) read(unit=lunit)this%anaattr%c
11907
11908if (associated(this%anavarattr%r)) read(unit=lunit)this%anavarattr%r
11909if (associated(this%anavarattr%i)) read(unit=lunit)this%anavarattr%i
11910if (associated(this%anavarattr%b)) read(unit=lunit)this%anavarattr%b
11911if (associated(this%anavarattr%d)) read(unit=lunit)this%anavarattr%d
11912if (associated(this%anavarattr%c)) read(unit=lunit)this%anavarattr%c
11913
11914if (associated(this%dativar%r)) read(unit=lunit)this%dativar%r
11915if (associated(this%dativar%i)) read(unit=lunit)this%dativar%i
11916if (associated(this%dativar%b)) read(unit=lunit)this%dativar%b
11917if (associated(this%dativar%d)) read(unit=lunit)this%dativar%d
11918if (associated(this%dativar%c)) read(unit=lunit)this%dativar%c
11919
11920if (associated(this%datiattr%r)) read(unit=lunit)this%datiattr%r
11921if (associated(this%datiattr%i)) read(unit=lunit)this%datiattr%i
11922if (associated(this%datiattr%b)) read(unit=lunit)this%datiattr%b
11923if (associated(this%datiattr%d)) read(unit=lunit)this%datiattr%d
11924if (associated(this%datiattr%c)) read(unit=lunit)this%datiattr%c
11925
11926if (associated(this%dativarattr%r)) read(unit=lunit)this%dativarattr%r
11927if (associated(this%dativarattr%i)) read(unit=lunit)this%dativarattr%i
11928if (associated(this%dativarattr%b)) read(unit=lunit)this%dativarattr%b
11929if (associated(this%dativarattr%d)) read(unit=lunit)this%dativarattr%d
11930if (associated(this%dativarattr%c)) read(unit=lunit)this%dativarattr%c
11931
11932call vol7d_alloc_vol (this)
11933
11934!! Volumi di valori e attributi per anagrafica e dati
11935
11936if (associated(this%volanar)) read(unit=lunit)this%volanar
11937if (associated(this%volanaattrr)) read(unit=lunit)this%volanaattrr
11938if (associated(this%voldatir)) read(unit=lunit)this%voldatir
11939if (associated(this%voldatiattrr)) read(unit=lunit)this%voldatiattrr
11940
11941if (associated(this%volanai)) read(unit=lunit)this%volanai
11942if (associated(this%volanaattri)) read(unit=lunit)this%volanaattri
11943if (associated(this%voldatii)) read(unit=lunit)this%voldatii
11944if (associated(this%voldatiattri)) read(unit=lunit)this%voldatiattri
11945
11946if (associated(this%volanab)) read(unit=lunit)this%volanab
11947if (associated(this%volanaattrb)) read(unit=lunit)this%volanaattrb
11948if (associated(this%voldatib)) read(unit=lunit)this%voldatib
11949if (associated(this%voldatiattrb)) read(unit=lunit)this%voldatiattrb
11950
11951if (associated(this%volanad)) read(unit=lunit)this%volanad
11952if (associated(this%volanaattrd)) read(unit=lunit)this%volanaattrd
11953if (associated(this%voldatid)) read(unit=lunit)this%voldatid
11954if (associated(this%voldatiattrd)) read(unit=lunit)this%voldatiattrd
11955
11956if (associated(this%volanac)) read(unit=lunit)this%volanac
11957if (associated(this%volanaattrc)) read(unit=lunit)this%volanaattrc
11958if (associated(this%voldatic)) read(unit=lunit)this%voldatic
11959if (associated(this%voldatiattrc)) read(unit=lunit)this%voldatiattrc
11960
11961if (.not. present(unit)) close(unit=lunit)
11962
11963end subroutine vol7d_read_from_file
11964
11965
11966! to double precision
11967elemental doubleprecision function doubledatd(voldat,var)
11968doubleprecision,intent(in) :: voldat
11969type(vol7d_var),intent(in) :: var
11970
11971doubledatd=voldat
11972
11973end function doubledatd
11974
11975
11976elemental doubleprecision function doubledatr(voldat,var)
11977real,intent(in) :: voldat
11978type(vol7d_var),intent(in) :: var
11979
11980if (c_e(voldat))then
11981 doubledatr=dble(voldat)
11982else
11983 doubledatr=dmiss
11984end if
11985
11986end function doubledatr
11987
11988
11989elemental doubleprecision function doubledati(voldat,var)
11990integer,intent(in) :: voldat
11991type(vol7d_var),intent(in) :: var
11992
11993if (c_e(voldat)) then
11994 if (c_e(var%scalefactor))then
11995 doubledati=dble(voldat)/10.d0**var%scalefactor
11996 else
11997 doubledati=dble(voldat)
11998 endif
11999else
12000 doubledati=dmiss
12001end if
12002
12003end function doubledati
12004
12005
12006elemental doubleprecision function doubledatb(voldat,var)
12007integer(kind=int_b),intent(in) :: voldat
12008type(vol7d_var),intent(in) :: var
12009
12010if (c_e(voldat)) then
12011 if (c_e(var%scalefactor))then
12012 doubledatb=dble(voldat)/10.d0**var%scalefactor
12013 else
12014 doubledatb=dble(voldat)
12015 endif
12016else
12017 doubledatb=dmiss
12018end if
12019
12020end function doubledatb
12021
12022
12023elemental doubleprecision function doubledatc(voldat,var)
12024CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
12025type(vol7d_var),intent(in) :: var
12026
12027doubledatc = c2d(voldat)
12028if (c_e(doubledatc) .and. c_e(var%scalefactor))then
12029 doubledatc=doubledatc/10.d0**var%scalefactor
12030end if
12031
12032end function doubledatc
12033
12034
12035! to integer
12036elemental integer function integerdatd(voldat,var)
12037doubleprecision,intent(in) :: voldat
12038type(vol7d_var),intent(in) :: var
12039
12040if (c_e(voldat))then
12041 if (c_e(var%scalefactor)) then
12042 integerdatd=nint(voldat*10d0**var%scalefactor)
12043 else
12044 integerdatd=nint(voldat)
12045 endif
12046else
12047 integerdatd=imiss
12048end if
12049
12050end function integerdatd
12051
12052
12053elemental integer function integerdatr(voldat,var)
12054real,intent(in) :: voldat
12055type(vol7d_var),intent(in) :: var
12056
12057if (c_e(voldat))then
12058 if (c_e(var%scalefactor)) then
12059 integerdatr=nint(voldat*10d0**var%scalefactor)
12060 else
12061 integerdatr=nint(voldat)
12062 endif
12063else
12064 integerdatr=imiss
12065end if
12066
12067end function integerdatr
12068
12069
12070elemental integer function integerdati(voldat,var)
12071integer,intent(in) :: voldat
12072type(vol7d_var),intent(in) :: var
12073
12074integerdati=voldat
12075
12076end function integerdati
12077
12078
12079elemental integer function integerdatb(voldat,var)
12080integer(kind=int_b),intent(in) :: voldat
12081type(vol7d_var),intent(in) :: var
12082
12083if (c_e(voldat))then
12084 integerdatb=voldat
12085else
12086 integerdatb=imiss
12087end if
12088
12089end function integerdatb
12090
12091
12092elemental integer function integerdatc(voldat,var)
12093CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
12094type(vol7d_var),intent(in) :: var
12095
12096integerdatc=c2i(voldat)
12097
12098end function integerdatc
12099
12100
12101! to real
12102elemental real function realdatd(voldat,var)
12103doubleprecision,intent(in) :: voldat
12104type(vol7d_var),intent(in) :: var
12105
12106if (c_e(voldat))then
12107 realdatd=real(voldat)
12108else
12109 realdatd=rmiss
12110end if
12111
12112end function realdatd
12113
12114
12115elemental real function realdatr(voldat,var)
12116real,intent(in) :: voldat
12117type(vol7d_var),intent(in) :: var
12118
12119realdatr=voldat
12120
12121end function realdatr
12122
12123
12124elemental real function realdati(voldat,var)
12125integer,intent(in) :: voldat
12126type(vol7d_var),intent(in) :: var
12127
12128if (c_e(voldat)) then
12129 if (c_e(var%scalefactor))then
12130 realdati=float(voldat)/10.**var%scalefactor
12131 else
12132 realdati=float(voldat)
12133 endif
12134else
12135 realdati=rmiss
12136end if
12137
12138end function realdati
12139
12140
12141elemental real function realdatb(voldat,var)
12142integer(kind=int_b),intent(in) :: voldat
12143type(vol7d_var),intent(in) :: var
12144
12145if (c_e(voldat)) then
12146 if (c_e(var%scalefactor))then
12147 realdatb=float(voldat)/10**var%scalefactor
12148 else
12149 realdatb=float(voldat)
12150 endif
12151else
12152 realdatb=rmiss
12153end if
12154
12155end function realdatb
12156
12157
12158elemental real function realdatc(voldat,var)
12159CHARACTER(len=vol7d_cdatalen),intent(in) :: voldat
12160type(vol7d_var),intent(in) :: var
12161
12162realdatc=c2r(voldat)
12163if (c_e(realdatc) .and. c_e(var%scalefactor))then
12164 realdatc=realdatc/10.**var%scalefactor
12165end if
12166
12167end function realdatc
12168
12169
12175FUNCTION realanavol(this, var) RESULT(vol)
12176TYPE(vol7d),INTENT(in) :: this
12177TYPE(vol7d_var),INTENT(in) :: var
12178REAL :: vol(SIZE(this%ana),size(this%network))
12179
12180CHARACTER(len=1) :: dtype
12181INTEGER :: indvar
12182
12183dtype = cmiss
12184indvar = index(this%anavar, var, type=dtype)
12185
12186IF (indvar > 0) THEN
12187 SELECT CASE (dtype)
12188 CASE("d")
12189 vol = realdat(this%volanad(:,indvar,:), var)
12190 CASE("r")
12191 vol = this%volanar(:,indvar,:)
12192 CASE("i")
12193 vol = realdat(this%volanai(:,indvar,:), var)
12194 CASE("b")
12195 vol = realdat(this%volanab(:,indvar,:), var)
12196 CASE("c")
12197 vol = realdat(this%volanac(:,indvar,:), var)
12198 CASE default
12199 vol = rmiss
12200 END SELECT
12201ELSE
12202 vol = rmiss
12203ENDIF
12204
12205END FUNCTION realanavol
12206
12207
12213FUNCTION integeranavol(this, var) RESULT(vol)
12214TYPE(vol7d),INTENT(in) :: this
12215TYPE(vol7d_var),INTENT(in) :: var
12216INTEGER :: vol(SIZE(this%ana),size(this%network))
12217
12218CHARACTER(len=1) :: dtype
12219INTEGER :: indvar
12220
12221dtype = cmiss
12222indvar = index(this%anavar, var, type=dtype)
12223
12224IF (indvar > 0) THEN
12225 SELECT CASE (dtype)
12226 CASE("d")
12227 vol = integerdat(this%volanad(:,indvar,:), var)
12228 CASE("r")
12229 vol = integerdat(this%volanar(:,indvar,:), var)
12230 CASE("i")
12231 vol = this%volanai(:,indvar,:)
12232 CASE("b")
12233 vol = integerdat(this%volanab(:,indvar,:), var)
12234 CASE("c")
12235 vol = integerdat(this%volanac(:,indvar,:), var)
12236 CASE default
12237 vol = imiss
12238 END SELECT
12239ELSE
12240 vol = imiss
12241ENDIF
12242
12243END FUNCTION integeranavol
12244
12245
12251subroutine move_datac (v7d,&
12252 indana,indtime,indlevel,indtimerange,indnetwork,&
12253 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
12254
12255TYPE(vol7d),intent(inout) :: v7d
12256
12257integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
12258integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
12259integer :: inddativar,inddativarattr
12260
12261
12262do inddativar=1,size(v7d%dativar%c)
12263
12264 if (c_e(v7d%voldatic(indana,indtime,indlevel,indtimerange,inddativar,indnetwork)) .and. &
12265 .not. c_e(v7d%voldatic(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
12266 ) then
12267
12268 ! dati
12269 v7d%voldatic &
12270 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
12271 v7d%voldatic &
12272 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
12273
12274
12275 ! attributi
12276 if (associated (v7d%dativarattr%i)) then
12277 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%c(inddativar))
12278 if (inddativarattr > 0 ) then
12279 v7d%voldatiattri &
12280 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12281 v7d%voldatiattri &
12282 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12283 end if
12284 end if
12285
12286 if (associated (v7d%dativarattr%r)) then
12287 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%c(inddativar))
12288 if (inddativarattr > 0 ) then
12289 v7d%voldatiattrr &
12290 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12291 v7d%voldatiattrr &
12292 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12293 end if
12294 end if
12295
12296 if (associated (v7d%dativarattr%d)) then
12297 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%c(inddativar))
12298 if (inddativarattr > 0 ) then
12299 v7d%voldatiattrd &
12300 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12301 v7d%voldatiattrd &
12302 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12303 end if
12304 end if
12305
12306 if (associated (v7d%dativarattr%b)) then
12307 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%c(inddativar))
12308 if (inddativarattr > 0 ) then
12309 v7d%voldatiattrb &
12310 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12311 v7d%voldatiattrb &
12312 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12313 end if
12314 end if
12315
12316 if (associated (v7d%dativarattr%c)) then
12317 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%c(inddativar))
12318 if (inddativarattr > 0 ) then
12319 v7d%voldatiattrc &
12320 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12321 v7d%voldatiattrc &
12322 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12323 end if
12324 end if
12325
12326 end if
12327
12328end do
12329
12330end subroutine move_datac
12331
12337subroutine move_datar (v7d,&
12338 indana,indtime,indlevel,indtimerange,indnetwork,&
12339 indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew)
12340
12341TYPE(vol7d),intent(inout) :: v7d
12342
12343integer,intent(in) :: indana,indtime,indlevel,indtimerange,indnetwork
12344integer,intent(in) :: indananew,indtimenew,indlevelnew,indtimerangenew,indnetworknew
12345integer :: inddativar,inddativarattr
12346
12347
12348do inddativar=1,size(v7d%dativar%r)
12349
12350 if (c_e(v7d%voldatir(indana,indtime,indlevel,indtimerange,inddativar,indnetwork)) .and. &
12351 .not. c_e(v7d%voldatir(indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew))&
12352 ) then
12353
12354 ! dati
12355 v7d%voldatir &
12356 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativar,indnetworknew) = &
12357 v7d%voldatir &
12358 (indana,indtime,indlevel,indtimerange,inddativar,indnetwork)
12359
12360
12361 ! attributi
12362 if (associated (v7d%dativarattr%i)) then
12363 inddativarattr = index(v7d%dativarattr%i,v7d%dativar%r(inddativar))
12364 if (inddativarattr > 0 ) then
12365 v7d%voldatiattri &
12366 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12367 v7d%voldatiattri &
12368 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12369 end if
12370 end if
12371
12372 if (associated (v7d%dativarattr%r)) then
12373 inddativarattr = index(v7d%dativarattr%r,v7d%dativar%r(inddativar))
12374 if (inddativarattr > 0 ) then
12375 v7d%voldatiattrr &
12376 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12377 v7d%voldatiattrr &
12378 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12379 end if
12380 end if
12381
12382 if (associated (v7d%dativarattr%d)) then
12383 inddativarattr = index(v7d%dativarattr%d,v7d%dativar%r(inddativar))
12384 if (inddativarattr > 0 ) then
12385 v7d%voldatiattrd &
12386 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12387 v7d%voldatiattrd &
12388 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12389 end if
12390 end if
12391
12392 if (associated (v7d%dativarattr%b)) then
12393 inddativarattr = index(v7d%dativarattr%b,v7d%dativar%r(inddativar))
12394 if (inddativarattr > 0 ) then
12395 v7d%voldatiattrb &
12396 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12397 v7d%voldatiattrb &
12398 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12399 end if
12400 end if
12401
12402 if (associated (v7d%dativarattr%c)) then
12403 inddativarattr = index(v7d%dativarattr%c,v7d%dativar%r(inddativar))
12404 if (inddativarattr > 0 ) then
12405 v7d%voldatiattrc &
12406 (indananew,indtimenew,indlevelnew,indtimerangenew,inddativarattr,indnetworknew,:) = &
12407 v7d%voldatiattrc &
12408 (indana,indtime,indlevel,indtimerange,inddativarattr,indnetwork,:)
12409 end if
12410 end if
12411
12412 end if
12413
12414end do
12415
12416end subroutine move_datar
12417
12418
12432subroutine v7d_rounding(v7din,v7dout,level,timerange,nostatproc)
12433type(vol7d),intent(inout) :: v7din
12434type(vol7d),intent(out) :: v7dout
12435type(vol7d_level),intent(in),optional :: level(:)
12436type(vol7d_timerange),intent(in),optional :: timerange(:)
12437!logical,intent(in),optional :: merge !< if there are data on more then one almost equal levels or timeranges
12438!! will be merged POINT BY POINT with priority for the fird data found ordered by icreasing var index
12439logical,intent(in),optional :: nostatproc
12440
12441integer :: nana,nlevel,ntime,ntimerange,nnetwork,nbin
12442integer :: iana,ilevel,itimerange,indl,indt,itime,inetwork
12443type(vol7d_level) :: roundlevel(size(v7din%level))
12444type(vol7d_timerange) :: roundtimerange(size(v7din%timerange))
12445type(vol7d) :: v7d_tmp
12446
12447
12448nbin=0
12449
12450if (associated(v7din%dativar%r)) nbin = nbin + size(v7din%dativar%r)
12451if (associated(v7din%dativar%i)) nbin = nbin + size(v7din%dativar%i)
12452if (associated(v7din%dativar%d)) nbin = nbin + size(v7din%dativar%d)
12453if (associated(v7din%dativar%b)) nbin = nbin + size(v7din%dativar%b)
12454
12455call init(v7d_tmp)
12456
12457roundlevel=v7din%level
12458
12459if (present(level))then
12460 do ilevel = 1, size(v7din%level)
12461 if ((any(v7din%level(ilevel) .almosteq. level))) then
12462 roundlevel(ilevel)=level(1)
12463 end if
12464 end do
12465end if
12466
12467roundtimerange=v7din%timerange
12468
12469if (present(timerange))then
12470 do itimerange = 1, size(v7din%timerange)
12471 if ((any(v7din%timerange(itimerange) .almosteq. timerange))) then
12472 roundtimerange(itimerange)=timerange(1)
12473 end if
12474 end do
12475end if
12476
12477!set istantaneous values everywere
12478!preserve p1 for forecast time
12479if (optio_log(nostatproc)) then
12480 roundtimerange(:)%timerange=254
12481 roundtimerange(:)%p2=0
12482end if
12483
12484
12485nana=size(v7din%ana)
12486nlevel=count_distinct(roundlevel,back=.true.)
12487ntime=size(v7din%time)
12488ntimerange=count_distinct(roundtimerange,back=.true.)
12489nnetwork=size(v7din%network)
12490
12491call init(v7d_tmp)
12492
12493if (nbin == 0) then
12494 call copy(v7din,v7d_tmp)
12495else
12496 call vol7d_convr(v7din,v7d_tmp)
12497end if
12498
12499v7d_tmp%level=roundlevel
12500v7d_tmp%timerange=roundtimerange
12501
12502do ilevel=1, size(v7d_tmp%level)
12503 indl=index(v7d_tmp%level,roundlevel(ilevel))
12504 do itimerange=1,size(v7d_tmp%timerange)
12505 indt=index(v7d_tmp%timerange,roundtimerange(itimerange))
12506
12507 if (indl /= ilevel .or. indt /= itimerange) then
12508
12509 do iana=1, nana
12510 do itime=1,ntime
12511 do inetwork=1,nnetwork
12512
12513 if (nbin > 0) then
12514 call move_datar (v7d_tmp,&
12515 iana,itime,ilevel,itimerange,inetwork,&
12516 iana,itime,indl,indt,inetwork)
12517 else
12518 call move_datac (v7d_tmp,&
12519 iana,itime,ilevel,itimerange,inetwork,&
12520 iana,itime,indl,indt,inetwork)
12521 end if
12522
12523 end do
12524 end do
12525 end do
12526
12527 end if
12528
12529 end do
12530end do
12531
12532! set to missing level and time > nlevel
12533do ilevel=nlevel+1,size(v7d_tmp%level)
12534 call init (v7d_tmp%level(ilevel))
12535end do
12536
12537do itimerange=ntimerange+1,size(v7d_tmp%timerange)
12538 call init (v7d_tmp%timerange(itimerange))
12539end do
12540
12541!copy with remove
12542CALL copy(v7d_tmp,v7dout,miss=.true.,lsort_timerange=.true.,lsort_level=.true.)
12543CALL delete(v7d_tmp)
12544
12545!call display(v7dout)
12546
12547end subroutine v7d_rounding
12548
12549
12550END MODULE vol7d_class
12551
12557
12558
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.