libsim Versione 7.2.6
volgrid6d_class.F90
1! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
2! authors:
3! Davide Cesari <dcesari@arpa.emr.it>
4! Paolo Patruno <ppatruno@arpa.emr.it>
5
6! This program is free software; you can redistribute it and/or
7! modify it under the terms of the GNU General Public License as
8! published by the Free Software Foundation; either version 2 of
9! the License, or (at your option) any later version.
10
11! This program is distributed in the hope that it will be useful,
12! but WITHOUT ANY WARRANTY; without even the implied warranty of
13! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14! GNU General Public License for more details.
15
16! You should have received a copy of the GNU General Public License
17! along with this program. If not, see <http://www.gnu.org/licenses/>.
18#include "config.h"
19!> \defgroup volgrid6d Libsim package, volgrid6d library.
20!! The libsim volgrid6d library contains classes for managing data on
21!! regular rectangular grids, tipically the output of a numerical
22!! weather prediction model, and for their import from a WMO GRIB file
23!! of from other geophysical data file formats. In order to compile
24!! and link programs using this library, you have to insert the
25!! required \c USE statements in the program units involved, specify
26!! the location of module files when compiling (tipically \c
27!! -I/usr/lib/gfortran/modules or \c -I/usr/lib64/gfortran/modules or
28!! \c -I/usr/include) and indicate the library name \c -lsim_volgrid6d
29!! when linking, assuming that the library has been installed in a
30!! default location.
31
32!> This module defines objects and methods for managing data volumes
33!! on rectangular georeferenced grids. The data are accomodated in a
34!! multi-dimensional array with 6 predefined dimensions. Different
35!! geographic coordinates and projections are supported, mainly
36!! inspired by grib coding standards. The \a volgrid6d object contains
37!! information and data on an homogeneous grid definition, while
38!! different grids are managed as arrays of \a volgrid6d objects.
39!! Every object contains also an identificator of the grid (\a grid_id
40!! object), carrying information about the driver used or which has to
41!! be used for import/export from/to file. With the help of \a
42!! gridinfo_def class, data can be imported and exported to the
43!! supported formats, mainly grib1 and grib2 through grib_api and many
44!! GIS-style formats through gdal.
45!!
46!! Simple example program \include example_vg6d_3.f90
47!! Example of transformation from volgrid6d to volgrid6d \include example_vg6d_5.f90
48!! Example of transformation from volgrid6d to vol7d \include example_vg6d_6.f90
49!! Example of transformation from vol7d to volgrid6d \include example_vg6d_7.f90
50!!
51!!\ingroup volgrid6d
52MODULE volgrid6d_class
57USE geo_proj_class
58USE grid_class
67!USE file_utilities
68#ifdef HAVE_DBALLE
70#endif
71IMPLICIT NONE
72
73character (len=255),parameter:: subcategory="volgrid6d_class"
74
75!> Object describing a rectangular, homogeneous gridded dataset
76type volgrid6d
77 type(griddim_def) :: griddim !< grid descriptor
78 TYPE(datetime),pointer :: time(:) !< time dimension descriptor
79 TYPE(vol7d_timerange),pointer :: timerange(:) !< timerange (forecast, analysis, statistically processed) dimension descriptor
80 TYPE(vol7d_level),pointer :: level(:) !< vertical level dimension descriptor
81 TYPE(volgrid6d_var),pointer :: var(:) !< physical variable dimension descriptor
82 TYPE(grid_id),POINTER :: gaid(:,:,:,:) !< array of grid identifiers, carrying information about the driver for import/export from/to file, indices are: (level,time,timerange,var)
83 REAL,POINTER :: voldati(:,:,:,:,:,:) !< array of data, indices are: (x,y,level,time,timerange,var)
84 integer :: time_definition !< time definition; 0=time is reference time ; 1=time is validity time
85 integer :: category = 0 !< log4fortran category
86end type volgrid6d
87
88!> Constructor, it creates a new instance of the object.
89INTERFACE init
90 MODULE PROCEDURE volgrid6d_init
91END INTERFACE
92
93!> Destructor, it releases every information and memory buffer
94!! associated with the object.
95INTERFACE delete
96 MODULE PROCEDURE volgrid6d_delete, volgrid6dv_delete
97END INTERFACE
98
99!> Import an object dirctly from a native file, from a \a gridinfo object
100!! or from a supported file format through a \a gridinfo object.
101INTERFACE import
102 MODULE PROCEDURE volgrid6d_read_from_file
103 MODULE PROCEDURE import_from_gridinfo, import_from_gridinfovv, &
104 volgrid6d_import_from_file
105END INTERFACE
106
107!> Export an object dirctly to a native file, to a \a gridinfo object
108!! or to a supported file format through a \a gridinfo object.
109INTERFACE export
110 MODULE PROCEDURE volgrid6d_write_on_file
111 MODULE PROCEDURE export_to_gridinfo, export_to_gridinfov, export_to_gridinfovv,&
112 volgrid6d_export_to_file
113END INTERFACE
114
115! methods for computing transformations through an initialised
116! grid_transform object, probably too low level to be interfaced
117INTERFACE compute
118 MODULE PROCEDURE volgrid6d_transform_compute, volgrid6d_v7d_transform_compute,&
119 v7d_volgrid6d_transform_compute, v7d_v7d_transform_compute
120END INTERFACE
121
122!> Transform between any combination of \a volgrid6d and \a vol7d objects
123!! by means of a \a transform_def object describing the transformation.
124INTERFACE transform
125 MODULE PROCEDURE volgrid6d_transform, volgrid6dv_transform,&
126 volgrid6d_v7d_transform, volgrid6dv_v7d_transform, v7d_volgrid6d_transform, &
127 v7d_v7d_transform
128END INTERFACE
129
130INTERFACE wind_rot
131 MODULE PROCEDURE vg6d_wind_rot
132END INTERFACE
133
134INTERFACE wind_unrot
135 MODULE PROCEDURE vg6d_wind_unrot
136END INTERFACE
137
138!> Display on standard output a description of the \a volgrid6d object
139!! provided.
140INTERFACE display
141 MODULE PROCEDURE display_volgrid6d,display_volgrid6dv
142END INTERFACE
143
144!> Reduce some dimensions (level and timerage) for semplification (rounding).
145!! You can use this for simplify and use variables in computation like alchimia
146!! where fields have to be on the same coordinate
147!! examples:
148!! means in time for short periods and istantaneous values
149!! 2 meter and surface levels
150!! If there are data on more then one almost equal levels or timeranges, the first var present (at least one point)
151!! will be taken (order is by icreasing var index).
152!! You can use predefined values for classic semplification
153!! almost_equal_levels and almost_equal_timeranges
154!! The level or timerange in output will be defined by the first element of level and timerange list
155INTERFACE rounding
156 MODULE PROCEDURE vg6d_rounding, vg6dv_rounding
157END INTERFACE
158
159private
160
161PUBLIC volgrid6d,init,delete,export,import,compute,transform, &
162 wind_rot,wind_unrot,vg6d_c2a,display,volgrid6d_alloc,volgrid6d_alloc_vol, &
163 volgrid_get_vol_2d, volgrid_set_vol_2d, volgrid_get_vol_3d, volgrid_set_vol_3d
164PUBLIC rounding, vg6d_reduce
165
166CONTAINS
167
168
169!> Constructor, it creates a new instance of the object.
170!! The constructor should be explicitly used only in rare cases,
171!! \a volgrid6d objects are usually created through the \a import
172!! interface.
173SUBROUTINE volgrid6d_init(this, griddim, time_definition, categoryappend)
174TYPE(volgrid6d) :: this !< object to be initialized
175TYPE(griddim_def),OPTIONAL :: griddim !< grid descriptor
176INTEGER,INTENT(IN),OPTIONAL :: time_definition !< 0=time is reference time; 1=time is validity time
177CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
178
179character(len=512) :: a_name
180
181if (present(categoryappend))then
182 call l4f_launcher(a_name,a_name_append=trim(subcategory)//"."//trim(categoryappend))
183else
184 call l4f_launcher(a_name,a_name_append=trim(subcategory))
185endif
186this%category=l4f_category_get(a_name)
187
188#ifdef DEBUG
189call l4f_category_log(this%category,l4f_debug,"init")
190#endif
191
192call init(this%griddim)
193
194if (present(griddim))then
195 call copy (griddim,this%griddim)
196end if
197
198CALL vol7d_var_features_init() ! initialise var features table once
199
200if(present(time_definition)) then
201 this%time_definition = time_definition
202else
203 this%time_definition = 0 !default to reference time
204end if
205
206nullify (this%time,this%timerange,this%level,this%var)
207nullify (this%gaid,this%voldati)
208
209END SUBROUTINE volgrid6d_init
210
211
212!> Allocate the dimension descriptors of the \a volgrid6d object.
213!! This method allocates the horizontal grid descriptor and the one
214!! dimensional arrays of the dimensions
215!! - time
216!! - vertical level
217!! - timerange
218!! - physical variable
219!!
220!! This method should be explicitly used only in rare cases, it is
221!! usually called implicitly through the \a import interface.
222SUBROUTINE volgrid6d_alloc(this, dim, ntime, nlevel, ntimerange, nvar, ini)
223TYPE(volgrid6d),INTENT(inout) :: this !< object whose decriptors should be allocated
224TYPE(grid_dim),INTENT(in),OPTIONAL :: dim !< horizontal grid size X, Y
225INTEGER,INTENT(in),OPTIONAL :: ntime !< number of time levels
226INTEGER,INTENT(in),OPTIONAL :: nlevel !< number of vertical levels
227INTEGER,INTENT(in),OPTIONAL :: ntimerange !< number of different timeranges
228INTEGER,INTENT(in),OPTIONAL :: nvar !< number of physical variables
229LOGICAL,INTENT(in),OPTIONAL :: ini !< if provided and \c .TRUE., for each allocated dimension descriptor the constructor is called without extra parameters, thus initializing everything as missing value
230
231INTEGER :: i, stallo
232LOGICAL :: linit
233
234#ifdef DEBUG
235call l4f_category_log(this%category,l4f_debug,"alloc")
236#endif
237
238IF (PRESENT(ini)) THEN
239 linit = ini
240ELSE
241 linit = .false.
242ENDIF
243
244
245if (present(dim)) call copy (dim,this%griddim%dim)
246
247
248IF (PRESENT(ntime)) THEN
249 IF (ntime >= 0) THEN
250 IF (ASSOCIATED(this%time)) DEALLOCATE(this%time)
251#ifdef DEBUG
252 call l4f_category_log(this%category,l4f_debug,"alloc ntime "//to_char(ntime))
253#endif
254 ALLOCATE(this%time(ntime),stat=stallo)
255 if (stallo /=0)then
256 call l4f_category_log(this%category,l4f_fatal,"allocating memory")
257 CALL raise_fatal_error()
258 end if
259 IF (linit) THEN
260 DO i = 1, ntime
261 this%time(i) = datetime_miss
262! CALL init(this%time(i)) ! senza argomento inizializza a zero non missing
263 ! baco di datetime?
264 ENDDO
265 ENDIF
266 ENDIF
267ENDIF
268IF (PRESENT(nlevel)) THEN
269 IF (nlevel >= 0) THEN
270 IF (ASSOCIATED(this%level)) DEALLOCATE(this%level)
271#ifdef DEBUG
272 call l4f_category_log(this%category,l4f_debug,"alloc nlevel "//to_char(nlevel))
273#endif
274 ALLOCATE(this%level(nlevel),stat=stallo)
275 if (stallo /=0)then
276 call l4f_category_log(this%category,l4f_fatal,"allocating memory")
277 CALL raise_fatal_error()
278 end if
279 IF (linit) THEN
280 DO i = 1, nlevel
281 CALL init(this%level(i))
282 ENDDO
283 ENDIF
284 ENDIF
285ENDIF
286IF (PRESENT(ntimerange)) THEN
287 IF (ntimerange >= 0) THEN
288 IF (ASSOCIATED(this%timerange)) DEALLOCATE(this%timerange)
289#ifdef DEBUG
290 call l4f_category_log(this%category,l4f_debug,"alloc ntimerange "//to_char(ntimerange))
291#endif
292 ALLOCATE(this%timerange(ntimerange),stat=stallo)
293 if (stallo /=0)then
294 call l4f_category_log(this%category,l4f_fatal,"allocating memory")
295 CALL raise_fatal_error()
296 end if
297 IF (linit) THEN
298 DO i = 1, ntimerange
299 CALL init(this%timerange(i))
300 ENDDO
301 ENDIF
302 ENDIF
303ENDIF
304IF (PRESENT(nvar)) THEN
305 IF (nvar >= 0) THEN
306 IF (ASSOCIATED(this%var)) DEALLOCATE(this%var)
307#ifdef DEBUG
308 call l4f_category_log(this%category,l4f_debug,"alloc nvar "//to_char(nvar))
309#endif
310 ALLOCATE(this%var(nvar),stat=stallo)
311 if (stallo /=0)then
312 call l4f_category_log(this%category,l4f_fatal,"allocating memory")
313 CALL raise_fatal_error()
314 end if
315 IF (linit) THEN
316 DO i = 1, nvar
317 CALL init(this%var(i))
318 ENDDO
319 ENDIF
320 ENDIF
321ENDIF
322
323end SUBROUTINE volgrid6d_alloc
324
325
326!> Allocate the data array of the \a volgrid6d object.
327!! This method allocates the main 6-dimensional data array
328!! \a this%voldati and the 4-dimensional \a grid_id array \a this%gaid
329!! with a shape dictated by the previous call(s) to \a vol7d_alloc().
330!! if any descriptor (except horizontal grid) has not been allocated
331!! yet, it is allocated here with a size of 1. This method should be
332!! explicitly used only in rare cases, it is usually called implicitly
333!! through the \a import interface.
334SUBROUTINE volgrid6d_alloc_vol(this, ini, inivol, decode)
335TYPE(volgrid6d),INTENT(inout) :: this !< object whose decriptors should be allocated
336LOGICAL,INTENT(in),OPTIONAL :: ini !< if provided and \c .TRUE., for each dimension descriptor not yet allocated and allocated here the constructor is called without extra parameters, thus initializing the element as missing value
337LOGICAL,INTENT(in),OPTIONAL :: inivol !< if provided and \c .FALSE., the allocated volumes will not be initialized to missing values
338LOGICAL,INTENT(in),OPTIONAL :: decode !< if provided and \c .TRUE., the \a this%voldati volume is allocated, otherwise only \a this%gaid will be allocated
339
340INTEGER :: stallo
341LOGICAL :: linivol
342
343#ifdef DEBUG
344call l4f_category_log(this%category,l4f_debug,"start alloc_vol")
345#endif
346
347IF (PRESENT(inivol)) THEN ! opposite condition, cannot use optio_log
348 linivol = inivol
349ELSE
350 linivol = .true.
351ENDIF
352
353IF (this%griddim%dim%nx > 0 .AND. this%griddim%dim%ny > 0) THEN
354! allocate dimension descriptors to a minimal size for having a
355! non-null volume
356 IF (.NOT. ASSOCIATED(this%var)) CALL volgrid6d_alloc(this, nvar=1, ini=ini)
357 IF (.NOT. ASSOCIATED(this%time)) CALL volgrid6d_alloc(this, ntime=1, ini=ini)
358 IF (.NOT. ASSOCIATED(this%level)) CALL volgrid6d_alloc(this, nlevel=1, ini=ini)
359 IF (.NOT. ASSOCIATED(this%timerange)) CALL volgrid6d_alloc(this, ntimerange=1, ini=ini)
360
361 IF (optio_log(decode)) THEN
362 IF (.NOT.ASSOCIATED(this%voldati)) THEN
363#ifdef DEBUG
364 CALL l4f_category_log(this%category,l4f_debug,"alloc voldati volume")
365#endif
366
367 ALLOCATE(this%voldati(this%griddim%dim%nx,this%griddim%dim%ny,&
368 SIZE(this%level), SIZE(this%time), &
369 SIZE(this%timerange), SIZE(this%var)),stat=stallo)
370 IF (stallo /= 0)THEN
371 CALL l4f_category_log(this%category,l4f_fatal,"allocating memory")
372 CALL raise_fatal_error()
373 ENDIF
374
375! this is not really needed if we can check other flags for a full
376! field missing values
377 IF (linivol) this%voldati = rmiss
378 this%voldati = rmiss
379 ENDIF
380 ENDIF
381
382 IF (.NOT.ASSOCIATED(this%gaid)) THEN
383#ifdef DEBUG
384 CALL l4f_category_log(this%category,l4f_debug,"alloc gaid volume")
385#endif
386 ALLOCATE(this%gaid(SIZE(this%level), SIZE(this%time), &
387 SIZE(this%timerange), SIZE(this%var)),stat=stallo)
388 IF (stallo /= 0)THEN
389 CALL l4f_category_log(this%category,l4f_fatal,"allocating memory")
390 CALL raise_fatal_error()
391 ENDIF
392
393 IF (linivol) THEN
394!!$ DO i=1 ,SIZE(this%gaid,1)
395!!$ DO ii=1 ,SIZE(this%gaid,2)
396!!$ DO iii=1 ,SIZE(this%gaid,3)
397!!$ DO iiii=1 ,SIZE(this%gaid,4)
398!!$ this%gaid(i,ii,iii,iiii) = grid_id_new() ! optimize?
399!!$ ENDDO
400!!$ ENDDO
401!!$ ENDDO
402!!$ ENDDO
403
404 this%gaid = grid_id_new()
405 ENDIF
406 ENDIF
407
408ELSE
409 CALL l4f_category_log(this%category,l4f_fatal,'volgrid6d_alloc_vol: &
410 &trying to allocate a volume with an invalid or unspecified horizontal grid')
411 CALL raise_fatal_error()
412ENDIF
413
414END SUBROUTINE volgrid6d_alloc_vol
415
416
417!> Return a 2-d pointer to a x-y slice of a volume. This method works
418!! both with volumes having allocated and non-allocated this%voldati
419!! array, and it returns a pointer to a 2-d slice either from the
420!! allocated this%voldati array or from the grid_id object on file or
421!! in memory. In the second case the pointer should be either
422!! ALLOCATE'd to the expected size or NULLIFY'ed, and if NULLIFY'ed,
423!! it is allocated within the method, thus it will have to be
424!! deallocated by the caller when not in use anymore. Since this
425!! method may be called many times by a program, it is optimized for
426!! speed and it does not make any check about the matching size of the
427!! pointer and the array or about the allocation status of \a this, so
428!! it should be called only when everything has been checked to be in
429!! good shape.
430SUBROUTINE volgrid_get_vol_2d(this, ilevel, itime, itimerange, ivar, voldati)
431TYPE(volgrid6d),INTENT(in) :: this !< object from which the slice has to be retrieved
432INTEGER,INTENT(in) :: ilevel !< index of vertical level of the slice
433INTEGER,INTENT(in) :: itime !< index of time level of the slice
434INTEGER,INTENT(in) :: itimerange !< index of timerange of the slice
435INTEGER,INTENT(in) :: ivar !< index of physical variable of the slice
436REAL,POINTER :: voldati(:,:) !< pointer to the data, if \a this%voldati is already allocated, it will just point to the requested slice, otherwise it will be allocated if and only if it is nullified on entry
437
438IF (ASSOCIATED(this%voldati)) THEN
439 voldati => this%voldati(:,:,ilevel,itime,itimerange,ivar)
440 RETURN
441ELSE
442 IF (.NOT.ASSOCIATED(voldati)) THEN
443 ALLOCATE(voldati(this%griddim%dim%nx,this%griddim%dim%ny))
444 ENDIF
445 CALL grid_id_decode_data(this%gaid(ilevel,itime,itimerange,ivar), voldati)
446ENDIF
447
448END SUBROUTINE volgrid_get_vol_2d
449
450
451!> Return a 3-d pointer to a x-y-z slice of a volume. This method works
452!! both with volumes having allocated and non-allocated this%voldati
453!! array, and it returns a pointer to a 3-d slice either from the
454!! allocated this%voldati array or from the grid_id object on file or
455!! in memory. In the second case the pointer should be either
456!! ALLOCATE'd to the expected size or NULLIFY'ed, and if NULLIFY'ed,
457!! it is allocated within the method, thus it will have to be
458!! deallocated by the caller when not in use anymore. Since this
459!! method may be called many times by a program, it is optimized for
460!! speed and it does not make any check about the matching size of the
461!! pointer and the array or about the allocation status of \a this, so
462!! it should be called only when everything has been checked to be in
463!! good shape.
464SUBROUTINE volgrid_get_vol_3d(this, itime, itimerange, ivar, voldati)
465TYPE(volgrid6d),INTENT(in) :: this !< object from which the slice has to be retrieved
466INTEGER,INTENT(in) :: itime !< index of time level of the slice
467INTEGER,INTENT(in) :: itimerange !< index of timerange of the slice
468INTEGER,INTENT(in) :: ivar !< index of physical variable of the slice
469REAL,POINTER :: voldati(:,:,:) !< pointer to the data, if \a this%voldati is already allocated, it will just point to the requested slice, otherwise it will be allocated if and only if it is nullified on entry
470
471INTEGER :: ilevel
472
473IF (ASSOCIATED(this%voldati)) THEN
474 voldati => this%voldati(:,:,:,itime,itimerange,ivar)
475 RETURN
476ELSE
477 IF (.NOT.ASSOCIATED(voldati)) THEN
478 ALLOCATE(voldati(this%griddim%dim%nx,this%griddim%dim%ny,SIZE(this%level)))
479 ENDIF
480!$OMP PARALLEL DEFAULT(SHARED)
481!$OMP MASTER
482 DO ilevel = 1, SIZE(this%level)
483!$OMP TASK FIRSTPRIVATE(ilevel)
484 CALL grid_id_decode_data(this%gaid(ilevel,itime,itimerange,ivar), &
485 voldati(:,:,ilevel))
486!$OMP END TASK
487 ENDDO
488!$OMP END MASTER
489!$OMP END PARALLEL
490ENDIF
491
492END SUBROUTINE volgrid_get_vol_3d
493
494
495!> Reset a 2-d x-y slice of a volume after the data have been modified.
496!! This method works both with volumes having allocated and
497!! non-allocated this%voldati array, and it updates the requested
498!! slice. In case \a this%voldati is already allocated, this is a
499!! no-operation while in the other case this method encodes the field
500!! provided into the grid_id object on file or in memory. Since this
501!! method may be called many times by a program, it is optimized for
502!! speed and it does not make any check about the matching size of the
503!! field and the array or about the allocation status of \a this, so
504!! it should be called only when everything has been checked to be in
505!! good shape.
506SUBROUTINE volgrid_set_vol_2d(this, ilevel, itime, itimerange, ivar, voldati)
507TYPE(volgrid6d),INTENT(inout) :: this !< object in which slice has to be updated
508INTEGER,INTENT(in) :: ilevel !< index of vertical level of the slice
509INTEGER,INTENT(in) :: itime !< index of time level of the slice
510INTEGER,INTENT(in) :: itimerange !< index of timerange of the slice
511INTEGER,INTENT(in) :: ivar !< index of physical variable of the slice
512REAL,INTENT(in) :: voldati(:,:) !< updated values of the slice
513
514IF (ASSOCIATED(this%voldati)) THEN
515 RETURN
516ELSE
517 CALL grid_id_encode_data(this%gaid(ilevel,itime,itimerange,ivar), voldati)
518ENDIF
519
520END SUBROUTINE volgrid_set_vol_2d
521
523!> Reset a 3-d x-y-z slice of a volume after the data have been modified.
524!! This method works both with volumes having allocated and
525!! non-allocated this%voldati array, and it updates the requested
526!! slice. In case \a this%voldati is already allocated, this is a
527!! no-operation while in the other case this method encodes the field
528!! provided into the grid_id object on file or in memory. Since this
529!! method may be called many times by a program, it is optimized for
530!! speed and it does not make any check about the matching size of the
531!! field and the array or about the allocation status of \a this, so
532!! it should be called only when everything has been checked to be in
533!! good shape.
534SUBROUTINE volgrid_set_vol_3d(this, itime, itimerange, ivar, voldati)
535TYPE(volgrid6d),INTENT(inout) :: this !< object in which slice has to be updated
536INTEGER,INTENT(in) :: itime !< index of time level of the slice
537INTEGER,INTENT(in) :: itimerange !< index of timerange of the slice
538INTEGER,INTENT(in) :: ivar !< index of physical variable of the slice
539REAL,INTENT(in) :: voldati(:,:,:) !< updated values of the slice
540
541INTEGER :: ilevel
542
543IF (ASSOCIATED(this%voldati)) THEN
544 RETURN
545ELSE
546!$OMP PARALLEL DEFAULT(SHARED)
547!$OMP MASTER
548 DO ilevel = 1, SIZE(this%level)
549!$OMP TASK FIRSTPRIVATE(ilevel)
550 CALL grid_id_encode_data(this%gaid(ilevel,itime,itimerange,ivar), &
551 voldati(:,:,ilevel))
552!$OMP END TASK
553 ENDDO
554!$OMP END MASTER
555!$OMP END PARALLEL
556ENDIF
557
558END SUBROUTINE volgrid_set_vol_3d
559
560
561!> Destructor, it releases every information and memory buffer
562!! associated with the object. It should be called also for objects
563!! crated through the \a import interface.
564SUBROUTINE volgrid6d_delete(this)
565TYPE(volgrid6d),INTENT(inout) :: this
566
567INTEGER :: i, ii, iii, iiii
568
569#ifdef DEBUG
570call l4f_category_log(this%category,l4f_debug,"delete")
571#endif
572
573if (associated(this%gaid))then
574
575 DO i=1 ,SIZE(this%gaid,1)
576 DO ii=1 ,SIZE(this%gaid,2)
577 DO iii=1 ,SIZE(this%gaid,3)
578 DO iiii=1 ,SIZE(this%gaid,4)
579 CALL delete(this%gaid(i,ii,iii,iiii))
580 ENDDO
581 ENDDO
582 ENDDO
583 ENDDO
584 DEALLOCATE(this%gaid)
585
586end if
587
588call delete(this%griddim)
589
590! call delete(this%time)
591! call delete(this%timerange)
592! call delete(this%level)
593! call delete(this%var)
594
595if (associated( this%time )) deallocate(this%time)
596if (associated( this%timerange )) deallocate(this%timerange)
597if (associated( this%level )) deallocate(this%level)
598if (associated( this%var )) deallocate(this%var)
599
600if (associated(this%voldati))deallocate(this%voldati)
601
602
603 !chiudo il logger
604call l4f_category_delete(this%category)
605
606END SUBROUTINE volgrid6d_delete
607
608
609!>\brief Scrittura su file di un volume Volgrid6d.
610!! Scrittura su file unformatted di un intero volume Volgrid6d.
611!! Il volume viene serializzato e scritto su file.
612!! Il file puo' essere aperto opzionalmente dall'utente. Si possono passare
613!! opzionalmente unità e nome del file altrimenti assegnati internamente a dei default; se assegnati internamente
614!! tali parametri saranno in output.
615!! Se non viene fornito il nome file viene utilizzato un file di default con nome pari al nome del programma in
616!! esecuzione con postfisso ".vg6d".
617!! Come parametro opzionale c'è la description che insieme alla data corrente viene inserita nell'header del file.
618subroutine volgrid6d_write_on_file (this,unit,description,filename,filename_auto)
619
620TYPE(volgrid6d),INTENT(IN) :: this !< volume volgrid6d da scrivere
621integer,optional,intent(inout) :: unit !< unità su cui scrivere; se passata =0 ritorna il valore rielaborato (default =rielaborato internamente con getlun )
622character(len=*),intent(in),optional :: filename !< nome del file su cui scrivere
623character(len=*),intent(out),optional :: filename_auto !< nome del file generato se "filename" è omesso
624character(len=*),INTENT(IN),optional :: description !< descrizione del volume
625
626integer :: lunit
627character(len=254) :: ldescription,arg,lfilename
628integer :: ntime, ntimerange, nlevel, nvar
629integer :: tarray(8)
630logical :: opened,exist
631
632#ifdef DEBUG
633call l4f_category_log(this%category,l4f_debug,"write on file")
634#endif
635
636ntime=0
637ntimerange=0
638nlevel=0
639nvar=0
640
641!call idate(im,id,iy)
642call date_and_time(values=tarray)
643call getarg(0,arg)
644
645if (present(description))then
646 ldescription=description
647else
648 ldescription="Volgrid6d generated by: "//trim(arg)
649end if
650
651if (.not. present(unit))then
652 lunit=getunit()
653else
654 if (unit==0)then
655 lunit=getunit()
656 unit=lunit
657 else
658 lunit=unit
659 end if
660end if
661
662lfilename=trim(arg)//".vg6d"
663if (index(arg,'/',back=.true.) > 0) lfilename=lfilename(index(arg,'/',back=.true.)+1 : )
664
665if (present(filename))then
666 if (filename /= "")then
667 lfilename=filename
668 end if
669end if
670
671if (present(filename_auto))filename_auto=lfilename
672
673
674inquire(unit=lunit,opened=opened)
675if (.not. opened) then
676 inquire(file=lfilename,exist=exist)
677 if (exist) CALL raise_error('file exist; cannot open new file')
678 if (.not.exist) open (unit=lunit,file=lfilename,form="UNFORMATTED")
679 !print *, "opened: ",lfilename
680end if
681
682if (associated(this%time)) ntime=size(this%time)
683if (associated(this%timerange)) ntimerange=size(this%timerange)
684if (associated(this%level)) nlevel=size(this%level)
685if (associated(this%var)) nvar=size(this%var)
686
687
688write(unit=lunit)ldescription
689write(unit=lunit)tarray
690
691call write_unit( this%griddim,lunit)
692write(unit=lunit) ntime, ntimerange, nlevel, nvar
693
694!! prime 4 dimensioni
695if (associated(this%time)) call write_unit(this%time, lunit)
696if (associated(this%level)) write(unit=lunit)this%level
697if (associated(this%timerange)) write(unit=lunit)this%timerange
698if (associated(this%var)) write(unit=lunit)this%var
699
700
701!! Volumi di valori dati
702
703if (associated(this%voldati)) write(unit=lunit)this%voldati
704
705if (.not. present(unit)) close(unit=lunit)
706
707end subroutine volgrid6d_write_on_file
708
709
710!>\brief Lettura da file di un volume Volgrid6d.
711!! Lettura da file unformatted di un intero volume Volgrid6d.
712!! Questa subroutine comprende volgrid6d_alloc e volgrid6d_alloc_vol.
713!! Il file puo' essere aperto opzionalmente dall'utente. Si possono passare
714!! opzionalmente unità e nome del file altrimenti assegnati internamente a dei default; se assegnati internamente
715!! tali parametri saranno in output.
716subroutine volgrid6d_read_from_file (this,unit,filename,description,tarray,filename_auto)
717
718TYPE(volgrid6d),INTENT(OUT) :: this !< Volume volgrid6d da leggere
719integer,intent(inout),optional :: unit !< unità su cui è stato aperto un file; se =0 rielaborato internamente (default = elaborato internamente con getunit)
720character(len=*),INTENT(in),optional :: filename !< nome del file eventualmente da aprire (default = nome dell'eseguibile.v7d)
721character(len=*),intent(out),optional :: filename_auto !< nome del file generato se "filename" è omesso
722character(len=*),INTENT(out),optional :: description !< descrizione del volume letto
723integer,intent(out),optional :: tarray(8) !< vettore come definito da "date_and_time" della data di scrittura del volume
724
725integer :: ntime, ntimerange, nlevel, nvar
726
727character(len=254) :: ldescription,lfilename,arg
728integer :: ltarray(8),lunit
729logical :: opened,exist
730
731#ifdef DEBUG
732call l4f_category_log(this%category,l4f_debug,"read from file")
733#endif
734
735call getarg(0,arg)
736
737if (.not. present(unit))then
738 lunit=getunit()
739else
740 if (unit==0)then
741 lunit=getunit()
742 unit=lunit
743 else
744 lunit=unit
745 end if
746end if
747
748lfilename=trim(arg)//".vg6d"
749if (index(arg,'/',back=.true.) > 0) lfilename=lfilename(index(arg,'/',back=.true.)+1 : )
750
751if (present(filename))then
752 if (filename /= "")then
753 lfilename=filename
754 end if
755end if
756
757if (present(filename_auto))filename_auto=lfilename
758
759
760inquire(unit=lunit,opened=opened)
761if (.not. opened) then
762 inquire(file=lfilename,exist=exist)
763 IF (.NOT. exist) CALL raise_fatal_error('file '//trim(lfilename)//' does not exist, cannot open')
764 open (unit=lunit,file=lfilename,form="UNFORMATTED")
765end if
766
767read(unit=lunit)ldescription
768read(unit=lunit)ltarray
769
770call l4f_log(l4f_info,"Info: reading volgrid6d from file: "//trim(lfilename))
771call l4f_log(l4f_info,"Info: description: "//trim(ldescription))
772!call l4f_log("Info: written on ",ltarray)
773
774if (present(description))description=ldescription
775if (present(tarray))tarray=ltarray
776
777
778call read_unit( this%griddim,lunit)
779read(unit=lunit) ntime, ntimerange, nlevel, nvar
780
781
782call volgrid6d_alloc (this, &
783 ntime=ntime, ntimerange=ntimerange, nlevel=nlevel, nvar=nvar)
784
785call volgrid6d_alloc_vol (this)
786
787if (associated(this%time)) call read_unit(this%time, lunit)
788if (associated(this%level)) read(unit=lunit)this%level
789if (associated(this%timerange)) read(unit=lunit)this%timerange
790if (associated(this%var)) read(unit=lunit)this%var
791
792
793!! Volumi di valori
794
795if (associated(this%voldati)) read(unit=lunit)this%voldati
796
797if (.not. present(unit)) close(unit=lunit)
798
799end subroutine volgrid6d_read_from_file
800
801
802!> Import a single \a gridinfo object into a \a volgrid6d object.
803!! This methods imports a single gridded field from a \a gridinfo
804!! object into a \a volgrid6d object, inserting it into the
805!! multidimensional structure of volgrid6d. The volgrid6d object must
806!! have been already initialized and the dimensions specified with
807!! volgrid6d_alloc(). If the \a force argument is missing or \a
808!! .FALSE. , the volgrid6d object dimension descriptors (time,
809!! timerange, vertical level, physical variable) must already have
810!! space for the corresponding values coming from gridinfo, otherwise
811!! the object will be rejected; this means that all the volgrid6d
812!! dimension descriptors should be correctly assigned. If \a force is
813!! \a .TRUE. , the gridinfo dimension descriptors that do not fit into
814!! available descriptors in the \a volgrid6d structure, will be
815!! accomodated in a empty (i.e. equal to missing value) descriptor, if
816!! available, otherwise the gridinfo will be rejected. The
817!! descriptor of the grid in the \a volgrid object is assigned to the
818!! descriptor contained in \a gridinfo if it is missing in \a volgrid,
819!! otherwise it is checked and the object is rejected if grids do not
820!! match.
821SUBROUTINE import_from_gridinfo(this, gridinfo, force, dup_mode, clone, &
822 isanavar)
823TYPE(volgrid6d),INTENT(inout) :: this !< object in which to import
824TYPE(gridinfo_def),INTENT(in) :: gridinfo !< gridinfo object to be imported
825LOGICAL,INTENT(in),OPTIONAL :: force !< if provided and \c .TRUE., the gridinfo is forced into an empty element of \a this, if required and possible
826INTEGER,INTENT(in),OPTIONAL :: dup_mode !< determines the behavior in case of duplicate metadata: if \a dup_mode is not provided or 0, a duplicate field overwrites, if \a dup_mode is 1, duplicate fields are merged with priority to the last
827LOGICAL , INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE. , clone the gaid's from \a gridinfo to \a this
828LOGICAL,INTENT(IN),OPTIONAL :: isanavar !< if provides and \a .TRUE., the gridinfo object is treated as time-independent and replicated for every time and timerange
829
830CHARACTER(len=255) :: type
831INTEGER :: itime0, itimerange0, itime1, itimerange1, itime, itimerange, &
832 ilevel, ivar, ldup_mode
833LOGICAL :: dup
834TYPE(datetime) :: correctedtime
835TYPE(vol7d_timerange) :: correctedtimerange
836REAL,ALLOCATABLE :: tmpgrid(:,:)
837
838IF (PRESENT(dup_mode)) THEN
839 ldup_mode = dup_mode
840ELSE
841 ldup_mode = 0
842ENDIF
843
844call get_val(this%griddim,proj_type=type)
845
846#ifdef DEBUG
847call l4f_category_log(this%category,l4f_debug,"import_from_gridinfo: "//trim(type))
848#endif
849
850if (.not. c_e(type))then
851 call copy(gridinfo%griddim, this%griddim)
852! ho gia` fatto init, altrimenti non potrei fare la get_val(this%griddim)
853! per cui meglio non ripetere
854! call init(this,gridinfo%griddim,categoryappend)
855 CALL volgrid6d_alloc_vol(this, ini=.true.) ! decode?
856
857else if (.not. (this%griddim == gridinfo%griddim ))then
858
859 CALL l4f_category_log(this%category,l4f_error, &
860 "volgrid and gridinfo grid type or size are different, gridinfo rejected")
861 CALL raise_error()
862 RETURN
863
864end if
865
866! Cerco gli indici del campo da inserire, se non trovo metto nel primo missing
867ilevel = index(this%level, gridinfo%level)
868IF (ilevel == 0 .AND. optio_log(force)) THEN
869 ilevel = index(this%level, vol7d_level_miss)
870 IF (ilevel /= 0) this%level(ilevel) = gridinfo%level
871ENDIF
872
873IF (ilevel == 0) THEN
874 CALL l4f_category_log(this%category,l4f_error, &
875 "volgrid6d: level not valid for volume, gridinfo rejected")
876 CALL raise_error()
877 RETURN
878ENDIF
879
880IF (optio_log(isanavar)) THEN ! assign to all times and timeranges
881 itime0 = 1
882 itime1 = SIZE(this%time)
883 itimerange0 = 1
884 itimerange1 = SIZE(this%timerange)
885ELSE ! usual case
886 correctedtime = gridinfo%time
887 IF (this%time_definition == 1 .OR. this%time_definition == 2) correctedtime = correctedtime + &
888 timedelta_new(sec=gridinfo%timerange%p1)
889 itime0 = index(this%time, correctedtime)
890 IF (itime0 == 0 .AND. optio_log(force)) THEN
891 itime0 = index(this%time, datetime_miss)
892 IF (itime0 /= 0) this%time(itime0) = correctedtime
893 ENDIF
894 IF (itime0 == 0) THEN
895 CALL l4f_category_log(this%category,l4f_error, &
896 "volgrid6d: time not valid for volume, gridinfo rejected")
897 CALL raise_error()
898 RETURN
899 ENDIF
900 itime1 = itime0
901
902 correctedtimerange = gridinfo%timerange
903 IF (this%time_definition == 2) correctedtimerange%p1 = 0
904 itimerange0 = index(this%timerange, correctedtimerange)
905 IF (itimerange0 == 0 .AND. optio_log(force)) THEN
906 itimerange0 = index(this%timerange, vol7d_timerange_miss)
907 IF (itimerange0 /= 0) this%timerange(itimerange0) = gridinfo%timerange
908 ENDIF
909 IF (itimerange0 == 0) THEN
910 CALL l4f_category_log(this%category,l4f_error, &
911 "volgrid6d: timerange not valid for volume, gridinfo rejected")
912 CALL raise_error()
913 RETURN
914 ENDIF
915 itimerange1 = itimerange0
916ENDIF
917
918ivar = index(this%var, gridinfo%var)
919IF (ivar == 0 .AND. optio_log(force)) THEN
920 ivar = index(this%var, volgrid6d_var_miss)
921 IF (ivar /= 0) this%var(ivar) = gridinfo%var
922ENDIF
923IF (ivar == 0) THEN
924 CALL l4f_category_log(this%category,l4f_error, &
925 "volgrid6d: var not valid for volume, gridinfo rejected")
926 CALL raise_error()
927 RETURN
928ENDIF
929
930DO itimerange = itimerange0, itimerange1
931 DO itime = itime0, itime1
932 IF (ASSOCIATED(this%gaid)) THEN
933 dup = .false.
934 IF (c_e(this%gaid(ilevel,itime,itimerange,ivar))) THEN
935 dup = .true.
936 CALL l4f_category_log(this%category,l4f_warn,"gaid exist: grib duplicated")
937! avoid memory leaks
938 IF (optio_log(clone)) CALL delete(this%gaid(ilevel,itime,itimerange,ivar))
939 ENDIF
940
941 IF (optio_log(clone)) THEN
942 CALL copy(gridinfo%gaid, this%gaid(ilevel,itime,itimerange,ivar))
943#ifdef DEBUG
944 CALL l4f_category_log(this%category,l4f_debug,"cloning to a new gaid")
945#endif
946 ELSE
947 this%gaid(ilevel,itime,itimerange,ivar) = gridinfo%gaid
948 ENDIF
949
950 IF (ASSOCIATED(this%voldati))THEN
951 IF (.NOT.dup .OR. ldup_mode == 0) THEN
952 this%voldati(:,:,ilevel,itime,itimerange,ivar) = decode_gridinfo(gridinfo)
953 ELSE IF (ldup_mode == 1) THEN
954 tmpgrid = decode_gridinfo(gridinfo) ! f2003 automatic allocation
955 WHERE(c_e(tmpgrid))
956 this%voldati(:,:,ilevel,itime,itimerange,ivar) = tmpgrid(:,:)
957 END WHERE
958 ELSE IF (ldup_mode == 2) THEN
959 WHERE(.NOT.c_e(this%voldati(:,:,ilevel,itime,itimerange,ivar)))
960 this%voldati(:,:,ilevel,itime,itimerange,ivar) = decode_gridinfo(gridinfo)
961 END WHERE
962 ENDIF
963 ENDIF
964
965 ELSE
966 CALL l4f_category_log(this%category,l4f_error, &
967 "gaid not allocated, you probably need to call volgrid6d_alloc_vol first")
968 CALL raise_error()
969 RETURN
970 ENDIF
971 ENDDO
972ENDDO
973
974
975END SUBROUTINE import_from_gridinfo
976
977
978!> Export a single grid of a \a volgrid6d object to a \a gridinfo_def object.
979!! A single 2d slice of a \a volgrid6d at a specified location is
980!! written into a \a gridinfo_def object, including the \a grid_id
981!! which can be used for the successive export to file.
982SUBROUTINE export_to_gridinfo(this, gridinfo, itime, itimerange, ilevel, ivar, &
983 gaid_template, clone)
984TYPE(volgrid6d),INTENT(in) :: this !< volume to be exported
985TYPE(gridinfo_def),INTENT(inout) :: gridinfo !< output gridinfo_def object
986INTEGER :: itime !< index within \a this of the element to be exported for the time dimension
987INTEGER :: itimerange !< index within \a this of the element to be exported for the timerange dimension
988INTEGER :: ilevel !< index within \a this of the element to be exported for the vertical level dimension
989INTEGER :: ivar !< index within \a this of the element to be exported for the variable dimension
990TYPE(grid_id),INTENT(in),OPTIONAL :: gaid_template !< \a grid_id template to be used for output data replacing the one contained in \a this
991LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE., clone the grid_id included in \a this rather than making a shallow copy
992
993TYPE(grid_id) :: gaid
994LOGICAL :: usetemplate
995REAL,POINTER :: voldati(:,:)
996TYPE(datetime) :: correctedtime
997
998#ifdef DEBUG
999CALL l4f_category_log(this%category,l4f_debug,"export_to_gridinfo")
1000#endif
1001
1002IF (.NOT.c_e(this%gaid(ilevel,itime,itimerange,ivar))) THEN
1003#ifdef DEBUG
1004 CALL l4f_category_log(this%category,l4f_debug,"empty gaid found, skipping export")
1005#endif
1006 RETURN
1007ENDIF
1008
1009usetemplate = .false.
1010IF (PRESENT(gaid_template)) THEN
1011 CALL copy(gaid_template, gaid)
1012#ifdef DEBUG
1013 CALL l4f_category_log(this%category,l4f_debug,"template cloned to a new gaid")
1014#endif
1015 usetemplate = c_e(gaid)
1016ENDIF
1017
1018IF (.NOT.usetemplate) THEN
1019 IF (optio_log(clone)) THEN
1020 CALL copy(this%gaid(ilevel,itime,itimerange,ivar), gaid)
1021#ifdef DEBUG
1022 CALL l4f_category_log(this%category,l4f_debug,"original gaid cloned to a new one")
1023#endif
1024 ELSE
1025 gaid = this%gaid(ilevel,itime,itimerange,ivar)
1026 ENDIF
1027ENDIF
1028
1029IF (this%time_definition == 1 .OR. this%time_definition == 2) THEN
1030 correctedtime = this%time(itime) - &
1031 timedelta_new(sec=this%timerange(itimerange)%p1)
1032ELSE
1033 correctedtime = this%time(itime)
1034ENDIF
1035
1036CALL init(gridinfo,gaid, this%griddim, correctedtime, this%timerange(itimerange), &
1037 this%level(ilevel), this%var(ivar))
1038
1039! reset the gridinfo, bad but necessary at this point for encoding the field
1040CALL export(gridinfo%griddim, gridinfo%gaid)
1041! encode the field
1042IF (ASSOCIATED(this%voldati)) THEN
1043 CALL encode_gridinfo(gridinfo, this%voldati(:,:,ilevel,itime,itimerange,ivar))
1044ELSE IF (usetemplate) THEN ! field must be forced into template in this case
1045 NULLIFY(voldati)
1046 CALL volgrid_get_vol_2d(this, ilevel, itime, itimerange, ivar, voldati)
1047 CALL encode_gridinfo(gridinfo, voldati)
1048 DEALLOCATE(voldati)
1049ENDIF
1050
1051END SUBROUTINE export_to_gridinfo
1052
1053
1054!> Import an array of \a gridinfo objects into an array of \a
1055!! volgrid6d objects.
1056!! This method imports an array of gridded fields from an \a
1057!! arrayof_gridinfo object into a suitable number of \a volgrid6d
1058!! objects. The number of \a volgrid6d allocated is determined by the
1059!! number of different grids encountered in \a arrayof_gridinfo.
1060!! Unlike the import for a single \a gridinfo, here the \a volgrid6d
1061!! object is a non-associated pointer to a 1-d array of uninitialized
1062!! objects, and all the dimension descriptors in each of the objects
1063!! are allocated and assigned within the method according to the data
1064!! contained in \a gridinfov.
1065!! If the \a anavar array argument is provided, all the input messages
1066!! whose variable maps to one of the B-table variables contained in \a
1067!! anavar are treated as time-independent (AKA anagraphic data,
1068!! station data, etc.), thus their time and timerange are ignored and
1069!! they are replicated for every time and timerange present in the
1070!! corresponding data volume.
1071SUBROUTINE import_from_gridinfovv(this, gridinfov, dup_mode, clone, decode, &
1072 time_definition, anavar, categoryappend)
1073TYPE(volgrid6d),POINTER :: this(:) !< object in which to import
1074TYPE(arrayof_gridinfo),INTENT(in) :: gridinfov !< array of gridinfo objects to be imported
1075INTEGER,INTENT(in),OPTIONAL :: dup_mode !< determines the behavior in case of duplicate metadata: if \a dup_mode is not provided or 0, a duplicate field overwrites, if \a dup_mode is 1, duplicate fields are merged with priority to the last
1076LOGICAL , INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE. , clone the gaid's from \a gridinfo to \a this
1077LOGICAL,INTENT(in),OPTIONAL :: decode !< if provided and \a .FALSE. the data volume in the elements of \a this is not allocated and successive work will be performed on grid_id's
1078INTEGER,INTENT(IN),OPTIONAL :: time_definition !< 0=time is reference time; 1=time is validity time
1079CHARACTER(len=*),INTENT(IN),OPTIONAL :: anavar(:) !< list of variables (B-table code equivalent) to be treated as time-independent data
1080CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
1081
1082INTEGER :: i, j, stallo
1083INTEGER :: ngrid, ntime, ntimerange, nlevel, nvar, ltime_definition
1084INTEGER :: category
1085CHARACTER(len=512) :: a_name
1086TYPE(datetime),ALLOCATABLE :: correctedtime(:)
1087LOGICAL,ALLOCATABLE :: isanavar(:)
1088TYPE(vol7d_var) :: lvar
1089TYPE(vol7d_timerange),ALLOCATABLE :: correctedtimerange(:)
1090
1091! category temporanea (altrimenti non possiamo loggare)
1092if (present(categoryappend))then
1093 call l4f_launcher(a_name,a_name_append=trim(subcategory)//"."//trim(categoryappend))
1094else
1095 call l4f_launcher(a_name,a_name_append=trim(subcategory))
1096endif
1097category=l4f_category_get(a_name)
1098
1099#ifdef DEBUG
1100call l4f_category_log(category,l4f_debug,"start import_from_gridinfovv")
1101#endif
1102
1103IF (PRESENT(time_definition)) THEN
1104 ltime_definition = max(min(time_definition, 2), 0)
1105ELSE
1106 ltime_definition = 0
1107ENDIF
1108
1109ngrid=count_distinct(gridinfov%array(1:gridinfov%arraysize)%griddim,back=.true.)
1110CALL l4f_category_log(category,l4f_info, t2c(ngrid)// &
1111 ' different grid definition(s) found in input data')
1112
1113ALLOCATE(this(ngrid),stat=stallo)
1114IF (stallo /= 0)THEN
1115 CALL l4f_category_log(category,l4f_fatal,"allocating memory")
1116 CALL raise_fatal_error()
1117ENDIF
1118DO i = 1, ngrid
1119 IF (PRESENT(categoryappend))THEN
1120 CALL init(this(i), time_definition=ltime_definition, categoryappend=trim(categoryappend)//"-vol"//t2c(i))
1121 ELSE
1122 CALL init(this(i), time_definition=ltime_definition, categoryappend="vol"//t2c(i))
1123 ENDIF
1124ENDDO
1125
1126this(:)%griddim=pack_distinct(gridinfov%array(1:gridinfov%arraysize)%griddim, &
1127 ngrid, back=.true.)
1128
1129! mark elements as ana variables (time-independent)
1130ALLOCATE(isanavar(gridinfov%arraysize))
1131isanavar(:) = .false.
1132IF (PRESENT(anavar)) THEN
1133 DO i = 1, gridinfov%arraysize
1134 DO j = 1, SIZE(anavar)
1135 lvar = convert(gridinfov%array(i)%var)
1136 IF (lvar%btable == anavar(j)) THEN
1137 isanavar(i) = .true.
1138 EXIT
1139 ENDIF
1140 ENDDO
1141 ENDDO
1142 CALL l4f_category_log(category,l4f_info,t2c(count(isanavar))//'/'// &
1143 t2c(gridinfov%arraysize)//' constant-data messages found in input data')
1144ENDIF
1145
1146IF (ltime_definition == 1 .OR. ltime_definition == 2) THEN ! verification time
1147 ALLOCATE(correctedtime(gridinfov%arraysize))
1148 correctedtime(:) = gridinfov%array(1:gridinfov%arraysize)%time
1149 DO i = 1, gridinfov%arraysize
1150 correctedtime(i) = correctedtime(i) + &
1151 timedelta_new(sec=gridinfov%array(i)%timerange%p1)
1152 ENDDO
1153ENDIF
1154IF (ltime_definition == 2) THEN ! set all to analysis
1155 ALLOCATE(correctedtimerange(gridinfov%arraysize))
1156 correctedtimerange(:) = gridinfov%array(1:gridinfov%arraysize)%timerange
1157 correctedtimerange(:)%p1 = 0
1158ENDIF
1159
1160DO i = 1, ngrid
1161 IF (PRESENT(anavar)) THEN
1162 j = count((this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1163 .AND. .NOT.isanavar(:))
1164 IF (j <= 0) THEN
1165 CALL l4f_category_log(category, l4f_fatal, 'grid n.'//t2c(i)// &
1166 ' has only constant data, this is not allowed')
1167 CALL l4f_category_log(category, l4f_fatal, 'please check anavar argument')
1168 CALL raise_fatal_error()
1169 ENDIF
1170 ENDIF
1171 IF (ltime_definition == 1 .OR. ltime_definition == 2) THEN ! verification time
1172 ntime = count_distinct(correctedtime, &
1173 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1174 .AND. .NOT.isanavar(:), back=.true.)
1175 ELSE
1176 ntime = count_distinct(gridinfov%array(1:gridinfov%arraysize)%time, &
1177 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1178 .AND. .NOT.isanavar(:), back=.true.)
1179 ENDIF
1180 IF (ltime_definition == 2) THEN ! set all to analysis
1181 ntimerange = count_distinct(correctedtimerange, &
1182 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1183 .AND. .NOT.isanavar(:), back=.true.)
1184 ELSE
1185 ntimerange = count_distinct(gridinfov%array(1:gridinfov%arraysize)%timerange, &
1186 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1187 .AND. .NOT.isanavar(:), back=.true.)
1188 ENDIF
1189 nlevel = count_distinct(gridinfov%array(1:gridinfov%arraysize)%level, &
1190 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1191 back=.true.)
1192 nvar = count_distinct(gridinfov%array(1:gridinfov%arraysize)%var, &
1193 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1194 back=.true.)
1195
1196#ifdef DEBUG
1197 CALL l4f_category_log(this(i)%category,l4f_debug,"alloc volgrid6d index: "//t2c(i))
1198#endif
1199
1200 CALL volgrid6d_alloc(this(i),this(i)%griddim%dim,ntime=ntime, &
1201 ntimerange=ntimerange,nlevel=nlevel,nvar=nvar)
1202
1203 IF (ltime_definition == 1 .OR. ltime_definition == 2) THEN ! verification time
1204 this(i)%time = pack_distinct(correctedtime, ntime, &
1205 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1206 .AND. .NOT.isanavar(:), back=.true.)
1207 ELSE
1208 this(i)%time = pack_distinct(gridinfov%array(1:gridinfov%arraysize)%time, ntime, &
1209 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1210 .AND. .NOT.isanavar(:), back=.true.)
1211 ENDIF
1212 CALL sort(this(i)%time)
1213
1214 IF (ltime_definition == 2) THEN ! set all to analysis
1215 this(i)%timerange = pack_distinct(correctedtimerange, ntimerange, &
1216 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1217 .AND. .NOT.isanavar(:), back=.true.)
1218 ELSE
1219 this(i)%timerange = pack_distinct(gridinfov%array(1:gridinfov%arraysize)%timerange, &
1220 ntimerange, mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim) &
1221 .AND. .NOT.isanavar(:), back=.true.)
1222 ENDIF
1223 CALL sort(this(i)%timerange)
1224
1225 this(i)%level=pack_distinct(gridinfov%array(1:gridinfov%arraysize)%level, &
1226 nlevel,mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1227 back=.true.)
1228 CALL sort(this(i)%level)
1229
1230 this(i)%var=pack_distinct(gridinfov%array(1:gridinfov%arraysize)%var, nvar, &
1231 mask=(this(i)%griddim == gridinfov%array(1:gridinfov%arraysize)%griddim), &
1232 back=.true.)
1233
1234#ifdef DEBUG
1235 CALL l4f_category_log(this(i)%category,l4f_debug,"alloc_vol volgrid6d index: "//t2c(i))
1236#endif
1237 CALL volgrid6d_alloc_vol(this(i), decode=decode)
1238
1239ENDDO
1240
1241IF (ltime_definition == 1 .OR. ltime_definition == 2) DEALLOCATE(correctedtime)
1242IF (ltime_definition == 2) DEALLOCATE(correctedtimerange)
1243
1244DO i = 1, gridinfov%arraysize
1245
1246#ifdef DEBUG
1247 CALL l4f_category_log(category,l4f_debug,"import from gridinfov index: "//t2c(i))
1248 CALL l4f_category_log(category,l4f_info, &
1249 "to volgrid6d index: "//t2c(index(this%griddim, gridinfov%array(i)%griddim)))
1250#endif
1251
1252 CALL import(this(index(this%griddim, gridinfov%array(i)%griddim)), &
1253 gridinfov%array(i), dup_mode=dup_mode, clone=clone, isanavar=isanavar(i))
1254
1255ENDDO
1256
1257!chiudo il logger temporaneo
1258CALL l4f_category_delete(category)
1260END SUBROUTINE import_from_gridinfovv
1261
1262
1263!> Export a \a volgrid6d object to an \a arrayof_gridinfo object.
1264!! The multidimensional \a volgrid6d structure is serialized into a
1265!! one-dimensional array of gridinfo_def objects, which is allocated
1266!! to the proper size if not already allocated, or it is extended
1267!! keeping the old data if any.
1268SUBROUTINE export_to_gridinfov(this, gridinfov, gaid_template, clone)
1269TYPE(volgrid6d),INTENT(inout) :: this !< volume to be exported
1270TYPE(arrayof_gridinfo),INTENT(inout) :: gridinfov !< output array of gridinfo_def objects
1271TYPE(grid_id),INTENT(in),OPTIONAL :: gaid_template !< \a grid_id template to be used for output data replacing the one contained in \a this
1272LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE., clone the grid_id included in \a this rather than making a shallow copy
1273
1274INTEGER :: i ,itime, itimerange, ilevel, ivar
1275INTEGER :: ntime, ntimerange, nlevel, nvar
1276TYPE(gridinfo_def) :: gridinfol
1277
1278#ifdef DEBUG
1279CALL l4f_category_log(this%category,l4f_debug,"start export_to_gridinfov")
1280#endif
1281
1282! this is necessary in order not to repeatedly and uselessly copy the
1283! same array of coordinates for each 2d grid during export, the
1284! side-effect is that the computed projection in this is lost
1285CALL dealloc(this%griddim%dim)
1286
1287i=0
1288ntime=size(this%time)
1289ntimerange=size(this%timerange)
1290nlevel=size(this%level)
1291nvar=size(this%var)
1292
1293DO itime=1,ntime
1294 DO itimerange=1,ntimerange
1295 DO ilevel=1,nlevel
1296 DO ivar=1,nvar
1297
1298 CALL init(gridinfol)
1299 CALL export(this, gridinfol, itime, itimerange, ilevel, ivar, &
1300 gaid_template=gaid_template, clone=clone)
1301 IF (c_e(gridinfol%gaid)) THEN
1302 CALL insert(gridinfov, gridinfol)
1303 ELSE
1304 CALL delete(gridinfol)
1305 ENDIF
1306
1307 ENDDO
1308 ENDDO
1309 ENDDO
1310ENDDO
1311
1312END SUBROUTINE export_to_gridinfov
1313
1314
1315!> Export an array of \a volgrid6d objects to an \a arrayof_gridinfo object.
1316!! The multidimensional \a volgrid6d structures are serialized into a
1317!! one-dimensional array of gridinfo_def objects, which is allocated
1318!! to the proper size if not already allocated, or it is extended
1319!! keeping the old data if any.
1320SUBROUTINE export_to_gridinfovv(this, gridinfov, gaid_template, clone)
1321!, &
1322! categoryappend)
1323TYPE(volgrid6d),INTENT(inout) :: this(:) !< volume array to be exported
1324TYPE(arrayof_gridinfo),INTENT(inout) :: gridinfov !< output array of gridinfo_def objects
1325TYPE(grid_id),INTENT(in),OPTIONAL :: gaid_template !< \a grid_id template to be used for output data replacing the one contained in \a this
1326LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE., clone the grid_id included in \a this rather than making a shallow copy
1327
1328INTEGER :: i
1329
1330DO i = 1, SIZE(this)
1331#ifdef DEBUG
1332 CALL l4f_category_log(this(i)%category,l4f_debug, &
1333 "export_to_gridinfovv grid index: "//t2c(i))
1334#endif
1335 CALL export(this(i), gridinfov, gaid_template=gaid_template, clone=clone)
1336ENDDO
1337
1338END SUBROUTINE export_to_gridinfovv
1339
1340
1341!> Import the content of a supported file (like grib or gdal-supported
1342!! format) into an array of \a volgrid6d objects. This method imports
1343!! a set of gridded fields from a file object into a suitable number
1344!! of \a volgrid6d objects. The data are imported by creating a
1345!! temporary \a gridinfo object, importing it into the \a volgrid6d
1346!! object cloning the gaid's and then destroying the gridinfo, so it
1347!! works similarly to volgrid6d_class::import_from_gridinfovv() method.
1348!! For a detailed explanation of the \a anavar argument, see the
1349!! documentation of volgrid6d_class::import_from_gridinfovv() method.
1350SUBROUTINE volgrid6d_import_from_file(this, filename, dup_mode, decode, &
1351 time_definition, anavar, categoryappend)
1352TYPE(volgrid6d),POINTER :: this(:) !< object in which to import
1353CHARACTER(len=*),INTENT(in) :: filename !< name of file from which to import
1354INTEGER,INTENT(in),OPTIONAL :: dup_mode !< determines the behavior in case of duplicate metadata: if \a dup_mode is not provided or 0, a duplicate field overwrites, if \a dup_mode is 1, duplicate fields are merged with priority to the last
1355LOGICAL,INTENT(in),OPTIONAL :: decode !< if provided and \a .FALSE. the data volume in the elements of \a this is not allocated and successive work will be performed on grid_id's
1356INTEGER,INTENT(IN),OPTIONAL :: time_definition !< 0=time is reference time; 1=time is validity time
1357CHARACTER(len=*),INTENT(IN),OPTIONAL :: anavar(:) !< list of variables (B-table code equivalent) to be treated as time-independent data
1358character(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
1359
1360TYPE(arrayof_gridinfo) :: gridinfo
1361INTEGER :: category
1362CHARACTER(len=512) :: a_name
1363
1364NULLIFY(this)
1365
1366IF (PRESENT(categoryappend))THEN
1367 CALL l4f_launcher(a_name,a_name_append= &
1368 trim(subcategory)//"."//trim(categoryappend))
1369ELSE
1370 CALL l4f_launcher(a_name,a_name_append=trim(subcategory))
1371ENDIF
1372category=l4f_category_get(a_name)
1373
1374CALL import(gridinfo, filename=filename, categoryappend=categoryappend)
1375
1376IF (gridinfo%arraysize > 0) THEN
1377
1378 CALL import(this, gridinfo, dup_mode=dup_mode, clone=.true., decode=decode, &
1379 time_definition=time_definition, anavar=anavar, &
1380 categoryappend=categoryappend)
1381
1382 CALL l4f_category_log(category,l4f_info,"deleting gridinfo")
1383 CALL delete(gridinfo)
1384
1385ELSE
1386 CALL l4f_category_log(category,l4f_info,"file does not contain gridded data")
1387ENDIF
1388
1389! close logger
1390CALL l4f_category_delete(category)
1391
1392END SUBROUTINE volgrid6d_import_from_file
1393
1394
1395!> High level method for exporting a volume array to file.
1396!! All the information contained into an array of \a volgrid6d
1397!! objects, i.e. dimension descriptors and data, is exported to a file
1398!! using the proper output driver (typically grib_api for grib
1399!! format). If a template is provided, it will determine the
1400!! characteristic of the output file, otherwise the \a grid_id
1401!! descriptors contained in the volgrid6d object will be used
1402SUBROUTINE volgrid6d_export_to_file(this, filename, gaid_template, categoryappend)
1403TYPE(volgrid6d) :: this(:) !< volume(s) to be exported
1404CHARACTER(len=*),INTENT(in) :: filename !< output file name
1405TYPE(grid_id),INTENT(in),OPTIONAL :: gaid_template !< template for the output file, if provided the grid_id information stored in the volgrid6d objects will be ignored
1406character(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
1407
1408TYPE(arrayof_gridinfo) :: gridinfo
1409INTEGER :: category
1410CHARACTER(len=512) :: a_name
1411
1412IF (PRESENT(categoryappend)) THEN
1413 CALL l4f_launcher(a_name,a_name_append=trim(subcategory)//"."//trim(categoryappend))
1414ELSE
1415 CALL l4f_launcher(a_name,a_name_append=trim(subcategory))
1416ENDIF
1417category=l4f_category_get(a_name)
1418
1419#ifdef DEBUG
1420CALL l4f_category_log(category,l4f_debug,"start export to file")
1421#endif
1422
1423CALL l4f_category_log(category,l4f_info,"writing volgrid6d to grib file: "//trim(filename))
1424
1425!IF (ASSOCIATED(this)) THEN
1426 CALL export(this, gridinfo, gaid_template=gaid_template, clone=.true.)
1427 IF (gridinfo%arraysize > 0) THEN
1428 CALL export(gridinfo, filename)
1429 CALL delete(gridinfo)
1430 ENDIF
1431!ELSE
1432! CALL l4f_category_log(category,L4F_INFO,"volume volgrid6d is not associated")
1433!ENDIF
1434
1435! close logger
1436CALL l4f_category_delete(category)
1437
1438END SUBROUTINE volgrid6d_export_to_file
1439
1440
1441!> Array destructor for \a volgrid6d class.
1442!! Delete an array of \a volgrid6d objects and deallocate the array
1443!! itself.
1444SUBROUTINE volgrid6dv_delete(this)
1445TYPE(volgrid6d),POINTER :: this(:) !< vector of volgrid6d object
1446
1447INTEGER :: i
1448
1449IF (ASSOCIATED(this)) THEN
1450 DO i = 1, SIZE(this)
1451#ifdef DEBUG
1452 CALL l4f_category_log(this(i)%category,l4f_debug, &
1453 "delete volgrid6d vector index: "//trim(to_char(i)))
1454#endif
1455 CALL delete(this(i))
1456 ENDDO
1457 DEALLOCATE(this)
1458ENDIF
1459
1460END SUBROUTINE volgrid6dv_delete
1461
1462
1463! Internal method for performing grid to grid computations
1464SUBROUTINE volgrid6d_transform_compute(this, volgrid6d_in, volgrid6d_out, &
1465 lev_out, var_coord_vol, clone)
1466TYPE(grid_transform),INTENT(in) :: this ! oggetto di trasformazione per il grigliato
1467type(volgrid6d), INTENT(in) :: volgrid6d_in ! oggetto da trasformare
1468type(volgrid6d), INTENT(inout) :: volgrid6d_out ! oggetto trasformato; deve essere completo (init, alloc, alloc_vol)
1469TYPE(vol7d_level),INTENT(in),OPTIONAL :: lev_out(:) ! vol7d_level object defining target vertical grid, for vertical interpolations
1470INTEGER,INTENT(in),OPTIONAL :: var_coord_vol ! index of variable defining vertical coordinate values in input volume
1471LOGICAL,INTENT(in),OPTIONAL :: clone ! se fornito e \c .TRUE., clona i gaid da volgrid6d_in a volgrid6d_out
1472
1473INTEGER :: ntime, ntimerange, inlevel, onlevel, nvar, &
1474 itime, itimerange, ilevel, ivar, levshift, levused, lvar_coord_vol, spos
1475REAL,POINTER :: voldatiin(:,:,:), voldatiout(:,:,:), coord_3d_in(:,:,:)
1476TYPE(vol7d_level) :: output_levtype
1477
1478
1479#ifdef DEBUG
1480call l4f_category_log(volgrid6d_in%category,l4f_debug,"start volgrid6d_transform_compute")
1481#endif
1482
1483ntime=0
1484ntimerange=0
1485inlevel=0
1486onlevel=0
1487nvar=0
1488lvar_coord_vol = optio_i(var_coord_vol)
1489
1490if (associated(volgrid6d_in%time))then
1491 ntime=size(volgrid6d_in%time)
1492 volgrid6d_out%time=volgrid6d_in%time
1493end if
1494
1495if (associated(volgrid6d_in%timerange))then
1496 ntimerange=size(volgrid6d_in%timerange)
1497 volgrid6d_out%timerange=volgrid6d_in%timerange
1498end if
1499
1500IF (ASSOCIATED(volgrid6d_in%level))THEN
1501 inlevel=SIZE(volgrid6d_in%level)
1502ENDIF
1503IF (PRESENT(lev_out)) THEN
1504 onlevel=SIZE(lev_out)
1505 volgrid6d_out%level=lev_out
1506ELSE IF (ASSOCIATED(volgrid6d_in%level))THEN
1507 onlevel=SIZE(volgrid6d_in%level)
1508 volgrid6d_out%level=volgrid6d_in%level
1509ENDIF
1510
1511if (associated(volgrid6d_in%var))then
1512 nvar=size(volgrid6d_in%var)
1513 volgrid6d_out%var=volgrid6d_in%var
1514end if
1515! allocate once for speed
1516IF (.NOT.ASSOCIATED(volgrid6d_in%voldati)) THEN
1517 ALLOCATE(voldatiin(volgrid6d_in%griddim%dim%nx, volgrid6d_in%griddim%dim%ny, &
1518 inlevel))
1519ENDIF
1520IF (.NOT.ASSOCIATED(volgrid6d_out%voldati)) THEN
1521 ALLOCATE(voldatiout(volgrid6d_out%griddim%dim%nx, volgrid6d_out%griddim%dim%ny, &
1522 onlevel))
1523ENDIF
1524
1525CALL get_val(this, levshift=levshift, levused=levused)
1526spos = imiss
1527IF (c_e(lvar_coord_vol)) THEN
1528 CALL get_val(this%trans, output_levtype=output_levtype)
1529 IF (output_levtype%level1 == 103 .OR. output_levtype%level1 == 108) THEN
1530 spos = firsttrue(volgrid6d_in%level(:) == vol7d_level_new(1))
1531 IF (spos == 0) THEN
1532 CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1533 'output level '//t2c(output_levtype%level1)// &
1534 ' requested, but height/press of surface not provided in volume')
1535 ENDIF
1536 IF (.NOT.c_e(levshift) .AND. .NOT.c_e(levused)) THEN
1537 CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1538 'internal inconsistence, levshift and levused undefined when they should be')
1539 ENDIF
1540 ENDIF
1541ENDIF
1542
1543DO ivar=1,nvar
1544! IF (c_e(var_coord_vol)) THEN
1545! IF (ivar == var_coord_vol) CYCLE ! skip coordinate variable in output
1546! ENDIF
1547 DO itimerange=1,ntimerange
1548 DO itime=1,ntime
1549! skip empty columns where possible, improve
1550 IF (c_e(levshift) .AND. c_e(levused)) THEN
1551 IF (.NOT.any(c_e( &
1552 volgrid6d_in%gaid(levshift+1:levshift+levused,itime,itimerange,ivar) &
1553 ))) cycle
1554 ENDIF
1555 DO ilevel = 1, min(inlevel,onlevel)
1556! if present gaid copy it
1557 IF (c_e(volgrid6d_in%gaid(ilevel,itime,itimerange,ivar)) .AND. .NOT. &
1558 c_e(volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))) THEN
1559
1560 IF (optio_log(clone)) THEN
1561 CALL copy(volgrid6d_in%gaid(ilevel,itime,itimerange,ivar),&
1562 volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))
1563#ifdef DEBUG
1564 CALL l4f_category_log(volgrid6d_in%category,l4f_debug, &
1565 "cloning gaid, level "//t2c(ilevel))
1566#endif
1567 ELSE
1568 volgrid6d_out%gaid(ilevel,itime,itimerange,ivar) = &
1569 volgrid6d_in%gaid(ilevel,itime,itimerange,ivar)
1570 ENDIF
1571 ENDIF
1572 ENDDO
1573! if out level are more, we have to clone anyway
1574 DO ilevel = min(inlevel,onlevel) + 1, onlevel
1575 IF (c_e(volgrid6d_in%gaid(inlevel,itime,itimerange,ivar)) .AND. .NOT. &
1576 c_e(volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))) then
1577
1578 CALL copy(volgrid6d_in%gaid(inlevel,itime,itimerange,ivar),&
1579 volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))
1580#ifdef DEBUG
1581 CALL l4f_category_log(volgrid6d_in%category,l4f_debug, &
1582 "forced cloning gaid, level "//t2c(inlevel)//"->"//t2c(ilevel))
1583#endif
1584 ENDIF
1585 ENDDO
1586
1587 IF (c_e(lvar_coord_vol)) THEN
1588 NULLIFY(coord_3d_in)
1589 CALL volgrid_get_vol_3d(volgrid6d_in, itime, itimerange, lvar_coord_vol, &
1590 coord_3d_in)
1591 IF (c_e(spos)) THEN ! compute difference wrt surface coordinate
1592 IF (spos == 0) THEN ! error condition, set all to missing and goodnight
1593 coord_3d_in(:,:,levshift+1:levshift+levused) = rmiss
1594 ELSE
1595 DO ilevel = levshift+1, levshift+levused
1596 WHERE(c_e(coord_3d_in(:,:,ilevel)) .AND. c_e(coord_3d_in(:,:,spos)))
1597 coord_3d_in(:,:,ilevel) = coord_3d_in(:,:,ilevel) - &
1598 coord_3d_in(:,:,spos)
1599 ELSEWHERE
1600 coord_3d_in(:,:,ilevel) = rmiss
1601 END WHERE
1602 ENDDO
1603 ENDIF
1604 ENDIF
1605 ENDIF
1606 CALL volgrid_get_vol_3d(volgrid6d_in, itime, itimerange, ivar, &
1607 voldatiin)
1608 IF (ASSOCIATED(volgrid6d_out%voldati)) & ! improve!!!!
1609 CALL volgrid_get_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
1610 voldatiout)
1611 IF (c_e(lvar_coord_vol)) THEN
1612 CALL compute(this, voldatiin, voldatiout, convert(volgrid6d_in%var(ivar)), &
1613 coord_3d_in(:,:,levshift+1:levshift+levused)) ! subset coord_3d_in
1614 ELSE
1615 CALL compute(this, voldatiin, voldatiout, convert(volgrid6d_in%var(ivar)))
1616 ENDIF
1617 CALL volgrid_set_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
1618 voldatiout)
1619 ENDDO
1620 ENDDO
1621ENDDO
1622
1623IF (c_e(lvar_coord_vol)) THEN
1624 DEALLOCATE(coord_3d_in)
1625ENDIF
1626IF (.NOT.ASSOCIATED(volgrid6d_in%voldati)) THEN
1627 DEALLOCATE(voldatiin)
1628ENDIF
1629IF (.NOT.ASSOCIATED(volgrid6d_out%voldati)) THEN
1630 DEALLOCATE(voldatiout)
1631ENDIF
1633
1634END SUBROUTINE volgrid6d_transform_compute
1635
1636
1637!> Performs the specified abstract transformation on the data provided.
1638!! The abstract transformation is specified by \a this parameter; the
1639!! corresponding specifical transformation (\a grid_transform object)
1640!! is created and destroyed internally. The output transformed object
1641!! is created internally and it does not require preliminary
1642!! initialisation.
1643SUBROUTINE volgrid6d_transform(this, griddim, volgrid6d_in, volgrid6d_out, &
1644 lev_out, volgrid6d_coord_in, maskgrid, maskbounds, clone, decode, categoryappend)
1645TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
1646TYPE(griddim_def),INTENT(in),OPTIONAL :: griddim !< griddim specifying the output grid (required by most transformation types)
1647TYPE(volgrid6d),INTENT(inout) :: volgrid6d_in !< object to be transformed, it is not modified, despite the INTENT(inout)
1648TYPE(volgrid6d),INTENT(out) :: volgrid6d_out !< transformed object, it does not require initialisation
1649TYPE(vol7d_level),INTENT(in),OPTIONAL,TARGET :: lev_out(:) !< vol7d_level object defining target vertical grid, for vertical interpolations
1650TYPE(volgrid6d),INTENT(in),OPTIONAL :: volgrid6d_coord_in !< object providing time constant input vertical coordinate for some kind of vertical interpolations
1651REAL,INTENT(in),OPTIONAL :: maskgrid(:,:) !< 2D field to be used for defining subareas according to its values, it must have the same shape as the field to be interpolated (for transformation subtype 'maskfill')
1652REAL,INTENT(in),OPTIONAL :: maskbounds(:) !< array of boundary values for defining a subset of valid points where the values of \a maskgrid are within the first and last value of \a maskbounds (for transformation type 'metamorphosis:maskfill')
1653LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \a .TRUE. , clone the \a gaid's from \a volgrid6d_in to \a volgrid6d_out
1654LOGICAL,INTENT(in),OPTIONAL :: decode !< determine whether the data in \a volgrid6d_out should be decoded or remain coded in gaid, if not provided, the decode status is taken from \a volgrid6d_in
1655CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
1656
1657TYPE(grid_transform) :: grid_trans
1658TYPE(vol7d_level),POINTER :: llev_out(:)
1659TYPE(vol7d_level) :: input_levtype, output_levtype
1660TYPE(vol7d_var) :: vcoord_var
1661INTEGER :: i, k, ntime, ntimerange, nlevel, nvar, var_coord_in, var_coord_vol, &
1662 cf_out, nxc, nyc, nxi, nyi, i3, i4, i5, i6, &
1663 ulstart, ulend, spos
1664REAL,ALLOCATABLE :: coord_3d_in(:,:,:)
1665TYPE(geo_proj) :: proj_in, proj_out
1666CHARACTER(len=80) :: trans_type
1667LOGICAL :: ldecode
1668LOGICAL,ALLOCATABLE :: mask_in(:)
1669
1670#ifdef DEBUG
1671call l4f_category_log(volgrid6d_in%category, l4f_debug, "start volgrid6d_transform")
1672#endif
1673
1674ntime=0
1675ntimerange=0
1676nlevel=0
1677nvar=0
1678
1679if (associated(volgrid6d_in%time)) ntime=size(volgrid6d_in%time)
1680if (associated(volgrid6d_in%timerange)) ntimerange=size(volgrid6d_in%timerange)
1681if (associated(volgrid6d_in%level)) nlevel=size(volgrid6d_in%level)
1682if (associated(volgrid6d_in%var)) nvar=size(volgrid6d_in%var)
1683
1684IF (ntime == 0 .OR. ntimerange == 0 .OR. nlevel == 0 .OR. nvar == 0) THEN
1685 CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1686 "trying to transform an incomplete volgrid6d object, ntime="//t2c(ntime)// &
1687 ' ntimerange='//t2c(ntimerange)//' nlevel='//t2c(nlevel)//' nvar='//t2c(nvar))
1688 CALL init(volgrid6d_out) ! initialize to empty
1689 CALL raise_error()
1690 RETURN
1691ENDIF
1692
1693CALL get_val(this, trans_type=trans_type)
1694
1695! store desired output component flag and unrotate if necessary
1696cf_out = imiss
1697IF (PRESENT(griddim) .AND. (trans_type == 'inter' .OR. trans_type == 'boxinter' &
1698 .OR. trans_type == 'stencilinter')) THEN ! improve condition!!
1699 CALL get_val(volgrid6d_in%griddim, proj=proj_in)
1700 CALL get_val(griddim, component_flag=cf_out, proj=proj_out)
1701! if different projections wind components must be referred to geographical system
1702 IF (proj_in /= proj_out) CALL vg6d_wind_unrot(volgrid6d_in)
1703ELSE IF (PRESENT(griddim)) THEN ! just get component_flag, the rest is rubbish
1704 CALL get_val(griddim, component_flag=cf_out)
1705ENDIF
1706
1707
1708var_coord_in = imiss
1709var_coord_vol = imiss
1710IF (trans_type == 'vertint') THEN
1711 IF (PRESENT(lev_out)) THEN
1712
1713! if volgrid6d_coord_in provided and allocated, check that it fits
1714 IF (PRESENT(volgrid6d_coord_in)) THEN
1715 IF (ASSOCIATED(volgrid6d_coord_in%voldati)) THEN
1716
1717! strictly 1 time and 1 timerange
1718 IF (SIZE(volgrid6d_coord_in%voldati,4) /= 1 .OR. &
1719 SIZE(volgrid6d_coord_in%voldati,5) /= 1) THEN
1720 CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1721 'volume providing constant input vertical coordinate must have &
1722 &only 1 time and 1 timerange')
1723 CALL init(volgrid6d_out) ! initialize to empty
1724 CALL raise_error()
1725 RETURN
1726 ENDIF
1727
1728! search for variable providing vertical coordinate
1729 CALL get_val(this, output_levtype=output_levtype)
1730 vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
1731 IF (.NOT.c_e(vcoord_var)) THEN
1732 CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1733 'requested output level type '//t2c(output_levtype%level1)// &
1734 ' does not correspond to any known physical variable for &
1735 &providing vertical coordinate')
1736 CALL init(volgrid6d_out) ! initialize to empty
1737 CALL raise_error()
1738 RETURN
1739 ENDIF
1740
1741 DO i = 1, SIZE(volgrid6d_coord_in%var)
1742 IF (convert(volgrid6d_coord_in%var(i)) == vcoord_var) THEN
1743 var_coord_in = i
1744 EXIT
1745 ENDIF
1746 ENDDO
1747
1748 IF (.NOT.c_e(var_coord_in)) THEN
1749 CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1750 'volume providing constant input vertical coordinate contains no &
1751 &variables matching output level type '//t2c(output_levtype%level1))
1752 CALL init(volgrid6d_out) ! initialize to empty
1753 CALL raise_error()
1754 RETURN
1755 ENDIF
1756 CALL l4f_category_log(volgrid6d_in%category, l4f_info, &
1757 'Coordinate for vertint found in coord volume at position '// &
1758 t2c(var_coord_in))
1759
1760! check horizontal grid
1761! this is too strict (component flag and so on)
1762! IF (volgrid6d_coord_in%griddim /= volgrid6d_in%griddim) THEN
1763 CALL get_val(volgrid6d_coord_in%griddim, nx=nxc, ny=nyc)
1764 CALL get_val(volgrid6d_in%griddim, nx=nxi, ny=nyi)
1765 IF (nxc /= nxi .OR. nyc /= nyi) THEN
1766 CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1767 'volume providing constant input vertical coordinate must have &
1768 &the same grid as the input')
1769 CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1770 'vertical coordinate: '//t2c(nxc)//'x'//t2c(nyc)// &
1771 ', input volume: '//t2c(nxi)//'x'//t2c(nyi))
1772 CALL init(volgrid6d_out) ! initialize to empty
1773 CALL raise_error()
1774 RETURN
1775 ENDIF
1776
1777! check vertical coordinate system
1778 CALL get_val(this, input_levtype=input_levtype)
1779 mask_in = & ! implicit allocation
1780 (volgrid6d_coord_in%level(:)%level1 == input_levtype%level1) .AND. &
1781 (volgrid6d_coord_in%level(:)%level2 == input_levtype%level2)
1782 ulstart = firsttrue(mask_in)
1783 ulend = lasttrue(mask_in)
1784 IF (ulstart == 0 .OR. ulend == 0) THEN
1785 CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1786 'coordinate file does not contain levels of type '// &
1787 t2c(input_levtype%level1)//'/'//t2c(input_levtype%level2)// &
1788 ' specified for input data')
1789 CALL init(volgrid6d_out) ! initialize to empty
1790 CALL raise_error()
1791 RETURN
1792 ENDIF
1793
1794 coord_3d_in = volgrid6d_coord_in%voldati(:,:,ulstart:ulend,1,1,var_coord_in) ! implicit allocation
1795! special case
1796 IF (output_levtype%level1 == 103 .OR. &
1797 output_levtype%level1 == 108) THEN ! surface coordinate needed
1798 spos = firsttrue(volgrid6d_coord_in%level(:) == vol7d_level_new(1))
1799 IF (spos == 0) THEN
1800 CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1801 'output level '//t2c(output_levtype%level1)// &
1802 ' requested, but height/press of surface not provided in coordinate file')
1803 CALL init(volgrid6d_out) ! initialize to empty
1804 CALL raise_error()
1805 RETURN
1806 ENDIF
1807 DO k = 1, SIZE(coord_3d_in,3)
1808 WHERE(c_e(coord_3d_in(:,:,k)) .AND. &
1809 c_e(volgrid6d_coord_in%voldati(:,:,spos,1,1,var_coord_in)))
1810 coord_3d_in(:,:,k) = coord_3d_in(:,:,k) - &
1811 volgrid6d_coord_in%voldati(:,:,spos,1,1,var_coord_in)
1812 ELSEWHERE
1813 coord_3d_in(:,:,k) = rmiss
1814 END WHERE
1815 ENDDO
1816 ENDIF
1817
1818 ENDIF
1819 ENDIF
1820
1821 IF (.NOT.c_e(var_coord_in)) THEN ! search for coordinate within volume
1822! search for variable providing vertical coordinate
1823 CALL get_val(this, output_levtype=output_levtype)
1824 vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
1825 IF (c_e(vcoord_var)) THEN
1826 DO i = 1, SIZE(volgrid6d_in%var)
1827 IF (convert(volgrid6d_in%var(i)) == vcoord_var) THEN
1828 var_coord_vol = i
1829 EXIT
1830 ENDIF
1831 ENDDO
1832
1833 IF (c_e(var_coord_vol)) THEN
1834 CALL l4f_category_log(volgrid6d_in%category, l4f_info, &
1835 'Coordinate for vertint found in input volume at position '// &
1836 t2c(var_coord_vol))
1837 ENDIF
1838
1839 ENDIF
1840 ENDIF
1841
1842 CALL init(volgrid6d_out, griddim=volgrid6d_in%griddim, &
1843 time_definition=volgrid6d_in%time_definition, categoryappend=categoryappend)
1844 IF (c_e(var_coord_in)) THEN
1845 CALL init(grid_trans, this, lev_in=volgrid6d_in%level, lev_out=lev_out, &
1846 coord_3d_in=coord_3d_in, categoryappend=categoryappend)
1847 ELSE
1848 CALL init(grid_trans, this, lev_in=volgrid6d_in%level, lev_out=lev_out, &
1849 categoryappend=categoryappend)
1850 ENDIF
1851
1852 CALL get_val(grid_trans, output_level_auto=llev_out) ! get levels if auto-generated
1853 IF (.NOT.ASSOCIATED(llev_out)) llev_out => lev_out
1854 nlevel = SIZE(llev_out)
1855 ELSE
1856 CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1857 'volgrid6d_transform: vertint requested but lev_out not provided')
1858 CALL init(volgrid6d_out) ! initialize to empty
1859 CALL raise_error()
1860 RETURN
1861 ENDIF
1862
1863ELSE
1864 CALL init(volgrid6d_out, griddim=griddim, &
1865 time_definition=volgrid6d_in%time_definition, categoryappend=categoryappend)
1866 CALL init(grid_trans, this, in=volgrid6d_in%griddim, out=volgrid6d_out%griddim, &
1867 maskgrid=maskgrid, maskbounds=maskbounds, categoryappend=categoryappend)
1868ENDIF
1869
1870
1871IF (c_e(grid_trans)) THEN ! transformation is valid
1872
1873 CALL volgrid6d_alloc(volgrid6d_out, ntime=ntime, nlevel=nlevel, &
1874 ntimerange=ntimerange, nvar=nvar)
1875
1876 IF (PRESENT(decode)) THEN ! explicitly set decode status
1877 ldecode = decode
1878 ELSE ! take it from input
1879 ldecode = ASSOCIATED(volgrid6d_in%voldati)
1880 ENDIF
1881! force decode if gaid is readonly
1882 decode_loop: DO i6 = 1,nvar
1883 DO i5 = 1, ntimerange
1884 DO i4 = 1, ntime
1885 DO i3 = 1, nlevel
1886 IF (c_e(volgrid6d_in%gaid(i3,i4,i5,i6))) THEN
1887 ldecode = ldecode .OR. grid_id_readonly(volgrid6d_in%gaid(i3,i4,i5,i6))
1888 EXIT decode_loop
1889 ENDIF
1890 ENDDO
1891 ENDDO
1892 ENDDO
1893 ENDDO decode_loop
1894
1895 IF (PRESENT(decode)) THEN
1896 IF (ldecode.NEQV.decode) THEN
1897 CALL l4f_category_log(volgrid6d_in%category, l4f_warn, &
1898 'volgrid6d_transform: decode status forced to .TRUE. because driver does not allow copy')
1899 ENDIF
1900 ENDIF
1901
1902 CALL volgrid6d_alloc_vol(volgrid6d_out, decode=ldecode)
1903
1904!ensure unproj was called
1905!call griddim_unproj(volgrid6d_out%griddim)
1906
1907 IF (trans_type == 'vertint') THEN
1908#ifdef DEBUG
1909 CALL l4f_category_log(volgrid6d_in%category, l4f_debug, &
1910 "volgrid6d_transform: vertint to "//t2c(nlevel)//" levels")
1911#endif
1912 CALL compute(grid_trans, volgrid6d_in, volgrid6d_out, lev_out=llev_out, &
1913 var_coord_vol=var_coord_vol, clone=clone)
1914 ELSE
1915 CALL compute(grid_trans, volgrid6d_in, volgrid6d_out, clone=clone)
1916 ENDIF
1917
1918 IF (cf_out == 0) THEN ! unrotated components are desired
1919 CALL wind_unrot(volgrid6d_out) ! unrotate if necessary
1920 ELSE IF (cf_out == 1) THEN ! rotated components are desired
1921 CALL wind_rot(volgrid6d_out) ! rotate if necessary
1922 ENDIF
1923
1924ELSE
1925! should log with grid_trans%category, but it is private
1926 CALL l4f_category_log(volgrid6d_in%category, l4f_error, &
1927 'volgrid6d_transform: transformation not valid')
1928 CALL raise_error()
1929ENDIF
1930
1931CALL delete (grid_trans)
1932
1933END SUBROUTINE volgrid6d_transform
1934
1935
1936!> Performs the specified abstract transformation on the arrays of
1937!! data provided. The abstract transformation is specified by \a this
1938!! parameter; the corresponding specifical transformation (\a
1939!! grid_transform object) is created and destroyed internally. The
1940!! output transformed object is created internally and it does not
1941!! require preliminary initialisation. According to the input data and
1942!! to the transformation type, the output array may have of one or
1943!! more \a volgrid6d elements on different grids.
1944SUBROUTINE volgrid6dv_transform(this, griddim, volgrid6d_in, volgrid6d_out, &
1945 lev_out, volgrid6d_coord_in, maskgrid, maskbounds, clone, decode, categoryappend)
1946TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
1947TYPE(griddim_def),INTENT(in),OPTIONAL :: griddim !< griddim specifying the output grid (required by most transformation types)
1948TYPE(volgrid6d),INTENT(inout) :: volgrid6d_in(:) !< object to be transformed, it is an array of volgrid6d objects, each of which will be transformed, it is not modified, despite the INTENT(inout)
1949TYPE(volgrid6d),POINTER :: volgrid6d_out(:) !< transformed object, it is a non associated pointer to an array of volgrid6d objects which will be allocated by the method
1950TYPE(vol7d_level),INTENT(in),OPTIONAL :: lev_out(:) !< vol7d_level object defining target vertical grid
1951TYPE(volgrid6d),INTENT(in),OPTIONAL :: volgrid6d_coord_in !< object providing time constant input vertical coordinate for some kind of vertical interpolations
1952REAL,INTENT(in),OPTIONAL :: maskgrid(:,:) !< 2D field to be used for defining subareas according to its values, it must have the same shape as the field to be interpolated (for transformation subtype 'maskfill')
1953REAL,INTENT(in),OPTIONAL :: maskbounds(:) !< array of boundary values for defining a subset of valid points where the values of \a maskgrid are within the first and last value of \a maskbounds (for transformation type 'metamorphosis:maskfill')
1954LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \a .TRUE. , clone the \a gaid's from \a volgrid6d_in to \a volgrid6d_out
1955LOGICAL,INTENT(in),OPTIONAL :: decode !< if provided and \a .FALSE. the data volume is not allocated, but work is performed on grid_id's
1956CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
1957
1958INTEGER :: i, stallo
1959
1960
1961allocate(volgrid6d_out(size(volgrid6d_in)),stat=stallo)
1962if (stallo /= 0)then
1963 call l4f_log(l4f_fatal,"allocating memory")
1964 call raise_fatal_error()
1965end if
1966
1967do i=1,size(volgrid6d_in)
1968 call transform(this, griddim, volgrid6d_in(i), volgrid6d_out(i), &
1969 lev_out=lev_out, volgrid6d_coord_in=volgrid6d_coord_in, &
1970 maskgrid=maskgrid, maskbounds=maskbounds, &
1971 clone=clone, decode=decode, categoryappend=categoryappend)
1972end do
1973
1974END SUBROUTINE volgrid6dv_transform
1975
1976
1977! Internal method for performing grid to sparse point computations
1978SUBROUTINE volgrid6d_v7d_transform_compute(this, volgrid6d_in, vol7d_out, &
1979 networkname, noconvert)
1980TYPE(grid_transform),INTENT(in) :: this ! oggetto di trasformazione per grigliato
1981type(volgrid6d), INTENT(in) :: volgrid6d_in ! oggetto da trasformare
1982type(vol7d), INTENT(inout) :: vol7d_out ! oggetto trasformato
1983CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname ! imposta il network in vol7d_out (default='generic')
1984LOGICAL,OPTIONAL,INTENT(in) :: noconvert !< do not try to match variable and convert values during transform
1985
1986INTEGER :: nntime, nana, ntime, ntimerange, nlevel, nvar, stallo
1987INTEGER :: itime, itimerange, ivar, inetwork
1988REAL,ALLOCATABLE :: voldatir_out(:,:,:)
1989TYPE(conv_func),POINTER :: c_func(:)
1990TYPE(datetime),ALLOCATABLE :: validitytime(:,:)
1991REAL,POINTER :: voldatiin(:,:,:)
1992
1993#ifdef DEBUG
1994call l4f_category_log(volgrid6d_in%category,l4f_debug,"start volgrid6d_v7d_transform_compute")
1995#endif
1996
1997ntime=0
1998ntimerange=0
1999nlevel=0
2000nvar=0
2001NULLIFY(c_func)
2002
2003if (present(networkname))then
2004 call init(vol7d_out%network(1),name=networkname)
2005else
2006 call init(vol7d_out%network(1),name='generic')
2007end if
2008
2009if (associated(volgrid6d_in%timerange))then
2010 ntimerange=size(volgrid6d_in%timerange)
2011 vol7d_out%timerange=volgrid6d_in%timerange
2012end if
2013
2014if (associated(volgrid6d_in%time))then
2015 ntime=size(volgrid6d_in%time)
2016
2017 if (vol7d_out%time_definition == volgrid6d_in%time_definition) then
2018
2019 ! i time sono definiti uguali: assegno
2020 vol7d_out%time=volgrid6d_in%time
2021
2022 else
2023 ! converto reference in validity
2024 allocate (validitytime(ntime,ntimerange),stat=stallo)
2025 if (stallo /=0)then
2026 call l4f_category_log(volgrid6d_in%category,l4f_fatal,"allocating memory")
2027 call raise_fatal_error()
2028 end if
2029
2030 do itime=1,ntime
2031 do itimerange=1,ntimerange
2032 if (vol7d_out%time_definition > volgrid6d_in%time_definition) then
2033 validitytime(itime,itimerange) = &
2034 volgrid6d_in%time(itime) + timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
2035 else
2036 validitytime(itime,itimerange) = &
2037 volgrid6d_in%time(itime) - timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
2038 end if
2039 end do
2040 end do
2041
2042 nntime = count_distinct(reshape(validitytime,(/ntime*ntimerange/)), back=.true.)
2043 vol7d_out%time=pack_distinct(reshape(validitytime,(/ntime*ntimerange/)), nntime,back=.true.)
2044
2045 end if
2046end if
2047
2048IF (ASSOCIATED(volgrid6d_in%level)) THEN
2049 nlevel = SIZE(volgrid6d_in%level)
2050 vol7d_out%level=volgrid6d_in%level
2051ENDIF
2052
2053IF (ASSOCIATED(volgrid6d_in%var)) THEN
2054 nvar = SIZE(volgrid6d_in%var)
2055 IF (.NOT. optio_log(noconvert)) THEN
2056 CALL vargrib2varbufr(volgrid6d_in%var, vol7d_out%dativar%r, c_func)
2057 ENDIF
2058ENDIF
2059
2060nana = SIZE(vol7d_out%ana)
2061
2062! allocate once for speed
2063IF (.NOT.ASSOCIATED(volgrid6d_in%voldati)) THEN
2064 ALLOCATE(voldatiin(volgrid6d_in%griddim%dim%nx, volgrid6d_in%griddim%dim%ny, &
2065 nlevel))
2066ENDIF
2067
2068ALLOCATE(voldatir_out(nana,1,nlevel),stat=stallo)
2069IF (stallo /= 0) THEN
2070 CALL l4f_category_log(volgrid6d_in%category,l4f_fatal,"allocating memory")
2071 CALL raise_fatal_error()
2072ENDIF
2073
2074inetwork=1
2075do itime=1,ntime
2076 do itimerange=1,ntimerange
2077! do ilevel=1,nlevel
2078 do ivar=1,nvar
2079
2080 !non è chiaro se questa sezione è utile o no
2081 !ossia il compute sotto sembra prevedere voldatir_out solo in out
2082!!$ if (vol7d_out%time_definition == volgrid6d_in%time_definition) then
2083!!$ voldatir_out=reshape(vol7d_out%voldatir(:,itime,ilevel,itimerange,ivar,inetwork),(/nana,1/))
2084!!$ else
2085!!$ voldatir_out=reshape(vol7d_out%voldatir(:,index(vol7d_out%time,validitytime(itime,ilevel)),ilevel,itimerange,ivar,inetwork),(/nana,1/))
2086!!$ end if
2087
2088 CALL volgrid_get_vol_3d(volgrid6d_in, itime, itimerange, ivar, &
2089 voldatiin)
2090
2091 CALL compute(this, voldatiin, voldatir_out, vol7d_out%dativar%r(ivar))
2092
2093 if (vol7d_out%time_definition == volgrid6d_in%time_definition) then
2094 vol7d_out%voldatir(:,itime,:,itimerange,ivar,inetwork) = &
2095 voldatir_out(:,1,:)
2096 else
2097 vol7d_out%voldatir(:,index(vol7d_out%time,validitytime(itime,itimerange)),:,itimerange,ivar,inetwork)=&
2098 reshape(voldatir_out,(/nana,nlevel/))
2099 end if
2100
2101! 1 indice della dimensione "anagrafica"
2102! 2 indice della dimensione "tempo"
2103! 3 indice della dimensione "livello verticale"
2104! 4 indice della dimensione "intervallo temporale"
2105! 5 indice della dimensione "variabile"
2106! 6 indice della dimensione "rete"
2107
2108 end do
2109! end do
2110 end do
2111end do
2112
2113deallocate(voldatir_out)
2114IF (.NOT.ASSOCIATED(volgrid6d_in%voldati)) THEN
2115 DEALLOCATE(voldatiin)
2116ENDIF
2117if (allocated(validitytime)) deallocate(validitytime)
2118
2119! Rescale valid data according to variable conversion table
2120IF (ASSOCIATED(c_func)) THEN
2121 DO ivar = 1, nvar
2122 CALL compute(c_func(ivar), vol7d_out%voldatir(:,:,:,:,ivar,:))
2123 ENDDO
2124 DEALLOCATE(c_func)
2125ENDIF
2126
2127end SUBROUTINE volgrid6d_v7d_transform_compute
2128
2129
2130!> Performs the specified abstract transformation on the data provided.
2131!! The abstract transformation is specified by \a this parameter; the
2132!! corresponding specifical transformation (\a grid_transform object)
2133!! is created and destroyed internally. The output transformed object
2134!! is created internally and it does not require preliminary
2135!! initialisation.
2136SUBROUTINE volgrid6d_v7d_transform(this, volgrid6d_in, vol7d_out, v7d, &
2137 maskgrid, maskbounds, networkname, noconvert, find_index, categoryappend)
2138TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
2139TYPE(volgrid6d),INTENT(inout) :: volgrid6d_in !< object to be transformed, it is not modified, despite the INTENT(inout)
2140TYPE(vol7d),INTENT(out) :: vol7d_out !< transformed object, it does not requires initialisation
2141TYPE(vol7d),INTENT(in),OPTIONAL :: v7d !< object containing a list of points over which transformation has to be done (required by some transformation types)
2142REAL,INTENT(in),OPTIONAL :: maskgrid(:,:) !< 2D field to be used for defining subareas according to its values, it must have the same shape as the field to be interpolated (for transformation type 'maskinter')
2143REAL,INTENT(in),OPTIONAL :: maskbounds(:) !< array of boundary values for defining subareas from the values of \a maskgrid, the number of subareas is SIZE(maskbounds) - 1, if not provided a default based on extreme values of \a makgrid is used
2144CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname !< set the output network name in vol7d_out (default='generic')
2145LOGICAL,OPTIONAL,INTENT(in) :: noconvert !< do not try to match variable and convert values during transform
2146PROCEDURE(basic_find_index),POINTER,OPTIONAL :: find_index
2147CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
2148
2149type(grid_transform) :: grid_trans
2150INTEGER :: ntime, ntimerange, nlevel, nvar, nana, time_definition, nnetwork, stallo
2151INTEGER :: itime, itimerange, inetwork
2152TYPE(datetime),ALLOCATABLE :: validitytime(:,:)
2153INTEGER,ALLOCATABLE :: point_index(:)
2154TYPE(vol7d) :: v7d_locana
2155
2156#ifdef DEBUG
2157call l4f_category_log(volgrid6d_in%category,l4f_debug,"start volgrid6d_v7d_transform")
2158#endif
2159
2160call vg6d_wind_unrot(volgrid6d_in)
2161
2162ntime=0
2163ntimerange=0
2164nlevel=0
2165nvar=0
2166nnetwork=1
2167
2168call get_val(this,time_definition=time_definition)
2169if (.not. c_e(time_definition)) then
2170 time_definition=1 ! default to validity time
2171endif
2172
2173IF (PRESENT(v7d)) THEN
2174 CALL vol7d_copy(v7d, v7d_locana)
2175ELSE
2176 CALL init(v7d_locana)
2177ENDIF
2178
2179if (associated(volgrid6d_in%timerange)) ntimerange=size(volgrid6d_in%timerange)
2180
2181if (associated(volgrid6d_in%time)) then
2182
2183 ntime=size(volgrid6d_in%time)
2184
2185 if (time_definition /= volgrid6d_in%time_definition) then
2186
2187 ! converto reference in validity
2188 allocate (validitytime(ntime,ntimerange),stat=stallo)
2189 if (stallo /=0)then
2190 call l4f_category_log(volgrid6d_in%category,l4f_fatal,"allocating memory")
2191 call raise_fatal_error()
2192 end if
2193
2194 do itime=1,ntime
2195 do itimerange=1,ntimerange
2196 if (time_definition > volgrid6d_in%time_definition) then
2197 validitytime(itime,itimerange) = &
2198 volgrid6d_in%time(itime) + timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
2199 else
2200 validitytime(itime,itimerange) = &
2201 volgrid6d_in%time(itime) - timedelta_new(sec=volgrid6d_in%timerange(itimerange)%p1)
2202 end if
2203 end do
2204 end do
2205
2206 ntime = count_distinct(reshape(validitytime,(/ntime*ntimerange/)), back=.true.)
2207 deallocate (validitytime)
2208
2209 end if
2210end if
2211
2212
2213if (associated(volgrid6d_in%level)) nlevel=size(volgrid6d_in%level)
2214if (associated(volgrid6d_in%var)) nvar=size(volgrid6d_in%var)
2215
2216CALL init(grid_trans, this, volgrid6d_in%griddim, v7d_locana, &
2217 maskgrid=maskgrid, maskbounds=maskbounds, find_index=find_index, &
2218 categoryappend=categoryappend)
2219CALL init (vol7d_out,time_definition=time_definition)
2220
2221IF (c_e(grid_trans)) THEN
2222
2223 nana=SIZE(v7d_locana%ana)
2224 CALL vol7d_alloc(vol7d_out, nana=nana, ntime=ntime, nlevel=nlevel, &
2225 ntimerange=ntimerange, ndativarr=nvar, nnetwork=nnetwork)
2226 vol7d_out%ana = v7d_locana%ana
2227
2228 CALL get_val(grid_trans, output_point_index=point_index)
2229 IF (ALLOCATED(point_index)) THEN
2230! check that size(point_index) == nana?
2231 CALL vol7d_alloc(vol7d_out, nanavari=1)
2232 CALL init(vol7d_out%anavar%i(1), 'B01192')
2233 ENDIF
2234
2235 CALL vol7d_alloc_vol(vol7d_out)
2236
2237 IF (ALLOCATED(point_index)) THEN
2238 DO inetwork = 1, nnetwork
2239 vol7d_out%volanai(:,1,inetwork) = point_index(:)
2240 ENDDO
2241 ENDIF
2242 CALL compute(grid_trans, volgrid6d_in, vol7d_out, networkname, noconvert)
2243ELSE
2244 CALL l4f_log(l4f_error, 'vg6d_v7d_transform: transformation not valid')
2245 CALL raise_error()
2246ENDIF
2247
2248CALL delete(grid_trans)
2249
2250#ifdef HAVE_DBALLE
2251! set variables to a conformal state
2252CALL vol7d_dballe_set_var_du(vol7d_out)
2253#endif
2254
2255CALL delete(v7d_locana)
2256
2257END SUBROUTINE volgrid6d_v7d_transform
2258
2259
2260!> Performs the specified abstract transformation on the arrays of
2261!! data provided. The abstract transformation is specified by \a this
2262!! parameter; the corresponding specifical transformation (\a
2263!! grid_transform object) is created and destroyed internally. The
2264!! output transformed object is created internally and it does not
2265!! require preliminary initialisation. The transformation performed on
2266!! each element of the input \a volgrid6d array object is merged into
2267!! a single \a vol7d output object.
2268SUBROUTINE volgrid6dv_v7d_transform(this, volgrid6d_in, vol7d_out, v7d, &
2269 maskgrid, maskbounds, networkname, noconvert, find_index, categoryappend)
2270TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
2271TYPE(volgrid6d),INTENT(inout) :: volgrid6d_in(:) !< object to be transformed, it is an array of volgrid6d objects, each of which will be transformed, it is not modified, despite the INTENT(inout)
2272TYPE(vol7d),INTENT(out) :: vol7d_out !< transformed object, it does not require initialisation
2273TYPE(vol7d),INTENT(in),OPTIONAL :: v7d !< object containing a list of points over which transformation has to be done (required by some transformation types)
2274REAL,INTENT(in),OPTIONAL :: maskgrid(:,:) !< 2D field to be used for defining subareas according to its values, it must have the same shape as the field to be interpolated (for transformation type 'maskinter')
2275REAL,INTENT(in),OPTIONAL :: maskbounds(:) !< array of boundary values for defining subareas from the values of \a maskgrid, the number of subareas is SIZE(maskbounds) - 1, if not provided a default based on extreme values of \a makgrid is used
2276CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname !< set the output network name in vol7d_out (default='generic')
2277LOGICAL,OPTIONAL,INTENT(in) :: noconvert !< do not try to match variable and convert values during transform
2278PROCEDURE(basic_find_index),POINTER,OPTIONAL :: find_index
2279CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
2280
2281integer :: i
2282TYPE(vol7d) :: v7dtmp
2283
2284
2285CALL init(v7dtmp)
2286CALL init(vol7d_out)
2287
2288DO i=1,SIZE(volgrid6d_in)
2289 CALL transform(this, volgrid6d_in(i), v7dtmp, v7d=v7d, &
2290 maskgrid=maskgrid, maskbounds=maskbounds, &
2291 networkname=networkname, noconvert=noconvert, find_index=find_index, &
2292 categoryappend=categoryappend)
2293 CALL vol7d_append(vol7d_out, v7dtmp)
2294ENDDO
2295
2296END SUBROUTINE volgrid6dv_v7d_transform
2297
2298
2299! Internal method for performing sparse point to grid computations
2300SUBROUTINE v7d_volgrid6d_transform_compute(this, vol7d_in, volgrid6d_out, networkname, gaid_template)
2301TYPE(grid_transform),INTENT(in) :: this ! object specifying the specific transformation
2302type(vol7d), INTENT(in) :: vol7d_in ! object to be transformed
2303type(volgrid6d), INTENT(inout) :: volgrid6d_out ! transformed object
2304CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname ! select the network to be processed from the \a vol7d input object, default the first network
2305TYPE(grid_id),OPTIONAL,INTENT(in) :: gaid_template ! the template (typically grib_api) to be associated with output data, it also helps in improving variable conversion
2306
2307integer :: nana, ntime, ntimerange, nlevel, nvar
2308INTEGER :: ilevel, itime, itimerange, ivar, inetwork
2309
2310REAL,POINTER :: voldatiout(:,:,:)
2311type(vol7d_network) :: network
2312TYPE(conv_func), pointer :: c_func(:)
2313!TODO category sarebbe da prendere da vol7d
2314#ifdef DEBUG
2315CALL l4f_category_log(volgrid6d_out%category, l4f_debug, &
2316 'start v7d_volgrid6d_transform_compute')
2317#endif
2318
2319ntime=0
2320ntimerange=0
2321nlevel=0
2322nvar=0
2323
2324IF (PRESENT(networkname)) THEN
2325 CALL init(network,name=networkname)
2326 inetwork = index(vol7d_in%network,network)
2327 IF (inetwork <= 0) THEN
2328 CALL l4f_category_log(volgrid6d_out%category,l4f_warn, &
2329 'network '//trim(networkname)//' not found, first network will be transformed')
2330 inetwork = 1
2331 ENDIF
2332ELSE
2333 inetwork = 1
2334ENDIF
2335
2336! no time_definition conversion implemented, output will be the same as input
2337if (associated(vol7d_in%time))then
2338 ntime=size(vol7d_in%time)
2339 volgrid6d_out%time=vol7d_in%time
2340end if
2341
2342if (associated(vol7d_in%timerange))then
2343 ntimerange=size(vol7d_in%timerange)
2344 volgrid6d_out%timerange=vol7d_in%timerange
2345end if
2346
2347if (associated(vol7d_in%level))then
2348 nlevel=size(vol7d_in%level)
2349 volgrid6d_out%level=vol7d_in%level
2350end if
2351
2352if (associated(vol7d_in%dativar%r))then
2353 nvar=size(vol7d_in%dativar%r)
2354 CALL varbufr2vargrib(vol7d_in%dativar%r, volgrid6d_out%var, c_func, gaid_template)
2355end if
2356
2357nana=SIZE(vol7d_in%voldatir, 1)
2358! allocate once for speed
2359IF (.NOT.ASSOCIATED(volgrid6d_out%voldati)) THEN
2360 ALLOCATE(voldatiout(volgrid6d_out%griddim%dim%nx, volgrid6d_out%griddim%dim%ny, &
2361 nlevel))
2362ENDIF
2363
2364DO ivar=1,nvar
2365 DO itimerange=1,ntimerange
2366 DO itime=1,ntime
2367
2368! clone the gaid template where I have data
2369 IF (PRESENT(gaid_template)) THEN
2370 DO ilevel = 1, nlevel
2371 IF (any(c_e(vol7d_in%voldatir(:,itime,ilevel,itimerange,ivar,inetwork)))) THEN
2372 CALL copy(gaid_template, volgrid6d_out%gaid(ilevel,itime,itimerange,ivar))
2373 ELSE
2374 volgrid6d_out%gaid(ilevel,itime,itimerange,ivar) = grid_id_new()
2375 ENDIF
2376 ENDDO
2377 ENDIF
2378
2379! get data
2380 IF (ASSOCIATED(volgrid6d_out%voldati)) & ! improve!!!!
2381 CALL volgrid_get_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
2382 voldatiout)
2383! do the interpolation
2384 CALL compute(this, &
2385 vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), voldatiout, &
2386 vol7d_in%dativar%r(ivar))
2387! rescale valid data according to variable conversion table
2388 IF (ASSOCIATED(c_func)) THEN
2389 CALL compute(c_func(ivar), voldatiout(:,:,:))
2390 ENDIF
2391! put data
2392 CALL volgrid_set_vol_3d(volgrid6d_out, itime, itimerange, ivar, &
2393 voldatiout)
2394
2395 ENDDO
2396 ENDDO
2397ENDDO
2398
2399IF (.NOT.ASSOCIATED(volgrid6d_out%voldati)) THEN
2400 DEALLOCATE(voldatiout)
2401ENDIF
2402IF (ASSOCIATED(c_func)) THEN
2403 DEALLOCATE(c_func)
2404ENDIF
2405
2406END SUBROUTINE v7d_volgrid6d_transform_compute
2407
2408
2409!> Performs the specified abstract transformation on the data provided.
2410!! The abstract transformation is specified by \a this parameter; the
2411!! corresponding specifical transformation (\a grid_transform object)
2412!! is created and destroyed internally. The output transformed object
2413!! is created internally and it does not require preliminary
2414!! initialisation.
2415SUBROUTINE v7d_volgrid6d_transform(this, griddim, vol7d_in, volgrid6d_out, &
2416 networkname, gaid_template, categoryappend)
2417TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
2418TYPE(griddim_def),INTENT(in),OPTIONAL :: griddim !< griddim specifying the output grid (required by most transformation types)
2419! TODO ripristinare intent(in) dopo le opportune modifiche in grid_class.F90
2420TYPE(vol7d),INTENT(inout) :: vol7d_in !< object to be transformed, it is not modified, despite the INTENT(inout)
2421TYPE(volgrid6d),INTENT(out) :: volgrid6d_out !< transformed object, it does not require initialisation
2422CHARACTER(len=*),OPTIONAL,INTENT(in) :: networkname !< select the network to be processed from the \a vol7d input object, default the first network
2423TYPE(grid_id),OPTIONAL,INTENT(in) :: gaid_template !< the template (typically grib_api) to be associated with output data, it also helps in improving variable conversion
2424CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
2425
2426type(grid_transform) :: grid_trans
2427integer :: ntime, ntimerange, nlevel, nvar
2428
2429
2430!TODO la category sarebbe da prendere da vol7d
2431!call l4f_category_log(vol7d_out%category,L4F_DEBUG,"start volgrid6d_transform")
2432
2433CALL vol7d_alloc_vol(vol7d_in) ! be safe
2434ntime=SIZE(vol7d_in%time)
2435ntimerange=SIZE(vol7d_in%timerange)
2436nlevel=SIZE(vol7d_in%level)
2437nvar=0
2438if (associated(vol7d_in%dativar%r)) nvar=size(vol7d_in%dativar%r)
2439
2440IF (nvar <= 0) THEN ! use vol7d category once it will be implemented
2441 CALL l4f_log(l4f_error, &
2442 "trying to transform a vol7d object incomplete or without real variables")
2443 CALL init(volgrid6d_out) ! initialize to empty
2444 CALL raise_error()
2445 RETURN
2446ENDIF
2447
2448CALL init(grid_trans, this, vol7d_in, griddim, categoryappend=categoryappend)
2449CALL init(volgrid6d_out, griddim, time_definition=vol7d_in%time_definition, &
2450 categoryappend=categoryappend)
2451
2452IF (c_e(grid_trans)) THEN
2453
2454 CALL volgrid6d_alloc(volgrid6d_out, griddim%dim, ntime=ntime, nlevel=nlevel, &
2455 ntimerange=ntimerange, nvar=nvar)
2456! can I avoid decode=.TRUE. here, using gaid_template?
2457 CALL volgrid6d_alloc_vol(volgrid6d_out, decode=.true.)
2458
2459 CALL compute(grid_trans, vol7d_in, volgrid6d_out, networkname, gaid_template)
2460
2461 CALL vg6d_wind_rot(volgrid6d_out)
2462ELSE
2463! should log with grid_trans%category, but it is private
2464 CALL l4f_log(l4f_error, 'v7d_vg6d_transform: transformation not valid')
2465 CALL raise_error()
2466ENDIF
2467
2468CALL delete(grid_trans)
2469
2470END SUBROUTINE v7d_volgrid6d_transform
2471
2472
2473! Internal method for performing sparse point to sparse point computations
2474SUBROUTINE v7d_v7d_transform_compute(this, vol7d_in, vol7d_out, lev_out, &
2475 var_coord_vol)
2476TYPE(grid_transform),INTENT(in) :: this ! oggetto di trasformazione per grigliato
2477type(vol7d), INTENT(in) :: vol7d_in ! oggetto da trasformare
2478type(vol7d), INTENT(inout) :: vol7d_out ! oggetto trasformato
2479TYPE(vol7d_level),INTENT(in),OPTIONAL :: lev_out(:) ! vol7d_level object defining target vertical grid, for vertical interpolations
2480INTEGER,INTENT(in),OPTIONAL :: var_coord_vol ! index of variable defining vertical coordinate values in input volume
2481
2482INTEGER :: itime, itimerange, ilevel, ivar, inetwork, &
2483 levshift, levused, lvar_coord_vol, spos
2484REAL,ALLOCATABLE :: coord_3d_in(:,:,:)
2485TYPE(vol7d_level) :: output_levtype
2486
2487lvar_coord_vol = optio_i(var_coord_vol)
2488vol7d_out%time(:) = vol7d_in%time(:)
2489vol7d_out%timerange(:) = vol7d_in%timerange(:)
2490IF (PRESENT(lev_out)) THEN
2491 vol7d_out%level(:) = lev_out(:)
2492ELSE
2493 vol7d_out%level(:) = vol7d_in%level(:)
2494ENDIF
2495vol7d_out%network(:) = vol7d_in%network(:)
2496IF (ASSOCIATED(vol7d_in%dativar%r)) THEN ! work only when real vars are available
2497 vol7d_out%dativar%r(:) = vol7d_in%dativar%r(:)
2498
2499 CALL get_val(this, levshift=levshift, levused=levused)
2500 spos = imiss
2501 IF (c_e(lvar_coord_vol)) THEN
2502 CALL get_val(this%trans, output_levtype=output_levtype)
2503 IF (output_levtype%level1 == 103 .OR. output_levtype%level1 == 108) THEN
2504 spos = firsttrue(vol7d_in%level(:) == vol7d_level_new(1))
2505 IF (spos == 0) THEN
2506 CALL l4f_log(l4f_error, &
2507 'output level '//t2c(output_levtype%level1)// &
2508 ' requested, but height/press of surface not provided in volume')
2509 ENDIF
2510 IF (.NOT.c_e(levshift) .AND. .NOT.c_e(levused)) THEN
2511 CALL l4f_log(l4f_error, &
2512 'internal inconsistence, levshift and levused undefined when they should be')
2513 ENDIF
2514 ALLOCATE(coord_3d_in(SIZE(vol7d_in%ana),1,SIZE(vol7d_in%level)))
2515 ENDIF
2516
2517 ENDIF
2518
2519 DO inetwork = 1, SIZE(vol7d_in%network)
2520 DO ivar = 1, SIZE(vol7d_in%dativar%r)
2521 DO itimerange = 1, SIZE(vol7d_in%timerange)
2522 DO itime = 1, SIZE(vol7d_in%time)
2523
2524! dirty trick to make voldatir look like a 2d-array of shape (nana,1)
2525 IF (c_e(lvar_coord_vol)) THEN
2526 IF (c_e(spos)) THEN ! compute difference wrt surface coordinate
2527 IF (spos == 0) THEN ! error condition, set all to missing and goodnight
2528 coord_3d_in(:,:,levshift+1:levshift+levused) = rmiss
2529 ELSE
2530 DO ilevel = levshift+1, levshift+levused
2531 WHERE(c_e(vol7d_in%voldatir(:,itime:itime,ilevel,itimerange,lvar_coord_vol,inetwork)) .AND. &
2532 c_e(vol7d_in%voldatir(:,itime:itime,spos,itimerange,lvar_coord_vol,inetwork)))
2533 coord_3d_in(:,:,ilevel) = vol7d_in%voldatir(:,itime:itime,ilevel,itimerange,lvar_coord_vol,inetwork) - &
2534 vol7d_in%voldatir(:,itime:itime,spos,itimerange,lvar_coord_vol,inetwork)
2535 ELSEWHERE
2536 coord_3d_in(:,:,ilevel) = rmiss
2537 END WHERE
2538 ENDDO
2539 ENDIF
2540 CALL compute(this, &
2541 vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), &
2542 vol7d_out%voldatir(:,itime:itime,:,itimerange,ivar,inetwork), &
2543 var=vol7d_in%dativar%r(ivar), &
2544 coord_3d_in=coord_3d_in)
2545 ELSE
2546 CALL compute(this, &
2547 vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), &
2548 vol7d_out%voldatir(:,itime:itime,:,itimerange,ivar,inetwork), &
2549 var=vol7d_in%dativar%r(ivar), &
2550 coord_3d_in=vol7d_in%voldatir(:,itime:itime,:,itimerange, &
2551 lvar_coord_vol,inetwork))
2552 ENDIF
2553 ELSE
2554 CALL compute(this, &
2555 vol7d_in%voldatir(:,itime,:,itimerange,ivar,inetwork), &
2556 vol7d_out%voldatir(:,itime:itime,:,itimerange,ivar,inetwork), &
2557 var=vol7d_in%dativar%r(ivar))
2558 ENDIF
2559 ENDDO
2560 ENDDO
2561 ENDDO
2562 ENDDO
2563
2564ENDIF
2565
2566END SUBROUTINE v7d_v7d_transform_compute
2567
2568
2569!> Performs the specified abstract transformation on the data provided.
2570!! The abstract transformation is specified by \a this parameter; the
2571!! corresponding specifical transformation (\a grid_transform object)
2572!! is created and destroyed internally. The output transformed object
2573!! is created internally and it does not require preliminary
2574!! initialisation. The success of the transformation can be checked
2575!! with the \a c_e method: c_e(vol7d_out).
2576SUBROUTINE v7d_v7d_transform(this, vol7d_in, vol7d_out, v7d, maskbounds, &
2577 lev_out, vol7d_coord_in, categoryappend)
2578TYPE(transform_def),INTENT(in) :: this !< object specifying the abstract transformation
2579TYPE(vol7d),INTENT(inout) :: vol7d_in !< object to be transformed, it is not modified, despite the INTENT(inout)
2580TYPE(vol7d),INTENT(out) :: vol7d_out !< transformed object, it does not require initialisation
2581TYPE(vol7d),INTENT(in),OPTIONAL :: v7d !< object containing a list of points over which transformation has to be done (required by some transformation types)
2582REAL,INTENT(in),OPTIONAL :: maskbounds(:) !< array of boundary values for defining a subset of valid points where the values of \a maskgrid are within the first and last value of \a maskbounds (for transformation type 'metamorphosis:maskfill')
2583TYPE(vol7d_level),INTENT(in),OPTIONAL,TARGET :: lev_out(:) !< vol7d_level object defining target vertical grid, for vertical interpolations
2584TYPE(vol7d),INTENT(in),OPTIONAL :: vol7d_coord_in !< object providing time constant input vertical coordinate for some kind of vertical interpolations
2585CHARACTER(len=*),INTENT(in),OPTIONAL :: categoryappend !< append this suffix to log4fortran namespace category
2586
2587INTEGER :: nvar, inetwork
2588TYPE(grid_transform) :: grid_trans
2589TYPE(vol7d_level),POINTER :: llev_out(:)
2590TYPE(vol7d_level) :: input_levtype, output_levtype
2591TYPE(vol7d_var) :: vcoord_var
2592REAL,ALLOCATABLE :: coord_3d_in(:,:,:)
2593INTEGER :: var_coord_in, var_coord_vol, i, k, ulstart, ulend, spos
2594INTEGER,ALLOCATABLE :: point_index(:)
2595TYPE(vol7d) :: v7d_locana, vol7d_tmpana
2596CHARACTER(len=80) :: trans_type
2597LOGICAL,ALLOCATABLE :: mask_in(:), point_mask(:)
2598
2599CALL vol7d_alloc_vol(vol7d_in) ! be safe
2600nvar=0
2601IF (ASSOCIATED(vol7d_in%dativar%r)) nvar=SIZE(vol7d_in%dativar%r)
2602
2603CALL init(v7d_locana)
2604IF (PRESENT(v7d)) v7d_locana = v7d
2605CALL init(vol7d_out, time_definition=vol7d_in%time_definition)
2606
2607CALL get_val(this, trans_type=trans_type)
2608
2609var_coord_vol = imiss
2610IF (trans_type == 'vertint') THEN
2611
2612 IF (PRESENT(lev_out)) THEN
2613
2614! if vol7d_coord_in provided and allocated, check that it fits
2615 var_coord_in = -1
2616 IF (PRESENT(vol7d_coord_in)) THEN
2617 IF (ASSOCIATED(vol7d_coord_in%voldatir) .AND. &
2618 ASSOCIATED(vol7d_coord_in%dativar%r)) THEN
2619
2620! strictly 1 time, 1 timerange and 1 network
2621 IF (SIZE(vol7d_coord_in%voldatir,2) /= 1 .OR. &
2622 SIZE(vol7d_coord_in%voldatir,4) /= 1 .OR. &
2623 SIZE(vol7d_coord_in%voldatir,6) /= 1) THEN
2624 CALL l4f_log(l4f_error, &
2625 'volume providing constant input vertical coordinate must have &
2626 &only 1 time, 1 timerange and 1 network')
2627 CALL raise_error()
2628 RETURN
2629 ENDIF
2630
2631! search for variable providing vertical coordinate
2632 CALL get_val(this, output_levtype=output_levtype)
2633 vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
2634 IF (.NOT.c_e(vcoord_var)) THEN
2635 CALL l4f_log(l4f_error, &
2636 'requested output level type '//t2c(output_levtype%level1)// &
2637 ' does not correspond to any known physical variable for &
2638 &providing vertical coordinate')
2639 CALL raise_error()
2640 RETURN
2641 ENDIF
2642
2643 var_coord_in = index(vol7d_coord_in%dativar%r, vcoord_var)
2644
2645 IF (var_coord_in <= 0) THEN
2646 CALL l4f_log(l4f_error, &
2647 'volume providing constant input vertical coordinate contains no &
2648 &real variables matching output level type '//t2c(output_levtype%level1))
2649 CALL raise_error()
2650 RETURN
2651 ENDIF
2652 CALL l4f_log(l4f_info, &
2653 'Coordinate for vertint found in coord volume at position '// &
2654 t2c(var_coord_in))
2655
2656! check vertical coordinate system
2657 CALL get_val(this, input_levtype=input_levtype)
2658 mask_in = & ! implicit allocation
2659 (vol7d_coord_in%level(:)%level1 == input_levtype%level1) .AND. &
2660 (vol7d_coord_in%level(:)%level2 == input_levtype%level2)
2661 ulstart = firsttrue(mask_in)
2662 ulend = lasttrue(mask_in)
2663 IF (ulstart == 0 .OR. ulend == 0) THEN
2664 CALL l4f_log(l4f_error, &
2665 'coordinate file does not contain levels of type '// &
2666 t2c(input_levtype%level1)//'/'//t2c(input_levtype%level2)// &
2667 ' specified for input data')
2668 CALL raise_error()
2669 RETURN
2670 ENDIF
2671
2672 coord_3d_in = vol7d_coord_in%voldatir(:,1:1,ulstart:ulend,1,var_coord_in,1) ! dirty 1:1, implicit allocation
2673! special case
2674 IF (output_levtype%level1 == 103 &
2675 .OR. output_levtype%level1 == 108) THEN ! surface coordinate needed
2676 spos = firsttrue(vol7d_coord_in%level(:) == vol7d_level_new(1))
2677 IF (spos == 0) THEN
2678 CALL l4f_log(l4f_error, &
2679 'output level '//t2c(output_levtype%level1)// &
2680 ' requested, but height/press of surface not provided in coordinate file')
2681 CALL raise_error()
2682 RETURN
2683 ENDIF
2684 DO k = 1, SIZE(coord_3d_in,3)
2685 WHERE(c_e(coord_3d_in(:,:,k)) .AND. &
2686 c_e(vol7d_coord_in%voldatir(:,1:1,spos,1,var_coord_in,1)))
2687 coord_3d_in(:,:,k) = coord_3d_in(:,:,k) - &
2688 vol7d_coord_in%voldatir(:,1:1,spos,1,var_coord_in,1)
2689 ELSEWHERE
2690 coord_3d_in(:,:,k) = rmiss
2691 END WHERE
2692 ENDDO
2693 ENDIF
2694
2695 ENDIF
2696 ENDIF
2697
2698 IF (var_coord_in <= 0) THEN ! search for coordinate within volume
2699! search for variable providing vertical coordinate
2700 CALL get_val(this, output_levtype=output_levtype)
2701 vcoord_var = vol7d_var_new(vol7d_level_to_var(output_levtype))
2702 IF (c_e(vcoord_var)) THEN
2703 DO i = 1, SIZE(vol7d_in%dativar%r)
2704 IF (vol7d_in%dativar%r(i) == vcoord_var) THEN
2705 var_coord_vol = i
2706 EXIT
2707 ENDIF
2708 ENDDO
2709
2710 IF (c_e(var_coord_vol)) THEN
2711 CALL l4f_log(l4f_info, &
2712 'Coordinate for vertint found in input volume at position '// &
2713 t2c(var_coord_vol))
2714 ENDIF
2715
2716 ENDIF
2717 ENDIF
2718
2719 IF (var_coord_in > 0) THEN
2720 CALL init(grid_trans, this, lev_in=vol7d_in%level, lev_out=lev_out, &
2721 coord_3d_in=coord_3d_in, categoryappend=categoryappend)
2722 ELSE
2723 CALL init(grid_trans, this, lev_in=vol7d_in%level, lev_out=lev_out, &
2724 categoryappend=categoryappend)
2725 ENDIF
2726
2727 CALL get_val(grid_trans, output_level_auto=llev_out) ! get levels if auto-generated
2728 IF (.NOT.associated(llev_out)) llev_out => lev_out
2729
2730 IF (c_e(grid_trans)) then! .AND. nvar > 0) THEN
2731
2732 CALL vol7d_alloc(vol7d_out, nana=SIZE(vol7d_in%ana), &
2733 ntime=SIZE(vol7d_in%time), ntimerange=SIZE(vol7d_in%timerange), &
2734 nlevel=SIZE(llev_out), nnetwork=SIZE(vol7d_in%network), ndativarr=nvar)
2735 vol7d_out%ana(:) = vol7d_in%ana(:)
2736
2737 CALL vol7d_alloc_vol(vol7d_out)
2738
2739! no need to check c_e(var_coord_vol) here since the presence of
2740! this%coord_3d_in (external) has precedence over coord_3d_in internal
2741! in grid_transform_compute
2742 CALL compute(grid_trans, vol7d_in, vol7d_out, llev_out, &
2743 var_coord_vol=var_coord_vol)
2744 ELSE
2745 CALL l4f_log(l4f_error, 'v7d_v7d_transform: transformation not valid')
2746 CALL raise_error()
2747 ENDIF
2748 ELSE
2749 CALL l4f_log(l4f_error, &
2750 'v7d_v7d_transform: vertint requested but lev_out not provided')
2751 CALL raise_error()
2752 ENDIF
2753
2754ELSE
2755
2756 CALL init(grid_trans, this, vol7d_in, v7d_locana, maskbounds=maskbounds, &
2757 categoryappend=categoryappend)
2758! if this init is successful, I am sure that v7d_locana%ana is associated
2759
2760 IF (c_e(grid_trans)) then! .AND. nvar > 0) THEN
2761
2762 CALL vol7d_alloc(vol7d_out, nana=SIZE(v7d_locana%ana), &
2763 ntime=SIZE(vol7d_in%time), ntimerange=SIZE(vol7d_in%timerange), &
2764 nlevel=SIZE(vol7d_in%level), nnetwork=SIZE(vol7d_in%network), ndativarr=nvar)
2765 vol7d_out%ana = v7d_locana%ana
2766
2767 CALL get_val(grid_trans, point_mask=point_mask, output_point_index=point_index)
2768
2769 IF (ALLOCATED(point_index)) THEN
2770 CALL vol7d_alloc(vol7d_out, nanavari=1)
2771 CALL init(vol7d_out%anavar%i(1), 'B01192')
2772 ENDIF
2773
2774 CALL vol7d_alloc_vol(vol7d_out)
2775
2776 IF (ALLOCATED(point_index)) THEN
2777 DO inetwork = 1, SIZE(vol7d_in%network)
2778 vol7d_out%volanai(:,1,inetwork) = point_index(:)
2779 ENDDO
2780 ENDIF
2781 CALL compute(grid_trans, vol7d_in, vol7d_out)
2782
2783 IF (ALLOCATED(point_mask)) THEN ! keep full ana
2784 IF (SIZE(point_mask) /= SIZE(vol7d_in%ana)) THEN
2785 CALL l4f_log(l4f_warn, &
2786 'v7d_v7d_transform: inconsistency in point size: '//t2c(SIZE(point_mask)) &
2787 //':'//t2c(SIZE(vol7d_in%ana)))
2788 ELSE
2789#ifdef DEBUG
2790 CALL l4f_log(l4f_debug, 'v7d_v7d_transform: merging ana from in to out')
2791#endif
2792 CALL vol7d_copy(vol7d_in, vol7d_tmpana, &
2793 lana=point_mask, lnetwork=(/.true./), &
2794 ltime=(/.false./), ltimerange=(/.false./), llevel=(/.false./))
2795 CALL vol7d_append(vol7d_out, vol7d_tmpana)
2796 ENDIF
2797 ENDIF
2798
2799 ELSE
2800 CALL l4f_log(l4f_error, 'v7d_v7d_transform: transformation not valid')
2801 CALL raise_error()
2802 ENDIF
2803
2804ENDIF
2805
2806CALL delete (grid_trans)
2807IF (.NOT. PRESENT(v7d)) CALL delete(v7d_locana)
2808
2809END SUBROUTINE v7d_v7d_transform
2810
2811
2812!> Unrotate the wind components.
2813!! It converts u and v components of vector quantities relative to the
2814!! defined grid in the direction of increasing x and y coordinates to
2815!! u and v components relative to easterly and notherly direction. The
2816!! original fields are overwritten.
2817!! \todo Check and correct wind component flag (to be moved in
2818!! griddim_def?)
2819subroutine vg6d_wind_unrot(this)
2820type(volgrid6d) :: this !< object containing wind to be unrotated
2821
2822integer :: component_flag
2823
2824call get_val(this%griddim,component_flag=component_flag)
2825
2826if (component_flag == 1) then
2827 call l4f_category_log(this%category,l4f_info, &
2828 "unrotating vector components")
2829 call vg6d_wind__un_rot(this,.false.)
2830 call set_val(this%griddim,component_flag=0)
2831else
2832 call l4f_category_log(this%category,l4f_info, &
2833 "no need to unrotate vector components")
2834end if
2835
2836end subroutine vg6d_wind_unrot
2837
2838
2839!> Rotate the wind components.
2840!! It converts u and v components of vector quantities
2841!! relative to easterly and notherly direction to
2842!! defined grid in the direction of increasing x and y coordinates.
2843!! The original fields are overwritten.
2844subroutine vg6d_wind_rot(this)
2845type(volgrid6d) :: this !< object containing wind to be rotated
2846
2847integer :: component_flag
2848
2849call get_val(this%griddim,component_flag=component_flag)
2850
2851if (component_flag == 0) then
2852 call l4f_category_log(this%category,l4f_info, &
2853 "rotating vector components")
2854 call vg6d_wind__un_rot(this,.true.)
2855 call set_val(this%griddim,component_flag=1)
2856else
2857 call l4f_category_log(this%category,l4f_info, &
2858 "no need to rotate vector components")
2859end if
2860
2861end subroutine vg6d_wind_rot
2862
2863
2864! Generic UnRotate the wind components.
2865SUBROUTINE vg6d_wind__un_rot(this,rot)
2866TYPE(volgrid6d) :: this ! object containing wind to be (un)rotated
2867LOGICAL :: rot ! if .true. rotate else unrotate
2868
2869INTEGER :: i, j, k, l, a11, a12, a21, a22, stallo
2870double precision,pointer :: rot_mat(:,:,:)
2871real,allocatable :: tmp_arr(:,:)
2872REAL,POINTER :: voldatiu(:,:), voldativ(:,:)
2873INTEGER,POINTER :: iu(:), iv(:)
2874
2875IF (.NOT.ASSOCIATED(this%var)) THEN
2876 CALL l4f_category_log(this%category, l4f_error, &
2877 "trying to unrotate an incomplete volgrid6d object")
2878 CALL raise_fatal_error()
2879! RETURN
2880ENDIF
2881
2882CALL volgrid6d_var_hor_comp_index(this%var, iu, iv)
2883IF (.NOT.ASSOCIATED(iu)) THEN
2884 CALL l4f_category_log(this%category,l4f_error, &
2885 "unrotation impossible")
2886 CALL raise_fatal_error()
2887! RETURN
2888ENDIF
2889
2890! Temporary workspace
2891ALLOCATE(tmp_arr(this%griddim%dim%nx, this%griddim%dim%ny),stat=stallo)
2892IF (stallo /= 0) THEN
2893 CALL l4f_category_log(this%category, l4f_fatal, "allocating memory")
2894 CALL raise_fatal_error()
2895ENDIF
2896! allocate once for speed
2897IF (.NOT.ASSOCIATED(this%voldati)) THEN
2898 ALLOCATE(voldatiu(this%griddim%dim%nx, this%griddim%dim%ny), &
2899 voldativ(this%griddim%dim%nx, this%griddim%dim%ny))
2900ENDIF
2901
2902CALL griddim_unproj(this%griddim)
2903CALL wind_unrot(this%griddim, rot_mat)
2904
2905a11=1
2906if (rot)then
2907 a12=2
2908 a21=3
2909else
2910 a12=3
2911 a21=2
2912end if
2913a22=4
2914
2915DO l = 1, SIZE(iu)
2916 DO k = 1, SIZE(this%timerange)
2917 DO j = 1, SIZE(this%time)
2918 DO i = 1, SIZE(this%level)
2919! get data
2920 CALL volgrid_get_vol_2d(this, i, j, k, iu(l), voldatiu)
2921 CALL volgrid_get_vol_2d(this, i, j, k, iv(l), voldativ)
2922! convert units forward
2923! CALL compute(conv_fwd(iu(l)), voldatiu)
2924! CALL compute(conv_fwd(iv(l)), voldativ)
2925
2926! multiply wind components by rotation matrix
2927 WHERE(voldatiu /= rmiss .AND. voldativ /= rmiss)
2928 tmp_arr(:,:) = real(voldatiu(:,:)*rot_mat(:,:,a11) + &
2929 voldativ(:,:)*rot_mat(:,:,a12))
2930 voldativ(:,:) = real(voldatiu(:,:)*rot_mat(:,:,a21) + &
2931 voldativ(:,:)*rot_mat(:,:,a22))
2932 voldatiu(:,:) = tmp_arr(:,:)
2933 END WHERE
2934! convert units backward
2935! CALL uncompute(conv_fwd(iu(l)), voldatiu)
2936! CALL uncompute(conv_fwd(iv(l)), voldativ)
2937! put data
2938 CALL volgrid_set_vol_2d(this, i, j, k, iu(l), voldatiu)
2939 CALL volgrid_set_vol_2d(this, i, j, k, iv(l), voldativ)
2940 ENDDO
2941 ENDDO
2942 ENDDO
2943ENDDO
2944
2945IF (.NOT.ASSOCIATED(this%voldati)) THEN
2946 DEALLOCATE(voldatiu, voldativ)
2947ENDIF
2948DEALLOCATE(rot_mat, tmp_arr, iu, iv)
2949
2950END SUBROUTINE vg6d_wind__un_rot
2951
2952
2953!!$ try to understand the problem:
2954!!$
2955!!$ case:
2956!!$
2957!!$ 1) we have only one volume: we have to provide the direction of shift
2958!!$ compute H and traslate on it
2959!!$ 2) we have two volumes:
2960!!$ 1) volume U and volume V: compute H and traslate on it
2961!!$ 2) volume U/V and volume H : translate U/V on H
2962!!$ 3) we have tree volumes: translate U and V on H
2963!!$
2964!!$ strange cases:
2965!!$ 1) do not have U in volume U
2966!!$ 2) do not have V in volume V
2967!!$ 3) we have others variables more than U and V in volumes U e V
2968!!$
2969!!$
2970!!$ so the steps are:
2971!!$ 1) find the volumes
2972!!$ 2) define or compute H grid
2973!!$ 3) trasform the volumes in H
2974
2975!!$ N.B.
2976!!$ case 1) for only one vg6d (U or V) is not managed, but
2977!!$ the not pubblic subroutines will work but you have to know what you want to do
2978
2979
2980!> Convert grids type C to type A.
2981!! Use this to interpolate data from grid type C to grid type A
2982!! Grids type are defined by Arakawa.
2983!!
2984!! We need to find these types of area in a vg6d array
2985!! T area of points with temterature etc.
2986!! U area of points with u components of winds
2987!! V area of points with v components of winds
2988!!
2989!! this method works if it finds
2990!! two volumes:
2991!! 1) volume U and volume V: compute H and traslate on it
2992!! 2) volume U/V and volume H : translate U/V on H
2993!! three volumes: translate U and V on H
2994!!
2995!! try to work well on more datasets at once
2996subroutine vg6d_c2a (this)
2997
2998TYPE(volgrid6d),INTENT(inout) :: this(:) !< vettor of volumes volgrid6d to elaborate
2999
3000integer :: ngrid,igrid,jgrid,ugrid,vgrid,tgrid
3001doubleprecision :: xmin, xmax, ymin, ymax
3002doubleprecision :: xmin_t, xmax_t, ymin_t, ymax_t
3003doubleprecision :: step_lon_t,step_lat_t
3004character(len=80) :: type_t,type
3005TYPE(griddim_def):: griddim_t
3006
3007ngrid=size(this)
3008
3009do igrid=1,ngrid
3010
3011 call init(griddim_t)
3012
3013 call get_val(this(igrid)%griddim,xmin=xmin_t, xmax=xmax_t, ymin=ymin_t, ymax=ymax_t,proj_type=type_t)
3014 step_lon_t=(xmax_t-xmin_t)/dble(this(igrid)%griddim%dim%nx-1)
3015 step_lat_t=(ymax_t-ymin_t)/dble(this(igrid)%griddim%dim%ny-1)
3016
3017 do jgrid=1,ngrid
3018
3019 ugrid=imiss
3020 vgrid=imiss
3021 tgrid=imiss
3022
3023#ifdef DEBUG
3024 call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: search U/V/T points:"//to_char(igrid)//to_char(jgrid))
3025 call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: test delta: "//to_char(step_lon_t)//to_char(step_lat_t))
3026#endif
3027
3028 if (this(igrid)%griddim == this(jgrid)%griddim ) cycle
3029
3030 if (this(igrid)%griddim%dim%nx == this(jgrid)%griddim%dim%nx .and. &
3031 this(igrid)%griddim%dim%ny == this(jgrid)%griddim%dim%ny ) then
3033 call get_val(this(jgrid)%griddim,xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,proj_type=type)
3034
3035 if (type_t /= type )cycle
3036
3037#ifdef DEBUG
3038 call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: test U "//&
3039 to_char(xmin)//to_char(xmax)//to_char(ymin)//to_char(ymax))
3040
3041 call l4f_category_log(this(igrid)%category,l4f_debug,"diff coordinate lon"//&
3042 to_char(abs(xmin - (xmin_t+step_lon_t/2.d0)))//&
3043 to_char(abs(xmax - (xmax_t+step_lon_t/2.d0))))
3044 call l4f_category_log(this(igrid)%category,l4f_debug,"diff coordinate lat"//&
3045 to_char(abs(ymin - (ymin_t+step_lat_t/2.d0)))//&
3046 to_char(abs(ymax - (ymax_t+step_lat_t/2.d0))))
3047#endif
3048
3049 if ( abs(xmin - (xmin_t+step_lon_t/2.d0)) < 1.d-3 .and. abs(xmax - (xmax_t+step_lon_t/2.d0)) < 1.d-3 ) then
3050 if ( abs(ymin - ymin_t) < 1.d-3 .and. abs(ymax - ymax_t) < 1.d-3 ) then
3051
3052#ifdef DEBUG
3053 call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: found U")
3054#endif
3055 ugrid=jgrid
3056 tgrid=igrid
3057
3058 end if
3059 end if
3060
3061#ifdef DEBUG
3062 call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: test V "//&
3063 to_char(xmin)//to_char(xmax)//to_char(ymin)//to_char(ymax))
3064#endif
3065
3066 if ( abs(ymin - (ymin_t+step_lat_t/2.d0)) < 1.d-3 .and. abs(ymax - (ymax_t+step_lat_t/2.d0)) < 1.d-3 ) then
3067 if ( abs(xmin - xmin_t) < 1.d-3 .and. abs(xmax - xmax_t) < 1.d-3 ) then
3068
3069#ifdef DEBUG
3070 call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: found V")
3071#endif
3072 vgrid=jgrid
3073 tgrid=igrid
3074
3075 end if
3076 end if
3077
3078
3079 ! test if we have U and V
3080
3081#ifdef DEBUG
3082 call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: test U and V"//&
3083 to_char(xmin_t)//to_char(xmax_t)//to_char(ymin_t)//to_char(ymax_t)//&
3084 to_char(xmin)//to_char(xmax)//to_char(ymin)//to_char(ymax))
3085
3086 call l4f_category_log(this(igrid)%category,l4f_debug,"UV diff coordinate lon"//&
3087 to_char(abs(xmin_t - xmin)-step_lon_t/2.d0)//&
3088 to_char(abs(xmax_t - xmax)-step_lon_t/2.d0))
3089 call l4f_category_log(this(igrid)%category,l4f_debug,"UV diff coordinate lat"//&
3090 to_char(abs(ymin_t - ymin) -step_lat_t/2.d0)//&
3091 to_char(abs(ymax_t - ymax)-step_lat_t/2.d0))
3092#endif
3093
3094 if ( abs(ymin - (ymin_t+step_lat_t/2.d0)) < 2.d-3 .and. abs(ymax - (ymax_t+step_lat_t/2.d0)) < 2.d-3 ) then
3095 if ( abs(xmin_t - (xmin+step_lon_t/2.d0)) < 2.d-3 .and. abs(xmax_t - (xmax+step_lon_t/2.d0)) < 2.d-3 ) then
3096
3097#ifdef DEBUG
3098 call l4f_category_log(this(igrid)%category,l4f_debug,"C grid: found U and V case up and right")
3099#endif
3100
3101 vgrid=jgrid
3102 ugrid=igrid
3103
3104 call init(griddim_t,xmin=xmin, xmax=xmax, ymin=ymin_t, ymax=ymax_t)
3105
3106 end if
3107 end if
3108 end if
3109
3110 ! abbiamo almeno due volumi: riportiamo U e/o V su T
3111 if (c_e(ugrid)) then
3112#ifdef DEBUG
3113 call l4f_category_log(this(igrid)%category,l4f_debug,"C grid U points: interpolate U on T "//to_char(tgrid)//to_char(ugrid))
3114#endif
3115 if (c_e(tgrid))then
3116 call vg6d_c2a_grid(this(ugrid),this(tgrid)%griddim,cgrid=1)
3117 else
3118 call vg6d_c2a_grid(this(ugrid),griddim_t,cgrid=1)
3119 end if
3120 call vg6d_c2a_mat(this(ugrid),cgrid=1)
3121 end if
3122
3123 if (c_e(vgrid)) then
3124#ifdef DEBUG
3125 call l4f_category_log(this(igrid)%category,l4f_debug,"C grid V points: interpolate V on T "//to_char(tgrid)//to_char(vgrid))
3126#endif
3127 if (c_e(tgrid))then
3128 call vg6d_c2a_grid(this(vgrid),this(tgrid)%griddim,cgrid=2)
3129 else
3130 call vg6d_c2a_grid(this(vgrid),griddim_t,cgrid=2)
3131 end if
3132 call vg6d_c2a_mat(this(vgrid),cgrid=2)
3133 end if
3134
3135 end do
3136
3137 call delete(griddim_t)
3138
3139end do
3140
3141
3142end subroutine vg6d_c2a
3143
3144
3145! Convert C grid to A grid
3146subroutine vg6d_c2a_grid(this,griddim_t,cgrid)
3147
3148type(volgrid6d),intent(inout) :: this ! object containing fields to be translated (U or V points)
3149type(griddim_def),intent(in),optional :: griddim_t ! object containing grid of T points
3150integer,intent(in) :: cgrid ! in C grid (Arakawa) we have 0=T,1=U,2=V points
3151
3152doubleprecision :: xmin, xmax, ymin, ymax
3153doubleprecision :: step_lon,step_lat
3154
3155
3156if (present(griddim_t)) then
3157
3158 call get_val(griddim_t,xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax)
3159 call set_val(this%griddim,xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax)
3160! improve grid definition precision
3161 CALL griddim_setsteps(this%griddim)
3162
3163else
3164
3165 select case (cgrid)
3166
3167 case(0)
3168
3169#ifdef DEBUG
3170 call l4f_category_log(this%category,l4f_debug,"C grid: T points, nothing to do")
3171#endif
3172 return
3173
3174 case (1)
3175
3176#ifdef DEBUG
3177 call l4f_category_log(this%category,l4f_debug,"C grid: U points, we need interpolation")
3178#endif
3179
3180 call get_val(this%griddim, xmin=xmin, xmax=xmax)
3181 step_lon=(xmax-xmin)/dble(this%griddim%dim%nx-1)
3182 xmin=xmin-step_lon/2.d0
3183 xmax=xmax-step_lon/2.d0
3184 call set_val(this%griddim, xmin=xmin, xmax=xmax)
3185! improve grid definition precision
3186 CALL griddim_setsteps(this%griddim)
3187
3188 case (2)
3189
3190#ifdef DEBUG
3191 call l4f_category_log(this%category,l4f_debug,"C grid: V points, we need interpolation")
3192#endif
3193
3194 call get_val(this%griddim, ymin=ymin, ymax=ymax)
3195 step_lat=(ymax-ymin)/dble(this%griddim%dim%ny-1)
3196 ymin=ymin-step_lat/2.d0
3197 ymax=ymax-step_lat/2.d0
3198 call set_val(this%griddim, ymin=ymin, ymax=ymax)
3199! improve grid definition precision
3200 CALL griddim_setsteps(this%griddim)
3201
3202 case default
3203
3204 call l4f_category_log(this%category,l4f_fatal,"C grid type not known")
3205 call raise_fatal_error ()
3206
3207 end select
3208
3209end if
3210
3211
3212call griddim_unproj(this%griddim)
3213
3214
3215end subroutine vg6d_c2a_grid
3216
3217! Convert C grid to A grid
3218subroutine vg6d_c2a_mat(this,cgrid)
3219
3220type(volgrid6d),intent(inout) :: this ! object containing fields to be translated
3221integer,intent(in) :: cgrid ! in C grid (Arakawa) we have 0=T,1=U,2=V points
3222
3223INTEGER :: i, j, k, iv, stallo
3224REAL,ALLOCATABLE :: tmp_arr(:,:)
3225REAL,POINTER :: voldatiuv(:,:)
3226
3227
3228IF (cgrid == 0) RETURN ! nothing to do
3229IF (cgrid /= 1 .AND. cgrid /= 2) THEN ! wrong cgrid
3230 CALL l4f_category_log(this%category,l4f_fatal,"C grid type "// &
3231 trim(to_char(cgrid))//" not known")
3232 CALL raise_error()
3233 RETURN
3234ENDIF
3235
3236! Temporary workspace
3237ALLOCATE(tmp_arr(this%griddim%dim%nx, this%griddim%dim%ny),stat=stallo)
3238if (stallo /=0)then
3239 call l4f_log(l4f_fatal,"allocating memory")
3240 call raise_fatal_error()
3241end if
3242
3243! allocate once for speed
3244IF (.NOT.ASSOCIATED(this%voldati)) THEN
3245 ALLOCATE(voldatiuv(this%griddim%dim%nx, this%griddim%dim%ny), stat=stallo)
3246 IF (stallo /= 0) THEN
3247 CALL l4f_log(l4f_fatal,"allocating memory")
3248 CALL raise_fatal_error()
3249 ENDIF
3250ENDIF
3251
3252IF (cgrid == 1) THEN ! U points to H points
3253 DO iv = 1, SIZE(this%var)
3254 DO k = 1, SIZE(this%timerange)
3255 DO j = 1, SIZE(this%time)
3256 DO i = 1, SIZE(this%level)
3257 tmp_arr(:,:) = rmiss
3258 CALL volgrid_get_vol_2d(this, i, j, k, iv, voldatiuv)
3259
3260! West boundary extrapolation
3261 WHERE(voldatiuv(1,:) /= rmiss .AND. voldatiuv(2,:) /= rmiss)
3262 tmp_arr(1,:) = voldatiuv(1,:) - (voldatiuv(2,:) - voldatiuv(1,:)) / 2.
3263 END WHERE
3264
3265! Rest of the field
3266 WHERE(voldatiuv(1:this%griddim%dim%nx-1,:) /= rmiss .AND. &
3267 voldatiuv(2:this%griddim%dim%nx,:) /= rmiss)
3268 tmp_arr(2:this%griddim%dim%nx,:) = &
3269 (voldatiuv(1:this%griddim%dim%nx-1,:) + &
3270 voldatiuv(2:this%griddim%dim%nx,:)) / 2.
3271 END WHERE
3272
3273 voldatiuv(:,:) = tmp_arr
3274 CALL volgrid_set_vol_2d(this, i, j, k, iv, voldatiuv)
3275 ENDDO
3276 ENDDO
3277 ENDDO
3278 ENDDO
3279
3280ELSE IF (cgrid == 2) THEN ! V points to H points
3281 DO iv = 1, SIZE(this%var)
3282 DO k = 1, SIZE(this%timerange)
3283 DO j = 1, SIZE(this%time)
3284 DO i = 1, SIZE(this%level)
3285 tmp_arr(:,:) = rmiss
3286 CALL volgrid_get_vol_2d(this, i, j, k, iv, voldatiuv)
3287
3288! South boundary extrapolation
3289 WHERE(voldatiuv(:,1) /= rmiss .AND. voldatiuv(:,2) /= rmiss)
3290 tmp_arr(:,1) = voldatiuv(:,1) - (voldatiuv(:,2) - voldatiuv(:,1)) / 2.
3291 END WHERE
3292
3293! Rest of the field
3294 WHERE(voldatiuv(:,1:this%griddim%dim%ny-1) /= rmiss .AND. &
3295 voldatiuv(:,2:this%griddim%dim%ny) /= rmiss)
3296 tmp_arr(:,2:this%griddim%dim%ny) = &
3297 (voldatiuv(:,1:this%griddim%dim%ny-1) + &
3298 voldatiuv(:,2:this%griddim%dim%ny)) / 2.
3299 END WHERE
3300
3301 voldatiuv(:,:) = tmp_arr
3302 CALL volgrid_set_vol_2d(this, i, j, k, iv, voldatiuv)
3303 ENDDO
3304 ENDDO
3305 ENDDO
3306 ENDDO
3307ENDIF ! tertium non datur
3308
3309IF (.NOT.ASSOCIATED(this%voldati)) THEN
3310 DEALLOCATE(voldatiuv)
3311ENDIF
3312DEALLOCATE (tmp_arr)
3313
3314end subroutine vg6d_c2a_mat
3315
3316
3317!> \brief Display object on screen
3318!!
3319!! Show brief content on screen.
3320subroutine display_volgrid6d (this)
3321
3322TYPE(volgrid6d),intent(in) :: this !< object to display
3323integer :: i
3324
3325#ifdef DEBUG
3326call l4f_category_log(this%category,l4f_debug,"ora mostro gridinfo " )
3327#endif
3328
3329print*,"----------------------- volgrid6d display ---------------------"
3330call display(this%griddim)
3331
3332IF (ASSOCIATED(this%time))then
3333 print*,"---- time vector ----"
3334 print*,"elements=",size(this%time)
3335 do i=1, size(this%time)
3336 call display(this%time(i))
3337 end do
3338end IF
3339
3340IF (ASSOCIATED(this%timerange))then
3341 print*,"---- timerange vector ----"
3342 print*,"elements=",size(this%timerange)
3343 do i=1, size(this%timerange)
3344 call display(this%timerange(i))
3345 end do
3346end IF
3347
3348IF (ASSOCIATED(this%level))then
3349 print*,"---- level vector ----"
3350 print*,"elements=",size(this%level)
3351 do i=1, size(this%level)
3352 call display(this%level(i))
3353 end do
3354end IF
3355
3356IF (ASSOCIATED(this%var))then
3357 print*,"---- var vector ----"
3358 print*,"elements=",size(this%var)
3359 do i=1, size(this%var)
3360 call display(this%var(i))
3361 end do
3362end IF
3363
3364IF (ASSOCIATED(this%gaid))then
3365 print*,"---- gaid vector (present mask only) ----"
3366 print*,"elements=",shape(this%gaid)
3367 print*,c_e(reshape(this%gaid,(/SIZE(this%gaid)/)))
3368end IF
3369
3370print*,"--------------------------------------------------------------"
3371
3372
3373end subroutine display_volgrid6d
3374
3375
3376!> \brief Display vector of object on screen
3377!!
3378!! Show brief content on screen.
3379subroutine display_volgrid6dv (this)
3380
3381TYPE(volgrid6d),intent(in) :: this(:) !< vector of object to display
3382integer :: i
3383
3384print*,"----------------------- volgrid6d vector ---------------------"
3385
3386print*,"elements=",size(this)
3387
3388do i=1, size(this)
3389
3390#ifdef DEBUG
3391 call l4f_category_log(this(i)%category,l4f_debug,"ora mostro il vettore volgrid6d" )
3392#endif
3393
3394 call display(this(i))
3395
3396end do
3397print*,"--------------------------------------------------------------"
3398
3399end subroutine display_volgrid6dv
3400
3401
3402!> Reduce some dimensions (level and timerage) for semplification (rounding).
3403!! Vector version of vg6dv_rounding
3404subroutine vg6dv_rounding(vg6din,vg6dout,level,timerange,nostatproc,merge)
3405type(volgrid6d),intent(in) :: vg6din(:) !< input volume
3406type(volgrid6d),intent(out),pointer :: vg6dout(:) !> output volume
3407type(vol7d_level),intent(in),optional :: level(:) !< almost equal level list
3408type(vol7d_timerange),intent(in),optional :: timerange(:) !< almost equal timerange list
3409logical,intent(in),optional :: merge !< if there are data on more then one almost equal levels or timeranges will be merged POINT BY POINT with priority for the fird data found ordered by icreasing var index
3410logical,intent(in),optional :: nostatproc !< do not take in account statistical processing code in timerange and P2
3411
3412integer :: i
3413
3414allocate(vg6dout(size(vg6din)))
3415
3416do i = 1, size(vg6din)
3417 call vg6d_rounding(vg6din(i),vg6dout(i),level,timerange,nostatproc,merge)
3418end do
3419
3420end subroutine vg6dv_rounding
3421
3422!> Reduce some dimensions (level and timerage) for semplification (rounding).
3423!! You can use this for simplify and use variables in computation like alchimia
3424!! where fields have to be on the same coordinate
3425!! examples:
3426!! means in time for short periods and istantaneous values
3427!! 2 meter and surface levels
3428!! If there are data on more then one almost equal levels or timeranges, the first var present (at least one point)
3429!! will be taken (order is by icreasing var index).
3430!! You can use predefined values for classic semplification
3431!! almost_equal_levels and almost_equal_timeranges
3432!! The level or timerange in output will be defined by the first element of level and timerange list
3433subroutine vg6d_rounding(vg6din,vg6dout,level,timerange,nostatproc,merge)
3434type(volgrid6d),intent(in) :: vg6din !< input volume
3435type(volgrid6d),intent(out) :: vg6dout !> output volume
3436type(vol7d_level),intent(in),optional :: level(:) !< almost equal level list
3437type(vol7d_timerange),intent(in),optional :: timerange(:) !< almost equal timerange list
3438logical,intent(in),optional :: merge !< if there are data on more then one almost equal levels or timeranges
3439!! will be merged POINT BY POINT with priority for the fird data found ordered by icreasing var index
3440logical,intent(in),optional :: nostatproc !< do not take in account statistical processing code in timerange and P2
3441
3442integer :: ilevel,itimerange
3443type(vol7d_level) :: roundlevel(size(vg6din%level))
3444type(vol7d_timerange) :: roundtimerange(size(vg6din%timerange))
3445
3446roundlevel=vg6din%level
3447
3448if (present(level))then
3449 do ilevel = 1, size(vg6din%level)
3450 if ((any(vg6din%level(ilevel) .almosteq. level))) then
3451 roundlevel(ilevel)=level(1)
3452 end if
3453 end do
3454end if
3455
3456roundtimerange=vg6din%timerange
3457
3458if (present(timerange))then
3459 do itimerange = 1, size(vg6din%timerange)
3460 if ((any(vg6din%timerange(itimerange) .almosteq. timerange))) then
3461 roundtimerange(itimerange)=timerange(1)
3462 end if
3463 end do
3464end if
3465
3466!set istantaneous values everywere
3467!preserve p1 for forecast time
3468if (optio_log(nostatproc)) then
3469 roundtimerange(:)%timerange=254
3470 roundtimerange(:)%p2=0
3471end if
3472
3473
3474call vg6d_reduce(vg6din,vg6dout,roundlevel,roundtimerange,merge)
3475
3476end subroutine vg6d_rounding
3477
3478!> Reduce some dimensions (level and timerage).
3479!! You can pass a volume and specify duplicated levels and timeranges in roundlevel and roundtimerange;
3480!! you get unique levels and timeranges in output.
3481!! If there are data on equal levels or timeranges, the first var present (at least one point)
3482!! will be taken (order is by icreasing var index).
3483!! you can specify merge and if there are data on equal levels or timeranges will be merged POINT BY POINT
3484!! with priority for the first data found ordered by icreasing var index (require to decode all the data)
3485!! Data are decoded only if needed so the output should be with or without voldata allocated
3486subroutine vg6d_reduce(vg6din,vg6dout,roundlevel,roundtimerange,merge)
3487type(volgrid6d),intent(in) :: vg6din !< input volume
3488type(volgrid6d),intent(out) :: vg6dout !< output volume
3489type(vol7d_level),intent(in) :: roundlevel(:) !< new level list to use for rounding
3490type(vol7d_timerange),intent(in) :: roundtimerange(:) !< new timerange list to use for rounding
3491logical,intent(in),optional :: merge !< if there are data on equal levels or timeranges will be merged POINT BY POINT with priority for the first data found ordered by icreasing var index (require to decode all the data)
3492
3493integer :: nlevel,ntime,ntimerange,nvar,ilevel,itimerange,ivar,indl,indt,itime,nx,ny
3494real,allocatable :: vol2d(:,:)
3495
3496nx=vg6din%griddim%dim%nx
3497ny=vg6din%griddim%dim%ny
3498nlevel=count_distinct(roundlevel,back=.true.)
3499ntime=size(vg6din%time)
3500ntimerange=count_distinct(roundtimerange,back=.true.)
3501nvar=size(vg6din%var)
3502
3503call init(vg6dout, vg6din%griddim, vg6din%time_definition, categoryappend="generated by vg6d_reduce")
3504call volgrid6d_alloc(vg6dout, vg6din%griddim%dim, ntime, nlevel, ntimerange, nvar)
3505
3506if ( ASSOCIATED(vg6din%voldati) .or. optio_log(merge)) then
3507 call volgrid6d_alloc_vol(vg6dout,inivol=.true.,decode=.true.)
3508 allocate(vol2d(nx,ny))
3509else
3510 call volgrid6d_alloc_vol(vg6dout,inivol=.true.,decode=.false.)
3511end if
3512
3513vg6dout%time=vg6din%time
3514vg6dout%var=vg6din%var
3515vg6dout%timerange=pack_distinct(roundtimerange,ntimerange,back=.true.)
3516vg6dout%level=pack_distinct(roundlevel,nlevel,back=.true.)
3517! sort modified dimensions
3518CALL sort(vg6dout%timerange)
3519CALL sort(vg6dout%level)
3520
3521do ilevel=1,size(vg6din%level)
3522 indl=index(vg6dout%level,roundlevel(ilevel))
3523 do itimerange=1,size(vg6din%timerange)
3524 indt=index(vg6dout%timerange,roundtimerange(itimerange))
3525 do ivar=1, nvar
3526 do itime=1,ntime
3527
3528 if ( ASSOCIATED(vg6din%voldati)) then
3529 vol2d=vg6din%voldati(:,:,ilevel,itime,itimerange,ivar)
3530 end if
3531
3532 if (optio_log(merge)) then
3533
3534 if ( .not. ASSOCIATED(vg6din%voldati)) then
3535 CALL grid_id_decode_data(vg6din%gaid(ilevel,itime,itimerange,ivar), vol2d)
3536 end if
3537
3538 !! merge present data point by point
3539 where (.not. c_e(vg6dout%voldati(:,:,indl,itime,indt,ivar)))
3540
3541 vg6dout%voldati(:,:,indl,itime,indt,ivar)=vol2d
3542
3543 end where
3544 else if ( ASSOCIATED(vg6din%voldati)) then
3545 if (.not. any(c_e(vg6dout%voldati(:,:,indl,itime,indt,ivar))))then
3546 vg6dout%voldati(:,:,indl,itime,indt,ivar)=vol2d
3547 end if
3548 end if
3549
3550 if (c_e(vg6din%gaid(ilevel,itime,itimerange,ivar)).and. .not. c_e(vg6dout%gaid(indl,itime,indt,ivar)))then
3551 call copy (vg6din%gaid(ilevel,itime,itimerange,ivar), vg6dout%gaid(indl,itime,indt,ivar))
3552 end if
3553 end do
3554 end do
3555 end do
3556end do
3557
3558if ( ASSOCIATED(vg6din%voldati) .or. optio_log(merge)) then
3559 deallocate(vol2d)
3560end if
3561
3562end subroutine vg6d_reduce
3563
3564
3565end module volgrid6d_class
3566
3567
3568
3569!>\example example_vg6d_3.f90
3570!!\brief Programma esempio semplice per gridinfo e volgrid6d.
3572!! Programma che importa da file un vettore di gridinfo poi lo importa in volgrid6d. Da volgrid6d viene di nuovo creato un vettore di gridinfo per poi exportare su file.
3573
3574!>\example example_vg6d_5.f90
3575!!\brief Programma trasformazione da volgrid6d a volgrid6d
3576!!
3577!! Legge grib da un file e li organizza in un vettore di strutture
3578!! volgrid6d mettendoli a disposizione per eventuali elaborazioni;
3579!! vengono poi riesportati a un file grib
3580
3581!>\example example_vg6d_8.f90
3582!! \brief Programma scrittura su file vettore di anagrafica
3583
3584!>\example example_vg6d_6.f90
3585!! \brief Programma trasformazione da volgrid6d a vol7d
3586
3587!>\example example_vg6d_7.f90
3588!! \brief Programma trasformazione da vol7d a volgrid6d
3589
3590!>\example example_vg6d_9.f90
3591!! \brief Example to create a grib editionNumber = 2 file from data generated in memory using a grib_api template
Method for inserting elements of the array at a desired position.
Operatore di valore assoluto di un intervallo.
Legge un oggetto datetime/timedelta o un vettore di oggetti datetime/timedelta da un file FORMATTED o...
Functions that return a trimmed CHARACTER representation of the input variable.
Restituiscono il valore dell'oggetto in forma di stringa stampabile.
Scrive un oggetto datetime/timedelta o un vettore di oggetti datetime/timedelta su un file FORMATTED ...
Copy an object, creating a fully new instance.
Method for returning the contents of the object.
Compute forward coordinate transformation from geographical system to projected system.
Method for setting the contents of the object.
Clone the object, creating a new independent instance of the object exactly equal to the starting one...
Decode and return the data array from a grid_id object associated to a gridinfo object.
Encode a data array into a grid_id object associated to a gridinfo object.
Index method.
Emit log message for a category with specific priority.
Convert a level type to a physical variable.
Destructor, it releases every information and memory buffer associated with the object.
Display on standard output a description of the volgrid6d object provided.
Export an object dirctly to a native file, to a gridinfo object or to a supported file format through...
Import an object dirctly from a native file, from a gridinfo object or from a supported file format t...
Constructor, it creates a new instance of the object.
Reduce some dimensions (level and timerage) for semplification (rounding).
Transform between any combination of volgrid6d and vol7d objects by means of a transform_def object d...
Apply the conversion function this to values.
This module defines usefull general purpose function and subroutine.
Classi per la gestione delle coordinate temporali.
Module for describing geographically referenced regular grids.
Module for defining the extension and coordinates of a rectangular georeferenced grid.
This module defines an abstract interface to different drivers for access to files containing gridded...
Module for defining transformations between rectangular georeferenced grids and between grids and spa...
Class for managing information about a single gridded georeferenced field, typically imported from an...
classe per la gestione del logging
Module for quickly interpreting the OPTIONAL parameters passed to a subprogram.
Classe per la gestione di un volume completo di dati osservati.
classe per import ed export di volumi da e in DB-All.e
Classe per la gestione dei livelli verticali in osservazioni meteo e affini.
Classe per la gestione degli intervalli temporali di osservazioni meteo e affini.
This module defines objects and methods for managing data volumes on rectangular georeferenced grids.
Class for managing physical variables in a grib 1/2 fashion.
Class for expressing an absolute time value.
This object completely describes a grid on a geographic projection.
Derived type describing the extension of a grid and the geographical coordinates of each point.
Derived type associated to a block/message/record/band of gridded data coming from a file-like object...
This object fully defines a transformation between a couple of particular griddim_def or vol7d object...
This object defines in an abstract way the type of transformation to be applied.
Derived type defining a dynamically extensible array of TYPE(gridinfo_def) elements.
Object describing a single gridded message/band.
Definisce un oggetto contenente i volumi anagrafica e dati e tutti i descrittori delle loro dimension...
Definisce il livello verticale di un'osservazione.
Definisce l'intervallo temporale di un'osservazione meteo.
Object describing a rectangular, homogeneous gridded dataset.
Class defining a real conversion function between units.
Definition of a physical variable in grib coding style.

Generated with Doxygen.