Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |