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

Last change on this file since 3680 was 3680, checked in by knoop, 3 years ago

Added cpp-option netcdf to inifor

  • Property svn:keywords set to Id
File size: 49.5 KB
Line 
1!> @file src/inifor_io.f90
2!------------------------------------------------------------------------------!
3! This file is part of the PALM model system.
4!
5! PALM is free software: you can redistribute it and/or modify it under the
6! terms of the GNU General Public License as published by the Free Software
7! Foundation, either version 3 of the License, or (at your option) any later
8! version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 2017-2018 Leibniz Universitaet Hannover
18! Copyright 2017-2018 Deutscher Wetterdienst Offenbach
19!------------------------------------------------------------------------------!
20!
21! Current revisions:
22! -----------------
23!
24!
25! Former revisions:
26! -----------------
27! $Id: inifor_io.f90 3680 2019-01-18 14:54:12Z knoop $
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
35! Prefixed all INIFOR modules with inifor_
36!
37!
38! 3615 2018-12-10 07:21:03Z raasch
39! bugfix: abort replaced by inifor_abort
40!
41! 3557 2018-11-22 16:01:22Z eckhard
42! Updated documentation, removed unused subroutine write_netcdf_variable_2d()
43!
44!
45! 3537 2018-11-20 10:53:14Z eckhard
46! New routine get_netcdf_dim_vector()
47!
48!
49! 3534 2018-11-19 15:35:16Z raasch
50! bugfix: INTENT attribute changed
51!
52! 3456 2018-10-30 14:29:54Z eckhard
53! NetCDf output of internal arrays only with --debug option
54!
55!
56! 3447 2018-10-29 15:52:54Z eckhard
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
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
71!
72! 3183 2018-07-27 14:25:55Z suehring
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
84!
85!
86! 3182 2018-07-27 13:36:03Z suehring
87! Initial revision
88!
89!
90!
91! Authors:
92! --------
93!> @author Eckhard Kadasch (Deutscher Wetterdienst, Offenbach)
94!
95! Description:
96! ------------
97!> The io module contains the functions needed to read and write netCDF data in
98!> INIFOR.
99!------------------------------------------------------------------------------!
100#if defined ( __netcdf )
101 MODULE inifor_io
102
103    USE inifor_control
104    USE inifor_defs,                                                           &
105        ONLY:  DATE, SNAME, PATH, PI, dp, hp, TO_RADIANS, TO_DEGREES, VERSION
106    USE inifor_types
107    USE inifor_util,                                                           &
108        ONLY:  add_hours_to, reverse, str, real_to_str
109    USE netcdf
110
111    IMPLICIT NONE
112
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!------------------------------------------------------------------------------!
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
128 CONTAINS
129
130!------------------------------------------------------------------------------!
131! Description:
132! ------------
133!> get_netcdf_variable_int() implements the integer variant for the
134!> get_netcdf_variable interface.
135!------------------------------------------------------------------------------!
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
142       INTEGER               ::  ncid
143       INTEGER, DIMENSION(3) ::  start, count
144
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) // "'."
164          CALL inifor_abort('get_netcdf_variable', message)
165
166       END IF
167
168       CALL check(nf90_close(ncid))
169 CALL run_control('time', 'read')
170
171    END SUBROUTINE get_netcdf_variable_int
172
173
174!------------------------------------------------------------------------------!
175! Description:
176! ------------
177!> get_netcdf_variable_real() implements the real variant for the
178!> get_netcdf_variable interface.
179!------------------------------------------------------------------------------!
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
186       INTEGER               ::  ncid
187       INTEGER, DIMENSION(3) ::  start, count
188
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) // "'."
208          CALL inifor_abort('get_netcdf_variable', message)
209
210       END IF
211
212       CALL check(nf90_close(ncid))
213 CALL run_control('time', 'read')
214
215    END SUBROUTINE get_netcdf_variable_real
216
217
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)
225
226       CHARACTER(LEN=*), INTENT(IN)         ::  filename
227       CHARACTER(LEN=*), INTENT(IN)         ::  coordname
228       REAL(dp), ALLOCATABLE, INTENT(INOUT) ::  coords(:)
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. &
234            nf90_inq_varid( ncid, coordname, varid ) .EQ. NF90_NOERR )  THEN
235
236          CALL check(nf90_inquire_variable( ncid, varid, dimids = dimids ))
237          CALL check(nf90_inquire_dimension( ncid, dimids(1), len = dimlen ))
238
239          ALLOCATE(coords(dimlen))
240          CALL check(nf90_get_var( ncid, varid, coords))
241
242       ELSE
243
244          message = "Failed to read '" // TRIM(coordname) // &
245             "' from file '" // TRIM(filename) // "'."
246          CALL inifor_abort('get_netcdf_dim_vector', message)
247
248       END IF
249
250    END SUBROUTINE get_netcdf_dim_vector
251
252
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!------------------------------------------------------------------------------!
259    SUBROUTINE get_input_dimensions(in_var, ncid)
260
261       TYPE(nc_var), INTENT(INOUT) ::  in_var
262       INTEGER, INTENT(IN)         ::  ncid
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
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!------------------------------------------------------------------------------!
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 " //                       &
302             TRIM(in_var % name) // " with " // TRIM(str(in_var % ndim)) //    &
303             " dimensions because only two- and and three-dimensional" //      &
304             " variables are supported."
305          CALL inifor_abort('get_netcdf_start_and_count', message)
306
307       END IF
308
309       start = (/ 1, 1, 1 /)
310       IF ( TRIM(in_var % name) .EQ. 'T_SO' )  THEN
311!
312!--       Skip depth = 0.0 for T_SO and reduce number of depths from 9 to 8
313          in_var % dimlen(3) = in_var % dimlen(3) - 1
314
315!
316!--       Start reading from second level, e.g. depth = 0.005 instead of 0.0
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
331!------------------------------------------------------------------------------!
332! Description:
333! ------------
334!> Routine for defining netCDF variables in the dynamic driver, INIFOR's netCDF
335!> output.
336!------------------------------------------------------------------------------!
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))
345        IF ( var % lod .GE. 0 )  THEN
346           CALL check(nf90_put_att(ncid, var % varid, "lod",           var % lod))
347        END IF
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
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!------------------------------------------------------------------------------!
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! ------------
380!> This routine parses and interpretes the command-line options and stores the
381!> resulting settings in the 'cfg' data structure.
382!------------------------------------------------------------------------------!
383    SUBROUTINE parse_command_line_arguments( cfg )
384
385       TYPE(inifor_config), INTENT(INOUT) ::  cfg
386
387       CHARACTER(LEN=PATH)                ::  option, arg
388       INTEGER                            ::  arg_count, i
389
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
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(' ') // &
405             "   - by setting the namelist parameters 'longitude' and 'latitude', or" // NEW_LINE(' ') // &
406             "   - by providing a static driver netCDF file via the -static command-line option."
407
408          i = 1
409          DO WHILE (i .LE. arg_count)
410
411             CALL GET_COMMAND_ARGUMENT( i, option )
412
413             SELECT CASE( TRIM(option) )
414
415             CASE( '--averaging-mode' )
416                CALL get_option_argument( i, arg )
417                cfg % averaging_mode = TRIM(arg)
418
419             CASE( '-date', '-d', '--date' )
420                CALL get_option_argument( i, arg )
421                cfg % start_date = TRIM(arg)
422
423             CASE( '--debug' )
424                cfg % debug = .TRUE.
425
426             CASE( '-z0', '-z', '--elevation' )
427                CALL get_option_argument( i, arg )
428                READ(arg, *) cfg % z0
429
430             CASE( '-p0', '-r', '--surface-pressure' )
431                cfg % p0_is_set = .TRUE.
432                CALL get_option_argument( i, arg )
433                READ(arg, *) cfg % p0
434
435             CASE( '-ug', '-u', '--geostrophic-u' )
436                cfg % ug_is_set = .TRUE.
437                CALL get_option_argument( i, arg )
438                READ(arg, *) cfg % ug
439
440             CASE( '-vg', '-v', '--geostrophic-v' )
441                cfg % vg_is_set = .TRUE.
442                CALL get_option_argument( i, arg )
443                READ(arg, *) cfg % vg
444
445             CASE( '-clon', '-clat' )
446                CALL inifor_abort('parse_command_line_arguments', message)
447
448             CASE( '-path', '-p', '--path' )
449                CALL get_option_argument( i, arg )
450                 cfg % input_path = TRIM(arg)
451
452             CASE( '-hhl', '-l', '--hhl-file' )
453                CALL get_option_argument( i, arg )
454                cfg % hhl_file = TRIM(arg)
455
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
465             CASE( '-static', '-t', '--static-driver' )
466                CALL get_option_argument( i, arg )
467                cfg % static_driver_file = TRIM(arg)
468
469             CASE( '-soil', '-s', '--soil-file')
470                CALL get_option_argument( i, arg )
471                cfg % soiltyp_file = TRIM(arg)
472
473             CASE( '--flow-prefix')
474                CALL get_option_argument( i, arg )
475                cfg % flow_prefix = TRIM(arg)
476                cfg % flow_prefix_is_set = .TRUE.
477   
478             CASE( '--radiation-prefix')
479                CALL get_option_argument( i, arg )
480                cfg % radiation_prefix = TRIM(arg)
481                cfg % radiation_prefix_is_set = .TRUE.
482   
483             CASE( '--soil-prefix')
484                CALL get_option_argument( i, arg )
485                cfg % soil_prefix = TRIM(arg)
486                cfg % soil_prefix_is_set = .TRUE.
487   
488             CASE( '--soilmoisture-prefix')
489                CALL get_option_argument( i, arg )
490                cfg % soilmoisture_prefix = TRIM(arg)
491                cfg % soilmoisture_prefix_is_set = .TRUE.
492
493             CASE( '-o', '--output' )
494                CALL get_option_argument( i, arg )
495                cfg % output_file = TRIM(arg)
496
497             CASE( '-n', '--namelist' )
498                CALL get_option_argument( i, arg )
499                cfg % namelist_file = TRIM(arg)
500
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
519             CASE DEFAULT
520                message = "unknown option '" // TRIM(option) // "'."
521                CALL inifor_abort('parse_command_line_arguments', message)
522
523             END SELECT
524
525             i = i + 1
526
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
538   
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
570!------------------------------------------------------------------------------!
571! Description:
572! ------------
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! ------------
640!> Get the argument of the i'th command line option, which is at the location
641!> i+1 of the argument list.
642!------------------------------------------------------------------------------!
643   SUBROUTINE get_option_argument(i, arg)
644      CHARACTER(LEN=PATH), INTENT(INOUT) ::  arg
645      INTEGER, INTENT(INOUT)             ::  i
646
647      i = i + 1
648      CALL GET_COMMAND_ARGUMENT(i, arg)
649
650   END SUBROUTINE
651
652
653!------------------------------------------------------------------------------!
654! Description:
655! ------------
656!> Checks the INIFOR configuration 'cfg' for plausibility.
657!------------------------------------------------------------------------------!
658   SUBROUTINE validate_config(cfg)
659      TYPE(inifor_config), INTENT(IN) ::  cfg
660
661      CALL verify_file(cfg % hhl_file, 'HHL')
662      CALL verify_file(cfg % namelist_file, 'NAMELIST')
663      CALL verify_file(cfg % soiltyp_file, 'SOILTYP')
664
665!
666!--   Only check optional static driver file name, if it has been given.
667      IF (TRIM(cfg % static_driver_file) .NE. '')  THEN
668         CALL verify_file(cfg % static_driver_file, 'static driver')
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."
679         CALL inifor_abort( 'validate_config', message )
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."
691         CALL inifor_abort( 'validate_config', message )
692      END SELECT
693
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."
702         CALL inifor_abort( 'validate_config', message )
703      END SELECT
704
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."
708         CALL inifor_abort( 'validate_config', message )
709      END IF
710
711   END SUBROUTINE validate_config
712
713
714!------------------------------------------------------------------------------!
715! Description:
716! ------------
717!> Check whether the given file is present on the filesystem.
718!------------------------------------------------------------------------------!
719   LOGICAL FUNCTION file_present(filename)
720      CHARACTER(LEN=PATH), INTENT(IN) ::  filename
721
722      INQUIRE(FILE=filename, EXIST=file_present)
723
724   END FUNCTION file_present
725
726
727!------------------------------------------------------------------------------!
728! Description:
729! ------------
730!> This routine initializes the dynamic driver file, i.e. INIFOR's netCDF output
731!> file.
732!>
733!> Besides writing metadata, such as global attributes, coordinates, variables,
734!> in the netCDF file, various netCDF IDs are saved for later, when INIFOR
735!> writes the actual data.
736!------------------------------------------------------------------------------!
737   SUBROUTINE setup_netcdf_dimensions(output_file, palm_grid,                  &
738                                      start_date_string, origin_lon, origin_lat)
739
740       TYPE(nc_file), INTENT(INOUT)      ::  output_file
741       TYPE(grid_definition), INTENT(IN) ::  palm_grid
742       CHARACTER (LEN=DATE), INTENT(IN)  ::  start_date_string
743       REAL(dp), INTENT(IN)              ::  origin_lon, origin_lat
744
745       CHARACTER (LEN=8)     ::  date_string
746       CHARACTER (LEN=10)    ::  time_string
747       CHARACTER (LEN=5)     ::  zone_string
748       CHARACTER (LEN=SNAME) ::  history_string
749       INTEGER               ::  ncid, nx, ny, nz, nt, dimids(3), dimvarids(3)
750       REAL(dp)              ::  z0
751
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
756!
757!--    Create the netCDF file as in netCDF-4/HDF5 format if __netcdf4 preprocessor flag is given
758#if defined( __netcdf4 )
759       CALL check(nf90_create(TRIM(output_file % name), OR(NF90_CLOBBER, NF90_HDF5), ncid))
760#else
761       CALL check(nf90_create(TRIM(output_file % name), NF90_CLOBBER, ncid))
762#endif
763
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
774!
775!------------------------------------------------------------------------------
776!- Section 2: Write global NetCDF attributes
777!------------------------------------------------------------------------------
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
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'))
787       CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'history',        TRIM(history_string)))
788       CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'references',     '--'))
789       CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'comment',        '--'))
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)'))))
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)
796       CALL check(nf90_put_att(ncid, NF90_GLOBAL, 'origin_z',       TRIM(real_to_str(z0, '(F18.13)'))))
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!------------------------------------------------------------------------------
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
814
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"))
821
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"))
825
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
832
833!
834!--    overwrite third dimid with the one of depth
835       CALL check(nf90_def_dim(ncid, "zsoil", SIZE(palm_grid % depths), dimids(3)) )
836!
837!--    save dimids for later
838       output_file % dimids_soil = dimids
839
840!
841!--    overwrite third dimvarid with the one of depth
842       CALL check(nf90_def_var(ncid, "zsoil", NF90_FLOAT, output_file % dimids_soil(3), dimvarids(3)))
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!
847!--    save dimvarids for later
848       output_file % dimvarids_soil = dimvarids
849!
850!------------------------------------------------------------------------------
851!- Section 2b: Define dimensions for cell faces/velocities
852!------------------------------------------------------------------------------
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
862
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"))
869
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"))
873
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
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"))
891       CALL check(nf90_put_att(ncid, output_file % dimvarid_time, "units",     &
892                               "seconds since " // start_date_string // " UTC"))
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
908!
909!--    TODO Read in soil depths from input file before this.
910       CALL check(nf90_put_var(ncid, output_file % dimvarids_soil(3), palm_grid%depths))
911
912!
913!--    Write time vector
914       CALL check(nf90_put_var(ncid, output_file % dimvarid_time, output_file % time))
915
916!
917!--    Close the file
918       CALL check(nf90_close(ncid))
919
920    END SUBROUTINE setup_netcdf_dimensions
921
922
923!------------------------------------------------------------------------------!
924! Description:
925! ------------
926!> Defines the netCDF variables to be written to the dynamic driver file
927!------------------------------------------------------------------------------!
928    SUBROUTINE setup_netcdf_variables(filename, output_variable_table, debug)
929
930       CHARACTER (LEN=*), INTENT(IN)        ::  filename
931       TYPE(nc_var), INTENT(INOUT), TARGET  ::  output_variable_table(:)
932       LOGICAL, INTENT(IN)                  ::  debug
933
934       TYPE(nc_var), POINTER                ::  var
935       INTEGER                              ::  i, ncid
936       LOGICAL                              ::  to_be_written
937
938       message = "Defining variables in dynamic driver '" // TRIM(filename) // "'."
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
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
952             message = "  variable #" // TRIM(str(i)) // " '" // TRIM(var%name) // "'."
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
964       message = "Dynamic driver '" // TRIM(filename) // "' initialized successfully."
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!------------------------------------------------------------------------------
1001!
1002!--    radiation budgets, precipitation
1003       IF (group % kind == 'running average' .OR.                              &
1004           group % kind == 'accumulated')  THEN
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)) ) // "."
1011             CALL inifor_abort('read_input_variables | accumulation', message)
1012          END IF
1013
1014!
1015!--       use two buffer arrays
1016          nbuffers = 2
1017          IF ( .NOT. ALLOCATED( buffer ) )  ALLOCATE( buffer(nbuffers) )
1018
1019!
1020!--       hour of the day
1021          hour = iter - 1
1022!
1023!--       chose correct buffer array
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)
1031          CALL get_netcdf_variable(input_file, input_var, buffer(buf_id) % array)
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
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)
1049          nbuffers = MAX( group % n_inputs, group % n_output_quantities )
1050          ALLOCATE( buffer(nbuffers) )
1051 CALL run_control('time', 'alloc')
1052         
1053!
1054!--       Read in all input variables, leave extra buffers-if any-untouched.
1055          DO ivar = 1, group % n_inputs
1056
1057             input_var => group % in_var_list(ivar)
1058
1059!
1060!            Check wheather P or PP is present in input file
1061             IF (input_var % name == 'P')  THEN
1062                input_var % name = TRIM( get_pressure_varname(input_file) )
1063 CALL run_control('time', 'read')
1064             END IF
1065
1066             CALL get_netcdf_variable(input_file, input_var, buffer(ivar) % array)
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
1077!------------------------------------------------------------------------------!
1078! Description:
1079! ------------
1080!> Select the appropriate buffer ID for accumulated COSMO input variables
1081!> depending on the current hour.
1082!------------------------------------------------------------------------------!
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))
1097           CALL inifor_abort('select_buffer', message)
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!------------------------------------------------------------------------------!
1108    CHARACTER(LEN=2) FUNCTION get_pressure_varname(input_file) RESULT(var)
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) // "'."
1126          CALL inifor_abort('get_pressure_var', message)
1127
1128       END IF
1129
1130       CALL check(nf90_close(ncid))
1131
1132    END FUNCTION get_pressure_varname
1133
1134
1135!------------------------------------------------------------------------------!
1136! Description:
1137! ------------
1138!> Read the given global attribute form the given netCDF file.
1139!------------------------------------------------------------------------------!
1140    FUNCTION get_netcdf_attribute(filename, attribute) RESULT(attribute_value)
1141
1142       CHARACTER(LEN=*), INTENT(IN) ::  filename, attribute
1143       REAL(dp)                     ::  attribute_value
1144
1145       INTEGER                      ::  ncid
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))
1150          CALL check(nf90_close(ncid))
1151
1152       ELSE
1153
1154          message = "Failed to read '" // TRIM(attribute) // &
1155                    "' from file '" // TRIM(filename) // "'."
1156          CALL inifor_abort('get_netcdf_attribute', message)
1157
1158       END IF
1159
1160    END FUNCTION get_netcdf_attribute
1161
1162
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!------------------------------------------------------------------------------!
1170    SUBROUTINE update_output(var, array, iter, output_file, cfg)
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
1175       TYPE(inifor_config)       ::  cfg
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
1184!
1185!--    Skip time dimension for output
1186       ndim = var % ndim
1187       IF ( var_is_time_dependent )  ndim = var % ndim - 1
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
1195!
1196!--    Reduce dimension of output array according to variable kind
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
1260       CASE ( 'constant scalar profile', 'geostrophic' )
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
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
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
1279       CASE DEFAULT
1280
1281           message = "Variable kind '" // TRIM(var % kind) //                  &
1282                    "' not recognized."
1283           CALL inifor_abort('update_output', message)
1284
1285       END SELECT
1286
1287       CALL check(nf90_close(ncid))
1288
1289    END SUBROUTINE update_output
1290
1291
1292!------------------------------------------------------------------------------!
1293! Description:
1294! ------------
1295!> Checks the status of a netCDF API call and aborts if an error occured
1296!------------------------------------------------------------------------------!
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) )
1304          CALL inifor_abort('io.check', message)
1305       END IF
1306
1307    END SUBROUTINE check
1308
1309 END MODULE inifor_io
1310#endif
Note: See TracBrowser for help on using the repository browser.