Home

Context Navigation

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

Last change on this file since 3731 was 3680, checked in by knoop, 6 years ago
Added cpp-option netcdf to inifor
Property svn:keywords set to `Id`
File size: 22.3 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |