FLEXPART v9.3.0
 All Classes Files Functions Variables
fp2nc4io_mod.F90
Go to the documentation of this file.
2 
3  !****************************************************************
4  ! *
5  ! Contains data and routines for dumping selected FLEXPART *
6  ! array variables to a NetCDF4 file. *
7  ! *
8  ! Don Morton (Boreal Scientific Computing LLC) *
9  ! *
10  ! May 2016 *
11  ! *
12  !****************************************************************
13 
14  USE par_mod
15  USE com_mod
16 
17  USE netcdf
18 
19  IMPLICIT NONE
20 
21  ! This variable should be in the range [1,9]. It has been suggested
22  ! that 2 offers reasonable compression in a reasonable time.
23  ! Higher values will offer more compression, but will take more time
24  INTEGER, PARAMETER :: DEFAULT_DEFLATE_LEVEL = 2
25  PRIVATE default_deflate_level
26 
27  ! These are valid variable names for the user of this module to reference
28  CHARACTER, DIMENSION(*), PARAMETER :: VALID_VARS = &
29 & (/ 't', 'u', 'v', 'w', 'q', &
30 & 'T', 'U', 'V', 'W', 'Q' /)
31  PRIVATE VALID_VARS
32 
33  ! Private routines in this module
34  PRIVATE private_dump_3dfield
35  PRIVATE private_read_3dfield
36  PRIVATE to_upper
37 
38 
39 CONTAINS
40 
42 
43  ! Prints the list of met variables that are considered valid in this
44  ! module
45 
46  IMPLICIT NONE
47  INTEGER :: i
48 
49  DO i=1,SIZE(valid_vars)
50  print *, valid_vars(i)
51  ENDDO
52 
53  END SUBROUTINE fp2nc4io_print_valid_vars
54 
55 
56 
57  LOGICAL FUNCTION fp2nc4io_vars_are_valid(num_vars, dump_vars)
58 
59  ! Returns True or False depending on whether all of the variables
60  ! in dump_vars are valid names (according to VALID_VARS)
61 
62  IMPLICIT NONE
63 
64  INTEGER, INTENT(IN) :: num_vars
65  CHARACTER, DIMENSION(num_vars), INTENT(IN) :: dump_vars ! var list
66 
67  LOGICAL :: all_good = .true.
68  INTEGER :: i
69 
70  DO i=1,num_vars
71  IF( .NOT. any(valid_vars == dump_vars(i)) ) THEN
72  all_good = .false.
73  ENDIF
74  ENDDO
75 
76  fp2nc4io_vars_are_valid = all_good
77 
78  END FUNCTION fp2nc4io_vars_are_valid
79 
80 
81 
82  SUBROUTINE fp2nc4io_dump(nc4_filepath, num_vars, dump_vars, deflate_level)
83 
84  ! Writes metadata plus variables in dump_vars to NetCDF4 file
85  ! All of the dumped variables come from FLEXPART modules
86  ! par_mod and com_mod
87 
88  IMPLICIT NONE
89 
90  CHARACTER(LEN=*), INTENT(IN) :: nc4_filepath ! Full path to dump file
91  INTEGER, INTENT(IN) :: num_vars ! Num variables in dump_vars
92  CHARACTER, DIMENSION(num_vars), INTENT(IN) :: dump_vars ! var list
93  INTEGER, OPTIONAL, INTENT(IN) :: deflate_level ! (should be 0-9)
94 
95 
96  INTEGER :: i, j, k
97  INTEGER :: ncfunc_retval ! NetCDF function call return values
98  INTEGER :: ncid ! NetCDF file id
99 
100 
101  ! Variables used by NetCDF routines
102  INTEGER :: x_dimid, y_dimid, z_dimid, dimids(3)
103  INTEGER :: varid
104  INTEGER :: deflevel ! Deflate level
105 
106 #ifdef TESTING
107  INTEGER :: nx_test, ny_test, nz_test
108  REAL, ALLOCATABLE, DIMENSION(:,:,:) :: testvar_array
109  CHARACTER(LEN=NF90_MAX_NAME) :: x_dimname, y_dimname, z_dimname
110  REAL :: orig_array_sum, test_array_sum
111 #endif
112 
113 #ifdef TESTING
114  print *,
115  print *, '*** Running in testing mode ***'
116  print *,
117 #endif
118 
119  ! Use default deflate level if it wasn't passed in, or if a bad
120  ! value was passed in.
121  IF (present(deflate_level)) THEN
122  IF (deflate_level < 0 .OR. deflate_level > 9) THEN
123  deflevel = default_deflate_level
124  ELSE
125  deflevel = deflate_level
126  ENDIF
127  ELSE
128  deflevel = default_deflate_level
129  ENDIF
130 
131  print *, 'Using deflate level: ', deflevel
132 
133  !!!!!!---------------------------------------------------
134  !!!!!! Now we are ready to dump to NetCDF4 file
135  !!!!!!---------------------------------------------------
136 
137  ncfunc_retval = nf90_create(nc4_filepath, &
138 & or(nf90_clobber, nf90_hdf5), ncid)
139  print *, 'Created file: ', trim(nc4_filepath)
140 
141  ! Define the dimensions, and get dimension ids passed back
142  ! The values nx, ny and nz come from FP par_mod
143  ncfunc_retval = nf90_def_dim(ncid, 'x', nx, x_dimid)
144  ncfunc_retval = nf90_def_dim(ncid, 'y', ny, y_dimid)
145  ncfunc_retval = nf90_def_dim(ncid, 'z', nz, z_dimid)
146  dimids = (/ x_dimid, y_dimid, z_dimid /)
147 
148  ! Write each of the 3d variables to the NetCDF file
149  DO i=1,num_vars
150  CALL private_dump_3dfield(ncid, dump_vars(i), dimids, deflevel)
151  print *, 'Dumped 3d field: ', dump_vars(i)
152  ENDDO
153 
154  ! Write the height field - variable 'height' is defined in com_mod
155  ncfunc_retval = nf90_def_var(ncid, 'height', nf90_double, &
156 & z_dimid, varid)
157 
158  ncfunc_retval = nf90_def_var_deflate(ncid, varid, &
159 & shuffle=0, &
160 & deflate=1, &
161 & deflate_level=deflevel)
162 
163  ncfunc_retval = nf90_put_var(ncid, varid, height(1:nz))
164 
165  ! Write some of the scalar metadata variables
166  ! dx, dy, xlon0, xlat0 are all defined in com_mod
167  ncfunc_retval = nf90_def_var(ncid, 'dx', nf90_double, varid)
168  ncfunc_retval = nf90_put_var(ncid, varid, dx)
169 
170  ncfunc_retval = nf90_def_var(ncid, 'dy', nf90_double, varid)
171  ncfunc_retval = nf90_put_var(ncid, varid, dy)
172 
173  ncfunc_retval = nf90_def_var(ncid, 'xlon0', nf90_double, varid)
174  ncfunc_retval = nf90_put_var(ncid, varid, xlon0)
175 
176  ncfunc_retval = nf90_def_var(ncid, 'ylat0', nf90_double, varid)
177  ncfunc_retval = nf90_put_var(ncid, varid, ylat0)
178 
179  ! All done, close the NetCDF file
180  ncfunc_retval = nf90_close(ncid)
181 
182 #ifdef TESTING
183  !!!!!!!!!!!!!! Reading !!!!!!!!!!!!!!!!!!!!
184  print *, "Opening nc file for reading"
185  ncfunc_retval = nf90_open(nc4_filepath, nf90_nowrite, ncid)
186 
187  ! Get dimensions
188  ncfunc_retval = nf90_inq_dimid(ncid, 'x', x_dimid)
189  ncfunc_retval = nf90_inquire_dimension(ncid, x_dimid, x_dimname, &
190 & nx_test)
191  print *, 'nx_test: ', nx_test
192 
193  ncfunc_retval = nf90_inq_dimid(ncid, 'y', y_dimid)
194  ncfunc_retval = nf90_inquire_dimension(ncid, y_dimid, y_dimname, &
195 & ny_test)
196  print *, 'ny_test: ', ny_test
197 
198  ncfunc_retval = nf90_inq_dimid(ncid, 'z', z_dimid)
199  ncfunc_retval = nf90_inquire_dimension(ncid, z_dimid, z_dimname, &
200 & nz_test)
201  print *, 'nz_test: ', nz_test
202 
203  ALLOCATE( testvar_array(0:nx_test-1, 0:ny_test-1, nz_test) )
204 
205  ! Read each variable and compare with original data
206  DO i=1,num_vars
207  CALL private_read_3dfield(ncid, dump_vars(i), &
208 & nx_test, ny_test, nz_test, &
209 & testvar_array)
210 
211 
212  IF (to_upper(dump_vars(i)) == 'U') THEN
213  orig_array_sum = sum( uu(0:nx_test-1, 0:ny_test-1, &
214 & 1:nz_test, 1) )
215  ELSEIF (to_upper(dump_vars(i)) == 'V') THEN
216  orig_array_sum = sum( vv(0:nx_test-1, 0:ny_test-1, &
217 & 1:nz_test, 1) )
218  ELSEIF (to_upper(dump_vars(i)) == 'T') THEN
219  orig_array_sum = sum( tt(0:nx_test-1, 0:ny_test-1, &
220 & 1:nz_test, 1) )
221  ELSEIF (to_upper(dump_vars(i)) == 'W') THEN
222  orig_array_sum = sum( ww(0:nx_test-1, 0:ny_test-1, &
223 & 1:nz_test, 1) )
224  ELSEIF (to_upper(dump_vars(i)) == 'Q') THEN
225  orig_array_sum = sum( qv(0:nx_test-1, 0:ny_test-1, &
226 & 1:nz_test, 1) )
227  ENDIF
228 
229  test_array_sum = sum( testvar_array(0:nx_test-1, 0:ny_test-1, &
230 & 1:nz_test) )
231 
232  print *, dump_vars(i), ': ', 'SUM of differences = ', &
233 & test_array_sum - orig_array_sum
234  IF ( abs(test_array_sum - orig_array_sum) .GT. 1.0e-3 ) THEN
235  print *, &
236 & 'WARNING WILL ROBINSON!: Sum of differences exceeds 1.0E-3'
237  ENDIF
238  ENDDO
239 
240  ncfunc_retval = nf90_close(ncid)
241  print *, 'Closed file: ', ncfunc_retval
242 
243 #endif
244 
245  END SUBROUTINE fp2nc4io_dump
246 
247 
248  SUBROUTINE private_dump_3dfield(ncid, varname, dimids, deflevel)
249 
250  ! Private routine meant to provide low level access for writing
251  ! specified varname to NetCDF4 file. It is assumed that the
252  ! NC4 file has already been opened, and that dimension id's have
253  ! already been obtained
254 
255  IMPLICIT NONE
256 
257  INTEGER, INTENT(IN) :: ncid ! NC4 file id
258  CHARACTER, INTENT(IN) :: varname
259  INTEGER, INTENT(IN) :: dimids(3) ! NC4 dimension ids
260  INTEGER, INTENT(IN) :: deflevel ! compression level
261 
262  ! NetCDF4 variables
263  CHARACTER :: nc_varname
264  INTEGER :: ncfunc_retval, varid
265 
266  ! Check that we have a valid varname. If not, buh-bye
267  IF( .NOT. any(valid_vars == varname) ) THEN
268  print *,
269  print *, 'fp2nc4io:private_dump_3d_field() bad var: ', varname
270  print *, ' ABORTING...'
271  print *,
272  stop
273  ENDIF
274 
275  ! Convert varname to upper case for use in NetCDF file
276  nc_varname = to_upper(varname)
277 
278  ! Create the variable in the NetCDF file
279  ncfunc_retval = nf90_def_var(ncid, nc_varname, nf90_double, &
280 & dimids, varid)
281 
282  ncfunc_retval = nf90_def_var_deflate(ncid, varid, &
283 & shuffle=0, &
284 & deflate=1, &
285 & deflate_level=deflevel)
286 
287  ! Write the data arrays
288  ! The values nx, ny and nz come from module com_mod
289  ! Likewise, the arrays uu, vv, tt, ww, qv are also from the
290  ! same module, and we assume they all have the same dimensions
291  ! (currently they do)
292  print *, 'Writing: ', nc_varname
293  IF (nc_varname == 'U') THEN
294  ncfunc_retval = nf90_put_var(ncid, varid, &
295 & uu(0:nx-1, 0:ny-1, 1:nz, 1))
296  ELSEIF (nc_varname == 'V') THEN
297  ncfunc_retval = nf90_put_var(ncid, varid, &
298 & vv(0:nx-1, 0:ny-1, 1:nz, 1))
299  ELSEIF (nc_varname == 'T') THEN
300  ncfunc_retval = nf90_put_var(ncid, varid, &
301 & tt(0:nx-1, 0:ny-1, 1:nz, 1))
302  ELSEIF (nc_varname == 'W') THEN
303  ncfunc_retval = nf90_put_var(ncid, varid, &
304 & ww(0:nx-1, 0:ny-1, 1:nz, 1))
305  ELSEIF (nc_varname == 'Q') THEN
306  ncfunc_retval = nf90_put_var(ncid, varid, &
307 & qv(0:nx-1, 0:ny-1, 1:nz, 1))
308  ELSE
309  print *,
310  print *, 'fp2nc4io:private_dump_3d_field() bad var: ', nc_varname
311  print *, ' ABORTING...'
312  print *,
313  ENDIF
314 
315  IF (ncfunc_retval /= 0) THEN
316  print *,
317  print *, '*** WARNING ***'
318  print *, ' fp2nc4io:private_dump_3d_field()'
319  print *, ' nf90_put_var returned error for var: ', nc_varname
320  print *,
321 
322  ENDIF
323 
324 
325  END SUBROUTINE private_dump_3dfield
326 
327 
328 
329  SUBROUTINE private_read_3dfield(ncid, varname, xdim, ydim, zdim, var_array)
330 
331  ! Private routine for reading full 3D array, specified by varname,
332  ! from NC4 file. Reads into preallocated array of size
333  ! xdim x ydim x zdim
334  IMPLICIT NONE
335 
336  INTEGER, INTENT(IN) :: ncid ! NC4 file id
337  CHARACTER, INTENT(IN) :: varname
338  INTEGER, INTENT(IN) :: xdim, ydim, zdim ! NC4 dimension ids
339  REAL, DIMENSION(xdim, ydim, zdim) :: var_array
340 
341  CHARACTER :: nc_varname
342  INTEGER :: ncfunc_retval, varid
343 
344  ! Check that we have a valid varname. If not, buh-bye
345  IF( .NOT. any(valid_vars == varname) ) THEN
346  print *,
347  print *, 'fp2nc4io:private_dump_3d_field() bad var: ', varname
348  print *, ' ABORTING...'
349  print *,
350  stop
351  ENDIF
352 
353  ! Convert varname to upper case for use in NetCDF file
354  nc_varname = to_upper(varname)
355 
356  ! Get the varid
357  ncfunc_retval = nf90_inq_varid(ncid, nc_varname, varid)
358 
359  ! Read the variable into var_array
360  ncfunc_retval = nf90_get_var(ncid, varid, var_array)
361 
362  END SUBROUTINE private_read_3dfield
363 
364 
365  CHARACTER FUNCTION to_upper(c)
366 
367  ! Utility function to convert lower case char to upper case
368 
369  IMPLICIT NONE
370 
371  CHARACTER, INTENT(IN) :: c
372 
373  INTEGER :: c_ascii_code
374 
375  c_ascii_code = iachar(c)
376  IF (c_ascii_code >= iachar("a") .AND. c_ascii_code <= iachar("z")) THEN
377  to_upper = achar(c_ascii_code - 32)
378  ELSE
379  to_upper = c
380  ENDIF
381 
382  END FUNCTION to_upper
383 
384 
385 
386 END MODULE fp2nc4io_mod
logical function fp2nc4io_vars_are_valid(num_vars, dump_vars)
subroutine fp2nc4io_print_valid_vars
subroutine, private private_dump_3dfield(ncid, varname, dimids, deflevel)
character function, private to_upper(c)
subroutine fp2nc4io_dump(nc4_filepath, num_vars, dump_vars, deflate_level)
subroutine, private private_read_3dfield(ncid, varname, xdim, ydim, zdim, var_array)