Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |