19 module volgrid6d_vapor_class
38 MODULE PROCEDURE volgrid6d_export_to_vapor
49 subroutine volgrid6d_export_to_vapor (this,normalize,rzscan,filename,filename_auto,reusevdf)
51 TYPE(volgrid6d),
INTENT(INOUT) :: this
52 logical,
intent(in) :: normalize
53 logical,
intent(in),
optional :: rzscan
54 character(len=*),
intent(in),
optional :: filename
55 character(len=*),
intent(out),
optional :: filename_auto
56 logical,
intent(in),
optional :: reusevdf
58 character(len=254) :: lfilename
59 integer :: ntime, ntimerange, ntimera, nlevel, nvar
61 integer :: ier, xyzdim(3),ivar,i,j
62 character(len=255),
allocatable :: varnames(:),vardescriptions(:),tsdescriptions(:)
63 doubleprecision :: extents(6),vdfmiss=rmiss
65 TYPE(conv_func),
pointer :: c_func(:)
66 TYPE(vol7d_var),
allocatable :: varbufr(:)
67 type(vol7d_var),
pointer :: dballevar(:)
68 CHARACTER(len=255) :: proj_type,mapprojection
70 integer :: zone, irzscan, indele
71 DOUBLE PRECISION :: xoff, yoff, ellips_smaj_axis, ellips_flatt
72 DOUBLE PRECISION :: longitude_south_pole, latitude_south_pole, angle_rotation
76 call vol7d_dballe_import_dballevar(dballevar)
78 if (optio_log(rzscan))
then 91 if (
present(filename))
then 92 if (filename /=
"")
then 97 if (
present(filename_auto))filename_auto=lfilename
100 inquire(file=lfilename,exist=exist)
102 if (optio_log(reusevdf))
then 103 if (.not. exist)
then 104 call l4f_category_log(this%category,l4f_error,
"file do not exist; cannot reuse file: "//trim(lfilename))
108 call l4f_category_log(this%category,l4f_error,
"file exist; cannot open new file: "//trim(lfilename))
112 call l4f_category_log(this%category,l4f_info,
"writing on file: "//trim(lfilename))
114 if (
associated(this%time)) ntime=
size(this%time)
115 if (
associated(this%timerange)) ntimerange=
size(this%timerange)
116 if (
associated(this%level)) nlevel=
size(this%level)
117 if (
associated(this%var)) nvar=
size(this%var)
119 if (
c_e(ntime) .and.
c_e(ntimerange) .and.
c_e(nlevel) .and.
c_e(nvar))
then 121 allocate(varnames(nvar),vardescriptions(nvar),varbufr(nvar))
123 call get_val (this%griddim, nx=xyzdim(1) , ny=xyzdim(2) )
129 if (
associated(this%voldati))
then 131 if (optio_log(normalize))
then 132 CALL vargrib2varbufr(this%var, varbufr, c_func)
135 IF (
ASSOCIATED(c_func))
THEN 137 call l4f_category_log(this%category,l4f_info,
"normalize is activated, so the volume data are changed in output")
139 this%voldati(:,:,:,:,:,ivar) =
convert(c_func(ivar),this%voldati(:,:,:,:,:,ivar))
148 j=firsttrue(varbufr(ivar)%btable == dballevar(:)%btable)
151 varbufr(ivar)%description = dballevar(j)%description
152 varbufr(ivar)%unit = dballevar(j)%unit
153 varbufr(ivar)%scalefactor = dballevar(j)%scalefactor
155 varnames(ivar) = varbufr(ivar)%btable
156 vardescriptions(ivar) = trim(varbufr(ivar)%description)//
"_"//trim(varbufr(ivar)%unit)
159 if (varnames(ivar) ==
"B10007")
then 165 varnames(ivar) =
"Vnotnormalized_"//
t2c(ivar)
166 vardescriptions(ivar) =
"None" 174 varnames(ivar) =
"V"//trim(
to_char(this%var(ivar)%number))
179 if (this%time_definition == 1)
then 181 allocate(tsdescriptions(ntimera))
184 tsdescriptions(i)=
to_char(this%time(i))
188 allocate(tsdescriptions(ntimera))
191 tsdescriptions(i)=
to_char(this%timerange(i))
196 call get_val (this%griddim, xmin=extents(1),ymin=extents(2), xmax=extents(4) , ymax=extents(5))
201 call get_val (this%griddim, proj_type=proj_type)
213 select case (proj_type)
216 call l4f_category_log(this%category,l4f_info,
"VDF: probably vapor do not support this projection ?: "//trim(proj_type))
217 mapprojection =
"+proj=latlon +ellps=sphere" 219 extents(1)=extents(1)*111177.d0
220 extents(2)=extents(2)*111177.d0
221 extents(3)=extents(3)*100000.d0
222 extents(4)=extents(4)*111177.d0
223 extents(5)=extents(5)*111177.d0
224 extents(6)=extents(6)*100000.d0
231 call l4f_category_log(this%category,l4f_info,
"VDF: vapor probably support this projection: "//trim(proj_type))
233 call get_val (this%griddim, longitude_south_pole=longitude_south_pole,&
234 latitude_south_pole=latitude_south_pole, angle_rotation=angle_rotation)
240 if (angle_rotation /= 0. )
then 241 call l4f_category_log(this%category,l4f_error,
"angle_rotation /= 0 not supported in vapor (proj)")
245 mapprojection =
"+proj=ob_tran +o_proj=latlong +o_lat_p="//
t2c(-latitude_south_pole)//&
246 "d +o_lon_p=0d +lon_0="//
t2c(longitude_south_pole)//
"d +ellps=sphere" 248 extents(1)=extents(1)*111177.d0
249 extents(2)=extents(2)*111177.d0
250 extents(3)=extents(3)*100000.d0
251 extents(4)=extents(4)*111177.d0
252 extents(5)=extents(5)*111177.d0
253 extents(6)=extents(6)*100000.d0
257 call l4f_category_log(this%category,l4f_warn,
"VDF: vapor do not support this projection: "//trim(proj_type))
258 call get_val (this%griddim, xmin=extents(1),ymin=extents(2), xmax=extents(4) , ymax=extents(5),&
259 zone=zone, xoff=xoff, yoff=yoff, ellips_smaj_axis=ellips_smaj_axis, ellips_flatt=ellips_flatt)
264 mapprojection =
"+proj=utm +zone="//
t2c(zone)
272 "VDF: proj or vdf (vapor) export do not support this projection: "//trim(proj_type))
273 mapprojection = cmiss
278 if(ier==0)
call l4f_category_log(this%category,l4f_info,
"VDF: projection parameter "//mapprojection)
281 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call create_metadata_from_file")
282 ier = vdf4f_create_metadata_from_file(lfilename)
284 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call create_metadata")
285 ier = vdf4f_create_metadata(xyzdim,vdctype=2)
287 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call set_missing_value")
288 if(ier==0) ier = vdf4f_set_missing_value(vdfmiss)
290 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call set_num_timesteps")
291 if(ier==0) ier = vdf4f_set_num_timesteps(ntimera)
293 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call vdf4f_set_comment")
294 if(ier==0) ier = vdf4f_set_comment(
"vogrid6d exported")
296 if(ier==0) ier = vdf4f_set_coord_system_type(coordsystemtype=
"cartesian")
298 if (
c_e(indele))
then 299 if(ier==0)
call l4f_category_log(this%category,l4f_info,
"VDF: ELEVATION (B10007) found: setting gridtype to layered")
300 extents(6) = maxval(this%voldati(:,:,:,:,:,indele),
c_e(this%voldati(:,:,:,:,:,indele)))
301 if(ier==0) ier = vdf4f_set_grid_type(gridtype=
"layered")
303 if(ier==0) ier = vdf4f_set_grid_type(gridtype=
"regular")
306 if(ier==0) ier = vdf4f_set_grid_extents(extents=extents)
307 if(ier==0) ier = vdf4f_set_map_projection(mapprojection=mapprojection)
315 if (
c_e(indele)) varnames(indele) =
"ELEVATION" 317 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call set_variables_names")
318 if(ier==0) ier = vdf4f_set_variables_names(nvar, varnames)
321 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call set_v_comment")
322 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call set_ts_comment")
325 if(ier==0) ier = vdf4f_set_ts_comment(i-1,tsdescriptions(i))
327 if(ier==0) ier = vdf4f_set_v_comment(i-1,varnames(j),vardescriptions(j))
334 varnames(ivar)=
"XY_"//
t2c(varnames(ivar))
337 if (
c_e(indele)) varnames(indele) =
"HGT" 339 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call vdf4f_set_variables_2d_xy")
340 if(ier==0) ier = vdf4f_set_variables_2d_xy(nvar, varnames)
344 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call write_metadata")
345 if(ier==0) ier = vdf4f_write_metadata(lfilename)
347 if(ier==0) ier = destroy_metadata_c()
352 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call vdf4f_create_writer")
353 if(ier==0) ier = vdf4f_create_writer(lfilename)
355 if (this%time_definition == 1)
then 357 if (ntimerange /= 1)
then 358 if(ier==0)
call l4f_category_log(this%category,l4f_warn,
"VDF: writing only first timerange, there are:"//
t2c(ntimerange))
361 if (.not.
c_e(indele))
call fill_underground_missing_values(this%voldati(:,:,:,:,1,:))
362 if(ier==0)
call l4f_category_log(this%category,l4f_info,
"scan VDF (vapor file) for times")
365 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call vdf4f_write")
366 if(ier==0) ier =
vdf4f_write(this%voldati(:,:,:,:,1,:), xyzdim, ntime, nvar, varnames, irzscan)
368 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call vdf4f_write_2d_xy")
369 if(ier==0) ier = vdf4f_write_2d_xy(this%voldati(:,:,1,:,1,:), xyzdim(:2), ntime, nvar ,varnames)
376 if(ier==0)
call l4f_category_log(this%category,l4f_warn,
"VDF: writing only fisth time, there are:"//
t2c(ntime))
379 if (.not.
c_e(indele))
call fill_underground_missing_values(this%voldati(:,:,:,1,:,:))
380 if(ier==0)
call l4f_category_log(this%category,l4f_info,
"scan VDF (vapor file) for timeranges")
383 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call vdf4f_write")
384 if(ier==0) ier =
vdf4f_write(this%voldati(:,:,:,1,:,:), xyzdim, ntimerange, nvar, varnames, irzscan)
386 if(ier==0)
call l4f_category_log(this%category,l4f_debug,
"VDF: call vdf4f_write_2d_xy")
387 if(ier==0) ier = vdf4f_write_2d_xy(this%voldati(:,:,1,1,:,:), xyzdim(:2), ntimerange, nvar, varnames)
393 if (ier /= 0)
call l4f_category_log(this%category,l4f_error,
"export to vdf: "//vdf4f_get_err_msg())
396 deallocate(varnames,vardescriptions,tsdescriptions,varbufr)
397 if (ier==0) ier = destroy_writer_c()
401 CALL raise_fatal_error(
"exporting to vdf")
406 call l4f_category_log(this%category,l4f_warn,
"volume with voldati not associated: not exported to vdf")
412 call l4f_category_log(this%category,l4f_warn,
"volume with some dimensions to 0: not exported to vdf")
429 subroutine fill_underground_missing_values(voldati)
430 real,
intent(inout) :: voldati(:,:,:,:,:)
432 integer :: x,y,z,tim,var,zz
434 do x=1,
size(voldati,1)
435 do y=1,
size(voldati,2)
436 do tim=1,
size(voldati,4)
437 do var=1,
size(voldati,5)
438 do z=1,
size(voldati,3)
440 if (.not.
c_e(voldati(x,y,z,tim,var)))
then 441 zz=firsttrue(
c_e(voldati(x,y,:,tim,var)))
443 voldati(x,y,z,tim,var)=voldati(x,y,firsttrue(
c_e(voldati(x,y,:,tim,var))),tim,var)
445 call l4f_log(l4f_warn,
"fill_underground_missing_values: there are only missing values in the full coloumn")
459 end subroutine fill_underground_missing_values
461 end subroutine volgrid6d_export_to_vapor
463 end module volgrid6d_vapor_class
Functions that return a trimmed CHARACTER representation of the input variable.
Represent level object in a pretty string.
Module for quickly interpreting the OPTIONAL parameters passed to a subprogram.
Apply the conversion function this to values.
Module for describing geographically referenced regular grids.
This module defines objects and methods for managing data volumes on rectangular georeferenced grids...
classe per import ed export di volumi da e in DB-All.e
Classe per la gestione di un volume completo di dati osservati.
Definitions of constants and functions for working with missing values.
interface to different architectures (cast some type)
classe per la gestione del logging
Class for managing physical variables in a grib 1/2 fashion.
Utilities for CHARACTER variables.
Emit log message for a category with specific priority.