Home

Context Navigation

source: palm/trunk/UTIL/inifor/src/inifor_io.f90 @ 3764

Last change on this file since 3764 was 3764, checked in by eckhard, 5 years ago
inifor: bugfix: removed dependency on radiation input files
Property svn:keywords set to `Id`
File size: 49.9 KB

Rev	Line
[3447]	1	!> @file src/inifor_io.f90
[2696]	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_io.f90 3764 2019-02-26 13:42:09Z eckhard $
[3764]	28	! Removed dependency on radiation input files
	29	!
	30	!
	31	! 3716 2019-02-05 17:02:38Z eckhard
[3716]	32	! Removed dependency on soilmoisture input files
	33	!
	34	!
	35	! 3680 2019-01-18 14:54:12Z knoop
[3678]	36	! Moved get_input_file_list() here from grid module, added check for presence of
	37	! input files
	38	!
	39	!
	40	!
	41	!
	42	! 3618 2018-12-10 13:25:22Z eckhard
[3618]	43	! Prefixed all INIFOR modules with inifor_
	44	!
	45	!
	46	! 3615 2018-12-10 07:21:03Z raasch
[3615]	47	! bugfix: abort replaced by inifor_abort
	48	!
	49	! 3557 2018-11-22 16:01:22Z eckhard
[3557]	50	! Updated documentation, removed unused subroutine write_netcdf_variable_2d()
	51	!
	52	!
	53	! 3537 2018-11-20 10:53:14Z eckhard
[3537]	54	! New routine get_netcdf_dim_vector()
	55	!
	56	!
	57	! 3534 2018-11-19 15:35:16Z raasch
[3534]	58	! bugfix: INTENT attribute changed
	59	!
	60	! 3456 2018-10-30 14:29:54Z eckhard
[3456]	61	! NetCDf output of internal arrays only with --debug option
	62	!
	63	!
	64	! 3447 2018-10-29 15:52:54Z eckhard
[3447]	65	! Removed INCLUDE statement for get_netcdf_variable()
	66	! Renamed source files for compatibilty with PALM build system
	67	!
	68	!
	69	! 3395 2018-10-22 17:32:49Z eckhard
[3395]	70	! Added command-line options for configuring the computation of geostrophic
	71	! winds (--averagin-mode, --averaging-angle)
	72	! Added command-line option --input-prefix for setting input file prefixes all
	73	! at once
	74	! Added --debug option for more verbose terminal output
	75	! Added option-specific *_is_set LOGICALs to indicate invocation from the
	76	! command-line
	77	! Improved error messages in case of empty file-name strings
	78	! Improved routine naming
[3262]	79	!
	80	! 3183 2018-07-27 14:25:55Z suehring
[3182]	81	! Introduced new PALM grid stretching
	82	! Updated variable names and metadata for PIDS v1.9 compatibility
	83	! Improved handling of the start date string
	84	! Better compatibility with older Intel compilers:
	85	! - avoiding implicit array allocation with new get_netcdf_variable()
	86	! subroutine instead of function
	87	! Improved command line interface:
	88	! - Added configuration validation
	89	! - New options to configure input file prefixes
	90	! - GNU-style short and long option names
	91	! - Added version and copyright output
[2696]	92	!
[3182]	93	!
[3183]	94	! 3182 2018-07-27 13:36:03Z suehring
[2696]	95	! Initial revision
	96	!
	97	!
	98	!
	99	! Authors:
	100	! --------
[3557]	101	!> @author Eckhard Kadasch (Deutscher Wetterdienst, Offenbach)
[2696]	102	!
	103	! Description:
	104	! ------------
	105	!> The io module contains the functions needed to read and write netCDF data in
	106	!> INIFOR.
	107	!------------------------------------------------------------------------------!
[3680]	108	#if defined ( __netcdf )
[3618]	109	MODULE inifor_io
[2696]	110
[3618]	111	USE inifor_control
	112	USE inifor_defs, &
[3182]	113	ONLY: DATE, SNAME, PATH, PI, dp, hp, TO_RADIANS, TO_DEGREES, VERSION
[3618]	114	USE inifor_types
	115	USE inifor_util, &
[3678]	116	ONLY: add_hours_to, reverse, str, real_to_str
[2696]	117	USE netcdf
	118
	119	IMPLICIT NONE
	120
[3557]	121	!------------------------------------------------------------------------------!
	122	! Description:
	123	! ------------
	124	!> get_netcdf_variable() reads the netCDF data and metadate for the netCDF
	125	!> variable 'in_var % name' from the file 'in_file'. The netCDF data array is
	126	!> stored in the 'buffer' array and metadata added to the respective members of
	127	!> the given 'in_var'.
	128	!------------------------------------------------------------------------------!
[3182]	129	INTERFACE get_netcdf_variable
	130	MODULE PROCEDURE get_netcdf_variable_int
	131	MODULE PROCEDURE get_netcdf_variable_real
	132	END INTERFACE get_netcdf_variable
	133
	134	PRIVATE :: get_netcdf_variable_int, get_netcdf_variable_real
	135
[2696]	136	CONTAINS
	137
[3557]	138	!------------------------------------------------------------------------------!
	139	! Description:
	140	! ------------
	141	!> get_netcdf_variable_int() implements the integer variant for the
	142	!> get_netcdf_variable interface.
	143	!------------------------------------------------------------------------------!
[3182]	144	SUBROUTINE get_netcdf_variable_int(in_file, in_var, buffer)
	145
	146	CHARACTER(LEN=PATH), INTENT(IN) :: in_file
	147	TYPE(nc_var), INTENT(INOUT) :: in_var
	148	INTEGER(hp), ALLOCATABLE, INTENT(INOUT) :: buffer(:,:,:)
	149
[3447]	150	INTEGER :: ncid
	151	INTEGER, DIMENSION(3) :: start, count
[3182]	152
[3447]	153	IF ( nf90_open( TRIM(in_file), NF90_NOWRITE, ncid ) .EQ. NF90_NOERR .AND. &
	154	nf90_inq_varid( ncid, in_var % name, in_var % varid ) .EQ. NF90_NOERR ) THEN
	155
	156	CALL get_input_dimensions(in_var, ncid)
	157
	158	CALL get_netcdf_start_and_count(in_var, start, count)
	159	CALL run_control('time', 'read')
	160
	161	ALLOCATE( buffer( count(1), count(2), count(3) ) )
	162	CALL run_control('time', 'alloc')
	163
	164	CALL check(nf90_get_var( ncid, in_var % varid, buffer, &
	165	start = start, &
	166	count = count ))
	167
	168	ELSE
	169
	170	message = "Failed to read '" // TRIM(in_var % name) // &
	171	"' from file '" // TRIM(in_file) // "'."
[3615]	172	CALL inifor_abort('get_netcdf_variable', message)
[3447]	173
	174	END IF
	175
	176	CALL check(nf90_close(ncid))
	177	CALL run_control('time', 'read')
	178
[3182]	179	END SUBROUTINE get_netcdf_variable_int
	180
	181
[3557]	182	!------------------------------------------------------------------------------!
	183	! Description:
	184	! ------------
	185	!> get_netcdf_variable_real() implements the real variant for the
	186	!> get_netcdf_variable interface.
	187	!------------------------------------------------------------------------------!
[3182]	188	SUBROUTINE get_netcdf_variable_real(in_file, in_var, buffer)
	189
	190	CHARACTER(LEN=PATH), INTENT(IN) :: in_file
	191	TYPE(nc_var), INTENT(INOUT) :: in_var
	192	REAL(dp), ALLOCATABLE, INTENT(INOUT) :: buffer(:,:,:)
	193
[3447]	194	INTEGER :: ncid
	195	INTEGER, DIMENSION(3) :: start, count
[3182]	196
[3447]	197	IF ( nf90_open( TRIM(in_file), NF90_NOWRITE, ncid ) .EQ. NF90_NOERR .AND. &
	198	nf90_inq_varid( ncid, in_var % name, in_var % varid ) .EQ. NF90_NOERR ) THEN
	199
	200	CALL get_input_dimensions(in_var, ncid)
	201
	202	CALL get_netcdf_start_and_count(in_var, start, count)
	203	CALL run_control('time', 'read')
	204
	205	ALLOCATE( buffer( count(1), count(2), count(3) ) )
	206	CALL run_control('time', 'alloc')
	207
	208	CALL check(nf90_get_var( ncid, in_var % varid, buffer, &
	209	start = start, &
	210	count = count ))
	211
	212	ELSE
	213
	214	message = "Failed to read '" // TRIM(in_var % name) // &
	215	"' from file '" // TRIM(in_file) // "'."
[3615]	216	CALL inifor_abort('get_netcdf_variable', message)
[3447]	217
	218	END IF
	219
	220	CALL check(nf90_close(ncid))
	221	CALL run_control('time', 'read')
	222
[3182]	223	END SUBROUTINE get_netcdf_variable_real
	224
	225
[3557]	226	!------------------------------------------------------------------------------!
	227	! Description:
	228	! ------------
	229	!> get_netcdf_dim_vector() reads the coordinate array 'coordname' from the
	230	!> netCDF file 'filename'.
	231	!------------------------------------------------------------------------------!
	232	SUBROUTINE get_netcdf_dim_vector(filename, coordname, coords)
[3537]	233
	234	CHARACTER(LEN=*), INTENT(IN) :: filename
[3557]	235	CHARACTER(LEN=*), INTENT(IN) :: coordname
	236	REAL(dp), ALLOCATABLE, INTENT(INOUT) :: coords(:)
[3537]	237
	238	INTEGER :: ncid, varid, dimlen
	239	INTEGER :: dimids(NF90_MAX_VAR_DIMS)
	240
	241	IF ( nf90_open( TRIM(filename), NF90_NOWRITE, ncid ) .EQ. NF90_NOERR .AND. &
[3557]	242	nf90_inq_varid( ncid, coordname, varid ) .EQ. NF90_NOERR ) THEN
[3537]	243
	244	CALL check(nf90_inquire_variable( ncid, varid, dimids = dimids ))
	245	CALL check(nf90_inquire_dimension( ncid, dimids(1), len = dimlen ))
	246
[3557]	247	ALLOCATE(coords(dimlen))
	248	CALL check(nf90_get_var( ncid, varid, coords))
[3537]	249
	250	ELSE
	251
[3557]	252	message = "Failed to read '" // TRIM(coordname) // &
[3537]	253	"' from file '" // TRIM(filename) // "'."
[3615]	254	CALL inifor_abort('get_netcdf_dim_vector', message)
[3537]	255
	256	END IF
	257
	258	END SUBROUTINE get_netcdf_dim_vector
	259
	260
[3557]	261	!------------------------------------------------------------------------------!
	262	! Description:
	263	! ------------
	264	!> get_input_dimensions() reads dimensions metadata of the netCDF variable given
	265	!> by 'in_var % name' into 'in_var' data structure.
	266	!------------------------------------------------------------------------------!
[3447]	267	SUBROUTINE get_input_dimensions(in_var, ncid)
	268
[3534]	269	TYPE(nc_var), INTENT(INOUT) :: in_var
	270	INTEGER, INTENT(IN) :: ncid
[3447]	271
	272	INTEGER :: i
	273
	274	CALL check(nf90_get_att( ncid, in_var % varid, "long_name", &
	275	in_var % long_name))
	276
	277	CALL check(nf90_get_att( ncid, in_var % varid, "units", &
	278	in_var % units ))
	279
	280	CALL check(nf90_inquire_variable( ncid, in_var % varid, &
	281	ndims = in_var % ndim, &
	282	dimids = in_var % dimids ))
	283
	284	DO i = 1, in_var % ndim
	285	CALL check(nf90_inquire_dimension( ncid, in_var % dimids(i), &
	286	name = in_var % dimname(i), &
	287	len = in_var % dimlen(i) ))
	288	END DO
	289
	290	END SUBROUTINE get_input_dimensions
	291
	292
[3557]	293	!------------------------------------------------------------------------------!
	294	! Description:
	295	! ------------
	296	!> get_netcdf_start_and_count() gets the start position and element counts for
	297	!> the given netCDF file. This information is used in get_netcdf_variable_int()
	298	!> and _real() for reading input variables..
	299	!------------------------------------------------------------------------------!
[3447]	300	SUBROUTINE get_netcdf_start_and_count(in_var, start, count)
	301
	302	TYPE(nc_var), INTENT(INOUT) :: in_var
	303	INTEGER, DIMENSION(3), INTENT(OUT) :: start, count
	304
	305	INTEGER :: ndim
	306
	307	IF ( in_var % ndim .LT. 2 .OR. in_var % ndim .GT. 4 ) THEN
	308
	309	message = "Failed reading NetCDF variable " // &
[3537]	310	TRIM(in_var % name) // " with " // TRIM(str(in_var % ndim)) // &
[3447]	311	" dimensions because only two- and and three-dimensional" // &
	312	" variables are supported."
[3615]	313	CALL inifor_abort('get_netcdf_start_and_count', message)
[3447]	314
	315	END IF
	316
	317	start = (/ 1, 1, 1 /)
	318	IF ( TRIM(in_var % name) .EQ. 'T_SO' ) THEN
[3557]	319	!
	320	!-- Skip depth = 0.0 for T_SO and reduce number of depths from 9 to 8
[3447]	321	in_var % dimlen(3) = in_var % dimlen(3) - 1
	322
[3557]	323	!
	324	!-- Start reading from second level, e.g. depth = 0.005 instead of 0.0
[3447]	325	start(3) = 2
	326	END IF
	327
	328	IF (in_var % ndim .EQ. 2) THEN
	329	in_var % dimlen(3) = 1
	330	ENDIF
	331
	332	ndim = MIN(in_var % ndim, 3)
	333	count = (/ 1, 1, 1 /)
	334	count(1:ndim) = in_var % dimlen(1:ndim)
	335
	336	END SUBROUTINE get_netcdf_start_and_count
	337
	338
[3557]	339	!------------------------------------------------------------------------------!
	340	! Description:
	341	! ------------
	342	!> Routine for defining netCDF variables in the dynamic driver, INIFOR's netCDF
	343	!> output.
	344	!------------------------------------------------------------------------------!
[2696]	345	SUBROUTINE netcdf_define_variable(var, ncid)
	346
	347	TYPE(nc_var), INTENT(INOUT) :: var
	348	INTEGER, INTENT(IN) :: ncid
	349
	350	CALL check(nf90_def_var(ncid, var % name, NF90_FLOAT, var % dimids(1:var % ndim), var % varid))
	351	CALL check(nf90_put_att(ncid, var % varid, "long_name", var % long_name))
	352	CALL check(nf90_put_att(ncid, var % varid, "units", var % units))
[3182]	353	IF ( var % lod .GE. 0 ) THEN
	354	CALL check(nf90_put_att(ncid, var % varid, "lod", var % lod))
	355	END IF
[2696]	356	CALL check(nf90_put_att(ncid, var % varid, "source", var % source))
	357	CALL check(nf90_put_att(ncid, var % varid, "_FillValue", NF90_FILL_REAL))
	358
	359	END SUBROUTINE netcdf_define_variable
	360
	361
[3557]	362	!------------------------------------------------------------------------------!
	363	! Description:
	364	! ------------
	365	!> netcdf_get_dimensions() reads in all dimensions and their lengths and stores
	366	!> them in the given the 'var' data structure. This information is used later
	367	!> for writing output variables in update_output().
	368	!------------------------------------------------------------------------------!
[2696]	369	SUBROUTINE netcdf_get_dimensions(var, ncid)
	370
	371	TYPE(nc_var), INTENT(INOUT) :: var
	372	INTEGER, INTENT(IN) :: ncid
	373	INTEGER :: i
	374	CHARACTER(SNAME) :: null
	375
	376	DO i = 1, var % ndim
	377	CALL check(nf90_inquire_dimension(ncid, var % dimids(i), &
	378	name = null, &
	379	len = var % dimlen(i) ) )
	380	END DO
	381
	382	END SUBROUTINE netcdf_get_dimensions
	383
	384
	385	!------------------------------------------------------------------------------!
	386	! Description:
	387	! ------------
[3557]	388	!> This routine parses and interpretes the command-line options and stores the
	389	!> resulting settings in the 'cfg' data structure.
[2696]	390	!------------------------------------------------------------------------------!
[3182]	391	SUBROUTINE parse_command_line_arguments( cfg )
[2696]	392
[3182]	393	TYPE(inifor_config), INTENT(INOUT) :: cfg
[2696]	394
[3182]	395	CHARACTER(LEN=PATH) :: option, arg
	396	INTEGER :: arg_count, i
[2696]	397
[3395]	398	cfg % p0_is_set = .FALSE.
	399	cfg % ug_is_set = .FALSE.
	400	cfg % vg_is_set = .FALSE.
	401	cfg % flow_prefix_is_set = .FALSE.
	402	cfg % input_prefix_is_set = .FALSE.
	403	cfg % radiation_prefix_is_set = .FALSE.
	404	cfg % soil_prefix_is_set = .FALSE.
	405	cfg % soilmoisture_prefix_is_set = .FALSE.
	406
[2696]	407	arg_count = COMMAND_ARGUMENT_COUNT()
	408	IF (arg_count .GT. 0) THEN
	409
	410	message = "The -clon and -clat command line options are depricated. " // &
	411	"Please remove them form your inifor command and specify the " // &
	412	"location of the PALM-4U origin either" // NEW_LINE(' ') // &
[3182]	413	" - by setting the namelist parameters 'longitude' and 'latitude', or" // NEW_LINE(' ') // &
[2696]	414	" - by providing a static driver netCDF file via the -static command-line option."
	415
[3182]	416	i = 1
	417	DO WHILE (i .LE. arg_count)
[2696]	418
	419	CALL GET_COMMAND_ARGUMENT( i, option )
	420
	421	SELECT CASE( TRIM(option) )
	422
[3395]	423	CASE( '--averaging-mode' )
	424	CALL get_option_argument( i, arg )
	425	cfg % averaging_mode = TRIM(arg)
	426
[3182]	427	CASE( '-date', '-d', '--date' )
	428	CALL get_option_argument( i, arg )
	429	cfg % start_date = TRIM(arg)
[2696]	430
[3395]	431	CASE( '--debug' )
	432	cfg % debug = .TRUE.
	433
[3182]	434	CASE( '-z0', '-z', '--elevation' )
	435	CALL get_option_argument( i, arg )
	436	READ(arg, *) cfg % z0
[2696]	437
[3182]	438	CASE( '-p0', '-r', '--surface-pressure' )
[3395]	439	cfg % p0_is_set = .TRUE.
[3182]	440	CALL get_option_argument( i, arg )
	441	READ(arg, *) cfg % p0
[2696]	442
[3182]	443	CASE( '-ug', '-u', '--geostrophic-u' )
[3395]	444	cfg % ug_is_set = .TRUE.
[3182]	445	CALL get_option_argument( i, arg )
	446	READ(arg, *) cfg % ug
[2696]	447
[3182]	448	CASE( '-vg', '-v', '--geostrophic-v' )
[3395]	449	cfg % vg_is_set = .TRUE.
[3182]	450	CALL get_option_argument( i, arg )
	451	READ(arg, *) cfg % vg
[2696]	452
[3182]	453	CASE( '-clon', '-clat' )
[3615]	454	CALL inifor_abort('parse_command_line_arguments', message)
[2696]	455
[3182]	456	CASE( '-path', '-p', '--path' )
	457	CALL get_option_argument( i, arg )
	458	cfg % input_path = TRIM(arg)
[2696]	459
[3182]	460	CASE( '-hhl', '-l', '--hhl-file' )
	461	CALL get_option_argument( i, arg )
[3395]	462	cfg % hhl_file = TRIM(arg)
[2696]	463
[3395]	464	CASE( '--input-prefix')
	465	CALL get_option_argument( i, arg )
	466	cfg % input_prefix = TRIM(arg)
	467	cfg % input_prefix_is_set = .TRUE.
	468
	469	CASE( '-a', '--averaging-angle' )
	470	CALL get_option_argument( i, arg )
	471	READ(arg, *) cfg % averaging_angle
	472
[3182]	473	CASE( '-static', '-t', '--static-driver' )
	474	CALL get_option_argument( i, arg )
[3395]	475	cfg % static_driver_file = TRIM(arg)
[2696]	476
[3182]	477	CASE( '-soil', '-s', '--soil-file')
	478	CALL get_option_argument( i, arg )
[3395]	479	cfg % soiltyp_file = TRIM(arg)
[2696]	480
[3182]	481	CASE( '--flow-prefix')
	482	CALL get_option_argument( i, arg )
[3395]	483	cfg % flow_prefix = TRIM(arg)
	484	cfg % flow_prefix_is_set = .TRUE.
	485
[3182]	486	CASE( '--radiation-prefix')
	487	CALL get_option_argument( i, arg )
[3395]	488	cfg % radiation_prefix = TRIM(arg)
	489	cfg % radiation_prefix_is_set = .TRUE.
	490
[3182]	491	CASE( '--soil-prefix')
	492	CALL get_option_argument( i, arg )
[3395]	493	cfg % soil_prefix = TRIM(arg)
	494	cfg % soil_prefix_is_set = .TRUE.
	495
[3182]	496	CASE( '--soilmoisture-prefix')
	497	CALL get_option_argument( i, arg )
[3395]	498	cfg % soilmoisture_prefix = TRIM(arg)
	499	cfg % soilmoisture_prefix_is_set = .TRUE.
[2696]	500
[3182]	501	CASE( '-o', '--output' )
	502	CALL get_option_argument( i, arg )
	503	cfg % output_file = TRIM(arg)
[2696]	504
[3182]	505	CASE( '-n', '--namelist' )
	506	CALL get_option_argument( i, arg )
	507	cfg % namelist_file = TRIM(arg)
[2696]	508
[3182]	509	CASE( '-mode', '-i', '--init-mode' )
	510	CALL get_option_argument( i, arg )
	511	cfg % ic_mode = TRIM(arg)
	512
	513	CASE( '-f', '--forcing-mode' )
	514	CALL get_option_argument( i, arg )
	515	cfg % bc_mode = TRIM(arg)
	516
	517	CASE( '--version' )
	518	CALL print_version()
	519	STOP
	520
	521	CASE( '--help' )
	522	CALL print_version()
	523	PRINT *, ""
	524	PRINT *, "For a list of command-line options have a look at the README file."
	525	STOP
	526
[2696]	527	CASE DEFAULT
[3182]	528	message = "unknown option '" // TRIM(option) // "'."
[3615]	529	CALL inifor_abort('parse_command_line_arguments', message)
[2696]	530
	531	END SELECT
	532
[3182]	533	i = i + 1
	534
[2696]	535	END DO
	536
	537	ELSE
	538
	539	message = "No arguments present, using default input and output files"
	540	CALL report('parse_command_line_arguments', message)
	541
	542	END IF
	543
	544	END SUBROUTINE parse_command_line_arguments
	545
[3182]	546
[3678]	547
	548	SUBROUTINE get_datetime_file_list( start_date_string, start_hour, end_hour, &
	549	step_hour, input_path, prefix, suffix, &
	550	file_list )
	551
	552	CHARACTER (LEN=DATE), INTENT(IN) :: start_date_string
	553	CHARACTER (LEN=*), INTENT(IN) :: prefix, suffix, input_path
	554	INTEGER, INTENT(IN) :: start_hour, end_hour, step_hour
	555	CHARACTER(LEN=*), ALLOCATABLE, INTENT(INOUT) :: file_list(:)
	556
	557	INTEGER :: number_of_intervals, hour, i
	558	CHARACTER(LEN=DATE) :: date_string
	559
	560	number_of_intervals = CEILING( REAL(end_hour - start_hour) / step_hour )
	561	ALLOCATE( file_list(number_of_intervals + 1) )
	562
	563	DO i = 0, number_of_intervals
	564
	565	hour = start_hour + i * step_hour
	566	date_string = add_hours_to(start_date_string, hour)
	567
	568	file_list(i+1) = TRIM(input_path) // TRIM(prefix) // &
	569	TRIM(date_string) // TRIM(suffix) // '.nc'
	570
	571	END DO
	572
	573	END SUBROUTINE get_datetime_file_list
	574
[3557]	575	!------------------------------------------------------------------------------!
	576	! Description:
	577	! ------------
[3678]	578	!> Establish a list of files based on the given start and end times and file
	579	!> prefixes and suffixes.
	580	!------------------------------------------------------------------------------!
	581	SUBROUTINE get_input_file_list( start_date_string, start_hour, end_hour, &
	582	step_hour, input_path, prefix, suffix, &
[3716]	583	file_list, nocheck )
[3678]	584
	585	CHARACTER (LEN=DATE), INTENT(IN) :: start_date_string
	586	CHARACTER (LEN=*), INTENT(IN) :: prefix, suffix, input_path
	587	INTEGER, INTENT(IN) :: start_hour, end_hour, step_hour
	588	CHARACTER(LEN=*), ALLOCATABLE, INTENT(INOUT) :: file_list(:)
[3716]	589	LOGICAL, OPTIONAL, INTENT(IN) :: nocheck
[3678]	590
	591	INTEGER :: number_of_intervals, hour, i
	592	CHARACTER(LEN=DATE) :: date_string
	593	CHARACTER(LEN=PATH) :: file_name
[3716]	594	LOGICAL :: check_files
[3678]	595
	596	CALL get_datetime_file_list( start_date_string, start_hour, end_hour, &
	597	step_hour, input_path, prefix, suffix, &
	598	file_list )
	599
[3716]	600	check_files = .TRUE.
	601	IF ( PRESENT ( nocheck ) ) THEN
	602	IF ( nocheck ) check_files = .FALSE.
	603	END IF
[3678]	604
[3716]	605	IF ( check_files ) THEN
[3678]	606
[3716]	607	tip = "Please check if you specified the correct file prefix " // &
	608	"using the options --input-prefix, --flow-prefix, etc."
	609
	610	DO i = 1, SIZE(file_list)
	611	CALL verify_file(file_list(i), 'input', tip)
	612	END DO
	613
	614	END IF
	615
[3678]	616	END SUBROUTINE get_input_file_list
	617
	618
	619	!------------------------------------------------------------------------------!
	620	! Description:
	621	! ------------
	622	!> Abort INIFOR if the given file is not present.
	623	!------------------------------------------------------------------------------!
	624	SUBROUTINE verify_file(file_name, file_kind, tip)
	625
	626	CHARACTER(LEN=*), INTENT(IN) :: file_name, file_kind
	627	CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: tip
	628
	629	IF (.NOT. file_present(file_name)) THEN
	630
	631	IF (LEN(TRIM(file_name)) == 0) THEN
	632
	633	message = "No name was given for the " // TRIM(file_kind) // " file."
	634
	635	ELSE
	636
	637	message = "The " // TRIM(file_kind) // " file '" // &
	638	TRIM(file_name) // "' was not found."
	639
	640	IF (PRESENT(tip)) THEN
	641	message = TRIM(message) // " " // TRIM(tip)
	642	END IF
	643
	644	END IF
	645
	646	CALL inifor_abort('verify_file', message)
	647
	648	END IF
	649
[3764]	650	message = "Set up input file name '" // TRIM(file_name) // "'"
	651	CALL report('verify_file', message)
	652
[3678]	653	END SUBROUTINE verify_file
	654
	655
	656	!------------------------------------------------------------------------------!
	657	! Description:
	658	! ------------
[3557]	659	!> Get the argument of the i'th command line option, which is at the location
	660	!> i+1 of the argument list.
	661	!------------------------------------------------------------------------------!
[3182]	662	SUBROUTINE get_option_argument(i, arg)
	663	CHARACTER(LEN=PATH), INTENT(INOUT) :: arg
	664	INTEGER, INTENT(INOUT) :: i
[2696]	665
[3182]	666	i = i + 1
	667	CALL GET_COMMAND_ARGUMENT(i, arg)
	668
	669	END SUBROUTINE
	670
	671
[3557]	672	!------------------------------------------------------------------------------!
	673	! Description:
	674	! ------------
	675	!> Checks the INIFOR configuration 'cfg' for plausibility.
	676	!------------------------------------------------------------------------------!
[3182]	677	SUBROUTINE validate_config(cfg)
	678	TYPE(inifor_config), INTENT(IN) :: cfg
	679
[3678]	680	CALL verify_file(cfg % hhl_file, 'HHL')
	681	CALL verify_file(cfg % namelist_file, 'NAMELIST')
	682	CALL verify_file(cfg % soiltyp_file, 'SOILTYP')
[3182]	683
[3557]	684	!
	685	!-- Only check optional static driver file name, if it has been given.
[3182]	686	IF (TRIM(cfg % static_driver_file) .NE. '') THEN
[3678]	687	CALL verify_file(cfg % static_driver_file, 'static driver')
[3182]	688	END IF
	689
	690	SELECT CASE( TRIM(cfg % ic_mode) )
	691	CASE( 'profile', 'volume')
	692	CASE DEFAULT
	693	message = "Initialization mode '" // TRIM(cfg % ic_mode) //&
	694	"' is not supported. " //&
	695	"Please select either 'profile' or 'volume', " //&
	696	"or omit the -i/--init-mode/-mode option entirely, which corresponds "//&
	697	"to the latter."
[3615]	698	CALL inifor_abort( 'validate_config', message )
[3182]	699	END SELECT
	700
	701
	702	SELECT CASE( TRIM(cfg % bc_mode) )
	703	CASE( 'real', 'ideal')
	704	CASE DEFAULT
	705	message = "Forcing mode '" // TRIM(cfg % bc_mode) //&
	706	"' is not supported. " //&
	707	"Please select either 'real' or 'ideal', " //&
	708	"or omit the -f/--forcing-mode option entirely, which corresponds "//&
	709	"to the latter."
[3615]	710	CALL inifor_abort( 'validate_config', message )
[3182]	711	END SELECT
	712
[3395]	713	SELECT CASE( TRIM(cfg % averaging_mode) )
	714	CASE( 'level', 'height')
	715	CASE DEFAULT
	716	message = "Averaging mode '" // TRIM(cfg % averaging_mode) //&
	717	"' is not supported. " //&
	718	"Please select either 'height' or 'level', " //&
	719	"or omit the --averaging-mode option entirely, which corresponds "//&
	720	"to the latter."
[3615]	721	CALL inifor_abort( 'validate_config', message )
[3395]	722	END SELECT
[3182]	723
[3395]	724	IF ( cfg % ug_is_set .NEQV. cfg % vg_is_set ) THEN
	725	message = "You specified only one component of the geostrophic " // &
	726	"wind. Please specify either both or none."
[3615]	727	CALL inifor_abort( 'validate_config', message )
[3395]	728	END IF
	729
[3182]	730	END SUBROUTINE validate_config
	731
	732
[3557]	733	!------------------------------------------------------------------------------!
	734	! Description:
	735	! ------------
	736	!> Check whether the given file is present on the filesystem.
	737	!------------------------------------------------------------------------------!
[3678]	738	LOGICAL FUNCTION file_present(filename)
[3182]	739	CHARACTER(LEN=PATH), INTENT(IN) :: filename
	740
[3678]	741	INQUIRE(FILE=filename, EXIST=file_present)
[3182]	742
	743	END FUNCTION file_present
	744
	745
[2696]	746	!------------------------------------------------------------------------------!
	747	! Description:
	748	! ------------
[3557]	749	!> This routine initializes the dynamic driver file, i.e. INIFOR's netCDF output
	750	!> file.
[2696]	751	!>
	752	!> Besides writing metadata, such as global attributes, coordinates, variables,
[3557]	753	!> in the netCDF file, various netCDF IDs are saved for later, when INIFOR
[2696]	754	!> writes the actual data.
	755	!------------------------------------------------------------------------------!
[3182]	756	SUBROUTINE setup_netcdf_dimensions(output_file, palm_grid, &
	757	start_date_string, origin_lon, origin_lat)
[2696]	758
	759	TYPE(nc_file), INTENT(INOUT) :: output_file
	760	TYPE(grid_definition), INTENT(IN) :: palm_grid
[3182]	761	CHARACTER (LEN=DATE), INTENT(IN) :: start_date_string
	762	REAL(dp), INTENT(IN) :: origin_lon, origin_lat
[2696]	763
[3182]	764	CHARACTER (LEN=8) :: date_string
	765	CHARACTER (LEN=10) :: time_string
	766	CHARACTER (LEN=5) :: zone_string
	767	CHARACTER (LEN=SNAME) :: history_string
[2696]	768	INTEGER :: ncid, nx, ny, nz, nt, dimids(3), dimvarids(3)
	769	REAL(dp) :: z0
	770
[3182]	771	message = "Initializing PALM-4U dynamic driver file '" // &
	772	TRIM(output_file % name) // "' and setting up dimensions."
	773	CALL report('setup_netcdf_dimensions', message)
	774
[3557]	775	!
	776	!-- Create the netCDF file as in netCDF-4/HDF5 format if __netcdf4 preprocessor flag is given
[3182]	777	#if defined( __netcdf4 )
[2696]	778	CALL check(nf90_create(TRIM(output_file % name), OR(NF90_CLOBBER, NF90_HDF5), ncid))
[3182]	779	#else
	780	CALL check(nf90_create(TRIM(output_file % name), NF90_CLOBBER, ncid))
	781	#endif
[2696]	782
[3395]	783	!------------------------------------------------------------------------------
	784	!- Section 1: Define NetCDF dimensions and coordinates
	785	!------------------------------------------------------------------------------
	786	nt = SIZE(output_file % time)
	787	nx = palm_grid % nx
	788	ny = palm_grid % ny
	789	nz = palm_grid % nz
	790	z0 = palm_grid % z0
	791
	792
[2696]	793	!
	794	!------------------------------------------------------------------------------
[3395]	795	!- Section 2: Write global NetCDF attributes
[2696]	796	!------------------------------------------------------------------------------
[3182]	797	CALL date_and_time(DATE=date_string, TIME=time_string, ZONE=zone_string)
	798	history_string = &
	799	'Created on '// date_string // &
	800	' at ' // time_string(1:2) // ':' // time_string(3:4) // &
	801	' (UTC' // zone_string // ')'
	802
[2696]	803	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'title', 'PALM input file for scenario ...'))
	804	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'institution', 'Deutscher Wetterdienst, Offenbach'))
	805	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'author', 'Eckhard Kadasch, eckhard.kadasch@dwd.de'))
[3182]	806	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'history', TRIM(history_string)))
[2696]	807	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'references', '--'))
	808	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'comment', '--'))
[3182]	809	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'origin_lat', TRIM(real_to_str(origin_lat*TO_DEGREES, '(F18.13)'))))
	810	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'origin_lon', TRIM(real_to_str(origin_lon*TO_DEGREES, '(F18.13)'))))
[3557]	811	!
	812	!-- FIXME: This is the elevation relative to COSMO-DE/D2 sea level and does
	813	!-- FIXME: not necessarily comply with DHHN2016 (c.f. PALM Input Data
	814	!-- FIXME: Standard v1.9., origin_z)
[3395]	815	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'origin_z', TRIM(real_to_str(z0, '(F18.13)'))))
[2696]	816	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'inifor_version', TRIM(VERSION)))
	817	CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'palm_version', '--'))
	818
	819	!
	820	!
	821	!------------------------------------------------------------------------------
	822	!- Section 2a: Define dimensions for cell centers (scalars in soil and atmosph.)
	823	!------------------------------------------------------------------------------
[3557]	824	!
	825	!-- reset dimids first
	826	dimids = (/0, 0, 0/)
	827	CALL check( nf90_def_dim(ncid, "x", nx+1, dimids(1)) )
	828	CALL check( nf90_def_dim(ncid, "y", ny+1, dimids(2)) )
	829	CALL check( nf90_def_dim(ncid, "z", nz, dimids(3)) )
	830	!
	831	!-- save dimids for later
	832	output_file % dimids_scl = dimids
[2696]	833
[3557]	834	!
	835	!-- reset dimvarids first
	836	dimvarids = (/0, 0, 0/)
	837	CALL check(nf90_def_var(ncid, "x", NF90_FLOAT, dimids(1), dimvarids(1)))
	838	CALL check(nf90_put_att(ncid, dimvarids(1), "standard_name", "x coordinate of cell centers"))
	839	CALL check(nf90_put_att(ncid, dimvarids(1), "units", "m"))
[2696]	840
[3557]	841	CALL check(nf90_def_var(ncid, "y", NF90_FLOAT, dimids(2), dimvarids(2)))
	842	CALL check(nf90_put_att(ncid, dimvarids(2), "standard_name", "y coordinate of cell centers"))
	843	CALL check(nf90_put_att(ncid, dimvarids(2), "units", "m"))
[2696]	844
[3557]	845	CALL check(nf90_def_var(ncid, "z", NF90_FLOAT, dimids(3), dimvarids(3)))
	846	CALL check(nf90_put_att(ncid, dimvarids(3), "standard_name", "z coordinate of cell centers"))
	847	CALL check(nf90_put_att(ncid, dimvarids(3), "units", "m"))
	848	!
	849	!-- save dimvarids for later
	850	output_file % dimvarids_scl = dimvarids
[2696]	851
[3557]	852	!
	853	!-- overwrite third dimid with the one of depth
[3182]	854	CALL check(nf90_def_dim(ncid, "zsoil", SIZE(palm_grid % depths), dimids(3)) )
[3557]	855	!
	856	!-- save dimids for later
	857	output_file % dimids_soil = dimids
[2696]	858
[3557]	859	!
	860	!-- overwrite third dimvarid with the one of depth
[3182]	861	CALL check(nf90_def_var(ncid, "zsoil", NF90_FLOAT, output_file % dimids_soil(3), dimvarids(3)))
[2696]	862	CALL check(nf90_put_att(ncid, dimvarids(3), "standard_name", "depth_below_land"))
	863	CALL check(nf90_put_att(ncid, dimvarids(3), "positive", "down"))
	864	CALL check(nf90_put_att(ncid, dimvarids(3), "units", "m"))
	865	!
[3557]	866	!-- save dimvarids for later
	867	output_file % dimvarids_soil = dimvarids
	868	!
[2696]	869	!------------------------------------------------------------------------------
	870	!- Section 2b: Define dimensions for cell faces/velocities
	871	!------------------------------------------------------------------------------
[3557]	872	!
	873	!-- reset dimids first
	874	dimids = (/0, 0, 0/)
	875	CALL check(nf90_def_dim(ncid, "xu", nx, dimids(1)) )
	876	CALL check(nf90_def_dim(ncid, "yv", ny, dimids(2)) )
	877	CALL check(nf90_def_dim(ncid, "zw", nz-1, dimids(3)) )
	878	!
	879	!-- save dimids for later
	880	output_file % dimids_vel = dimids
[2696]	881
[3557]	882	!
	883	!-- reset dimvarids first
	884	dimvarids = (/0, 0, 0/)
	885	CALL check(nf90_def_var(ncid, "xu", NF90_FLOAT, dimids(1), dimvarids(1)))
	886	CALL check(nf90_put_att(ncid, dimvarids(1), "standard_name", "x coordinate of cell faces"))
	887	CALL check(nf90_put_att(ncid, dimvarids(1), "units", "m"))
[2696]	888
[3557]	889	CALL check(nf90_def_var(ncid, "yv", NF90_FLOAT, dimids(2), dimvarids(2)))
	890	CALL check(nf90_put_att(ncid, dimvarids(2), "standard_name", "y coordinate of cell faces"))
	891	CALL check(nf90_put_att(ncid, dimvarids(2), "units", "m"))
[2696]	892
[3557]	893	CALL check(nf90_def_var(ncid, "zw", NF90_FLOAT, dimids(3), dimvarids(3)))
	894	CALL check(nf90_put_att(ncid, dimvarids(3), "standard_name", "z coordinate of cell faces"))
	895	CALL check(nf90_put_att(ncid, dimvarids(3), "units", "m"))
	896	!
	897	!-- save dimvarids for later
	898	output_file % dimvarids_vel = dimvarids
[2696]	899
	900	!
	901	!------------------------------------------------------------------------------
	902	!- Section 2c: Define time dimension
	903	!------------------------------------------------------------------------------
	904	CALL check(nf90_def_dim(ncid, "time", nt, output_file % dimid_time) )
	905	CALL check(nf90_def_var(ncid, "time", NF90_FLOAT, &
	906	output_file % dimid_time, &
	907	output_file % dimvarid_time))
	908	CALL check(nf90_put_att(ncid, output_file % dimvarid_time, "standard_name", "time"))
	909	CALL check(nf90_put_att(ncid, output_file % dimvarid_time, "long_name", "time"))
[3182]	910	CALL check(nf90_put_att(ncid, output_file % dimvarid_time, "units", &
	911	"seconds since " // start_date_string // " UTC"))
[2696]	912
	913	CALL check(nf90_enddef(ncid))
	914
	915	!
	916	!------------------------------------------------------------------------------
	917	!- Section 3: Write grid coordinates
	918	!------------------------------------------------------------------------------
	919	CALL check(nf90_put_var(ncid, output_file % dimvarids_scl(1), palm_grid%x))
	920	CALL check(nf90_put_var(ncid, output_file % dimvarids_scl(2), palm_grid%y))
	921	CALL check(nf90_put_var(ncid, output_file % dimvarids_scl(3), palm_grid%z))
	922
	923	CALL check(nf90_put_var(ncid, output_file % dimvarids_vel(1), palm_grid%xu))
	924	CALL check(nf90_put_var(ncid, output_file % dimvarids_vel(2), palm_grid%yv))
	925	CALL check(nf90_put_var(ncid, output_file % dimvarids_vel(3), palm_grid%zw))
	926
[3557]	927	!
	928	!-- TODO Read in soil depths from input file before this.
[2696]	929	CALL check(nf90_put_var(ncid, output_file % dimvarids_soil(3), palm_grid%depths))
	930
[3557]	931	!
	932	!-- Write time vector
[2696]	933	CALL check(nf90_put_var(ncid, output_file % dimvarid_time, output_file % time))
	934
[3557]	935	!
	936	!-- Close the file
[2696]	937	CALL check(nf90_close(ncid))
	938
	939	END SUBROUTINE setup_netcdf_dimensions
	940
	941
[3557]	942	!------------------------------------------------------------------------------!
	943	! Description:
	944	! ------------
	945	!> Defines the netCDF variables to be written to the dynamic driver file
	946	!------------------------------------------------------------------------------!
[3456]	947	SUBROUTINE setup_netcdf_variables(filename, output_variable_table, debug)
[2696]	948
	949	CHARACTER (LEN=*), INTENT(IN) :: filename
	950	TYPE(nc_var), INTENT(INOUT), TARGET :: output_variable_table(:)
[3456]	951	LOGICAL, INTENT(IN) :: debug
	952
[2696]	953	TYPE(nc_var), POINTER :: var
	954	INTEGER :: i, ncid
[3456]	955	LOGICAL :: to_be_written
[2696]	956
[3182]	957	message = "Defining variables in dynamic driver '" // TRIM(filename) // "'."
[2696]	958	CALL report('setup_netcdf_variables', message)
	959
	960	CALL check(nf90_open(TRIM(filename), NF90_WRITE, ncid))
	961	CALL check(nf90_redef(ncid))
	962
	963	DO i = 1, SIZE(output_variable_table)
	964
	965	var => output_variable_table(i)
	966
[3456]	967	to_be_written = ( var % to_be_processed .AND. .NOT. var % is_internal) .OR. &
	968	( var % is_internal .AND. debug )
	969
	970	IF ( to_be_written ) THEN
[3182]	971	message = " variable #" // TRIM(str(i)) // " '" // TRIM(var%name) // "'."
[2696]	972	CALL report('setup_netcdf_variables', message)
	973
	974	CALL netcdf_define_variable(var, ncid)
	975	CALL netcdf_get_dimensions(var, ncid)
	976	END IF
	977
	978	END DO
	979
	980	CALL check(nf90_enddef(ncid))
	981	CALL check(nf90_close(ncid))
	982
[3182]	983	message = "Dynamic driver '" // TRIM(filename) // "' initialized successfully."
[2696]	984	CALL report('setup_netcdf_variables', message)
	985
	986	END SUBROUTINE setup_netcdf_variables
	987
	988
	989	!------------------------------------------------------------------------------!
	990	! Description:
	991	! ------------
	992	!> This routine reads and returns all input variables of the given IO group
	993	!> It accomodates the data by allocating a container variable that represents a
	994	!> list of arrays of the same length as the groups variable list. (This list
	995	!> will typically contain one or two items.) After the container, its members
	996	!> are allocated one by one with the appropriate, possibly different,
	997	!> dimensions.
	998	!>
	999	!> The 'group' is an INTENT(INOUT) variable so that 'get_netcdf_variable()' can
	1000	!> record netCDF IDs in the 'in_var_list()' member variable.
	1001	!------------------------------------------------------------------------------!
	1002	SUBROUTINE read_input_variables(group, iter, buffer)
	1003	TYPE(io_group), INTENT(INOUT), TARGET :: group
	1004	INTEGER, INTENT(IN) :: iter
	1005	TYPE(container), ALLOCATABLE, INTENT(INOUT) :: buffer(:)
	1006	INTEGER :: hour, buf_id
	1007	TYPE(nc_var), POINTER :: input_var
	1008	CHARACTER(LEN=PATH), POINTER :: input_file
	1009	INTEGER :: ivar, nbuffers
	1010
	1011	message = "Reading data for I/O group '" // TRIM(group % in_var_list(1) % name) // "'."
	1012	CALL report('read_input_variables', message)
	1013
	1014	input_file => group % in_files(iter)
	1015
	1016	!
	1017	!------------------------------------------------------------------------------
	1018	!- Section 1: Load input buffers for accumulated variables
	1019	!------------------------------------------------------------------------------
[3557]	1020	!
	1021	!-- radiation budgets, precipitation
[2696]	1022	IF (group % kind == 'running average' .OR. &
[3557]	1023	group % kind == 'accumulated') THEN
[2696]	1024
	1025	IF (SIZE(group % in_var_list) .GT. 1 ) THEN
	1026	message = "I/O groups may not contain more than one " // &
	1027	"accumulated variable. Group '" // TRIM(group % kind) //&
	1028	"' contains " // &
	1029	TRIM( str(SIZE(group % in_var_list)) ) // "."
[3615]	1030	CALL inifor_abort('read_input_variables \| accumulation', message)
[2696]	1031	END IF
	1032
[3557]	1033	!
	1034	!-- use two buffer arrays
[2696]	1035	nbuffers = 2
	1036	IF ( .NOT. ALLOCATED( buffer ) ) ALLOCATE( buffer(nbuffers) )
	1037
[3557]	1038	!
	1039	!-- hour of the day
	1040	hour = iter - 1
	1041	!
	1042	!-- chose correct buffer array
[2696]	1043	buf_id = select_buffer(hour)
	1044
	1045	CALL run_control('time', 'read')
	1046	IF ( ALLOCATED(buffer(buf_id) % array) ) DEALLOCATE(buffer(buf_id) % array)
	1047	CALL run_control('time', 'alloc')
	1048
	1049	input_var => group % in_var_list(1)
[3182]	1050	CALL get_netcdf_variable(input_file, input_var, buffer(buf_id) % array)
[2696]	1051	CALL report('read_input_variables', "Read accumulated " // TRIM(group % in_var_list(1) % name))
	1052
	1053	IF ( input_var % is_upside_down ) CALL reverse(buffer(buf_id) % array)
	1054	CALL run_control('time', 'comp')
	1055
	1056	!------------------------------------------------------------------------------
	1057	!- Section 2: Load input buffers for normal I/O groups
	1058	!------------------------------------------------------------------------------
	1059	ELSE
	1060
[3557]	1061	!
	1062	!-- Allocate one input buffer per input_variable. If more quantities
	1063	!-- have to be computed than input variables exist in this group,
	1064	!-- allocate more buffers. For instance, in the thermodynamics group,
	1065	!-- there are three input variabels (PP, T, Qv) and four quantities
	1066	!-- necessart (P, Theta, Rho, qv) for the corresponding output fields
	1067	!-- (p0, Theta, qv, ug, and vg)
[3395]	1068	nbuffers = MAX( group % n_inputs, group % n_output_quantities )
[2696]	1069	ALLOCATE( buffer(nbuffers) )
	1070	CALL run_control('time', 'alloc')
	1071
[3557]	1072	!
	1073	!-- Read in all input variables, leave extra buffers-if any-untouched.
[3395]	1074	DO ivar = 1, group % n_inputs
[2696]	1075
	1076	input_var => group % in_var_list(ivar)
	1077
[3557]	1078	!
	1079	! Check wheather P or PP is present in input file
[2696]	1080	IF (input_var % name == 'P') THEN
[3395]	1081	input_var % name = TRIM( get_pressure_varname(input_file) )
[2696]	1082	CALL run_control('time', 'read')
	1083	END IF
	1084
[3182]	1085	CALL get_netcdf_variable(input_file, input_var, buffer(ivar) % array)
[2696]	1086
	1087	IF ( input_var % is_upside_down ) CALL reverse(buffer(ivar) % array)
	1088	CALL run_control('time', 'comp')
	1089
	1090	END DO
	1091	END IF
	1092
	1093	END SUBROUTINE read_input_variables
	1094
	1095
[3557]	1096	!------------------------------------------------------------------------------!
	1097	! Description:
	1098	! ------------
	1099	!> Select the appropriate buffer ID for accumulated COSMO input variables
	1100	!> depending on the current hour.
	1101	!------------------------------------------------------------------------------!
[2696]	1102	INTEGER FUNCTION select_buffer(hour)
	1103	INTEGER, INTENT(IN) :: hour
	1104	INTEGER :: step
	1105
	1106	select_buffer = 0
	1107	step = MODULO(hour, 3) + 1
	1108
	1109	SELECT CASE(step)
	1110	CASE(1, 3)
	1111	select_buffer = 1
	1112	CASE(2)
	1113	select_buffer = 2
	1114	CASE DEFAULT
	1115	message = "Invalid step '" // TRIM(str(step))
[3615]	1116	CALL inifor_abort('select_buffer', message)
[2696]	1117	END SELECT
	1118	END FUNCTION select_buffer
	1119
	1120
	1121	!------------------------------------------------------------------------------!
	1122	! Description:
	1123	! ------------
	1124	!> Checks if the input_file contains the absolute pressure, 'P', or the pressure
	1125	!> perturbation, 'PP', and returns the appropriate string.
	1126	!------------------------------------------------------------------------------!
[3395]	1127	CHARACTER(LEN=2) FUNCTION get_pressure_varname(input_file) RESULT(var)
[2696]	1128	CHARACTER(LEN=*) :: input_file
	1129	INTEGER :: ncid, varid
	1130
	1131	CALL check(nf90_open( TRIM(input_file), NF90_NOWRITE, ncid ))
	1132	IF ( nf90_inq_varid( ncid, 'P', varid ) .EQ. NF90_NOERR ) THEN
	1133
	1134	var = 'P'
	1135
	1136	ELSE IF ( nf90_inq_varid( ncid, 'PP', varid ) .EQ. NF90_NOERR ) THEN
	1137
	1138	var = 'PP'
	1139	CALL report('get_pressure_var', 'Using PP instead of P')
	1140
	1141	ELSE
	1142
	1143	message = "Failed to read '" // TRIM(var) // &
	1144	"' from file '" // TRIM(input_file) // "'."
[3615]	1145	CALL inifor_abort('get_pressure_var', message)
[2696]	1146
	1147	END IF
	1148
	1149	CALL check(nf90_close(ncid))
	1150
[3395]	1151	END FUNCTION get_pressure_varname
[2696]	1152
	1153
[3557]	1154	!------------------------------------------------------------------------------!
	1155	! Description:
	1156	! ------------
	1157	!> Read the given global attribute form the given netCDF file.
	1158	!------------------------------------------------------------------------------!
[2696]	1159	FUNCTION get_netcdf_attribute(filename, attribute) RESULT(attribute_value)
	1160
	1161	CHARACTER(LEN=*), INTENT(IN) :: filename, attribute
	1162	REAL(dp) :: attribute_value
	1163
[3557]	1164	INTEGER :: ncid
[2696]	1165
	1166	IF ( nf90_open( TRIM(filename), NF90_NOWRITE, ncid ) == NF90_NOERR ) THEN
	1167
	1168	CALL check(nf90_get_att(ncid, NF90_GLOBAL, TRIM(attribute), attribute_value))
[3182]	1169	CALL check(nf90_close(ncid))
[2696]	1170
	1171	ELSE
	1172
	1173	message = "Failed to read '" // TRIM(attribute) // &
	1174	"' from file '" // TRIM(filename) // "'."
[3615]	1175	CALL inifor_abort('get_netcdf_attribute', message)
[2696]	1176
	1177	END IF
	1178
	1179	END FUNCTION get_netcdf_attribute
	1180
	1181
[3557]	1182
	1183	!------------------------------------------------------------------------------!
	1184	! Description:
	1185	! ------------
	1186	!> Updates the dynamic driver with the interpolated field of the current
	1187	!> variable at the current time step.
	1188	!------------------------------------------------------------------------------!
[3456]	1189	SUBROUTINE update_output(var, array, iter, output_file, cfg)
[2696]	1190	TYPE(nc_var), INTENT(IN) :: var
	1191	REAL(dp), INTENT(IN) :: array(:,:,:)
	1192	INTEGER, INTENT(IN) :: iter
	1193	TYPE(nc_file), INTENT(IN) :: output_file
[3456]	1194	TYPE(inifor_config) :: cfg
[2696]	1195
	1196	INTEGER :: ncid, ndim, start(4), count(4)
	1197	LOGICAL :: var_is_time_dependent
	1198
	1199	var_is_time_dependent = ( &
	1200	var % dimids( var % ndim ) == output_file % dimid_time &
	1201	)
	1202
[3557]	1203	!
	1204	!-- Skip time dimension for output
[3182]	1205	ndim = var % ndim
	1206	IF ( var_is_time_dependent ) ndim = var % ndim - 1
[2696]	1207
	1208	start(:) = (/1,1,1,1/)
	1209	start(ndim+1) = iter
	1210	count(1:ndim) = var%dimlen(1:ndim)
	1211
	1212	CALL check(nf90_open(output_file % name, NF90_WRITE, ncid))
	1213
[3557]	1214	!
	1215	!-- Reduce dimension of output array according to variable kind
[2696]	1216	SELECT CASE (TRIM(var % kind))
	1217
	1218	CASE ( 'init scalar profile', 'init u profile', 'init v profile', &
	1219	'init w profile' )
	1220
	1221	CALL check(nf90_put_var( ncid, var%varid, array(1,1,:) ) )
	1222
	1223	CASE ( 'init soil', 'init scalar', 'init u', 'init v', 'init w' )
	1224
	1225	CALL check(nf90_put_var( ncid, var%varid, array(:,:,:) ) )
	1226
	1227	CASE( 'left scalar', 'right scalar', 'left w', 'right w' )
	1228
	1229	CALL check(nf90_put_var( ncid, var%varid, array(1,:,:), &
	1230	start=start(1:ndim+1), &
	1231	count=count(1:ndim) ) )
	1232
	1233
	1234	IF (.NOT. SIZE(array, 2) .EQ. var % dimlen(1)) THEN
	1235	PRINT *, "inifor: update_output: Dimension ", 1, " of variable ", &
	1236	TRIM(var % name), " (", var % dimlen(1), &
	1237	") does not match the dimension of the output array (", &
	1238	SIZE(array, 2), ")."
	1239	STOP
	1240	END IF
	1241
	1242
	1243	CASE( 'north scalar', 'south scalar', 'north w', 'south w' )
	1244
	1245	CALL check(nf90_put_var( ncid, var%varid, array(:,1,:), &
	1246	start=start(1:ndim+1), &
	1247	count=count(1:ndim) ) )
	1248
	1249
	1250	CASE( 'surface forcing', 'top scalar', 'top w' )
	1251
	1252	CALL check(nf90_put_var( ncid, var%varid, array(:,:,1), &
	1253	start=start(1:ndim+1), &
	1254	count=count(1:ndim) ) )
	1255
	1256	CASE ( 'left u', 'right u', 'left v', 'right v' )
	1257
	1258	CALL check(nf90_put_var( ncid, var%varid, array(1,:,:), &
	1259	start=start(1:ndim+1), &
	1260	count=count(1:ndim) ) )
	1261
	1262	CASE ( 'north u', 'south u', 'north v', 'south v' )
	1263
	1264	CALL check(nf90_put_var( ncid, var%varid, array(:,1,:), &
	1265	start=start(1:ndim+1), &
	1266	count=count(1:ndim) ) )
	1267
	1268	CASE ( 'top u', 'top v' )
	1269
	1270	CALL check(nf90_put_var( ncid, var%varid, array(:,:,1), &
	1271	start=start(1:ndim+1), &
	1272	count=count(1:ndim) ) )
	1273
	1274	CASE ( 'time series' )
	1275
	1276	CALL check(nf90_put_var( ncid, var%varid, array(1,1,1), &
	1277	start=start(1:ndim+1) ) )
	1278
[3456]	1279	CASE ( 'constant scalar profile', 'geostrophic' )
[2696]	1280
	1281	CALL check(nf90_put_var( ncid, var%varid, array(1,1,:), &
	1282	start=start(1:ndim+1), &
	1283	count=count(1:ndim) ) )
	1284
[3456]	1285	CASE ( 'internal profile' )
	1286
	1287	IF ( cfg % debug ) THEN
	1288	CALL check(nf90_put_var( ncid, var%varid, array(1,1,:), &
	1289	start=start(1:ndim+1), &
	1290	count=count(1:ndim) ) )
	1291	END IF
	1292
[3182]	1293	CASE ( 'large-scale scalar forcing', 'large-scale w forcing' )
	1294
	1295	message = "Doing nothing in terms of writing large-scale forings."
	1296	CALL report('update_output', message)
	1297
[2696]	1298	CASE DEFAULT
	1299
	1300	message = "Variable kind '" // TRIM(var % kind) // &
	1301	"' not recognized."
[3615]	1302	CALL inifor_abort('update_output', message)
[2696]	1303
	1304	END SELECT
	1305
	1306	CALL check(nf90_close(ncid))
	1307
	1308	END SUBROUTINE update_output
	1309
	1310
[3557]	1311	!------------------------------------------------------------------------------!
	1312	! Description:
	1313	! ------------
	1314	!> Checks the status of a netCDF API call and aborts if an error occured
	1315	!------------------------------------------------------------------------------!
[2696]	1316	SUBROUTINE check(status)
	1317
	1318	INTEGER, INTENT(IN) :: status
	1319
	1320	IF (status /= nf90_noerr) THEN
	1321	message = "NetCDF API call failed with error: " // &
	1322	TRIM( nf90_strerror(status) )
[3615]	1323	CALL inifor_abort('io.check', message)
[2696]	1324	END IF
	1325
	1326	END SUBROUTINE check
	1327
[3618]	1328	END MODULE inifor_io
[3680]	1329	#endif

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |