Home

Context Navigation

source: palm/trunk/UTIL/inifor/src/io.f90 @ 3254

Last change on this file since 3254 was 3183, checked in by suehring, 6 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 34.1 KB

Line
1	!> @file src/io.f90
2	!------------------------------------------------------------------------------!
3	! This file is part of the PALM model system.
4	!
5	! PALM is free software: you can redistribute it and/or modify it under the
6	! terms of the GNU General Public License as published by the Free Software
7	! Foundation, either version 3 of the License, or (at your option) any later
8	! version.
9	!
10	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
13	!
14	! You should have received a copy of the GNU General Public License along with
15	! PALM. If not, see <http://www.gnu.org/licenses/>.
16	!
17	! Copyright 2017-2018 Leibniz Universitaet Hannover
18	! Copyright 2017-2018 Deutscher Wetterdienst Offenbach
19	!------------------------------------------------------------------------------!
20	!
21	! Current revisions:
22	! -----------------
23	!
24	!
25	! Former revisions:
26	! -----------------
27	! $Id: io.f90 3183 2018-07-27 14:25:55Z suehring $
28	! Introduced new PALM grid stretching
29	! Updated variable names and metadata for PIDS v1.9 compatibility
30	! Improved handling of the start date string
31	! Better compatibility with older Intel compilers:
32	! - avoiding implicit array allocation with new get_netcdf_variable()
33	! subroutine instead of function
34	! Improved command line interface:
35	! - Added configuration validation
36	! - New options to configure input file prefixes
37	! - GNU-style short and long option names
38	! - Added version and copyright output
39	!
40	!
41	! 3182 2018-07-27 13:36:03Z suehring
42	! Initial revision
43	!
44	!
45	!
46	! Authors:
47	! --------
48	! @author Eckhard Kadasch
49	!
50	! Description:
51	! ------------
52	!> The io module contains the functions needed to read and write netCDF data in
53	!> INIFOR.
54	!------------------------------------------------------------------------------!
55	MODULE io
56
57	USE control
58	USE defs, &
59	ONLY: DATE, SNAME, PATH, PI, dp, hp, TO_RADIANS, TO_DEGREES, VERSION
60	USE netcdf
61	USE types
62	USE util, &
63	ONLY: reverse, str, real_to_str
64
65	IMPLICIT NONE
66
67	INTERFACE get_netcdf_variable
68	MODULE PROCEDURE get_netcdf_variable_int
69	MODULE PROCEDURE get_netcdf_variable_real
70	END INTERFACE get_netcdf_variable
71
72	PRIVATE :: get_netcdf_variable_int, get_netcdf_variable_real
73
74	CONTAINS
75
76	SUBROUTINE get_netcdf_variable_int(in_file, in_var, buffer)
77
78	CHARACTER(LEN=PATH), INTENT(IN) :: in_file
79	TYPE(nc_var), INTENT(INOUT) :: in_var
80	INTEGER(hp), ALLOCATABLE, INTENT(INOUT) :: buffer(:,:,:)
81
82	INCLUDE 'get_netcdf_variable.inc'
83
84	END SUBROUTINE get_netcdf_variable_int
85
86
87	SUBROUTINE get_netcdf_variable_real(in_file, in_var, buffer)
88
89	CHARACTER(LEN=PATH), INTENT(IN) :: in_file
90	TYPE(nc_var), INTENT(INOUT) :: in_var
91	REAL(dp), ALLOCATABLE, INTENT(INOUT) :: buffer(:,:,:)
92
93	INCLUDE 'get_netcdf_variable.inc'
94
95	END SUBROUTINE get_netcdf_variable_real
96
97
98	SUBROUTINE netcdf_define_variable(var, ncid)
99
100	TYPE(nc_var), INTENT(INOUT) :: var
101	INTEGER, INTENT(IN) :: ncid
102
103	CALL check(nf90_def_var(ncid, var % name, NF90_FLOAT, var % dimids(1:var % ndim), var % varid))
104	CALL check(nf90_put_att(ncid, var % varid, "long_name", var % long_name))
105	CALL check(nf90_put_att(ncid, var % varid, "units", var % units))
106	IF ( var % lod .GE. 0 ) THEN
107	CALL check(nf90_put_att(ncid, var % varid, "lod", var % lod))
108	END IF
109	CALL check(nf90_put_att(ncid, var % varid, "source", var % source))
110	CALL check(nf90_put_att(ncid, var % varid, "_FillValue", NF90_FILL_REAL))
111
112	END SUBROUTINE netcdf_define_variable
113
114
115	SUBROUTINE netcdf_get_dimensions(var, ncid)
116
117	TYPE(nc_var), INTENT(INOUT) :: var
118	INTEGER, INTENT(IN) :: ncid
119	INTEGER :: i
120	CHARACTER(SNAME) :: null
121
122	! Remember dimension lenghts for NetCDF writing routine
123	DO i = 1, var % ndim
124	CALL check(nf90_inquire_dimension(ncid, var % dimids(i), &
125	name = null, &
126	len = var % dimlen(i) ) )
127	END DO
128
129	END SUBROUTINE netcdf_get_dimensions
130
131
132	!------------------------------------------------------------------------------!
133	! Description:
134	! ------------
135	!> This routine initializes Inifor. This includes parsing command-line
136	!> arguments, setting the names of the input and output file names as well as
137	!> the name of the input namelist and, subsequently, reading in and setting grid
138	!> parameters for the PALM-4U computational grid.
139	!------------------------------------------------------------------------------!
140	SUBROUTINE parse_command_line_arguments( cfg )
141
142	TYPE(inifor_config), INTENT(INOUT) :: cfg
143
144	CHARACTER(LEN=PATH) :: option, arg
145	INTEGER :: arg_count, i
146
147	arg_count = COMMAND_ARGUMENT_COUNT()
148	IF (arg_count .GT. 0) THEN
149
150	! Every option should have an argument.
151	!IF ( MOD(arg_count, 2) .NE. 0 ) THEN
152	! message = "Syntax error in command line."
153	! CALL abort('parse_command_line_arguments', message)
154	!END IF
155
156	message = "The -clon and -clat command line options are depricated. " // &
157	"Please remove them form your inifor command and specify the " // &
158	"location of the PALM-4U origin either" // NEW_LINE(' ') // &
159	" - by setting the namelist parameters 'longitude' and 'latitude', or" // NEW_LINE(' ') // &
160	" - by providing a static driver netCDF file via the -static command-line option."
161
162	i = 1
163	DO WHILE (i .LE. arg_count)
164
165	CALL GET_COMMAND_ARGUMENT( i, option )
166
167	SELECT CASE( TRIM(option) )
168
169	CASE( '-date', '-d', '--date' )
170	CALL get_option_argument( i, arg )
171	cfg % start_date = TRIM(arg)
172
173	! Elevation of the PALM-4U domain above sea level
174	CASE( '-z0', '-z', '--elevation' )
175	CALL get_option_argument( i, arg )
176	READ(arg, *) cfg % z0
177
178	! surface pressure, at z0
179	CASE( '-p0', '-r', '--surface-pressure' )
180	CALL get_option_argument( i, arg )
181	READ(arg, *) cfg % p0
182
183	! geostrophic wind in x direction
184	CASE( '-ug', '-u', '--geostrophic-u' )
185	CALL get_option_argument( i, arg )
186	READ(arg, *) cfg % ug
187
188	! geostrophic wind in y direction
189	CASE( '-vg', '-v', '--geostrophic-v' )
190	CALL get_option_argument( i, arg )
191	READ(arg, *) cfg % vg
192
193	! domain centre geographical longitude and latitude
194	CASE( '-clon', '-clat' )
195	CALL abort('parse_command_line_arguments', message)
196	!READ(arg, *) lambda_cg
197	!lambda_cg = lambda_cg * TO_RADIANS
198	!READ(arg, *) phi_cg
199	!phi_cg = phi_cg * TO_RADIANS
200
201	CASE( '-path', '-p', '--path' )
202	CALL get_option_argument( i, arg )
203	cfg % input_path = TRIM(arg)
204
205	CASE( '-hhl', '-l', '--hhl-file' )
206	CALL get_option_argument( i, arg )
207	cfg % hhl_file = TRIM(arg)
208
209	CASE( '-static', '-t', '--static-driver' )
210	CALL get_option_argument( i, arg )
211	cfg % static_driver_file = TRIM(arg)
212
213	CASE( '-soil', '-s', '--soil-file')
214	CALL get_option_argument( i, arg )
215	cfg % soiltyp_file = TRIM(arg)
216
217	CASE( '--flow-prefix')
218	CALL get_option_argument( i, arg )
219	cfg % flow_prefix = TRIM(arg)
220
221	CASE( '--radiation-prefix')
222	CALL get_option_argument( i, arg )
223	cfg % radiation_prefix = TRIM(arg)
224
225	CASE( '--soil-prefix')
226	CALL get_option_argument( i, arg )
227	cfg % soil_prefix = TRIM(arg)
228
229	CASE( '--soilmoisture-prefix')
230	CALL get_option_argument( i, arg )
231	cfg % soilmoisture_prefix = TRIM(arg)
232
233	CASE( '-o', '--output' )
234	CALL get_option_argument( i, arg )
235	cfg % output_file = TRIM(arg)
236
237	CASE( '-n', '--namelist' )
238	CALL get_option_argument( i, arg )
239	cfg % namelist_file = TRIM(arg)
240
241	! initial condition mode: 'profile' / 'volume'
242	CASE( '-mode', '-i', '--init-mode' )
243	CALL get_option_argument( i, arg )
244	cfg % ic_mode = TRIM(arg)
245
246	! boundary conditions / forcing mode: 'ideal' / 'real'
247	CASE( '-f', '--forcing-mode' )
248	CALL get_option_argument( i, arg )
249	cfg % bc_mode = TRIM(arg)
250
251	CASE( '--version' )
252	CALL print_version()
253	STOP
254
255	CASE( '--help' )
256	CALL print_version()
257	PRINT *, ""
258	PRINT *, "For a list of command-line options have a look at the README file."
259	STOP
260
261	CASE DEFAULT
262	message = "unknown option '" // TRIM(option) // "'."
263	CALL abort('parse_command_line_arguments', message)
264
265	END SELECT
266
267	i = i + 1
268
269	END DO
270
271	ELSE
272
273	message = "No arguments present, using default input and output files"
274	CALL report('parse_command_line_arguments', message)
275
276	END IF
277
278	END SUBROUTINE parse_command_line_arguments
279
280
281	SUBROUTINE get_option_argument(i, arg)
282	CHARACTER(LEN=PATH), INTENT(INOUT) :: arg
283	INTEGER, INTENT(INOUT) :: i
284
285	i = i + 1
286	CALL GET_COMMAND_ARGUMENT(i, arg)
287
288	END SUBROUTINE
289
290
291	SUBROUTINE validate_config(cfg)
292	TYPE(inifor_config), INTENT(IN) :: cfg
293	LOGICAL :: all_files_present
294
295	all_files_present = .TRUE.
296	all_files_present = all_files_present .AND. file_present(cfg % hhl_file)
297	all_files_present = all_files_present .AND. file_present(cfg % namelist_file)
298	all_files_present = all_files_present .AND. file_present(cfg % output_file)
299	all_files_present = all_files_present .AND. file_present(cfg % soiltyp_file)
300
301	! Only check optional static driver file name, if it has been given.
302	IF (TRIM(cfg % static_driver_file) .NE. '') THEN
303	all_files_present = all_files_present .AND. file_present(cfg % static_driver_file)
304	END IF
305
306	IF (.NOT. all_files_present) THEN
307	message = "INIFOR configuration invalid; some input files are missing."
308	CALL abort( 'validate_config', message )
309	END IF
310
311
312	SELECT CASE( TRIM(cfg % ic_mode) )
313	CASE( 'profile', 'volume')
314	CASE DEFAULT
315	message = "Initialization mode '" // TRIM(cfg % ic_mode) //&
316	"' is not supported. " //&
317	"Please select either 'profile' or 'volume', " //&
318	"or omit the -i/--init-mode/-mode option entirely, which corresponds "//&
319	"to the latter."
320	CALL abort( 'validate_config', message )
321	END SELECT
322
323
324	SELECT CASE( TRIM(cfg % bc_mode) )
325	CASE( 'real', 'ideal')
326	CASE DEFAULT
327	message = "Forcing mode '" // TRIM(cfg % bc_mode) //&
328	"' is not supported. " //&
329	"Please select either 'real' or 'ideal', " //&
330	"or omit the -f/--forcing-mode option entirely, which corresponds "//&
331	"to the latter."
332	CALL abort( 'validate_config', message )
333	END SELECT
334
335
336	END SUBROUTINE validate_config
337
338
339	LOGICAL FUNCTION file_present(filename)
340	CHARACTER(LEN=PATH), INTENT(IN) :: filename
341
342	INQUIRE(FILE=filename, EXIST=file_present)
343
344	IF (.NOT. file_present) THEN
345	message = "The given file '" // "' does not exist."
346	CALL report('file_present', message)
347	END IF
348
349	END FUNCTION file_present
350
351
352	!------------------------------------------------------------------------------!
353	! Description:
354	! ------------
355	!> This routine initializes the Inifor output file, i.e. the PALM-4U
356	!> initializing and forcing data as a NetCDF file.
357	!>
358	!> Besides writing metadata, such as global attributes, coordinates, variables,
359	!> in the NetCDF file, various NetCDF IDs are saved for later, when Inifor
360	!> writes the actual data.
361	!------------------------------------------------------------------------------!
362	SUBROUTINE setup_netcdf_dimensions(output_file, palm_grid, &
363	start_date_string, origin_lon, origin_lat)
364
365	TYPE(nc_file), INTENT(INOUT) :: output_file
366	TYPE(grid_definition), INTENT(IN) :: palm_grid
367	CHARACTER (LEN=DATE), INTENT(IN) :: start_date_string
368	REAL(dp), INTENT(IN) :: origin_lon, origin_lat
369
370	CHARACTER (LEN=8) :: date_string
371	CHARACTER (LEN=10) :: time_string
372	CHARACTER (LEN=5) :: zone_string
373	CHARACTER (LEN=SNAME) :: history_string
374	INTEGER :: ncid, nx, ny, nz, nt, dimids(3), dimvarids(3)
375	REAL(dp) :: z0
376
377	message = "Initializing PALM-4U dynamic driver file '" // &
378	TRIM(output_file % name) // "' and setting up dimensions."
379	CALL report('setup_netcdf_dimensions', message)
380
381	! Create the NetCDF file. NF90_CLOBBER selects overwrite mode.
382	#if defined( __netcdf4 )
383	CALL check(nf90_create(TRIM(output_file % name), OR(NF90_CLOBBER, NF90_HDF5), ncid))
384	#else
385	CALL check(nf90_create(TRIM(output_file % name), NF90_CLOBBER, ncid))
386	#endif
387
388	!
389	!------------------------------------------------------------------------------
390	!- Section 1: Write global NetCDF attributes
391	!------------------------------------------------------------------------------
392	CALL date_and_time(DATE=date_string, TIME=time_string, ZONE=zone_string)
393	history_string = &
394	'Created on '// date_string // &
395	' at ' // time_string(1:2) // ':' // time_string(3:4) // &
396	' (UTC' // zone_string // ')'
397
398	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'title', 'PALM input file for scenario ...'))
399	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'institution', 'Deutscher Wetterdienst, Offenbach'))
400	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'author', 'Eckhard Kadasch, eckhard.kadasch@dwd.de'))
401	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'history', TRIM(history_string)))
402	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'references', '--'))
403	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'comment', '--'))
404	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'origin_lat', TRIM(real_to_str(origin_lat*TO_DEGREES, '(F18.13)'))))
405	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'origin_lon', TRIM(real_to_str(origin_lon*TO_DEGREES, '(F18.13)'))))
406	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'inifor_version', TRIM(VERSION)))
407	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'palm_version', '--'))
408
409	!
410	!------------------------------------------------------------------------------
411	!- Section 2: Define NetCDF dimensions and coordinates
412	!------------------------------------------------------------------------------
413	nt = SIZE(output_file % time)
414	nx = palm_grid % nx
415	ny = palm_grid % ny
416	nz = palm_grid % nz
417	z0 = palm_grid % z0
418
419	!
420	!------------------------------------------------------------------------------
421	!- Section 2a: Define dimensions for cell centers (scalars in soil and atmosph.)
422	!------------------------------------------------------------------------------
423	dimids = (/0, 0, 0/) ! reset dimids
424	CALL check( nf90_def_dim(ncid, "x", nx+1, dimids(1)) )
425	CALL check( nf90_def_dim(ncid, "y", ny+1, dimids(2)) )
426	CALL check( nf90_def_dim(ncid, "z", nz, dimids(3)) )
427	output_file % dimids_scl = dimids ! save dimids for later
428
429	dimvarids = (/0, 0, 0/) ! reset dimvarids
430	CALL check(nf90_def_var(ncid, "x", NF90_FLOAT, dimids(1), dimvarids(1)))
431	CALL check(nf90_put_att(ncid, dimvarids(1), "standard_name", "x coordinate of cell centers"))
432	CALL check(nf90_put_att(ncid, dimvarids(1), "units", "m"))
433
434	CALL check(nf90_def_var(ncid, "y", NF90_FLOAT, dimids(2), dimvarids(2)))
435	CALL check(nf90_put_att(ncid, dimvarids(2), "standard_name", "y coordinate of cell centers"))
436	CALL check(nf90_put_att(ncid, dimvarids(2), "units", "m"))
437
438	CALL check(nf90_def_var(ncid, "z", NF90_FLOAT, dimids(3), dimvarids(3)))
439	CALL check(nf90_put_att(ncid, dimvarids(3), "standard_name", "z coordinate of cell centers"))
440	CALL check(nf90_put_att(ncid, dimvarids(3), "units", "m"))
441	output_file % dimvarids_scl = dimvarids ! save dimvarids for later
442
443	! overwrite third dimid with the one of depth
444	CALL check(nf90_def_dim(ncid, "zsoil", SIZE(palm_grid % depths), dimids(3)) )
445	output_file % dimids_soil = dimids ! save dimids for later
446
447	! overwrite third dimvarid with the one of depth
448	CALL check(nf90_def_var(ncid, "zsoil", NF90_FLOAT, output_file % dimids_soil(3), dimvarids(3)))
449	CALL check(nf90_put_att(ncid, dimvarids(3), "standard_name", "depth_below_land"))
450	CALL check(nf90_put_att(ncid, dimvarids(3), "positive", "down"))
451	CALL check(nf90_put_att(ncid, dimvarids(3), "units", "m"))
452	output_file % dimvarids_soil = dimvarids ! save dimvarids for later
453	!
454	!------------------------------------------------------------------------------
455	!- Section 2b: Define dimensions for cell faces/velocities
456	!------------------------------------------------------------------------------
457	dimids = (/0, 0, 0/) ! reset dimids
458	CALL check(nf90_def_dim(ncid, "xu", nx, dimids(1)) )
459	CALL check(nf90_def_dim(ncid, "yv", ny, dimids(2)) )
460	CALL check(nf90_def_dim(ncid, "zw", nz-1, dimids(3)) )
461	output_file % dimids_vel = dimids ! save dimids for later
462
463	dimvarids = (/0, 0, 0/) ! reset dimvarids
464	CALL check(nf90_def_var(ncid, "xu", NF90_FLOAT, dimids(1), dimvarids(1)))
465	CALL check(nf90_put_att(ncid, dimvarids(1), "standard_name", "x coordinate of cell faces"))
466	CALL check(nf90_put_att(ncid, dimvarids(1), "units", "m"))
467
468	CALL check(nf90_def_var(ncid, "yv", NF90_FLOAT, dimids(2), dimvarids(2)))
469	CALL check(nf90_put_att(ncid, dimvarids(2), "standard_name", "y coordinate of cell faces"))
470	CALL check(nf90_put_att(ncid, dimvarids(2), "units", "m"))
471
472	CALL check(nf90_def_var(ncid, "zw", NF90_FLOAT, dimids(3), dimvarids(3)))
473	CALL check(nf90_put_att(ncid, dimvarids(3), "standard_name", "z coordinate of cell faces"))
474	CALL check(nf90_put_att(ncid, dimvarids(3), "units", "m"))
475	output_file % dimvarids_vel = dimvarids ! save dimvarids for later
476
477	!
478	!------------------------------------------------------------------------------
479	!- Section 2c: Define time dimension
480	!------------------------------------------------------------------------------
481	CALL check(nf90_def_dim(ncid, "time", nt, output_file % dimid_time) )
482	CALL check(nf90_def_var(ncid, "time", NF90_FLOAT, &
483	output_file % dimid_time, &
484	output_file % dimvarid_time))
485	CALL check(nf90_put_att(ncid, output_file % dimvarid_time, "standard_name", "time"))
486	CALL check(nf90_put_att(ncid, output_file % dimvarid_time, "long_name", "time"))
487	CALL check(nf90_put_att(ncid, output_file % dimvarid_time, "units", &
488	"seconds since " // start_date_string // " UTC"))
489
490	CALL check(nf90_enddef(ncid))
491
492	!
493	!------------------------------------------------------------------------------
494	!- Section 3: Write grid coordinates
495	!------------------------------------------------------------------------------
496	CALL check(nf90_put_var(ncid, output_file % dimvarids_scl(1), palm_grid%x))
497	CALL check(nf90_put_var(ncid, output_file % dimvarids_scl(2), palm_grid%y))
498	CALL check(nf90_put_var(ncid, output_file % dimvarids_scl(3), palm_grid%z))
499
500	CALL check(nf90_put_var(ncid, output_file % dimvarids_vel(1), palm_grid%xu))
501	CALL check(nf90_put_var(ncid, output_file % dimvarids_vel(2), palm_grid%yv))
502	CALL check(nf90_put_var(ncid, output_file % dimvarids_vel(3), palm_grid%zw))
503
504	! TODO Read in soil depths from input file before this.
505	CALL check(nf90_put_var(ncid, output_file % dimvarids_soil(3), palm_grid%depths))
506
507	! Write time vector
508	CALL check(nf90_put_var(ncid, output_file % dimvarid_time, output_file % time))
509
510	! Close the file
511	CALL check(nf90_close(ncid))
512
513	END SUBROUTINE setup_netcdf_dimensions
514
515
516	SUBROUTINE setup_netcdf_variables(filename, output_variable_table)
517
518	CHARACTER (LEN=*), INTENT(IN) :: filename
519	TYPE(nc_var), INTENT(INOUT), TARGET :: output_variable_table(:)
520	TYPE(nc_var), POINTER :: var
521	INTEGER :: i, ncid
522
523	message = "Defining variables in dynamic driver '" // TRIM(filename) // "'."
524	CALL report('setup_netcdf_variables', message)
525
526	CALL check(nf90_open(TRIM(filename), NF90_WRITE, ncid))
527	CALL check(nf90_redef(ncid))
528
529	DO i = 1, SIZE(output_variable_table)
530
531	var => output_variable_table(i)
532
533	IF ( var % to_be_processed ) THEN
534	message = " variable #" // TRIM(str(i)) // " '" // TRIM(var%name) // "'."
535	CALL report('setup_netcdf_variables', message)
536
537	CALL netcdf_define_variable(var, ncid)
538	CALL netcdf_get_dimensions(var, ncid)
539	END IF
540
541	END DO
542
543	CALL check(nf90_enddef(ncid))
544	CALL check(nf90_close(ncid))
545
546	message = "Dynamic driver '" // TRIM(filename) // "' initialized successfully."
547	CALL report('setup_netcdf_variables', message)
548
549	END SUBROUTINE setup_netcdf_variables
550
551
552	!------------------------------------------------------------------------------!
553	! Description:
554	! ------------
555	!> This routine reads and returns all input variables of the given IO group
556	!> It accomodates the data by allocating a container variable that represents a
557	!> list of arrays of the same length as the groups variable list. (This list
558	!> will typically contain one or two items.) After the container, its members
559	!> are allocated one by one with the appropriate, possibly different,
560	!> dimensions.
561	!>
562	!> The 'group' is an INTENT(INOUT) variable so that 'get_netcdf_variable()' can
563	!> record netCDF IDs in the 'in_var_list()' member variable.
564	!------------------------------------------------------------------------------!
565	SUBROUTINE read_input_variables(group, iter, buffer)
566	TYPE(io_group), INTENT(INOUT), TARGET :: group
567	INTEGER, INTENT(IN) :: iter
568	TYPE(container), ALLOCATABLE, INTENT(INOUT) :: buffer(:)
569	INTEGER :: hour, buf_id
570	TYPE(nc_var), POINTER :: input_var
571	CHARACTER(LEN=PATH), POINTER :: input_file
572	INTEGER :: ivar, nbuffers
573
574	message = "Reading data for I/O group '" // TRIM(group % in_var_list(1) % name) // "'."
575	CALL report('read_input_variables', message)
576
577	input_file => group % in_files(iter)
578
579	!
580	!------------------------------------------------------------------------------
581	!- Section 1: Load input buffers for accumulated variables
582	!------------------------------------------------------------------------------
583	IF (group % kind == 'running average' .OR. &
584	group % kind == 'accumulated') THEN ! radiation budgets, precipitation
585
586	IF (SIZE(group % in_var_list) .GT. 1 ) THEN
587	message = "I/O groups may not contain more than one " // &
588	"accumulated variable. Group '" // TRIM(group % kind) //&
589	"' contains " // &
590	TRIM( str(SIZE(group % in_var_list)) ) // "."
591	CALL abort('read_input_variables \| accumulation', message)
592	END IF
593
594	! use two buffer arrays
595	nbuffers = 2
596	IF ( .NOT. ALLOCATED( buffer ) ) ALLOCATE( buffer(nbuffers) )
597
598	! chose correct buffer array
599	hour = iter - 1! hour of the day
600	buf_id = select_buffer(hour)
601
602	CALL run_control('time', 'read')
603	IF ( ALLOCATED(buffer(buf_id) % array) ) DEALLOCATE(buffer(buf_id) % array)
604	CALL run_control('time', 'alloc')
605
606	input_var => group % in_var_list(1)
607	CALL get_netcdf_variable(input_file, input_var, buffer(buf_id) % array)
608	CALL report('read_input_variables', "Read accumulated " // TRIM(group % in_var_list(1) % name))
609
610	IF ( input_var % is_upside_down ) CALL reverse(buffer(buf_id) % array)
611	CALL run_control('time', 'comp')
612
613	!------------------------------------------------------------------------------
614	!- Section 2: Load input buffers for normal I/O groups
615	!------------------------------------------------------------------------------
616	ELSE
617
618	nbuffers = SIZE( group % in_var_list )
619	ALLOCATE( buffer(nbuffers) )
620	CALL run_control('time', 'alloc')
621
622	DO ivar = 1, nbuffers
623
624	input_var => group % in_var_list(ivar)
625
626	! Check wheather P or PP is present in input file
627	IF (input_var % name == 'P') THEN
628	input_var % name = TRIM( get_pressure_var(input_file) )
629	CALL run_control('time', 'read')
630	END IF
631
632	CALL get_netcdf_variable(input_file, input_var, buffer(ivar) % array)
633
634	IF ( input_var % is_upside_down ) CALL reverse(buffer(ivar) % array)
635	CALL run_control('time', 'comp')
636
637	END DO
638	END IF
639
640	END SUBROUTINE read_input_variables
641
642
643	INTEGER FUNCTION select_buffer(hour)
644	INTEGER, INTENT(IN) :: hour
645	INTEGER :: step
646
647	select_buffer = 0
648	step = MODULO(hour, 3) + 1
649
650	SELECT CASE(step)
651	CASE(1, 3)
652	select_buffer = 1
653	CASE(2)
654	select_buffer = 2
655	CASE DEFAULT
656	message = "Invalid step '" // TRIM(str(step))
657	CALL abort('select_buffer', message)
658	END SELECT
659	END FUNCTION select_buffer
660
661
662	!------------------------------------------------------------------------------!
663	! Description:
664	! ------------
665	!> Checks if the input_file contains the absolute pressure, 'P', or the pressure
666	!> perturbation, 'PP', and returns the appropriate string.
667	!------------------------------------------------------------------------------!
668	CHARACTER(LEN=2) FUNCTION get_pressure_var(input_file) RESULT(var)
669	CHARACTER(LEN=*) :: input_file
670	INTEGER :: ncid, varid
671
672	CALL check(nf90_open( TRIM(input_file), NF90_NOWRITE, ncid ))
673	IF ( nf90_inq_varid( ncid, 'P', varid ) .EQ. NF90_NOERR ) THEN
674
675	var = 'P'
676
677	ELSE IF ( nf90_inq_varid( ncid, 'PP', varid ) .EQ. NF90_NOERR ) THEN
678
679	var = 'PP'
680	CALL report('get_pressure_var', 'Using PP instead of P')
681
682	ELSE
683
684	message = "Failed to read '" // TRIM(var) // &
685	"' from file '" // TRIM(input_file) // "'."
686	CALL abort('get_pressure_var', message)
687
688	END IF
689
690	CALL check(nf90_close(ncid))
691
692	END FUNCTION get_pressure_var
693
694
695	FUNCTION get_netcdf_attribute(filename, attribute) RESULT(attribute_value)
696
697	CHARACTER(LEN=*), INTENT(IN) :: filename, attribute
698	REAL(dp) :: attribute_value
699
700	INTEGER :: ncid
701
702	IF ( nf90_open( TRIM(filename), NF90_NOWRITE, ncid ) == NF90_NOERR ) THEN
703
704	CALL check(nf90_get_att(ncid, NF90_GLOBAL, TRIM(attribute), attribute_value))
705	CALL check(nf90_close(ncid))
706
707	ELSE
708
709	message = "Failed to read '" // TRIM(attribute) // &
710	"' from file '" // TRIM(filename) // "'."
711	CALL abort('get_netcdf_attribute', message)
712
713	END IF
714
715	END FUNCTION get_netcdf_attribute
716
717
718	SUBROUTINE update_output(var, array, iter, output_file)
719	TYPE(nc_var), INTENT(IN) :: var
720	REAL(dp), INTENT(IN) :: array(:,:,:)
721	INTEGER, INTENT(IN) :: iter
722	TYPE(nc_file), INTENT(IN) :: output_file
723
724	INTEGER :: ncid, ndim, start(4), count(4)
725	LOGICAL :: var_is_time_dependent
726
727	var_is_time_dependent = ( &
728	var % dimids( var % ndim ) == output_file % dimid_time &
729	)
730
731	! Skip time dimension for output
732	ndim = var % ndim
733	IF ( var_is_time_dependent ) ndim = var % ndim - 1
734
735	start(:) = (/1,1,1,1/)
736	start(ndim+1) = iter
737	count(1:ndim) = var%dimlen(1:ndim)
738
739	CALL check(nf90_open(output_file % name, NF90_WRITE, ncid))
740
741	! Reduce dimension of output array according to variable kind
742	SELECT CASE (TRIM(var % kind))
743
744	CASE ( 'init scalar profile', 'init u profile', 'init v profile', &
745	'init w profile' )
746
747	CALL check(nf90_put_var( ncid, var%varid, array(1,1,:) ) )
748
749	CASE ( 'init soil', 'init scalar', 'init u', 'init v', 'init w' )
750
751	CALL check(nf90_put_var( ncid, var%varid, array(:,:,:) ) )
752
753	CASE( 'left scalar', 'right scalar', 'left w', 'right w' )
754
755	CALL check(nf90_put_var( ncid, var%varid, array(1,:,:), &
756	start=start(1:ndim+1), &
757	count=count(1:ndim) ) )
758
759
760	IF (.NOT. SIZE(array, 2) .EQ. var % dimlen(1)) THEN
761	PRINT *, "inifor: update_output: Dimension ", 1, " of variable ", &
762	TRIM(var % name), " (", var % dimlen(1), &
763	") does not match the dimension of the output array (", &
764	SIZE(array, 2), ")."
765	STOP
766	END IF
767
768
769	CASE( 'north scalar', 'south scalar', 'north w', 'south w' )
770
771	CALL check(nf90_put_var( ncid, var%varid, array(:,1,:), &
772	start=start(1:ndim+1), &
773	count=count(1:ndim) ) )
774
775
776	CASE( 'surface forcing', 'top scalar', 'top w' )
777
778	CALL check(nf90_put_var( ncid, var%varid, array(:,:,1), &
779	start=start(1:ndim+1), &
780	count=count(1:ndim) ) )
781
782	CASE ( 'left u', 'right u', 'left v', 'right v' )
783
784	CALL check(nf90_put_var( ncid, var%varid, array(1,:,:), &
785	start=start(1:ndim+1), &
786	count=count(1:ndim) ) )
787
788	CASE ( 'north u', 'south u', 'north v', 'south v' )
789
790	CALL check(nf90_put_var( ncid, var%varid, array(:,1,:), &
791	start=start(1:ndim+1), &
792	count=count(1:ndim) ) )
793
794	CASE ( 'top u', 'top v' )
795
796	CALL check(nf90_put_var( ncid, var%varid, array(:,:,1), &
797	start=start(1:ndim+1), &
798	count=count(1:ndim) ) )
799
800	CASE ( 'time series' )
801
802	CALL check(nf90_put_var( ncid, var%varid, array(1,1,1), &
803	start=start(1:ndim+1) ) )
804
805	CASE ( 'constant scalar profile' )
806
807	CALL check(nf90_put_var( ncid, var%varid, array(1,1,:), &
808	start=start(1:ndim+1), &
809	count=count(1:ndim) ) )
810
811	CASE ( 'large-scale scalar forcing', 'large-scale w forcing' )
812
813	message = "Doing nothing in terms of writing large-scale forings."
814	CALL report('update_output', message)
815
816	CASE DEFAULT
817
818	message = "Variable kind '" // TRIM(var % kind) // &
819	"' not recognized."
820	CALL abort('update_output', message)
821
822	END SELECT
823
824	CALL check(nf90_close(ncid))
825
826	END SUBROUTINE update_output
827
828
829	SUBROUTINE write_netcdf_variable_2d(var, iter, output_file, buffer)
830	TYPE(nc_var), INTENT(IN) :: var
831	INTEGER, INTENT(IN) :: iter
832	TYPE(nc_file), INTENT(IN) :: output_file
833	REAL(dp), INTENT(IN) :: buffer(:,:,:)
834
835	INTEGER :: ncid, ndim_out, start(4), count(4)
836	LOGICAL :: last_dimension_is_time
837
838	ndim_out = var % ndim
839
840	last_dimension_is_time = var % dimids( var % ndim ) == output_file % dimid_time
841	IF ( last_dimension_is_time ) THEN
842	ndim_out = ndim_out - 1
843	END IF
844
845	start(:) = (/1,1,1,iter/)
846	count(1:ndim_out) = var%dimlen(1:ndim_out)
847
848	CALL check(nf90_open(output_file % name, NF90_WRITE, ncid))
849
850	IF (TRIM(var % kind) .EQ. 'left/right scalar') THEN
851
852	CALL check(nf90_put_var( ncid, var%varid, buffer(1,:,:) ) )
853
854	ELSE IF (TRIM(var % kind) .EQ. 'north/south scalar') THEN
855
856	CALL check(nf90_put_var( ncid, var%varid, buffer(:,1,:) ) )
857
858	ELSE IF (TRIM(var % kind) .EQ. 'top scalar') THEN
859
860	CALL check(nf90_put_var( ncid, var%varid, buffer(:,:,1) ) )
861	ELSE
862
863	CALL check(nf90_put_var( ncid, var%varid, buffer ) )
864
865	END IF
866	CALL check(nf90_close(ncid))
867
868	END SUBROUTINE write_netcdf_variable_2d
869
870
871	SUBROUTINE check(status)
872
873	INTEGER, INTENT(IN) :: status
874
875	IF (status /= nf90_noerr) THEN
876	message = "NetCDF API call failed with error: " // &
877	TRIM( nf90_strerror(status) )
878	CALL abort('io.check', message)
879	END IF
880
881	END SUBROUTINE check
882
883	END MODULE io

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |