Home

Context Navigation

source: palm/trunk/UTIL/inifor/src/inifor.f90 @ 3525

Last change on this file since 3525 was 3456, checked in by eckhard, 6 years ago
inifor: Removed surface forcing and internal arrays from netCDF output
Property svn:keywords set to `Id`
File size: 21.0 KB

Line
1	!> @file src/inifor.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: inifor.f90 3456 2018-10-30 14:29:54Z kanani $
28	! NetCDf output of internal arrays only with --debug option
29	!
30	!
31	! 3401 2018-10-23 09:33:31Z eckhard
32	! Re-compute geostrophic winds every time step
33	!
34	! 3395 2018-10-22 17:32:49Z eckhard
35	! Added main loop support for computation of geostrophic winds and surface
36	! pressure
37	! Cleaned up terminal output, show some meVssages only with --debug option
38	!
39	! 3183 2018-07-27 14:25:55Z suehring
40	! Introduced new PALM grid stretching
41	! Renamend initial-condition mode variable 'mode' to 'ic_mode'
42	! Improved log messages
43	!
44	!
45	! 3182 2018-07-27 13:36:03Z suehring
46	! Initial revision
47	!
48	!
49	!
50	! Authors:
51	! --------
52	! @author Eckhard Kadasch
53	!
54	! Description:
55	! ------------
56	!> INIFOR is an interpolation tool for generating meteorological initialization
57	!> and forcing data for the urban climate model PALM-4U. The required
58	!> meteorological fields are interpolated from output data of the mesoscale
59	!> model COSMO-DE. This is the main program file.
60	!------------------------------------------------------------------------------!
61	PROGRAM inifor
62
63	USE control
64	USE defs
65	USE grid, &
66	ONLY: setup_parameters, setup_grids, setup_variable_tables, &
67	setup_io_groups, fini_grids, fini_variables, fini_io_groups, &
68	fini_file_lists, preprocess, origin_lon, origin_lat, &
69	output_file, io_group_list, output_var_table, &
70	cosmo_grid, palm_grid, nx, ny, nz, p0, cfg, f3, &
71	averaging_width_ns, averaging_width_ew, phi_n, lambda_n, &
72	lam_centre, phi_centre
73	USE io
74	USE transform, &
75	ONLY: average_profile, interpolate_2d, interpolate_3d, &
76	geostrophic_winds, extrapolate_density, extrapolate_pressure, &
77	get_surface_pressure
78	USE types
79
80	IMPLICIT NONE
81
82	INTEGER :: igroup
83	INTEGER :: ivar
84	INTEGER :: iter
85	REAL(dp), ALLOCATABLE, DIMENSION(:,:,:) :: output_arr
86	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_centre
87	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: ug_arr
88	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: vg_arr
89	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_north
90	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_south
91	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_east
92	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: rho_west
93	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: p_north
94	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: p_south
95	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: p_east
96	REAL(dp), ALLOCATABLE, DIMENSION(:), TARGET :: p_west
97	REAL(dp), POINTER, DIMENSION(:) :: internal_arr
98	REAL(dp), POINTER, DIMENSION(:) :: ug_vg_arr
99	TYPE(nc_var), POINTER :: output_var
100	TYPE(io_group), POINTER :: group
101	TYPE(container), ALLOCATABLE :: input_buffer(:)
102	LOGICAL, SAVE :: ug_vg_have_been_computed = .FALSE.
103	LOGICAL, SAVE :: debugging_output = .TRUE.
104
105	!> \mainpage About INIFOR
106	!> ...
107	!
108	!------------------------------------------------------------------------------
109	!- Section 1: Initialization
110	!------------------------------------------------------------------------------
111	CALL run_control('init', 'void')
112
113	! Initialize INIFOR's parameters from command-line interface and namelists
114	CALL setup_parameters()
115
116	! Initialize all grids, including interpolation neighbours and weights
117	CALL setup_grids()
118	CALL run_control('time', 'init')
119
120	! Initialize the netCDF output file and define dimensions
121	CALL setup_netcdf_dimensions(output_file, palm_grid, cfg % start_date, &
122	origin_lon, origin_lat)
123	CALL run_control('time', 'write')
124
125	! Set up the tables containing the input and output variables and set
126	! the corresponding netCDF dimensions for each output variable
127	CALL setup_variable_tables(cfg % ic_mode)
128	CALL run_control('time', 'write')
129
130	! Add the output variables to the netCDF output file
131	CALL setup_netcdf_variables(output_file % name, output_var_table, &
132	cfg % debug)
133
134	CALL setup_io_groups()
135	CALL run_control('time', 'init')
136
137	!------------------------------------------------------------------------------
138	!- Section 2: Main loop
139	!------------------------------------------------------------------------------
140
141	DO igroup = 1, SIZE(io_group_list)
142
143	group => io_group_list(igroup)
144	IF ( group % to_be_processed ) THEN
145
146	DO iter = 1, group % nt
147
148	!------------------------------------------------------------------------------
149	!- Section 2.1: Read and preprocess input data
150	!------------------------------------------------------------------------------
151	CALL read_input_variables(group, iter, input_buffer)
152	CALL run_control('time', 'read')
153
154	CALL preprocess(group, input_buffer, cosmo_grid, iter)
155	CALL run_control('time', 'comp')
156
157	!TODO: move this assertion into 'preprocess'.
158	IF ( .NOT. ALL(input_buffer(:) % is_preprocessed .AND. .TRUE.) ) THEN
159	message = "Input buffers for group '" // TRIM(group % kind) // &
160	"' could not be preprocessed sucessfully."
161	CALL abort('main loop', message)
162	END IF
163
164	!------------------------------------------------------------------------------
165	!- Section 2.2: Interpolate each output variable of the group
166	!------------------------------------------------------------------------------
167	DO ivar = 1, group % nv
168
169	output_var => group % out_vars( ivar )
170
171	IF ( output_var % to_be_processed .AND. &
172	iter .LE. output_var % nt ) THEN
173
174	message = "Processing '" // TRIM(output_var % name) // &
175	"' (" // TRIM(output_var % kind) // &
176	"), iteration " // TRIM(str(iter)) //" of " // &
177	TRIM(str(output_var % nt))
178	CALL report('main loop', message)
179
180	SELECT CASE( TRIM(output_var % task) )
181
182	CASE( 'interpolate_2d' )
183
184	SELECT CASE( TRIM(output_var % kind) )
185
186	CASE( 'init soil' )
187
188	ALLOCATE( output_arr( 0:output_var % grid % nx, &
189	0:output_var % grid % ny, &
190	SIZE(output_var % grid % depths) ) )
191
192	CASE ( 'surface forcing' )
193
194	ALLOCATE( output_arr( 0:output_var % grid % nx, &
195	0:output_var % grid % ny, 1 ) )
196
197	CASE DEFAULT
198
199	message = "'" // TRIM(output_var % kind) // "' is not a soil variable"
200	CALL abort("main loop", message)
201
202	END SELECT
203	CALL run_control('time', 'alloc')
204
205	CALL interpolate_2d(input_buffer(output_var % input_id) % array(:,:,:), &
206	output_arr(:,:,:), output_var % intermediate_grid, output_var)
207	CALL run_control('time', 'comp')
208
209
210	CASE ( 'interpolate_3d' )
211
212	ALLOCATE( output_arr( 0:output_var % grid % nx, &
213	0:output_var % grid % ny, &
214	1:output_var % grid % nz ) )
215
216	CALL run_control('time', 'alloc')
217	CALL interpolate_3d( &
218	input_buffer(output_var % input_id) % array(:,:,:), &
219	output_arr(:,:,:), &
220	output_var % intermediate_grid, &
221	output_var % grid)
222	CALL run_control('time', 'comp')
223
224	CASE ( 'average profile' )
225
226	ALLOCATE( output_arr( 0:output_var % grid % nx, &
227	0:output_var % grid % ny, &
228	1:output_var % grid % nz ) )
229	CALL run_control('time', 'alloc')
230
231	CALL average_profile( &
232	input_buffer(output_var % input_id) % array(:,:,:), &
233	output_arr(0,0,:), &
234	output_var % averaging_grid)
235
236	IF ( TRIM(output_var % name) == &
237	'surface_forcing_surface_pressure' ) THEN
238
239	IF ( cfg % p0_is_set ) THEN
240	output_arr(0,0,1) = p0
241	ELSE
242	CALL get_surface_pressure( &
243	output_arr(0,0,:), rho_centre, &
244	output_var % averaging_grid)
245	END IF
246
247	END IF
248	CALL run_control('time', 'comp')
249
250	CASE ( 'internal profile' )
251
252	message = "Averaging of internal profile for variable '" //&
253	TRIM(output_var % name) // "' is not supported."
254
255	SELECT CASE (TRIM(output_var % name))
256
257	CASE('internal_density_centre')
258	ALLOCATE( rho_centre( 1:output_var % grid % nz ) )
259	internal_arr => rho_centre
260
261	CASE('internal_density_north')
262	ALLOCATE( rho_north( 1:output_var % grid % nz ) )
263	internal_arr => rho_north
264
265	CASE('internal_density_south')
266	ALLOCATE( rho_south( 1:output_var % grid % nz ) )
267	internal_arr => rho_south
268
269	CASE('internal_density_east')
270	ALLOCATE( rho_east( 1:output_var % grid % nz) )
271	internal_arr => rho_east
272
273	CASE('internal_density_west')
274	ALLOCATE( rho_west( 1:output_var % grid % nz ) )
275	internal_arr => rho_west
276
277	CASE('internal_pressure_north')
278	ALLOCATE( p_north( 1:output_var % grid % nz ) )
279	internal_arr => p_north
280
281	CASE('internal_pressure_south')
282	ALLOCATE( p_south( 1:output_var % grid % nz ) )
283	internal_arr => p_south
284
285	CASE('internal_pressure_east')
286	ALLOCATE( p_east( 1:output_var % grid % nz) )
287	internal_arr => p_east
288
289	CASE('internal_pressure_west')
290	ALLOCATE( p_west( 1:output_var % grid % nz ) )
291	internal_arr => p_west
292
293	CASE DEFAULT
294	CALL abort('main loop', message)
295
296	END SELECT
297	CALL run_control('time', 'alloc')
298
299
300	CALL average_profile( &
301	input_buffer(output_var % input_id) % array(:,:,:),&
302	internal_arr(:), &
303	output_var % averaging_grid)
304
305	SELECT CASE (TRIM(output_var % name))
306
307	CASE('internal_density_centre', &
308	'internal_density_north', &
309	'internal_density_south', &
310	'internal_density_east', &
311	'internal_density_west')
312	CALL extrapolate_density(internal_arr, &
313	output_var % averaging_grid)
314
315	CASE('internal_pressure_north')
316	CALL extrapolate_pressure(internal_arr, rho_north, &
317	output_var % averaging_grid)
318
319	CASE('internal_pressure_south')
320	CALL extrapolate_pressure(internal_arr, rho_south, &
321	output_var % averaging_grid)
322
323	CASE('internal_pressure_east')
324	CALL extrapolate_pressure(internal_arr, rho_east, &
325	output_var % averaging_grid)
326
327	CASE('internal_pressure_west')
328	CALL extrapolate_pressure(internal_arr, rho_west, &
329	output_var % averaging_grid)
330
331	CASE DEFAULT
332	CALL abort('main loop', message)
333
334	END SELECT
335
336	IF (.TRUE.) THEN
337	ALLOCATE( output_arr(1,1,1:output_var % grid % nz) )
338	output_arr(1,1,:) = internal_arr(:)
339	END IF
340	CALL run_control('time', 'comp')
341
342
343	! This case gets called twice, the first time for ug, the
344	! second time for vg. We compute ug and vg at the first call
345	! and keep vg (and ug for that matter) around for the second
346	! call.
347	CASE ( 'geostrophic winds' )
348
349	IF (.NOT. ug_vg_have_been_computed ) THEN
350	ALLOCATE( ug_arr(output_var % grid % nz) )
351	ALLOCATE( vg_arr(output_var % grid % nz) )
352
353	IF ( cfg % ug_is_set ) THEN
354	ug_arr = cfg % ug
355	vg_arr = cfg % vg
356	ELSE
357	CALL geostrophic_winds( p_north, p_south, p_east, &
358	p_west, rho_centre, f3, &
359	averaging_width_ew, &
360	averaging_width_ns, &
361	phi_n, lambda_n, phi_centre, lam_centre, &
362	ug_arr, vg_arr )
363	END IF
364
365	ug_vg_have_been_computed = .TRUE.
366
367	END IF
368
369	! Prepare output of geostrophic winds
370	SELECT CASE(TRIM(output_var % name))
371	CASE ('ls_forcing_ug')
372	ug_vg_arr => ug_arr
373	CASE ('ls_forcing_vg')
374	ug_vg_arr => vg_arr
375	END SELECT
376
377	ALLOCATE( output_arr(1,1,output_var % grid % nz) )
378	output_arr(1,1,:) = ug_vg_arr(:)
379
380	CASE ( 'average scalar' )
381
382	ALLOCATE( output_arr(1,1,1) )
383	CALL run_control('time', 'alloc')
384	output_arr(1,1,1) = p0
385	CALL run_control('time', 'comp')
386
387	CASE ( 'set profile' )
388
389	ALLOCATE( output_arr( 1, 1, 1:nz ) )
390	CALL run_control('time', 'alloc')
391
392	SELECT CASE (TRIM(output_var % name))
393
394	!CASE('ls_forcing_ug')
395	! output_arr(1, 1, :) = ug
396
397	!CASE('ls_forcing_vg')
398	! output_arr(1, 1, :) = vg
399
400	CASE('nudging_tau')
401	output_arr(1, 1, :) = NUDGING_TAU
402
403	CASE DEFAULT
404	message = "'" // TRIM(output_var % name) // &
405	"' is not a valid '" // TRIM(output_var % kind) //&
406	"' variable kind."
407	CALL abort('main loop', message)
408	END SELECT
409	CALL run_control('time', 'comp')
410
411	CASE('average large-scale profile')
412	message = "Averaging of large-scale forcing profiles " //&
413	"has not been implemented, yet."
414	CALL abort('main loop', message)
415	!ALLOCATE( output_arr( 1, 1, 1:nz ) )
416
417	CASE DEFAULT
418	message = "Processing task '" // TRIM(output_var % task) //&
419	"' not recognized."
420	CALL abort('', message)
421
422	END SELECT
423	CALL run_control('time', 'comp')
424
425	!------------------------------------------------------------------------------
426	!- Section 2.3: Write current time step of current variable
427	!------------------------------------------------------------------------------
428	IF (.NOT. output_var % is_internal .OR. debugging_output) THEN
429	message = "Writing variable '" // TRIM(output_var%name) // "'."
430	CALL report('main loop', message)
431	CALL update_output(output_var, output_arr, iter, &
432	output_file, cfg)
433	CALL run_control('time', 'write')
434	END IF
435
436	IF (ALLOCATED(output_arr)) DEALLOCATE(output_arr)
437	CALL run_control('time', 'alloc')
438
439	END IF
440
441	END DO ! ouput variables
442
443	ug_vg_have_been_computed = .FALSE.
444	IF ( group % kind == 'thermodynamics' ) THEN
445	DEALLOCATE( rho_centre )
446	DEALLOCATE( ug_arr )
447	DEALLOCATE( vg_arr )
448	IF ( .NOT. cfg % ug_is_set ) THEN
449	DEALLOCATE( rho_north )
450	DEALLOCATE( rho_south )
451	DEALLOCATE( rho_east )
452	DEALLOCATE( rho_west )
453	DEALLOCATE( p_north )
454	DEALLOCATE( p_south )
455	DEALLOCATE( p_east )
456	DEALLOCATE( p_west )
457	END IF
458	END IF
459
460	IF ( group % kind == 'running average' .OR. &
461	group % kind == 'accumulated' ) THEN
462	! Keep input buffer around for averaged (radiation) and
463	! accumulated COSMO-DE quantities (precipitation).
464	ELSE
465	CALL report('main loop', 'Deallocating input buffer', cfg % debug)
466	DEALLOCATE(input_buffer)
467	END IF
468	CALL run_control('time', 'alloc')
469
470	END DO ! time steps / input files
471
472	IF (ALLOCATED(input_buffer)) THEN
473	CALL report('main loop', 'Deallocating input buffer', cfg % debug)
474	DEALLOCATE(input_buffer)
475	END IF
476	CALL run_control('time', 'alloc')
477
478	ELSE
479
480	message = "Skipping IO group " // TRIM(str(igroup)) // " '" // TRIM(group % kind) // "'"
481	IF ( ALLOCATED(group % in_var_list) ) THEN
482	message = TRIM(message) // " with input variable '" // &
483	TRIM(group % in_var_list(1) % name) // "'."
484	END IF
485
486	CALL report('main loop', message, cfg % debug)
487
488	END IF ! IO group % to_be_processed
489
490	END DO ! IO groups
491
492	!------------------------------------------------------------------------------
493	!- Section 3: Clean up.
494	!------------------------------------------------------------------------------
495	CALL fini_file_lists()
496	CALL fini_io_groups()
497	CALL fini_variables()
498	!CALL fini_grids()
499	CALL run_control('time', 'alloc')
500	CALL run_control('report', 'void')
501
502	message = "Finished writing dynamic driver '" // TRIM(output_file % name) // &
503	"' successfully."
504	CALL report('main loop', message)
505
506
507	END PROGRAM inifor

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |