Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |