source: palm/trunk/SOURCE/diagnostic_output_quantities_mod.f90 @ 4869

Last change on this file since 4869 was 4866, checked in by suehring, 4 years ago

Implemented vertical passive scalar flux

  • Property svn:keywords set to Id
File size: 97.2 KB
RevLine 
[3994]1!> @file diagnostic_output_quantities_mod.f90
[4583]2!--------------------------------------------------------------------------------------------------!
[3994]3! This file is part of the PALM model system.
4!
[4583]5! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
6! Public License as published by the Free Software Foundation, either version 3 of the License, or
7! (at your option) any later version.
[3994]8!
[4583]9! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
10! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
11! Public License for more details.
[3994]12!
[4583]13! You should have received a copy of the GNU General Public License along with PALM. If not, see
14! <http://www.gnu.org/licenses/>.
[3994]15!
[4828]16! Copyright 1997-2021 Leibniz Universitaet Hannover
[4583]17!--------------------------------------------------------------------------------------------------!
[3994]18!
19! Current revisions:
20! ------------------
[4671]21!
[4866]22!
[3994]23! Former revisions:
24! -----------------
25! $Id: diagnostic_output_quantities_mod.f90 4866 2021-02-05 17:50:14Z schwenkel $
[4866]26! Implemented vertical passive scalar flux
27! Implemented air temperature as diagnostic output quantity (previous commit)
28!
29! 4861 2021-02-01 10:51:45Z raasch
[4768]30! Enable 3D data output also with 64-bit precision
31!
32! 4757 2020-10-26 10:23:38Z schwenkel
[4757]33! Implement relative humidity as diagnostic output quantity
34!
35! 4671 2020-09-09 20:27:58Z pavelkrc
[4671]36! Implementation of downward facing USM and LSM surfaces
[4861]37!
[4671]38! 4583 2020-06-29 12:36:47Z raasch
[4583]39! file re-formatted to follow the PALM coding standard
40!
41! 4560 2020-06-11 12:19:47Z suehring
[4560]42! - Bugfix in calculation of vertical momentum and scalar fluxes
43! - remove averaged output variables from PUBLIC list
[4583]44!
[4560]45! 4535 2020-05-15 12:07:23Z raasch
[4535]46! bugfix for restart data format query
[4583]47!
[4535]48! 4518 2020-05-04 15:44:28Z suehring
[4583]49! * Define arrays over ghost points in order to allow for standard mpi-io treatment. By this
50!   modularization of restart-data input is possible with the module interface.
[4518]51! * Bugfix: add missing restart input of wtheta_av, wq_av, wu_av, and wv_av.
[4583]52!
[4518]53! 4517 2020-05-03 14:29:30Z raasch
[4457]54! use statement for exchange horiz added,
55! bugfix for call of exchange horiz 2d
[4583]56!
[4457]57! 4431 2020-02-27 23:23:01Z gronemeier
[4431]58! added wspeed and wdir output; bugfix: set fill_value in case of masked output
59!
60! 4360 2020-01-07 11:25:50Z suehring
[4583]61! added output of wu, wv, wtheta and wq to enable covariance calculation according to temporal EC
62! method
[4431]63!
[4351]64! 4346 2019-12-18 11:55:56Z motisi
[4583]65! Introduction of wall_flags_total_0, which currently sets bits based on static topography
66! information used in wall_flags_static_0
[4431]67!
[4346]68! 4331 2019-12-10 18:25:02Z suehring
[4331]69! - Modularize 2-m potential temperature output
70! - New output for 10-m wind speed
[4431]71!
[4331]72! 4329 2019-12-10 15:46:36Z motisi
[4329]73! Renamed wall_flags_0 to wall_flags_static_0
[4431]74!
[4329]75! 4182 2019-08-22 15:20:23Z scharf
[4182]76! Corrected "Former revisions" section
[4431]77!
[4182]78! 4167 2019-08-16 11:01:48Z suehring
[4583]79! Changed behaviour of masked output over surface to follow terrain and ignore buildings
80! (J.Resler, T.Gronemeier)
[4431]81!
[4167]82! 4157 2019-08-14 09:19:12Z suehring
[4583]83! Initialization restructured, in order to work also when data output during spin-up is enabled.
[4431]84!
[4157]85! 4132 2019-08-02 12:34:17Z suehring
[4132]86! Bugfix in masked data output
[4431]87!
[4132]88! 4069 2019-07-01 14:05:51Z Giersch
[4583]89! Masked output running index mid has been introduced as a local variable to avoid runtime error
90! (Loop variable has been modified) in time_integration
[4431]91!
[4069]92! 4039 2019-06-18 10:32:41Z suehring
[4583]93! - Add output of uu, vv, ww to enable variance calculation according temporal EC method
[4039]94! - Allocate arrays only when they are required
95! - Formatting adjustment
96! - Rename subroutines
97! - Further modularization
[4431]98!
[4039]99! 3998 2019-05-23 13:38:11Z suehring
[4431]100! Bugfix in gathering all output strings
101!
[3998]102! 3995 2019-05-22 18:59:54Z suehring
[4583]103! Avoid compiler warnings about unused variable and fix string operation which is not allowed with
104! PGI compiler
[4431]105!
[3995]106! 3994 2019-05-22 18:08:09Z suehring
[3994]107! Initial revision
108!
[4182]109! Authors:
110! --------
[3994]111! @author Farah Kanani-Suehring
[4182]112!
[4431]113!
[3994]114! Description:
115! ------------
116!> ...
[4583]117!--------------------------------------------------------------------------------------------------!
[3994]118 MODULE diagnostic_output_quantities_mod
[4431]119
[4583]120    USE arrays_3d,                                                                                 &
121        ONLY:  ddzu,                                                                               &
[4757]122               exner,                                                                              &
123               hyp,                                                                                &
[4583]124               pt,                                                                                 &
125               q,                                                                                  &
[4866]126               s,                                                                                  &
[4583]127               u,                                                                                  &
128               v,                                                                                  &
129               w,                                                                                  &
130               zu,                                                                                 &
[4331]131               zw
132
[4583]133    USE basic_constants_and_equations_mod,                                                         &
[4861]134        ONLY:  degc_to_k,                                                                          &
135               kappa,                                                                              &
[4757]136               pi,                                                                                 &
137               magnus,                                                                             &
138               rd_d_rv
[4331]139
[4583]140    USE control_parameters,                                                                        &
141        ONLY:  current_timestep_number,                                                            &
142               data_output,                                                                        &
[4757]143               data_output_pr,                                                                     &
144               humidity,                                                                           &
[4583]145               message_string,                                                                     &
[4866]146               passive_scalar,                                                                     &
[4583]147               restart_data_format_output,                                                         &
[4331]148               varnamelength
[4431]149!
[4583]150!     USE cpulog,                                                                                   &
[3994]151!         ONLY:  cpu_log, log_point
[4039]152
[4583]153    USE exchange_horiz_mod,                                                                        &
[4457]154        ONLY:  exchange_horiz_2d
155
[4583]156   USE grid_variables,                                                                             &
[4039]157        ONLY:  ddx, ddy
[4331]158
[4583]159    USE indices,                                                                                   &
160        ONLY:  nbgp,                                                                               &
161               nxl,                                                                                &
162               nxlg,                                                                               &
163               nxr,                                                                                &
164               nxrg,                                                                               &
165               nyn,                                                                                &
166               nyng,                                                                               &
167               nys,                                                                                &
168               nysg,                                                                               &
169               nzb,                                                                                &
170               nzt,                                                                                &
[4346]171               wall_flags_total_0
[4331]172
[3994]173    USE kinds
174
[4757]175    USE profil_parameter,                                                                          &
176        ONLY:  dopr_index
177
[4518]178    USE restart_data_mpi_io_mod,                                                                   &
179        ONLY:  rd_mpi_io_check_array,                                                              &
180               rrd_mpi_io,                                                                         &
181               wrd_mpi_io
182
[4757]183    USE statistics,                                                                                &
184        ONLY:  hom,                                                                                &
185               pr_palm,                                                                            &
186               rmask,                                                                              &
187               statistic_regions,                                                                  &
188               sums_l
189
[4583]190    USE surface_mod,                                                                               &
191        ONLY:  surf_def_h,                                                                         &
192               surf_lsm_h,                                                                         &
193               surf_type,                                                                          &
[4331]194               surf_usm_h
[3994]195
196    IMPLICIT NONE
197
198    CHARACTER(LEN=varnamelength), DIMENSION(500) ::  do_all = ' '
199
200    INTEGER(iwp) ::  timestep_number_at_prev_calc = 0  !< ...at previous diagnostic output calculation
[4431]201
[4039]202    LOGICAL ::  initialized_diagnostic_output_quantities = .FALSE. !< flag indicating whether output is initialized
203    LOGICAL ::  prepared_diagnostic_output_quantities = .FALSE.    !< flag indicating whether output is p
[3994]204
[4331]205    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  pt_2m     !< 2-m air potential temperature
206    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  pt_2m_av  !< averaged 2-m air potential temperature
207    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  uv_10m    !< horizontal wind speed at 10m
208    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  uv_10m_av !< averaged horizontal wind speed at 10m
209
[4757]210    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  rh     !< relative humidity
211    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  rh_av  !< avg. relative humidity
[4861]212    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ta     !< air temperature
213    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ta_av  !< avg. air temperature
[3994]214    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ti     !< rotation(u,v,w) aka turbulence intensity
[4431]215    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ti_av  !< avg. rotation(u,v,w) aka turbulence intensity
216    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  u_center     !< u at center of grid box
217    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  u_center_av  !< mean of u_center
218    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  uu           !< uu
219    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  uu_av        !< mean of uu
220    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wspeed       !< horizontal wind speed
221    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wspeed_av    !< mean of horizotal wind speed
222    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  v_center     !< v at center of grid box
223    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  v_center_av  !< mean of v_center
224    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  vv           !< vv
225    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  vv_av        !< mean of vv
226    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wdir         !< wind direction
227    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wdir_av      !< mean wind direction
228    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ww           !< ww
229    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ww_av        !< mean of ww
230    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wu           !< wu
231    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wu_av        !< mean of wu
232    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wv           !< wv
233    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wv_av        !< mean of wv
234    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wtheta       !< wtheta
235    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wtheta_av    !< mean of wtheta
236    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wq           !< wq
237    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  wq_av        !< mean of wq
[4866]238    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ws           !< ws
239    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ws_av        !< mean of ws
[3994]240
241
242    SAVE
243
244    PRIVATE
245
246!
247!-- Public variables
[4583]248    PUBLIC do_all,                                                                                 &
249           initialized_diagnostic_output_quantities,                                               &
250           prepared_diagnostic_output_quantities,                                                  &
[4560]251           timestep_number_at_prev_calc
[4431]252!
253!-- Public routines
[4583]254    PUBLIC doq_3d_data_averaging,                                                                  &
255           doq_calculate,                                                                          &
256           doq_check_data_output,                                                                  &
[4757]257           doq_check_data_output_pr,                                                               &
[4583]258           doq_define_netcdf_grid,                                                                 &
259           doq_init,                                                                               &
260           doq_output_2d,                                                                          &
261           doq_output_3d,                                                                          &
262           doq_output_mask,                                                                        &
[4757]263           doq_statistics,                                                                         &
[4583]264           doq_rrd_local,                                                                          &
[4039]265           doq_wrd_local
[3994]266
267
[4039]268    INTERFACE doq_3d_data_averaging
269       MODULE PROCEDURE doq_3d_data_averaging
[4431]270    END INTERFACE doq_3d_data_averaging
[3994]271
[4039]272    INTERFACE doq_calculate
273       MODULE PROCEDURE doq_calculate
274    END INTERFACE doq_calculate
[3994]275
[4039]276    INTERFACE doq_check_data_output
277       MODULE PROCEDURE doq_check_data_output
278    END INTERFACE doq_check_data_output
[4431]279
[4757]280    INTERFACE doq_check_data_output_pr
281       MODULE PROCEDURE doq_check_data_output_pr
282    END INTERFACE doq_check_data_output_pr
283
[4039]284    INTERFACE doq_define_netcdf_grid
285       MODULE PROCEDURE doq_define_netcdf_grid
286    END INTERFACE doq_define_netcdf_grid
[4431]287
[4039]288    INTERFACE doq_output_2d
289       MODULE PROCEDURE doq_output_2d
290    END INTERFACE doq_output_2d
[4431]291
[4039]292    INTERFACE doq_output_3d
293       MODULE PROCEDURE doq_output_3d
294    END INTERFACE doq_output_3d
[4431]295
[4039]296    INTERFACE doq_output_mask
297       MODULE PROCEDURE doq_output_mask
298    END INTERFACE doq_output_mask
[4431]299
[4039]300    INTERFACE doq_init
301       MODULE PROCEDURE doq_init
302    END INTERFACE doq_init
[3994]303
[4757]304    INTERFACE doq_statistics
305       MODULE PROCEDURE doq_statistics
306    END INTERFACE doq_statistics
307
[4039]308    INTERFACE doq_prepare
309       MODULE PROCEDURE doq_prepare
310    END INTERFACE doq_prepare
[4431]311
[4518]312   INTERFACE doq_rrd_local
313       MODULE PROCEDURE doq_rrd_local_ftn
314       MODULE PROCEDURE doq_rrd_local_mpi
315    END INTERFACE doq_rrd_local
[4431]316
[4039]317    INTERFACE doq_wrd_local
318       MODULE PROCEDURE doq_wrd_local
319    END INTERFACE doq_wrd_local
[3994]320
321 CONTAINS
[4431]322
[4583]323!--------------------------------------------------------------------------------------------------!
[4039]324! Description:
325! ------------
[4583]326!> Sum up and time-average diagnostic output quantities as well as allocate the array necessary for
327!> storing the average.
328!--------------------------------------------------------------------------------------------------!
[4039]329 SUBROUTINE doq_3d_data_averaging( mode, variable )
[3994]330
[4583]331    USE control_parameters,                                                                        &
[4039]332        ONLY:  average_count_3d
333
[4431]334    CHARACTER (LEN=*) ::  mode     !<
335    CHARACTER (LEN=*) ::  variable !<
[4039]336
337    INTEGER(iwp) ::  i !<
338    INTEGER(iwp) ::  j !<
339    INTEGER(iwp) ::  k !<
340
[4583]341
[4039]342    IF ( mode == 'allocate' )  THEN
343
344       SELECT CASE ( TRIM( variable ) )
345
[4757]346          CASE ( 'rh' )
347             IF ( .NOT. ALLOCATED( rh_av ) )  THEN
348                ALLOCATE( rh_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
349             ENDIF
350             rh_av = 0.0_wp
351
[4861]352          CASE ( 'ta' )
353             IF ( .NOT. ALLOCATED( ta_av ) )  THEN
354                ALLOCATE( ta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
355             ENDIF
356             ta_av = 0.0_wp
357
[4039]358          CASE ( 'ti' )
359             IF ( .NOT. ALLOCATED( ti_av ) )  THEN
[4518]360                ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]361             ENDIF
362             ti_av = 0.0_wp
[4431]363
[4039]364          CASE ( 'uu' )
365             IF ( .NOT. ALLOCATED( uu_av ) )  THEN
[4518]366                ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]367             ENDIF
368             uu_av = 0.0_wp
[4431]369
[4039]370          CASE ( 'vv' )
371             IF ( .NOT. ALLOCATED( vv_av ) )  THEN
[4518]372                ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]373             ENDIF
374             vv_av = 0.0_wp
[4431]375
[4039]376          CASE ( 'ww' )
377             IF ( .NOT. ALLOCATED( ww_av ) )  THEN
[4518]378                ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]379             ENDIF
380             ww_av = 0.0_wp
[4431]381
[4351]382           CASE ( 'wu' )
383             IF ( .NOT. ALLOCATED( wu_av ) )  THEN
[4518]384                ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]385             ENDIF
386             wu_av = 0.0_wp
[4431]387
[4351]388           CASE ( 'wv' )
389             IF ( .NOT. ALLOCATED( wv_av ) )  THEN
[4518]390                ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]391             ENDIF
392             wv_av = 0.0_wp
[4431]393
[4351]394           CASE ( 'wtheta' )
395             IF ( .NOT. ALLOCATED( wtheta_av ) )  THEN
[4518]396                ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]397             ENDIF
398             wtheta_av = 0.0_wp
[4431]399
[4351]400           CASE ( 'wq' )
401             IF ( .NOT. ALLOCATED( wq_av ) )  THEN
[4518]402                ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]403             ENDIF
404             wq_av = 0.0_wp
[4431]405
[4866]406          CASE ( 'ws' )
407             IF ( .NOT. ALLOCATED( ws_av ) )  THEN
408                ALLOCATE( ws_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
409             ENDIF
410             ws_av = 0.0_wp
411
[4331]412          CASE ( 'theta_2m*' )
413             IF ( .NOT. ALLOCATED( pt_2m_av ) )  THEN
414                ALLOCATE( pt_2m_av(nysg:nyng,nxlg:nxrg) )
415             ENDIF
416             pt_2m_av = 0.0_wp
[4431]417
[4331]418          CASE ( 'wspeed_10m*' )
419             IF ( .NOT. ALLOCATED( uv_10m_av ) )  THEN
420                ALLOCATE( uv_10m_av(nysg:nyng,nxlg:nxrg) )
421             ENDIF
422             uv_10m_av = 0.0_wp
[4039]423
[4431]424          CASE ( 'wspeed' )
425             IF ( .NOT. ALLOCATED( wspeed_av ) )  THEN
[4518]426                ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]427             ENDIF
428             wspeed_av = 0.0_wp
429
430          CASE ( 'wdir' )
431             IF ( .NOT. ALLOCATED( u_center_av ) )  THEN
[4518]432                ALLOCATE( u_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]433             ENDIF
434             IF ( .NOT. ALLOCATED( v_center_av ) )  THEN
[4518]435                ALLOCATE( v_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]436             ENDIF
437             u_center_av = 0.0_wp
438             v_center_av = 0.0_wp
439
[4039]440          CASE DEFAULT
441             CONTINUE
442
443       END SELECT
444
445    ELSEIF ( mode == 'sum' )  THEN
446
447       SELECT CASE ( TRIM( variable ) )
[4431]448
[4757]449          CASE ( 'rh' )
450             IF ( ALLOCATED( rh_av ) )  THEN
451                DO  i = nxl, nxr
452                   DO  j = nys, nyn
453                      DO  k = nzb, nzt+1
454                         rh_av(k,j,i) = rh_av(k,j,i) + rh(k,j,i)
455                      ENDDO
456                   ENDDO
457                ENDDO
458             ENDIF
459
[4861]460          CASE ( 'ta' )
461             IF ( ALLOCATED( ta_av ) )  THEN
462                DO  i = nxl, nxr
463                   DO  j = nys, nyn
464                      DO  k = nzb, nzt+1
465                         ta_av(k,j,i) = ta_av(k,j,i) + ta(k,j,i)
466                      ENDDO
467                   ENDDO
468                ENDDO
469             ENDIF
470
[4039]471          CASE ( 'ti' )
[4583]472             IF ( ALLOCATED( ti_av ) )  THEN
[4039]473                DO  i = nxl, nxr
474                   DO  j = nys, nyn
475                      DO  k = nzb, nzt+1
476                         ti_av(k,j,i) = ti_av(k,j,i) + ti(k,j,i)
477                      ENDDO
478                   ENDDO
479                ENDDO
480             ENDIF
[4431]481
[4039]482          CASE ( 'uu' )
[4583]483             IF ( ALLOCATED( uu_av ) )  THEN
[4039]484                DO  i = nxl, nxr
485                   DO  j = nys, nyn
486                      DO  k = nzb, nzt+1
487                         uu_av(k,j,i) = uu_av(k,j,i) + uu(k,j,i)
488                      ENDDO
489                   ENDDO
490                ENDDO
491             ENDIF
[4431]492
[4039]493          CASE ( 'vv' )
[4583]494             IF ( ALLOCATED( vv_av ) )  THEN
[4039]495                DO  i = nxl, nxr
496                   DO  j = nys, nyn
497                      DO  k = nzb, nzt+1
498                         vv_av(k,j,i) = vv_av(k,j,i) + vv(k,j,i)
499                      ENDDO
500                   ENDDO
501                ENDDO
502             ENDIF
[4431]503
[4039]504          CASE ( 'ww' )
[4583]505             IF ( ALLOCATED( ww_av ) )  THEN
[4039]506                DO  i = nxl, nxr
507                   DO  j = nys, nyn
508                      DO  k = nzb, nzt+1
509                         ww_av(k,j,i) = ww_av(k,j,i) + ww(k,j,i)
510                      ENDDO
511                   ENDDO
512                ENDDO
513             ENDIF
[4431]514
[4351]515          CASE ( 'wu' )
[4583]516             IF ( ALLOCATED( wu_av ) )  THEN
[4351]517                DO  i = nxl, nxr
518                   DO  j = nys, nyn
519                      DO  k = nzb, nzt+1
520                         wu_av(k,j,i) = wu_av(k,j,i) + wu(k,j,i)
521                      ENDDO
522                   ENDDO
523                ENDDO
524             ENDIF
[4431]525
[4351]526          CASE ( 'wv' )
[4583]527             IF ( ALLOCATED( wv_av ) )  THEN
[4351]528                DO  i = nxl, nxr
529                   DO  j = nys, nyn
530                      DO  k = nzb, nzt+1
531                         wv_av(k,j,i) = wv_av(k,j,i) + wv(k,j,i)
532                      ENDDO
533                   ENDDO
534                ENDDO
535             ENDIF
[4431]536
[4351]537          CASE ( 'wtheta' )
[4583]538             IF ( ALLOCATED( wtheta_av ) )  THEN
[4351]539                DO  i = nxl, nxr
540                   DO  j = nys, nyn
541                      DO  k = nzb, nzt+1
542                         wtheta_av(k,j,i) = wtheta_av(k,j,i) + wtheta(k,j,i)
543                      ENDDO
544                   ENDDO
545                ENDDO
546             ENDIF
[4431]547
[4351]548          CASE ( 'wq' )
[4583]549             IF ( ALLOCATED( wq_av ) )  THEN
[4351]550                DO  i = nxl, nxr
551                   DO  j = nys, nyn
552                      DO  k = nzb, nzt+1
553                         wq_av(k,j,i) = wq_av(k,j,i) + wq(k,j,i)
554                      ENDDO
555                   ENDDO
556                ENDDO
557             ENDIF
[4431]558
[4866]559          CASE ( 'ws' )
560             IF ( ALLOCATED( ws_av ) )  THEN
561                DO  i = nxl, nxr
562                   DO  j = nys, nyn
563                      DO  k = nzb, nzt+1
564                         ws_av(k,j,i) = ws_av(k,j,i) + ws(k,j,i)
565                      ENDDO
566                   ENDDO
567                ENDDO
568             ENDIF
569
[4331]570          CASE ( 'theta_2m*' )
[4583]571             IF ( ALLOCATED( pt_2m_av ) )  THEN
[4331]572                DO  i = nxl, nxr
573                   DO  j = nys, nyn
574                      pt_2m_av(j,i) = pt_2m_av(j,i) + pt_2m(j,i)
575                   ENDDO
576                ENDDO
577             ENDIF
[4039]578
[4331]579          CASE ( 'wspeed_10m*' )
[4583]580             IF ( ALLOCATED( uv_10m_av ) )  THEN
[4331]581                DO  i = nxl, nxr
582                   DO  j = nys, nyn
583                      uv_10m_av(j,i) = uv_10m_av(j,i) + uv_10m(j,i)
584                   ENDDO
585                ENDDO
586             ENDIF
587
[4431]588          CASE ( 'wspeed' )
[4583]589            IF ( ALLOCATED( wspeed_av ) )  THEN
[4431]590               DO  i = nxl, nxr
591                  DO  j = nys, nyn
592                     DO  k = nzb, nzt+1
593                         wspeed_av(k,j,i) = wspeed_av(k,j,i) + wspeed(k,j,i)
594                     ENDDO
595                  ENDDO
596               ENDDO
597            ENDIF
598
599          CASE ( 'wdir' )
[4583]600             IF ( ALLOCATED( u_center_av )  .AND.  ALLOCATED( v_center_av ) )  THEN
[4431]601                DO  i = nxl, nxr
602                   DO  j = nys, nyn
603                      DO  k = nzb, nzt+1
604                        u_center_av(k,j,i) = u_center_av(k,j,i) + u_center(k,j,i)
605                        v_center_av(k,j,i) = v_center_av(k,j,i) + v_center(k,j,i)
606                      ENDDO
607                   ENDDO
608                ENDDO
609             ENDIF
610
[4039]611          CASE DEFAULT
612             CONTINUE
613
614       END SELECT
615
616    ELSEIF ( mode == 'average' )  THEN
617
618       SELECT CASE ( TRIM( variable ) )
619
[4757]620          CASE ( 'rh' )
621             IF ( ALLOCATED( rh_av ) )  THEN
622                DO  i = nxl, nxr
623                   DO  j = nys, nyn
624                      DO  k = nzb, nzt+1
625                         rh_av(k,j,i) = rh_av(k,j,i) / REAL( average_count_3d, KIND=wp )
626                      ENDDO
627                   ENDDO
628                ENDDO
629             ENDIF
630
[4861]631          CASE ( 'ta' )
632             IF ( ALLOCATED( ta_av ) )  THEN
633                DO  i = nxl, nxr
634                   DO  j = nys, nyn
635                      DO  k = nzb, nzt+1
636                         ta_av(k,j,i) = ta_av(k,j,i) / REAL( average_count_3d, KIND=wp )
637                      ENDDO
638                   ENDDO
639                ENDDO
640             ENDIF
641
[4039]642          CASE ( 'ti' )
[4583]643             IF ( ALLOCATED( ti_av ) )  THEN
[4039]644                DO  i = nxl, nxr
645                   DO  j = nys, nyn
646                      DO  k = nzb, nzt+1
647                         ti_av(k,j,i) = ti_av(k,j,i) / REAL( average_count_3d, KIND=wp )
648                      ENDDO
649                   ENDDO
650                ENDDO
651             ENDIF
[4431]652
[4039]653          CASE ( 'uu' )
[4583]654             IF ( ALLOCATED( uu_av ) )  THEN
[4039]655                DO  i = nxl, nxr
656                   DO  j = nys, nyn
657                      DO  k = nzb, nzt+1
658                         uu_av(k,j,i) = uu_av(k,j,i) / REAL( average_count_3d, KIND=wp )
659                      ENDDO
660                   ENDDO
661                ENDDO
662             ENDIF
[4431]663
[4039]664          CASE ( 'vv' )
[4583]665             IF ( ALLOCATED( vv_av ) )  THEN
[4039]666                DO  i = nxl, nxr
667                   DO  j = nys, nyn
668                      DO  k = nzb, nzt+1
669                         vv_av(k,j,i) = vv_av(k,j,i) / REAL( average_count_3d, KIND=wp )
670                      ENDDO
671                   ENDDO
672                ENDDO
673             ENDIF
[4431]674
[4039]675          CASE ( 'ww' )
[4583]676             IF ( ALLOCATED( ww_av ) )  THEN
[4039]677                DO  i = nxl, nxr
678                   DO  j = nys, nyn
679                      DO  k = nzb, nzt+1
680                         ww_av(k,j,i) = ww_av(k,j,i) / REAL( average_count_3d, KIND=wp )
681                      ENDDO
682                   ENDDO
683                ENDDO
684             ENDIF
685
[4351]686          CASE ( 'wu' )
[4583]687             IF ( ALLOCATED( wu_av ) )  THEN
[4351]688                DO  i = nxl, nxr
689                   DO  j = nys, nyn
690                      DO  k = nzb, nzt+1
691                         wu_av(k,j,i) = wu_av(k,j,i) / REAL( average_count_3d, KIND=wp )
692                      ENDDO
693                   ENDDO
694                ENDDO
695             ENDIF
[4431]696
[4351]697          CASE ( 'wv' )
[4583]698             IF ( ALLOCATED( wv_av ) )  THEN
[4351]699                DO  i = nxl, nxr
700                   DO  j = nys, nyn
701                      DO  k = nzb, nzt+1
702                         wv_av(k,j,i) = wv_av(k,j,i) / REAL( average_count_3d, KIND=wp )
703                      ENDDO
704                   ENDDO
705                ENDDO
706             ENDIF
[4431]707
[4351]708          CASE ( 'wtheta' )
[4583]709             IF ( ALLOCATED( wtheta_av ) )  THEN
[4351]710                DO  i = nxl, nxr
711                   DO  j = nys, nyn
712                      DO  k = nzb, nzt+1
713                         wtheta_av(k,j,i) = wtheta_av(k,j,i) / REAL( average_count_3d, KIND=wp )
714                      ENDDO
715                   ENDDO
716                ENDDO
717             ENDIF
[4431]718
[4351]719          CASE ( 'wq' )
[4583]720             IF ( ALLOCATED( wq_av ) )  THEN
[4351]721                DO  i = nxl, nxr
722                   DO  j = nys, nyn
723                      DO  k = nzb, nzt+1
724                         wq_av(k,j,i) = wq_av(k,j,i) / REAL( average_count_3d, KIND=wp )
725                      ENDDO
726                   ENDDO
727                ENDDO
728             ENDIF
729
[4866]730          CASE ( 'ws' )
731             IF ( ALLOCATED( ws_av ) )  THEN
732                DO  i = nxl, nxr
733                   DO  j = nys, nyn
734                      DO  k = nzb, nzt+1
735                         ws_av(k,j,i) = ws_av(k,j,i) / REAL( average_count_3d, KIND=wp )
736                      ENDDO
737                   ENDDO
738                ENDDO
739             ENDIF
740
[4331]741         CASE ( 'theta_2m*' )
[4583]742            IF ( ALLOCATED( pt_2m_av ) )  THEN
[4331]743               DO  i = nxlg, nxrg
744                  DO  j = nysg, nyng
745                     pt_2m_av(j,i) = pt_2m_av(j,i) / REAL( average_count_3d, KIND=wp )
746                  ENDDO
747               ENDDO
[4457]748               CALL exchange_horiz_2d( pt_2m_av )
[4331]749            ENDIF
750
751         CASE ( 'wspeed_10m*' )
[4583]752            IF ( ALLOCATED( uv_10m_av ) )  THEN
[4331]753               DO  i = nxlg, nxrg
754                  DO  j = nysg, nyng
755                     uv_10m_av(j,i) = uv_10m_av(j,i) / REAL( average_count_3d, KIND=wp )
756                  ENDDO
757               ENDDO
[4457]758               CALL exchange_horiz_2d( uv_10m_av )
[4331]759            ENDIF
760
[4431]761         CASE ( 'wspeed' )
[4583]762             IF ( ALLOCATED( wspeed_av ) )  THEN
[4431]763                DO  i = nxl, nxr
764                   DO  j = nys, nyn
765                      DO  k = nzb, nzt+1
766                         wspeed_av(k,j,i) = wspeed_av(k,j,i) / REAL( average_count_3d, KIND=wp )
767                      ENDDO
768                   ENDDO
769                ENDDO
770             ENDIF
771
772          CASE ( 'wdir' )
[4583]773             IF ( ALLOCATED( u_center_av )  .AND.  ALLOCATED( v_center_av ) )  THEN
[4431]774
775                IF ( .NOT. ALLOCATED( wdir_av ) )  THEN
[4518]776                   ALLOCATE( wdir_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]777                ENDIF
778                wdir_av = 0.0_wp
779
780                DO  i = nxl, nxr
781                   DO  j = nys, nyn
782                      DO  k = nzb, nzt+1
783                         u_center_av(k,j,i) = u_center_av(k,j,i) / REAL( average_count_3d, KIND=wp )
784                         v_center_av(k,j,i) = v_center_av(k,j,i) / REAL( average_count_3d, KIND=wp )
[4583]785                         wdir_av(k,j,i) = ATAN2( u_center_av(k,j,i), v_center_av(k,j,i) )          &
786                                          / pi * 180.0_wp + 180.0_wp
[4431]787                      ENDDO
788                   ENDDO
789                ENDDO
790             ENDIF
791
[4039]792       END SELECT
793
794    ENDIF
795
796
[4431]797 END SUBROUTINE doq_3d_data_averaging
798
[4583]799!--------------------------------------------------------------------------------------------------!
[3994]800! Description:
801! ------------
[4039]802!> Check data output for diagnostic output
[4583]803!--------------------------------------------------------------------------------------------------!
[4331]804 SUBROUTINE doq_check_data_output( var, unit, i, ilen, k )
[3994]805
[4039]806    IMPLICIT NONE
[3994]807
[4431]808    CHARACTER (LEN=*) ::  unit  !<
[4039]809    CHARACTER (LEN=*) ::  var   !<
[3994]810
[4331]811    INTEGER(iwp), OPTIONAL, INTENT(IN) ::  i     !< Current element of data_output
812    INTEGER(iwp), OPTIONAL, INTENT(IN) ::  ilen  !< Length of current entry in data_output
813    INTEGER(iwp), OPTIONAL, INTENT(IN) ::  k     !< Output is xy mode? 0 = no, 1 = yes
814
[4583]815
[4039]816    SELECT CASE ( TRIM( var ) )
817
[4757]818       CASE ( 'rh' )
819          unit = '%'
820
[4861]821       CASE ( 'ta' )
822          unit = 'degree_C'
823
[4039]824       CASE ( 'ti' )
825          unit = '1/s'
[4431]826
[4039]827       CASE ( 'uu' )
828          unit = 'm2/s2'
[4431]829
[4039]830       CASE ( 'vv' )
831          unit = 'm2/s2'
[4431]832
[4039]833       CASE ( 'ww' )
834          unit = 'm2/s2'
[4431]835
[4351]836       CASE ( 'wu' )
837          unit = 'm2/s2'
[4431]838
[4351]839       CASE ( 'wv' )
840          unit = 'm2/s2'
[4431]841
[4351]842       CASE ( 'wtheta' )
843          unit = 'Km/s'
[4431]844
[4351]845       CASE ( 'wq' )
846          unit = 'm/s'
[4431]847
[4866]848       CASE ( 'ws' )
849          IF ( .NOT.  passive_scalar )  THEN
850             message_string = 'data_output = ' // TRIM( var ) //                                   &
851                              ' is not implemented for passive_scalar = .FALSE.'
852             CALL message( 'diagnostic_output', 'PA0093', 1, 2, 0, 6, 0 )
853          ENDIF
854          unit = 'm/s'
855
[4431]856       CASE ( 'wspeed' )
857          unit = 'm/s'
858
859       CASE ( 'wdir' )
860          unit = 'degree'
[4331]861!
862!--    Treat horizotal cross-section output quanatities
863       CASE ( 'theta_2m*', 'wspeed_10m*' )
864!
865!--       Check if output quantity is _xy only.
866          IF ( k == 0  .OR.  data_output(i)(ilen-2:ilen) /= '_xy' )  THEN
[4583]867             message_string = 'illegal value for data_output: "' //                                &
868                              TRIM( var ) // '" & only 2d-horizontal ' //                          &
[4331]869                              'cross sections are allowed for this value'
870             CALL message( 'diagnostic_output', 'PA0111', 1, 2, 0, 6, 0 )
871          ENDIF
872
873          IF ( TRIM( var ) == 'theta_2m*'   )  unit = 'K'
874          IF ( TRIM( var ) == 'wspeed_10m*' )  unit = 'm/s'
875
[4039]876       CASE DEFAULT
877          unit = 'illegal'
878
879    END SELECT
880
881
882 END SUBROUTINE doq_check_data_output
[4431]883
[4757]884
[4583]885!--------------------------------------------------------------------------------------------------!
[4757]886! Description:
887! ------------
888!> Set the unit of user defined profile output quantities. For those variables not recognized by the
889!> user, the parameter unit is set to "illegal", which tells the calling routine that the
890!> output variable is not defined and leads to a program abort.
891!--------------------------------------------------------------------------------------------------!
892 SUBROUTINE doq_check_data_output_pr( variable, var_count, unit, dopr_unit )
893
894    CHARACTER (LEN=*) ::  unit      !<
895    CHARACTER (LEN=*) ::  variable  !<
896    CHARACTER (LEN=*) ::  dopr_unit !< local value of dopr_unit
897
898    INTEGER(iwp) ::  pr_index       !<
899    INTEGER(iwp) ::  var_count      !<
900
901    SELECT CASE ( TRIM( variable ) )
902
903       CASE ( 'rh' )
904          IF ( .NOT.  humidity )  THEN
905             message_string = 'data_output_pr = ' // TRIM( data_output_pr(var_count) ) //       &
906                              ' requires humidity'
907             CALL message( 'check_parameters', 'PA0358', 1, 2, 0, 6, 0 )
908          ENDIF
909          pr_index = 130
910          dopr_index(var_count) = pr_index
911          dopr_unit     = '%'
912          unit = dopr_unit
913          hom(:,2,pr_index,:)  = SPREAD( zu, 2, statistic_regions+1 )
914
915       CASE DEFAULT
916          unit = 'illegal'
917
918    END SELECT
919
920
921 END SUBROUTINE doq_check_data_output_pr
922
923
924!--------------------------------------------------------------------------------------------------!
[4039]925!
926! Description:
927! ------------
928!> Subroutine defining appropriate grid for netcdf variables.
[4583]929!--------------------------------------------------------------------------------------------------!
[4039]930 SUBROUTINE doq_define_netcdf_grid( variable, found, grid_x, grid_y, grid_z )
[4431]931
[3994]932    IMPLICIT NONE
[4039]933
934    CHARACTER (LEN=*), INTENT(IN)  ::  variable    !<
[4431]935    CHARACTER (LEN=*), INTENT(OUT) ::  grid_x      !<
936    CHARACTER (LEN=*), INTENT(OUT) ::  grid_y      !<
937    CHARACTER (LEN=*), INTENT(OUT) ::  grid_z      !<
[4039]938
[4583]939    LOGICAL, INTENT(OUT)           ::  found       !<
940
941
[4039]942    found  = .TRUE.
[4431]943
[4039]944    SELECT CASE ( TRIM( variable ) )
[3995]945!
[4039]946!--    s grid
[4757]947       CASE ( 'rh', 'rh_xy', 'rh_xz', 'rh_yz',                                                     &
[4861]948              'ta', 'ta_xy', 'ta_xz', 'ta_yz',                                                     &
[4757]949              'ti', 'ti_xy', 'ti_xz', 'ti_yz',                                                     &
[4583]950              'wspeed', 'wspeed_xy', 'wspeed_xz', 'wspeed_yz',                                     &
[4560]951              'wdir', 'wdir_xy', 'wdir_xz', 'wdir_yz' )
[3994]952
[4039]953          grid_x = 'x'
954          grid_y = 'y'
[4331]955          grid_z = 'zu'
[4039]956!
[4331]957!--    s grid surface variables
958       CASE ( 'theta_2m*_xy', 'wspeed_10m*' )
959
960          grid_x = 'x'
961          grid_y = 'y'
962          grid_z = 'zu'
963!
[4039]964!--    u grid
965       CASE ( 'uu', 'uu_xy', 'uu_xz', 'uu_yz' )
[3994]966
[4039]967          grid_x = 'xu'
968          grid_y = 'y'
969          grid_z = 'zu'
970!
971!--    v grid
972       CASE ( 'vv', 'vv_xy', 'vv_xz', 'vv_yz'  )
973
974          grid_x = 'x'
975          grid_y = 'yv'
976          grid_z = 'zu'
977!
978!--    w grid
[4583]979       CASE ( 'ww', 'ww_xy', 'ww_xz', 'ww_yz',                                                     &
980              'wu', 'wu_xy', 'wu_xz', 'wu_yz',                                                     &
981              'wv', 'wv_xy', 'wv_xz', 'wv_yz',                                                     &
982              'wtheta', 'wtheta_xy', 'wtheta_xz', 'wtheta_yz',                                     &
[4866]983              'wq', 'wq_xy', 'wq_xz', 'wq_yz',                                                     &
984              'ws', 'ws_xy', 'ws_xz', 'ws_yz' )
[4039]985
986          grid_x = 'x'
987          grid_y = 'y'
988          grid_z = 'zw'
989
990       CASE DEFAULT
991          found  = .FALSE.
992          grid_x = 'none'
993          grid_y = 'none'
994          grid_z = 'none'
995
996    END SELECT
997
998
999 END SUBROUTINE doq_define_netcdf_grid
[4431]1000
[4583]1001!--------------------------------------------------------------------------------------------------!
[4039]1002!
1003! Description:
1004! ------------
1005!> Subroutine defining 2D output variables
[4583]1006!--------------------------------------------------------------------------------------------------!
1007 SUBROUTINE doq_output_2d( av, variable, found, grid, mode, local_pf, two_d, nzb_do, nzt_do,       &
1008                           fill_value )
[4039]1009
1010
1011    IMPLICIT NONE
1012
[4431]1013    CHARACTER (LEN=*) ::  grid     !<
1014    CHARACTER (LEN=*) ::  mode     !<
1015    CHARACTER (LEN=*) ::  variable !<
[4039]1016
1017    INTEGER(iwp) ::  av       !< value indicating averaged or non-averaged output
1018    INTEGER(iwp) ::  flag_nr  !< number of the topography flag (0: scalar, 1: u, 2: v, 3: w)
1019    INTEGER(iwp) ::  i        !< grid index x-direction
1020    INTEGER(iwp) ::  j        !< grid index y-direction
1021    INTEGER(iwp) ::  k        !< grid index z-direction
[4431]1022    INTEGER(iwp) ::  nzb_do   !<
1023    INTEGER(iwp) ::  nzt_do   !<
[4039]1024
1025    LOGICAL ::  found             !< true if variable is in list
1026    LOGICAL ::  resorted          !< true if array is resorted
1027    LOGICAL ::  two_d             !< flag parameter that indicates 2D variables (horizontal cross sections)
1028
1029    REAL(wp) ::  fill_value       !< value for the _FillValue attribute
[4431]1030
1031    REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) ::  local_pf !<
[4039]1032    REAL(wp), DIMENSION(:,:,:), POINTER ::                 to_be_resorted  !< points to array which needs to be resorted for output
[4431]1033
[4583]1034
[4039]1035    flag_nr  = 0
1036    found    = .TRUE.
1037    resorted = .FALSE.
1038    two_d    = .FALSE.
1039
1040    SELECT CASE ( TRIM( variable ) )
1041
[4757]1042       CASE ( 'rh_xy', 'rh_xz', 'rh_yz' )
1043           IF ( av == 0 )  THEN
1044              to_be_resorted => rh
1045           ELSE
1046              IF ( .NOT. ALLOCATED( rh_av ) )  THEN
1047                 ALLOCATE( rh_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1048                 rh_av = REAL( fill_value, KIND = wp )
1049              ENDIF
1050              to_be_resorted => rh_av
1051           ENDIF
1052           flag_nr = 0
1053
1054           IF ( mode == 'xy' )  grid = 'zu'
1055
[4861]1056       CASE ( 'ta_xy', 'ta_xz', 'ta_yz' )
1057           IF ( av == 0 )  THEN
1058              to_be_resorted => ta
1059           ELSE
1060              IF ( .NOT. ALLOCATED( ta_av ) )  THEN
1061                 ALLOCATE( ta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1062                 ta_av = REAL( fill_value, KIND = wp )
1063              ENDIF
1064              to_be_resorted => ta_av
1065           ENDIF
1066           flag_nr = 0
1067
1068           IF ( mode == 'xy' )  grid = 'zu'
1069
[4039]1070       CASE ( 'ti_xy', 'ti_xz', 'ti_yz' )
1071           IF ( av == 0 )  THEN
1072              to_be_resorted => ti
1073           ELSE
[4583]1074              IF ( .NOT. ALLOCATED( ti_av ) )  THEN
[4518]1075                 ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]1076                 ti_av = REAL( fill_value, KIND = wp )
1077              ENDIF
1078              to_be_resorted => ti_av
1079           ENDIF
1080           flag_nr = 0
[4431]1081
[4039]1082           IF ( mode == 'xy' )  grid = 'zu'
[4431]1083
[4039]1084       CASE ( 'uu_xy', 'uu_xz', 'uu_yz' )
1085          IF ( av == 0 )  THEN
1086             to_be_resorted => uu
1087          ELSE
[4583]1088             IF ( .NOT. ALLOCATED( uu_av ) )  THEN
[4518]1089                ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]1090                uu_av = REAL( fill_value, KIND = wp )
1091             ENDIF
1092             to_be_resorted => uu_av
1093          ENDIF
1094          flag_nr = 1
[4431]1095
[4039]1096          IF ( mode == 'xy' )  grid = 'zu'
[4431]1097
[4039]1098       CASE ( 'vv_xy', 'vv_xz', 'vv_yz' )
1099          IF ( av == 0 )  THEN
1100             to_be_resorted => vv
1101          ELSE
[4583]1102             IF ( .NOT. ALLOCATED( vv_av ) )  THEN
[4518]1103                ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]1104                vv_av = REAL( fill_value, KIND = wp )
1105             ENDIF
1106             to_be_resorted => vv_av
1107          ENDIF
1108          flag_nr = 2
[4431]1109
[4039]1110          IF ( mode == 'xy' )  grid = 'zu'
[4431]1111
[4039]1112       CASE ( 'ww_xy', 'ww_xz', 'ww_yz' )
1113          IF ( av == 0 )  THEN
1114             to_be_resorted => ww
1115          ELSE
[4583]1116             IF ( .NOT. ALLOCATED( ww_av ) )  THEN
[4518]1117                ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]1118                ww_av = REAL( fill_value, KIND = wp )
1119             ENDIF
1120             to_be_resorted => ww_av
1121          ENDIF
1122          flag_nr = 3
[4431]1123
[4039]1124          IF ( mode == 'xy' )  grid = 'zw'
[4431]1125
[4351]1126       CASE ( 'wu_xy', 'wu_xz', 'wu_yz' )
1127          IF ( av == 0 )  THEN
1128             to_be_resorted => wu
1129          ELSE
[4583]1130             IF ( .NOT. ALLOCATED( wu_av ) )  THEN
[4518]1131                ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]1132                wu_av = REAL( fill_value, KIND = wp )
1133             ENDIF
1134             to_be_resorted => wu_av
1135          ENDIF
1136          flag_nr = 0
[4431]1137
[4351]1138          IF ( mode == 'xy' )  grid = 'zw'
[4431]1139
[4351]1140       CASE ( 'wv_xy', 'wv_xz', 'wv_yz' )
1141          IF ( av == 0 )  THEN
1142             to_be_resorted => wv
1143          ELSE
[4583]1144             IF ( .NOT. ALLOCATED( wv_av ) )  THEN
[4518]1145                ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]1146                wv_av = REAL( fill_value, KIND = wp )
1147             ENDIF
1148             to_be_resorted => wv_av
1149          ENDIF
1150          flag_nr = 0
[4431]1151
[4351]1152          IF ( mode == 'xy' )  grid = 'zw'
[4431]1153
[4351]1154       CASE ( 'wtheta_xy', 'wtheta_xz', 'wtheta_yz' )
1155          IF ( av == 0 )  THEN
1156             to_be_resorted => wtheta
1157          ELSE
[4583]1158             IF ( .NOT. ALLOCATED( wtheta_av ) )  THEN
[4518]1159                ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]1160                wtheta_av = REAL( fill_value, KIND = wp )
1161             ENDIF
1162             to_be_resorted => wtheta_av
1163          ENDIF
1164          flag_nr = 0
[4431]1165
[4351]1166          IF ( mode == 'xy' )  grid = 'zw'
[4431]1167
[4351]1168       CASE ( 'wq_xy', 'wq_xz', 'wq_yz' )
1169          IF ( av == 0 )  THEN
1170             to_be_resorted => wq
1171          ELSE
[4583]1172             IF ( .NOT. ALLOCATED( wq_av ) )  THEN
[4518]1173                ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]1174                wq_av = REAL( fill_value, KIND = wp )
1175             ENDIF
1176             to_be_resorted => wq_av
1177          ENDIF
1178          flag_nr = 0
[4431]1179
[4351]1180          IF ( mode == 'xy' )  grid = 'zw'
[4431]1181
[4866]1182       CASE ( 'ws_xy', 'ws_xz', 'ws_yz' )
1183          IF ( av == 0 )  THEN
1184             to_be_resorted => ws
1185          ELSE
1186             IF ( .NOT. ALLOCATED( ws_av ) )  THEN
1187                ALLOCATE( ws_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1188                ws_av = REAL( fill_value, KIND = wp )
1189             ENDIF
1190             to_be_resorted => ws_av
1191          ENDIF
1192          flag_nr = 0
1193
1194          IF ( mode == 'xy' )  grid = 'zw'
1195
[4331]1196       CASE ( 'theta_2m*_xy' )        ! 2d-array
1197          IF ( av == 0 )  THEN
1198             DO  i = nxl, nxr
1199                DO  j = nys, nyn
1200                   local_pf(i,j,nzb+1) = pt_2m(j,i)
1201                ENDDO
1202             ENDDO
1203          ELSE
[4583]1204             IF ( .NOT. ALLOCATED( pt_2m_av ) )  THEN
[4331]1205                ALLOCATE( pt_2m_av(nysg:nyng,nxlg:nxrg) )
1206                pt_2m_av = REAL( fill_value, KIND = wp )
1207             ENDIF
1208             DO  i = nxl, nxr
1209                DO  j = nys, nyn
1210                   local_pf(i,j,nzb+1) = pt_2m_av(j,i)
1211                ENDDO
1212             ENDDO
1213          ENDIF
1214          resorted = .TRUE.
1215          two_d    = .TRUE.
1216          grid     = 'zu1'
[4431]1217
[4331]1218       CASE ( 'wspeed_10m*_xy' )        ! 2d-array
1219          IF ( av == 0 )  THEN
1220             DO  i = nxl, nxr
1221                DO  j = nys, nyn
1222                   local_pf(i,j,nzb+1) = uv_10m(j,i)
1223                ENDDO
1224             ENDDO
1225          ELSE
[4583]1226             IF ( .NOT. ALLOCATED( uv_10m_av ) )  THEN
[4331]1227                ALLOCATE( uv_10m_av(nysg:nyng,nxlg:nxrg) )
1228                uv_10m_av = REAL( fill_value, KIND = wp )
1229             ENDIF
1230             DO  i = nxl, nxr
1231                DO  j = nys, nyn
1232                   local_pf(i,j,nzb+1) = uv_10m_av(j,i)
1233                ENDDO
1234             ENDDO
1235          ENDIF
1236          resorted = .TRUE.
1237          two_d    = .TRUE.
1238          grid     = 'zu1'
[4039]1239
[4431]1240       CASE ( 'wspeed_xy', 'wspeed_xz', 'wspeed_yz' )
1241          IF ( av == 0 )  THEN
1242             to_be_resorted => wspeed
1243          ELSE
[4583]1244             IF ( .NOT. ALLOCATED( wspeed_av ) )  THEN
[4518]1245                ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]1246                wspeed_av = REAL( fill_value, KIND = wp )
1247             ENDIF
1248             to_be_resorted => wspeed_av
1249          ENDIF
1250          flag_nr = 0
1251
1252          IF ( mode == 'xy' )  grid = 'zu'
1253
1254       CASE ( 'wdir_xy', 'wdir_xz', 'wdir_yz' )
1255          IF ( av == 0 )  THEN
1256             to_be_resorted => wdir
1257          ELSE
[4583]1258             IF ( .NOT. ALLOCATED( wdir_av ) )  THEN
[4518]1259                ALLOCATE( wdir_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]1260                wdir_av = REAL( fill_value, KIND = wp )
1261             ENDIF
1262             to_be_resorted => wdir_av
1263          ENDIF
1264          flag_nr = 0
1265
1266          IF ( mode == 'xy' )  grid = 'zu'
1267
[4039]1268       CASE DEFAULT
1269          found = .FALSE.
1270          grid  = 'none'
1271
1272    END SELECT
[4431]1273
1274    IF ( found  .AND.  .NOT. resorted )  THEN
[4039]1275       DO  i = nxl, nxr
1276          DO  j = nys, nyn
1277             DO  k = nzb_do, nzt_do
[4583]1278                local_pf(i,j,k) = MERGE( to_be_resorted(k,j,i),                                    &
1279                                         REAL( fill_value, KIND = wp ),                            &
1280                                         BTEST( wall_flags_total_0(k,j,i), flag_nr ) )
[4039]1281             ENDDO
1282          ENDDO
1283       ENDDO
[3994]1284    ENDIF
[4431]1285
[4039]1286 END SUBROUTINE doq_output_2d
[4431]1287
1288
[4583]1289!--------------------------------------------------------------------------------------------------!
[4039]1290!
1291! Description:
1292! ------------
1293!> Subroutine defining 3D output variables
[4583]1294!--------------------------------------------------------------------------------------------------!
1295 SUBROUTINE doq_output_3d( av, variable, found, local_pf, fill_value, nzb_do, nzt_do )
[4431]1296
[4039]1297    IMPLICIT NONE
[3994]1298
[4431]1299    CHARACTER (LEN=*) ::  variable !<
[3994]1300
[4039]1301    INTEGER(iwp) ::  av       !< index indicating averaged or instantaneous output
1302    INTEGER(iwp) ::  flag_nr  !< number of the topography flag (0: scalar, 1: u, 2: v, 3: w)
1303    INTEGER(iwp) ::  i        !< index variable along x-direction
1304    INTEGER(iwp) ::  j        !< index variable along y-direction
1305    INTEGER(iwp) ::  k        !< index variable along z-direction
1306    INTEGER(iwp) ::  nzb_do   !< lower limit of the data output (usually 0)
1307    INTEGER(iwp) ::  nzt_do   !< vertical upper limit of the data output (usually nz_do3d)
[3994]1308
[4039]1309    LOGICAL ::  found             !< true if variable is in list
1310    LOGICAL ::  resorted          !< true if array is resorted
[3994]1311
[4039]1312    REAL(wp) ::  fill_value       !< value for the _FillValue attribute
1313
[4768]1314    REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) ::  local_pf        !<
[4039]1315    REAL(wp), DIMENSION(:,:,:), POINTER ::                 to_be_resorted  !< points to array which needs to be resorted for output
1316
[4583]1317
[4039]1318    flag_nr  = 0
1319    found    = .TRUE.
1320    resorted = .FALSE.
[4431]1321
[4039]1322    SELECT CASE ( TRIM( variable ) )
1323
[4757]1324       CASE ( 'rh' )
1325          IF ( av == 0 )  THEN
1326             to_be_resorted => rh
1327          ELSE
1328             IF ( .NOT. ALLOCATED( rh ) )  THEN
1329                ALLOCATE( rh(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1330                rh = REAL( fill_value, KIND = wp )
1331             ENDIF
1332             to_be_resorted => rh
1333          ENDIF
1334          flag_nr = 0
1335
[4861]1336       CASE ( 'ta' )
1337          IF ( av == 0 )  THEN
1338             to_be_resorted => ta
1339          ELSE
1340             IF ( .NOT. ALLOCATED( ta ) )  THEN
1341                ALLOCATE( ta(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1342                ta = REAL( fill_value, KIND = wp )
1343             ENDIF
1344             to_be_resorted => ta
1345          ENDIF
1346          flag_nr = 0
1347
[4039]1348       CASE ( 'ti' )
1349          IF ( av == 0 )  THEN
1350             to_be_resorted => ti
1351          ELSE
[4583]1352             IF ( .NOT. ALLOCATED( ti_av ) )  THEN
[4518]1353                ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]1354                ti_av = REAL( fill_value, KIND = wp )
1355             ENDIF
1356             to_be_resorted => ti_av
1357          ENDIF
1358          flag_nr = 0
[4431]1359
[4039]1360       CASE ( 'uu' )
1361          IF ( av == 0 )  THEN
1362             to_be_resorted => uu
1363          ELSE
[4583]1364             IF ( .NOT. ALLOCATED( uu_av ) )  THEN
[4518]1365                ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]1366                uu_av = REAL( fill_value, KIND = wp )
1367             ENDIF
1368             to_be_resorted => uu_av
1369          ENDIF
1370          flag_nr = 1
[4431]1371
[4039]1372       CASE ( 'vv' )
1373          IF ( av == 0 )  THEN
1374             to_be_resorted => vv
1375          ELSE
[4583]1376             IF ( .NOT. ALLOCATED( vv_av ) )  THEN
[4518]1377                ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]1378                vv_av = REAL( fill_value, KIND = wp )
1379             ENDIF
1380             to_be_resorted => vv_av
1381          ENDIF
1382          flag_nr = 2
[4431]1383
[4039]1384       CASE ( 'ww' )
1385          IF ( av == 0 )  THEN
1386             to_be_resorted => ww
1387          ELSE
[4583]1388             IF ( .NOT. ALLOCATED( ww_av ) )  THEN
[4518]1389                ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]1390                ww_av = REAL( fill_value, KIND = wp )
1391             ENDIF
1392             to_be_resorted => ww_av
1393          ENDIF
1394          flag_nr = 3
1395
[4351]1396       CASE ( 'wu' )
1397          IF ( av == 0 )  THEN
1398             to_be_resorted => wu
1399          ELSE
[4583]1400             IF ( .NOT. ALLOCATED( wu_av ) )  THEN
[4518]1401                ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]1402                wu_av = REAL( fill_value, KIND = wp )
1403             ENDIF
1404             to_be_resorted => wu_av
1405          ENDIF
1406          flag_nr = 0
1407
1408       CASE ( 'wv' )
1409          IF ( av == 0 )  THEN
1410             to_be_resorted => wv
1411          ELSE
[4583]1412             IF ( .NOT. ALLOCATED( wv_av ) )  THEN
[4518]1413                ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]1414                wv_av = REAL( fill_value, KIND = wp )
1415             ENDIF
1416             to_be_resorted => wv_av
1417          ENDIF
1418          flag_nr = 0
[4431]1419
[4351]1420       CASE ( 'wtheta' )
1421          IF ( av == 0 )  THEN
1422             to_be_resorted => wtheta
1423          ELSE
[4583]1424             IF ( .NOT. ALLOCATED( wtheta_av ) )  THEN
[4518]1425                ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]1426                wtheta_av = REAL( fill_value, KIND = wp )
1427             ENDIF
1428             to_be_resorted => wtheta_av
1429          ENDIF
1430          flag_nr = 0
[4431]1431
[4351]1432       CASE ( 'wq' )
1433          IF ( av == 0 )  THEN
1434             to_be_resorted => wq
1435          ELSE
[4583]1436             IF ( .NOT. ALLOCATED( wq_av ) )  THEN
[4518]1437                ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]1438                wq_av = REAL( fill_value, KIND = wp )
1439             ENDIF
1440             to_be_resorted => wq_av
1441          ENDIF
1442          flag_nr = 0
1443
[4866]1444       CASE ( 'ws' )
1445          IF ( av == 0 )  THEN
1446             to_be_resorted => ws
1447          ELSE
1448             IF ( .NOT. ALLOCATED( ws_av ) )  THEN
1449                ALLOCATE( ws_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1450                ws_av = REAL( fill_value, KIND = wp )
1451             ENDIF
1452             to_be_resorted => ws_av
1453          ENDIF
1454          flag_nr = 0
1455
[4431]1456       CASE ( 'wspeed' )
1457          IF ( av == 0 )  THEN
1458             to_be_resorted => wspeed
1459          ELSE
[4583]1460             IF ( .NOT. ALLOCATED( wspeed_av ) )  THEN
[4518]1461                ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]1462                wspeed_av = REAL( fill_value, KIND = wp )
1463             ENDIF
1464             to_be_resorted => wspeed_av
1465          ENDIF
1466          flag_nr = 0
1467
1468       CASE ( 'wdir' )
1469          IF ( av == 0 )  THEN
1470             to_be_resorted => wdir
1471          ELSE
[4583]1472             IF ( .NOT. ALLOCATED( wdir_av ) )  THEN
[4518]1473                ALLOCATE( wdir_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]1474                wdir_av = REAL( fill_value, KIND = wp )
1475             ENDIF
1476             to_be_resorted => wdir_av
1477          ENDIF
1478          flag_nr = 0
1479
[4039]1480       CASE DEFAULT
1481          found = .FALSE.
1482
1483    END SELECT
[4431]1484
1485    IF ( found  .AND.  .NOT. resorted )  THEN
[4039]1486       DO  i = nxl, nxr
1487          DO  j = nys, nyn
1488             DO  k = nzb_do, nzt_do
[4583]1489                local_pf(i,j,k) = MERGE( to_be_resorted(k,j,i),                                    &
1490                                         REAL( fill_value, KIND = wp ),                            &
1491                                         BTEST( wall_flags_total_0(k,j,i), flag_nr ) )
[4039]1492             ENDDO
1493          ENDDO
1494       ENDDO
1495    ENDIF
1496
1497 END SUBROUTINE doq_output_3d
[4431]1498
[4583]1499
1500!--------------------------------------------------------------------------------------------------!
1501!
[3994]1502! Description:
1503! ------------
[4583]1504!> Resorts the user-defined output quantity with indices (k,j,i) to a temporary array with indices
1505!> (i,j,k) for masked data output.
1506!--------------------------------------------------------------------------------------------------!
[4069]1507 SUBROUTINE doq_output_mask( av, variable, found, local_pf, mid )
[4431]1508
[4039]1509    USE control_parameters
[4431]1510
[4039]1511    USE indices
[4167]1512
[4039]1513    IMPLICIT NONE
[3994]1514
[4583]1515    REAL(wp), PARAMETER   ::  fill_value = -9999.0_wp       !< value for the _FillValue attribute
1516
[4167]1517    CHARACTER (LEN=*) ::  variable   !<
1518    CHARACTER (LEN=5) ::  grid       !< flag to distinquish between staggered grids
[3994]1519
[4167]1520    INTEGER(iwp) ::  av              !< index indicating averaged or instantaneous output
1521    INTEGER(iwp) ::  flag_nr         !< number of the topography flag (0: scalar, 1: u, 2: v, 3: w)
1522    INTEGER(iwp) ::  i               !< index variable along x-direction
1523    INTEGER(iwp) ::  j               !< index variable along y-direction
1524    INTEGER(iwp) ::  k               !< index variable along z-direction
1525    INTEGER(iwp) ::  im              !< loop index for masked variables
1526    INTEGER(iwp) ::  jm              !< loop index for masked variables
1527    INTEGER(iwp) ::  kk              !< masked output index variable along z-direction
1528    INTEGER(iwp) ::  mid             !< masked output running index
1529    INTEGER(iwp) ::  ktt             !< k index of highest horizontal surface
[3994]1530
[4167]1531    LOGICAL      ::  found           !< true if variable is in list
1532    LOGICAL      ::  resorted        !< true if array is resorted
[3994]1533
[4583]1534    REAL(wp), DIMENSION(mask_size_l(mid,1),mask_size_l(mid,2),mask_size_l(mid,3)) ::  local_pf   !<
[3994]1535
[4583]1536    REAL(wp), DIMENSION(:,:,:), POINTER  ::  to_be_resorted  !< points to array which needs to be resorted for output
[4167]1537
[4583]1538
[4039]1539    flag_nr  = 0
1540    found    = .TRUE.
1541    resorted = .FALSE.
1542    grid     = 's'
1543
1544    SELECT CASE ( TRIM( variable ) )
1545
[4757]1546       CASE ( 'rh' )
1547          IF ( av == 0 )  THEN
1548             to_be_resorted => rh
1549          ELSE
1550             to_be_resorted => rh_av
1551          ENDIF
1552          grid = 's'
1553          flag_nr = 0
1554
[4861]1555       CASE ( 'ta' )
1556          IF ( av == 0 )  THEN
1557             to_be_resorted => ta
1558          ELSE
1559             to_be_resorted => ta_av
1560          ENDIF
1561          grid = 's'
1562          flag_nr = 0
1563
[4039]1564       CASE ( 'ti' )
1565          IF ( av == 0 )  THEN
1566             to_be_resorted => ti
1567          ELSE
1568             to_be_resorted => ti_av
1569          ENDIF
1570          grid = 's'
1571          flag_nr = 0
[4431]1572
[4039]1573       CASE ( 'uu' )
1574          IF ( av == 0 )  THEN
1575             to_be_resorted => uu
1576          ELSE
1577             to_be_resorted => uu_av
1578          ENDIF
1579          grid = 'u'
1580          flag_nr = 1
[4431]1581
[4039]1582       CASE ( 'vv' )
1583          IF ( av == 0 )  THEN
1584             to_be_resorted => vv
1585          ELSE
1586             to_be_resorted => vv_av
1587          ENDIF
1588          grid = 'v'
1589          flag_nr = 2
[4431]1590
[4039]1591       CASE ( 'ww' )
1592          IF ( av == 0 )  THEN
1593             to_be_resorted => ww
1594          ELSE
1595             to_be_resorted => ww_av
1596          ENDIF
1597          grid = 'w'
1598          flag_nr = 3
[4431]1599
[4351]1600       CASE ( 'wu' )
1601          IF ( av == 0 )  THEN
1602             to_be_resorted => wu
1603          ELSE
1604             to_be_resorted => wu_av
1605          ENDIF
1606          grid = 's'
1607          flag_nr = 0
[4431]1608
[4351]1609       CASE ( 'wv' )
1610          IF ( av == 0 )  THEN
1611             to_be_resorted => wv
1612          ELSE
1613             to_be_resorted => wv_av
1614          ENDIF
1615          grid = 's'
1616          flag_nr = 0
[4431]1617
[4351]1618       CASE ( 'wtheta' )
1619          IF ( av == 0 )  THEN
1620             to_be_resorted => wtheta
1621          ELSE
1622             to_be_resorted => wtheta_av
1623          ENDIF
1624          grid = 's'
1625          flag_nr = 0
[4431]1626
[4351]1627       CASE ( 'wq' )
1628          IF ( av == 0 )  THEN
1629             to_be_resorted => wq
1630          ELSE
1631             to_be_resorted => wq_av
1632          ENDIF
1633          grid = 's'
1634          flag_nr = 0
[4039]1635
[4866]1636       CASE ( 'ws' )
1637          IF ( av == 0 )  THEN
1638             to_be_resorted => ws
1639          ELSE
1640             to_be_resorted => ws_av
1641          ENDIF
1642          grid = 's'
1643          flag_nr = 0
1644
[4431]1645       CASE ( 'wspeed' )
1646          IF ( av == 0 )  THEN
1647             to_be_resorted => wspeed
1648          ELSE
1649             to_be_resorted => wspeed_av
1650          ENDIF
1651          grid = 's'
1652          flag_nr = 0
1653
1654       CASE ( 'wdir' )
1655          IF ( av == 0 )  THEN
1656             to_be_resorted => wdir
1657          ELSE
1658             to_be_resorted => wdir_av
1659          ENDIF
1660          grid = 's'
1661          flag_nr = 0
1662
[4039]1663       CASE DEFAULT
1664          found = .FALSE.
1665
1666    END SELECT
[4431]1667
[4039]1668    IF ( found  .AND.  .NOT. resorted )  THEN
1669       IF ( .NOT. mask_surface(mid) )  THEN
1670!
1671!--       Default masked output
1672          DO  i = 1, mask_size_l(mid,1)
1673             DO  j = 1, mask_size_l(mid,2)
1674                DO  k = 1, mask_size_l(mid,3)
[4583]1675                   local_pf(i,j,k) = MERGE( to_be_resorted(mask_k(mid,k),                          &
1676                                                           mask_j(mid,j),                          &
1677                                                           mask_i(mid,i)),                         &
1678                                            REAL( fill_value, KIND = wp ),                         &
1679                                            BTEST( wall_flags_total_0(mask_k(mid,k),               &
1680                                                                      mask_j(mid,j),               &
1681                                                                      mask_i(mid,i)),              &
[4431]1682                                                   flag_nr ) )
[4039]1683                ENDDO
1684             ENDDO
1685          ENDDO
1686
1687       ELSE
1688!
1689!--       Terrain-following masked output
1690          DO  i = 1, mask_size_l(mid,1)
1691             DO  j = 1, mask_size_l(mid,2)
1692!
[4167]1693!--             Get k index of the highest terraing surface
1694                im = mask_i(mid,i)
1695                jm = mask_j(mid,j)
[4583]1696                ktt = MINLOC( MERGE( 1, 0, BTEST( wall_flags_total_0(:,jm,im), 5 ) ), DIM=1 ) - 1
[4167]1697                DO  k = 1, mask_size_l(mid,3)
1698                   kk = MIN( ktt+mask_k(mid,k), nzt+1 )
[4039]1699!
[4167]1700!--                Set value if not in building
[4346]1701                   IF ( BTEST( wall_flags_total_0(kk,jm,im), 6 ) )  THEN
[4167]1702                      local_pf(i,j,k) = fill_value
1703                   ELSE
1704                      local_pf(i,j,k) = to_be_resorted(kk,jm,im)
1705                   ENDIF
[4039]1706                ENDDO
1707             ENDDO
1708          ENDDO
1709
1710       ENDIF
1711    ENDIF
[4431]1712
[4039]1713 END SUBROUTINE doq_output_mask
1714
[4583]1715!--------------------------------------------------------------------------------------------------!
[4039]1716! Description:
1717! ------------
1718!> Allocate required arrays
[4583]1719!--------------------------------------------------------------------------------------------------!
[4039]1720 SUBROUTINE doq_init
1721
[3994]1722    IMPLICIT NONE
[4431]1723
[4039]1724    INTEGER(iwp) ::  ivar   !< loop index over all 2d/3d/mask output quantities
[4157]1725
[4583]1726
[4039]1727!
1728!-- Next line is to avoid compiler warnings about unused variables
1729    IF ( timestep_number_at_prev_calc == 0 )  CONTINUE
[4157]1730!
1731!-- Preparatory steps and initialization of output arrays
1732    IF ( .NOT.  prepared_diagnostic_output_quantities )  CALL doq_prepare
[3994]1733
[4039]1734    initialized_diagnostic_output_quantities = .FALSE.
[4431]1735
[4039]1736    ivar = 1
1737
[4431]1738    DO  WHILE ( ivar <= SIZE( do_all ) )
1739
[4039]1740       SELECT CASE ( TRIM( do_all(ivar) ) )
1741!
[4757]1742!--       Allocate array for 'relative humidity'
1743          CASE ( 'rh' )
1744             IF ( .NOT. ALLOCATED( rh ) )  THEN
1745                ALLOCATE( rh(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1746                rh = 0.0_wp
1747             ENDIF
1748!
[4861]1749!--       Allocate array for 'air temperature'
1750          CASE ( 'ta' )
1751             IF ( .NOT. ALLOCATED( ta ) )  THEN
1752                ALLOCATE( ta(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1753                ta = 0.0_wp
1754             ENDIF
1755!
[4039]1756!--       Allocate array for 'turbulence intensity'
1757          CASE ( 'ti' )
1758             IF ( .NOT. ALLOCATED( ti ) )  THEN
[4518]1759                ALLOCATE( ti(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]1760                ti = 0.0_wp
1761             ENDIF
1762!
1763!--       Allocate array for uu
1764          CASE ( 'uu' )
1765             IF ( .NOT. ALLOCATED( uu ) )  THEN
[4518]1766                ALLOCATE( uu(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]1767                uu = 0.0_wp
1768             ENDIF
1769!
1770!--       Allocate array for vv
1771          CASE ( 'vv' )
1772             IF ( .NOT. ALLOCATED( vv ) )  THEN
[4518]1773                ALLOCATE( vv(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]1774                vv = 0.0_wp
1775             ENDIF
1776!
1777!--       Allocate array for ww
1778          CASE ( 'ww' )
1779             IF ( .NOT. ALLOCATED( ww ) )  THEN
[4518]1780                ALLOCATE( ww(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4039]1781                ww = 0.0_wp
1782             ENDIF
[4331]1783!
[4351]1784!--       Allocate array for wu
1785          CASE ( 'wu' )
1786             IF ( .NOT. ALLOCATED( wu ) )  THEN
[4518]1787                ALLOCATE( wu(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]1788                wu = 0.0_wp
1789             ENDIF
1790!
1791!--       Allocate array for wv
1792          CASE ( 'wv' )
1793             IF ( .NOT. ALLOCATED( wv ) )  THEN
[4518]1794                ALLOCATE( wv(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]1795                wv = 0.0_wp
1796             ENDIF
1797!
1798!--       Allocate array for wtheta
1799          CASE ( 'wtheta' )
1800             IF ( .NOT. ALLOCATED( wtheta ) )  THEN
[4518]1801                ALLOCATE( wtheta(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]1802                wtheta = 0.0_wp
1803             ENDIF
1804!
1805!--       Allocate array for wq
1806          CASE ( 'wq' )
1807             IF ( .NOT. ALLOCATED( wq ) )  THEN
[4518]1808                ALLOCATE( wq(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4351]1809                wq = 0.0_wp
1810             ENDIF
1811!
[4866]1812!--       Allocate array for ws
1813          CASE ( 'ws' )
1814             IF ( .NOT. ALLOCATED( ws ) )  THEN
1815                ALLOCATE( ws(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
1816                ws = 0.0_wp
1817             ENDIF
1818!
[4331]1819!--       Allocate array for 2-m potential temperature
1820          CASE ( 'theta_2m*' )
1821             IF ( .NOT. ALLOCATED( pt_2m ) )  THEN
[4518]1822                ALLOCATE( pt_2m(nysg:nyng,nxlg:nxrg) )
[4331]1823                pt_2m = 0.0_wp
1824             ENDIF
1825!
1826!--       Allocate array for 10-m wind speed
1827          CASE ( 'wspeed_10m*' )
1828             IF ( .NOT. ALLOCATED( uv_10m ) )  THEN
[4518]1829                ALLOCATE( uv_10m(nysg:nyng,nxlg:nxrg) )
[4331]1830                uv_10m = 0.0_wp
1831             ENDIF
[4431]1832!
1833!--       Allocate array for wspeed
1834          CASE ( 'wspeed' )
1835             IF ( .NOT. ALLOCATED( wspeed ) )  THEN
[4518]1836                ALLOCATE( wspeed(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]1837                wspeed = 0.0_wp
1838             ENDIF
[4039]1839
[4431]1840!
1841!--       Allocate array for wdir
1842          CASE ( 'wdir' )
1843             IF ( .NOT. ALLOCATED( u_center ) )  THEN
[4518]1844                ALLOCATE( u_center(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]1845                u_center = 0.0_wp
1846             ENDIF
1847             IF ( .NOT. ALLOCATED( v_center ) )  THEN
[4518]1848                ALLOCATE( v_center(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]1849                v_center = 0.0_wp
1850             ENDIF
1851             IF ( .NOT. ALLOCATED( wdir ) )  THEN
[4518]1852                ALLOCATE( wdir(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[4431]1853                wdir = 0.0_wp
1854             ENDIF
1855
[4039]1856       END SELECT
1857
1858       ivar = ivar + 1
1859    ENDDO
1860
1861    initialized_diagnostic_output_quantities = .TRUE.
[4431]1862
[4039]1863 END SUBROUTINE doq_init
1864
1865
1866!--------------------------------------------------------------------------------------------------!
1867! Description:
1868! ------------
[4757]1869!> Calculation of module-specific statistics, i.e. horizontally averaged profiles and time series.
1870!> This is called for every statistic region sr, but at least for the region "total domain" (sr=0).
1871!--------------------------------------------------------------------------------------------------!
1872 SUBROUTINE doq_statistics( mode, sr, tn )
1873
1874
1875    CHARACTER (LEN=*) ::  mode   !<
1876
1877!     INTEGER(iwp) ::  i    !<
1878!     INTEGER(iwp) ::  j    !<
1879!     INTEGER(iwp) ::  k    !<
1880    INTEGER(iwp) ::  sr   !<
1881    INTEGER(iwp) ::  tn   !<
1882!
1883!-- Next line is to avoid compiler warning about unused variables. Please remove.
1884    IF ( sr == 0  .OR.  tn == 0 )  CONTINUE
1885
1886    IF ( mode == 'profiles' )  THEN
1887
1888    ELSEIF ( mode == 'time_series' )  THEN
1889
1890    ENDIF
1891
1892 END SUBROUTINE doq_statistics
1893
1894
1895!--------------------------------------------------------------------------------------------------!
1896! Description:
1897! ------------
[4039]1898!> Calculate diagnostic quantities
1899!--------------------------------------------------------------------------------------------------!
1900 SUBROUTINE doq_calculate
1901
1902    IMPLICIT NONE
1903
[4331]1904    INTEGER(iwp) ::  i          !< grid index x-dimension
[4431]1905    INTEGER(iwp) ::  j          !< grid index y-dimension
[4331]1906    INTEGER(iwp) ::  k          !< grid index z-dimension
[4039]1907    INTEGER(iwp) ::  ivar       !< loop index over all 2d/3d/mask output quantities
[4431]1908
[4757]1909    REAL(wp) ::  temp           !< temperature
1910    REAL(wp) ::  e_s            !< saturation vapor pressure
1911    REAL(wp) ::  q_s            !< saturation mixing ratio
1912
[4331]1913    TYPE(surf_type), POINTER ::  surf     !< surf-type array, used to generalize subroutines
[3994]1914
1915
1916!     CALL cpu_log( log_point(41), 'calculate_quantities', 'start' )
[4157]1917
[3994]1918!
[4583]1919!-- Save timestep number to check in time_integration if doq_calculate has been called already,
1920!-- since the CALL occurs at two locations, but the calculations need to be done only once per
1921!-- timestep.
[3994]1922    timestep_number_at_prev_calc = current_timestep_number
1923
[4039]1924    ivar = 1
[3994]1925
[4431]1926    DO  WHILE ( ivar <= SIZE( do_all ) )
[3994]1927
[4039]1928       SELECT CASE ( TRIM( do_all(ivar) ) )
[3994]1929!
[4757]1930!--       rh (relative humidity)
1931          CASE ( 'rh' )
1932             DO  i = nxl, nxr
1933                DO  j = nys, nyn
1934                   DO  k = nzb+1, nzt
1935                      IF ( humidity )  THEN
1936                         temp = exner(k) * pt(k,j,i)
1937                         e_s = magnus( temp )
1938                         q_s = rd_d_rv * e_s / ( hyp(k) - e_s )
1939                         rh(k,j,i) = q(k,j,i) / q_s * 100.0_wp                                     &
1940                                     * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 1 ) )
1941                      ENDIF
1942                   ENDDO
1943                ENDDO
1944             ENDDO
1945!
[4861]1946!--       ta (air temperature)
1947          CASE ( 'ta' )
1948             DO  i = nxl, nxr
1949                DO  j = nys, nyn
1950                   DO  k = nzb+1, nzt
1951                      ta(k,j,i) = exner(k) * pt(k,j,i) - degc_to_k
1952                   ENDDO
1953                ENDDO
1954             ENDDO
1955!
[3994]1956!--       Calculate 'turbulence intensity' from rot[(u,v,w)] at scalar grid point
1957          CASE ( 'ti' )
1958             DO  i = nxl, nxr
1959                DO  j = nys, nyn
1960                   DO  k = nzb+1, nzt
[4583]1961                      ti(k,j,i) = 0.25_wp * SQRT(                                                  &
1962                                       (   (   w(k,j+1,i) + w(k-1,j+1,i)                           &
1963                                             - w(k,j-1,i) - w(k-1,j-1,i) ) * ddy                   &
1964                                         - (   v(k+1,j,i) + v(k+1,j+1,i)                           &
1965                                             - v(k-1,j,i) - v(k-1,j+1,i) ) * ddzu(k) )**2          &
1966                                     + (   (   u(k+1,j,i) + u(k+1,j,i+1)                           &
1967                                             - u(k-1,j,i) - u(k-1,j,i+1) ) * ddzu(k)               &
1968                                         - (   w(k,j,i+1) + w(k-1,j,i+1)                           &
1969                                             - w(k,j,i-1) - w(k-1,j,i-1) ) * ddx     )**2          &
1970                                     + (   (   v(k,j,i+1) + v(k,j+1,i+1)                           &
1971                                             - v(k,j,i-1) - v(k,j+1,i-1) ) * ddx                   &
1972                                         - (   u(k,j+1,i) + u(k,j+1,i+1)                           &
1973                                             - u(k,j-1,i) - u(k,j-1,i+1) ) * ddy     )**2  )       &
1974                                  * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
[3994]1975                   ENDDO
1976                ENDDO
[4431]1977             ENDDO
[4039]1978!
1979!--       uu
1980          CASE ( 'uu' )
1981             DO  i = nxl, nxr
1982                DO  j = nys, nyn
1983                   DO  k = nzb+1, nzt
[4583]1984                      uu(k,j,i) = u(k,j,i) * u(k,j,i)                                              &
1985                                  * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 1 ) )
[4039]1986                   ENDDO
1987                ENDDO
[3994]1988             ENDDO
[4039]1989!
1990!--       vv
1991          CASE ( 'vv' )
1992             DO  i = nxl, nxr
1993                DO  j = nys, nyn
1994                   DO  k = nzb+1, nzt
[4583]1995                      vv(k,j,i) = v(k,j,i) * v(k,j,i)                                              &
1996                                  * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 2 ) )
[4039]1997                   ENDDO
1998                ENDDO
1999             ENDDO
2000!
2001!--       ww
2002          CASE ( 'ww' )
2003             DO  i = nxl, nxr
2004                DO  j = nys, nyn
2005                   DO  k = nzb+1, nzt-1
[4583]2006                      ww(k,j,i) = w(k,j,i) * w(k,j,i)                                              &
2007                                  * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 3 ) )
[4039]2008                   ENDDO
2009                ENDDO
2010             ENDDO
[4331]2011!
[4351]2012!--       wu
2013          CASE ( 'wu' )
2014             DO  i = nxl, nxr
2015                DO  j = nys, nyn
2016                   DO  k = nzb+1, nzt-1
[4583]2017                      wu(k,j,i) = w(k,j,i)                                                         &
2018                                  * 0.25_wp * ( u(k,j,i) + u(k,j,i+1) + u(k+1,j,i) + u(k+1,j,i+1) )&
2019                                  * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
[4351]2020                   ENDDO
2021                ENDDO
2022             ENDDO
2023!
2024!--       wv
2025          CASE ( 'wv' )
2026             DO  i = nxl, nxr
2027                DO  j = nys, nyn
2028                   DO  k = nzb+1, nzt-1
[4583]2029                      wv(k,j,i) = w(k,j,i)                                                         &
2030                                  * 0.25_wp * ( v(k,j,i) + v(k,j+1,i) + v(k+1,j,i) + v(k+1,j+1,i) )&
2031                                  * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) )
[4351]2032                   ENDDO
2033                ENDDO
2034             ENDDO
2035!
2036!--       wtheta
2037          CASE ( 'wtheta' )
2038             DO  i = nxl, nxr
2039                DO  j = nys, nyn
2040                   DO  k = nzb+1, nzt-1
[4583]2041                      wtheta(k,j,i) = w(k,j,i)                                                     &
2042                                     *  0.5_wp  * ( pt(k,j,i) + pt(k+1,j,i) )                      &
2043                                     * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 3 ))
[4351]2044                   ENDDO
2045                ENDDO
2046             ENDDO
2047!
2048!--       wq
2049          CASE ( 'wq' )
2050             DO  i = nxl, nxr
2051                DO  j = nys, nyn
2052                   DO  k = nzb+1, nzt-1
[4583]2053                      wq(k,j,i) = w(k,j,i) * 0.5_wp * ( q(k,j,i) + q(k+1,j,i) )                    &
2054                                  * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 3 ) )
[4351]2055                   ENDDO
2056                ENDDO
2057             ENDDO
2058!
[4866]2059!--       ws
2060          CASE ( 'ws' )
2061             DO  i = nxl, nxr
2062                DO  j = nys, nyn
2063                   DO  k = nzb+1, nzt-1
2064                      ws(k,j,i) = w(k,j,i) * 0.5_wp * ( s(k,j,i) + s(k+1,j,i) )                    &
2065                                  * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 3 ) )
2066                   ENDDO
2067                ENDDO
2068             ENDDO
2069!
[4331]2070!--       2-m potential temperature
2071          CASE ( 'theta_2m*' )
2072!
[4583]2073!--          2-m potential temperature is caluclated from surface arrays. In case the 2m level is
2074!--          below the first grid point, MOST is applied, else, linear interpolation between two
2075!--          vertical grid levels is applied. To access all surfaces, iterate over all horizontally-
[4331]2076!--          upward facing surface types.
2077             surf => surf_def_h(0)
2078             CALL calc_pt_2m
[4671]2079             surf => surf_lsm_h(0)
[4331]2080             CALL calc_pt_2m
[4671]2081             surf => surf_usm_h(0)
[4331]2082             CALL calc_pt_2m
2083!
2084!--       10-m wind speed
2085          CASE ( 'wspeed_10m*' )
2086!
[4583]2087!--          10-m wind speed is caluclated from surface arrays. In case the 10m level is below the
2088!--          first grid point, MOST is applied, else, linear interpolation between two vertical grid
2089!--          levels is applied. To access all surfaces, iterate over all horizontally-upward facing
2090!--          surface types.
[4331]2091             surf => surf_def_h(0)
2092             CALL calc_wind_10m
[4671]2093             surf => surf_lsm_h(0)
[4331]2094             CALL calc_wind_10m
[4671]2095             surf => surf_usm_h(0)
[4331]2096             CALL calc_wind_10m
[4431]2097!
[4583]2098!--       Horizontal wind speed
[4431]2099          CASE ( 'wspeed' )
2100             DO  i = nxl, nxr
2101                DO  j = nys, nyn
2102                   DO  k = nzb, nzt+1
[4583]2103                      wspeed(k,j,i) = SQRT( ( 0.5_wp * ( u(k,j,i) + u(k,j,i+1) ) )**2              &
2104                                          + ( 0.5_wp * ( v(k,j,i) + v(k,j+1,i) ) )**2 )            &
2105                                     * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ))
[4431]2106                   ENDDO
2107                ENDDO
2108             ENDDO
[3994]2109
[4431]2110!
[4583]2111!--       Horizontal wind direction
[4431]2112          CASE ( 'wdir' )
2113             DO  i = nxl, nxr
2114                DO  j = nys, nyn
2115                   DO  k = nzb, nzt+1
2116                      u_center(k,j,i) = 0.5_wp * ( u(k,j,i) + u(k,j,i+1) )
2117                      v_center(k,j,i) = 0.5_wp * ( v(k,j,i) + v(k,j+1,i) )
2118
[4583]2119                      wdir(k,j,i) = ATAN2( u_center(k,j,i), v_center(k,j,i) )                      &
2120                                    / pi * 180.0_wp + 180.0_wp
[4431]2121                   ENDDO
2122                ENDDO
2123             ENDDO
2124
[4039]2125       END SELECT
[3994]2126
[4039]2127       ivar = ivar + 1
[3994]2128    ENDDO
2129
2130!     CALL cpu_log( log_point(41), 'calculate_quantities', 'stop' )
2131
[4331]2132!
[4583]2133!-- The following block contains subroutines to calculate diagnostic surface quantities.
[4331]2134    CONTAINS
[4583]2135!--------------------------------------------------------------------------------------------------!
[4331]2136! Description:
2137! ------------
2138!> Calculation of 2-m potential temperature.
[4583]2139!--------------------------------------------------------------------------------------------------!
[4331]2140       SUBROUTINE calc_pt_2m
2141
[4583]2142          USE surface_layer_fluxes_mod,                                                            &
[4331]2143              ONLY:  psi_h
2144
2145          IMPLICIT NONE
2146
2147          INTEGER(iwp) ::  kk     !< running index along the z-dimension
2148          INTEGER(iwp) ::  m      !< running index for surface elements
[4431]2149
[4331]2150          DO  m = 1, surf%ns
2151
2152             i = surf%i(m)
2153             j = surf%j(m)
2154             k = surf%k(m)
2155!
[4583]2156!--          If 2-m level is below the first grid level, MOST is used for calculation of
2157!--          2-m temperature.
[4331]2158             IF ( surf%z_mo(m) > 2.0_wp )  THEN
[4583]2159                pt_2m(j,i) = surf%pt_surface(m) + surf%ts(m) / kappa                               &
2160                                * ( LOG( 2.0_wp / surf%z0h(m) )                                    &
2161                                    - psi_h( 2.0_wp      / surf%ol(m) )                            &
2162                                    + psi_h( surf%z0h(m) / surf%ol(m) ) )
[4331]2163!
[4583]2164!--          If 2-m level is above the first grid level, 2-m temperature is linearly interpolated
2165!--          between the two nearest vertical grid levels. Note, since 2-m temperature is only
2166!--          computed for horizontal upward-facing surfaces, only a vertical interpolation is
2167!--          necessary.
[4331]2168             ELSE
2169!
[4431]2170!--             zw(k-1) defines the height of the surface.
[4331]2171                kk = k
2172                DO WHILE ( zu(kk) - zw(k-1) < 2.0_wp  .AND.  kk <= nzt )
[4431]2173                   kk = kk + 1
[4331]2174                ENDDO
2175!
[4431]2176!--             kk defines the index of the first grid level >= 2m.
[4583]2177                pt_2m(j,i) = pt(kk-1,j,i) +                                                        &
2178                              ( zw(k-1) + 2.0_wp - zu(kk-1)     ) *                                &
2179                              ( pt(kk,j,i)       - pt(kk-1,j,i) ) /                                &
[4331]2180                              ( zu(kk)           - zu(kk-1)     )
2181             ENDIF
2182
2183          ENDDO
2184
2185       END SUBROUTINE calc_pt_2m
[4431]2186
[4583]2187!--------------------------------------------------------------------------------------------------!
[4331]2188! Description:
2189! ------------
2190!> Calculation of 10-m wind speed.
[4583]2191!--------------------------------------------------------------------------------------------------!
[4331]2192       SUBROUTINE calc_wind_10m
2193
[4583]2194          USE surface_layer_fluxes_mod,                                                            &
[4331]2195              ONLY:  psi_m
2196
2197          IMPLICIT NONE
2198
2199          INTEGER(iwp) ::  kk     !< running index along the z-dimension
2200          INTEGER(iwp) ::  m      !< running index for surface elements
2201
2202          REAL(wp) ::  uv_l !< wind speed at lower grid point
2203          REAL(wp) ::  uv_u !< wind speed at upper grid point
[4431]2204
[4583]2205
[4331]2206          DO  m = 1, surf%ns
2207
2208             i = surf%i(m)
2209             j = surf%j(m)
2210             k = surf%k(m)
2211!
[4583]2212!--          If 10-m level is below the first grid level, MOST is used for calculation of 10-m
2213!--          temperature.
[4331]2214             IF ( surf%z_mo(m) > 10.0_wp )  THEN
[4583]2215                uv_10m(j,i) = surf%us(m) / kappa                                                   &
2216                              * ( LOG( 10.0_wp /  surf%z0(m) )                                     &
2217                                  - psi_m( 10.0_wp    / surf%ol(m) )                               &
2218                                  + psi_m( surf%z0(m) / surf%ol(m) ) )
[4331]2219!
[4583]2220!--          If 10-m level is above the first grid level, 10-m wind speed is linearly interpolated
2221!--          between the two nearest vertical grid levels. Note, since 10-m temperature is only
2222!--          computed for horizontal upward-facing surfaces, only a vertical interpolation is
2223!--          necessary.
[4331]2224             ELSE
2225!
[4431]2226!--             zw(k-1) defines the height of the surface.
[4331]2227                kk = k
2228                DO WHILE ( zu(kk) - zw(k-1) < 10.0_wp  .AND.  kk <= nzt )
[4431]2229                   kk = kk + 1
[4331]2230                ENDDO
2231!
2232!--             kk defines the index of the first grid level >= 10m.
[4583]2233                uv_l = SQRT( ( 0.5_wp * ( u(kk-1,j,i) + u(kk-1,j,i+1) ) )**2                       &
[4331]2234                           + ( 0.5_wp * ( v(kk-1,j,i) + v(kk-1,j+1,i) ) )**2 )
2235
[4583]2236                uv_u = SQRT( ( 0.5_wp * ( u(kk,j,i)   + u(kk,j,i+1)   ) )**2                       &
[4331]2237                           + ( 0.5_wp * ( v(kk,j,i)   + v(kk,j+1,i)   ) )**2 )
2238
[4583]2239                uv_10m(j,i) = uv_l + ( zw(k-1) + 10.0_wp - zu(kk-1) ) *                            &
2240                                     ( uv_u              - uv_l     ) /                            &
[4331]2241                                     ( zu(kk)            - zu(kk-1) )
2242
2243             ENDIF
2244
2245          ENDDO
2246
2247       END SUBROUTINE calc_wind_10m
2248
[4039]2249 END SUBROUTINE doq_calculate
[3994]2250
2251
[4583]2252!--------------------------------------------------------------------------------------------------!
[3994]2253! Description:
2254! ------------
[4583]2255!> Preparation of the diagnostic output, counting of the module-specific output quantities and
2256!> gathering of the output names.
2257!--------------------------------------------------------------------------------------------------!
[4039]2258 SUBROUTINE doq_prepare
[3994]2259
[4583]2260    USE control_parameters,                                                                        &
[4069]2261        ONLY:  do2d, do3d, domask, masks
[3994]2262
2263    IMPLICIT NONE
2264
[4583]2265    CHARACTER (LEN=varnamelength), DIMENSION(0:1,500) ::  do2d_var = ' '  !< label array for 2d output quantities
[3994]2266
2267    INTEGER(iwp) ::  av         !< index defining type of output, av=0 instantaneous, av=1 averaged
2268    INTEGER(iwp) ::  ivar       !< loop index
2269    INTEGER(iwp) ::  ivar_all   !< loop index
2270    INTEGER(iwp) ::  l          !< index for cutting string
[4431]2271    INTEGER(iwp) ::  mid          !< masked output running index
[3994]2272
[4583]2273
[3994]2274    prepared_diagnostic_output_quantities = .FALSE.
2275
2276    ivar     = 1
2277    ivar_all = 1
2278
2279    DO  av = 0, 1
2280!
2281!--    Remove _xy, _xz, or _yz from string
2282       l = MAX( 3, LEN_TRIM( do2d(av,ivar) ) )
[3998]2283       do2d_var(av,ivar)(1:l-3) = do2d(av,ivar)(1:l-3)
[3994]2284!
[4039]2285!--    Gather 2d output quantity names.
[4583]2286!--    Check for double occurrence of output quantity, e.g. by _xy, _yz, _xz.
[3994]2287       DO  WHILE ( do2d_var(av,ivar)(1:1) /= ' ' )
2288          IF ( .NOT.  ANY( do_all == do2d_var(av,ivar) ) )  THEN
2289             do_all(ivar_all) = do2d_var(av,ivar)
2290          ENDIF
2291          ivar = ivar + 1
2292          ivar_all = ivar_all + 1
2293          l = MAX( 3, LEN_TRIM( do2d(av,ivar) ) )
[3998]2294          do2d_var(av,ivar)(1:l-3) = do2d(av,ivar)(1:l-3)
[3994]2295       ENDDO
2296
2297       ivar = 1
2298!
[4039]2299!--    Gather 3d output quantity names
[3994]2300       DO  WHILE ( do3d(av,ivar)(1:1) /= ' ' )
2301          do_all(ivar_all) = do3d(av,ivar)
2302          ivar = ivar + 1
2303          ivar_all = ivar_all + 1
2304       ENDDO
2305
2306       ivar = 1
2307!
[4583]2308!--    Gather masked output quantity names. Also check for double output e.g. by different masks.
[3994]2309       DO  mid = 1, masks
2310          DO  WHILE ( domask(mid,av,ivar)(1:1) /= ' ' )
[4039]2311             IF ( .NOT.  ANY( do_all == domask(mid,av,ivar) ) )  THEN
[4132]2312                do_all(ivar_all) = domask(mid,av,ivar)
[4039]2313             ENDIF
[3994]2314
2315             ivar = ivar + 1
2316             ivar_all = ivar_all + 1
2317          ENDDO
2318          ivar = 1
2319       ENDDO
2320
2321    ENDDO
2322
2323    prepared_diagnostic_output_quantities = .TRUE.
2324
[4039]2325 END SUBROUTINE doq_prepare
[4431]2326
[4583]2327!--------------------------------------------------------------------------------------------------!
[4039]2328! Description:
2329! ------------
2330!> Subroutine reads local (subdomain) restart data
[4583]2331!--------------------------------------------------------------------------------------------------!
2332 SUBROUTINE doq_rrd_local_ftn( k, nxlf, nxlc, nxl_on_file, nxrf, nxrc, nxr_on_file, nynf, nync,    &
2333                               nyn_on_file, nysf, nysc, nys_on_file, tmp_2d, tmp_3d, found )
[4431]2334
[4518]2335    USE control_parameters
2336
2337    USE indices
2338
2339    USE kinds
2340
2341    USE pegrid
2342
2343
2344    IMPLICIT NONE
2345
2346    INTEGER(iwp) ::  k               !<
2347    INTEGER(iwp) ::  nxlc            !<
2348    INTEGER(iwp) ::  nxlf            !<
2349    INTEGER(iwp) ::  nxl_on_file     !<
2350    INTEGER(iwp) ::  nxrc            !<
2351    INTEGER(iwp) ::  nxrf            !<
2352    INTEGER(iwp) ::  nxr_on_file     !<
2353    INTEGER(iwp) ::  nync            !<
2354    INTEGER(iwp) ::  nynf            !<
2355    INTEGER(iwp) ::  nyn_on_file     !<
2356    INTEGER(iwp) ::  nysc            !<
2357    INTEGER(iwp) ::  nysf            !<
2358    INTEGER(iwp) ::  nys_on_file     !<
2359
2360    LOGICAL, INTENT(OUT)  :: found
2361
2362    REAL(wp), DIMENSION(nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_2d   !<
2363
2364    REAL(wp), DIMENSION(nzb:nzt+1,nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_3d   !<
2365
2366
2367    found = .TRUE.
2368
2369    SELECT CASE ( restart_string(1:length) )
2370
2371       CASE ( 'pt_2m_av' )
2372          IF ( .NOT. ALLOCATED( pt_2m_av ) )  THEN
2373             ALLOCATE( pt_2m_av(nysg:nyng,nxlg:nxrg) )
2374          ENDIF
2375          IF ( k == 1 )  READ ( 13 )  tmp_2d
2376          pt_2m_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  =                                     &
2377                      tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
[4583]2378
[4757]2379       CASE ( 'rh_av' )
2380          IF ( .NOT. ALLOCATED( rh_av ) )  THEN
2381             ALLOCATE( rh_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2382          ENDIF
2383          IF ( k == 1 )  READ ( 13 )  tmp_3d
2384          rh_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                       &
2385                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2386
[4861]2387       CASE ( 'ta_av' )
2388          IF ( .NOT. ALLOCATED( ta_av ) )  THEN
2389             ALLOCATE( ta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2390          ENDIF
2391          IF ( k == 1 )  READ ( 13 )  tmp_3d
2392          ta_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                       &
2393                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2394
[4518]2395       CASE ( 'ti_av' )
2396          IF ( .NOT. ALLOCATED( ti_av ) )  THEN
2397             ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2398          ENDIF
2399          IF ( k == 1 )  READ ( 13 )  tmp_3d
2400          ti_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                       &
2401                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2402
2403       CASE ( 'u_center_av' )
2404          IF ( .NOT. ALLOCATED( u_center_av ) )  THEN
2405             ALLOCATE( u_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2406          ENDIF
2407          IF ( k == 1 )  READ ( 13 )  tmp_3d
2408          u_center_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                 &
2409                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2410
2411       CASE ( 'uu_av' )
2412          IF ( .NOT. ALLOCATED( uu_av ) )  THEN
2413             ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2414          ENDIF
2415          IF ( k == 1 )  READ ( 13 )  tmp_3d
2416          uu_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                       &
2417                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2418
2419       CASE ( 'uv_10m_av' )
2420          IF ( .NOT. ALLOCATED( uv_10m_av ) )  THEN
2421             ALLOCATE( uv_10m_av(nysg:nyng,nxlg:nxrg) )
2422          ENDIF
2423          IF ( k == 1 )  READ ( 13 )  tmp_2d
2424          uv_10m_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  =                                    &
2425                       tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2426
2427       CASE ( 'v_center_av' )
2428          IF ( .NOT. ALLOCATED( v_center_av ) )  THEN
2429             ALLOCATE( v_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2430          ENDIF
2431          IF ( k == 1 )  READ ( 13 )  tmp_3d
2432          v_center_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                 &
2433                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2434
2435       CASE ( 'vv_av' )
2436          IF ( .NOT. ALLOCATED( vv_av ) )  THEN
2437             ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2438          ENDIF
2439          IF ( k == 1 )  READ ( 13 )  tmp_3d
2440          vv_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                       &
2441                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2442
2443       CASE ( 'wtheta_av' )
2444          IF ( .NOT. ALLOCATED( wtheta_av ) )  THEN
2445             ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2446          ENDIF
2447          IF ( k == 1 )  READ ( 13 )  tmp_3d
2448          wtheta_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                   &
2449                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2450
2451       CASE ( 'wq_av' )
2452          IF ( .NOT. ALLOCATED( wq_av ) )  THEN
2453             ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2454          ENDIF
2455          IF ( k == 1 )  READ ( 13 )  tmp_3d
2456          wq_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                       &
2457                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2458
[4866]2459       CASE ( 'ws_av' )
2460          IF ( .NOT. ALLOCATED( ws_av ) )  THEN
2461             ALLOCATE( ws_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2462          ENDIF
2463          IF ( k == 1 )  READ ( 13 )  tmp_3d
2464          ws_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                       &
2465                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2466
[4518]2467       CASE ( 'wspeed_av' )
2468          IF ( .NOT. ALLOCATED( wspeed_av ) )  THEN
2469             ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2470          ENDIF
2471          IF ( k == 1 )  READ ( 13 )  tmp_3d
2472          wspeed_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                   &
2473                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2474
2475       CASE ( 'wu_av' )
2476          IF ( .NOT. ALLOCATED( wu_av ) )  THEN
2477             ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2478          ENDIF
2479          IF ( k == 1 )  READ ( 13 )  tmp_3d
2480          wu_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                       &
2481                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2482
2483       CASE ( 'wv_av' )
2484          IF ( .NOT. ALLOCATED( wv_av ) )  THEN
2485             ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2486          ENDIF
2487          IF ( k == 1 )  READ ( 13 )  tmp_3d
2488          wv_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                       &
2489                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2490
2491       CASE ( 'ww_av' )
2492          IF ( .NOT. ALLOCATED( ww_av ) )  THEN
2493             ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2494          ENDIF
2495          IF ( k == 1 )  READ ( 13 )  tmp_3d
2496          ww_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =                                       &
2497                      tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
2498
2499
2500       CASE DEFAULT
2501
2502          found = .FALSE.
2503
2504    END SELECT
2505
2506 END SUBROUTINE doq_rrd_local_ftn
2507
[4583]2508!--------------------------------------------------------------------------------------------------!
[4039]2509! Description:
2510! ------------
[4518]2511!> Read module-specific local restart data arrays (MPI-IO).
[4583]2512!--------------------------------------------------------------------------------------------------!
[4518]2513 SUBROUTINE doq_rrd_local_mpi
2514
2515    LOGICAL ::  array_found  !< control flad indicating whether the array is found in the file or not
2516
2517
2518    CALL rd_mpi_io_check_array( 'pt_2m_av' , found = array_found )
2519    IF ( array_found )  THEN
2520       IF ( .NOT. ALLOCATED( pt_2m_av ) )  ALLOCATE( pt_2m_av(nysg:nyng,nxlg:nxrg) )
2521       CALL rrd_mpi_io( 'pt_2m_av', pt_2m_av )
2522    ENDIF
2523
[4861]2524    CALL rd_mpi_io_check_array( 'ta_av' , found = array_found )
2525    IF ( array_found )  THEN
2526       IF ( .NOT. ALLOCATED( ta_av ) )  ALLOCATE( ta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2527       CALL rrd_mpi_io( 'ta_av', ta_av )
2528    ENDIF
2529
[4518]2530    CALL rd_mpi_io_check_array( 'ti_av' , found = array_found )
2531    IF ( array_found )  THEN
2532       IF ( .NOT. ALLOCATED( ti_av ) )  ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2533       CALL rrd_mpi_io( 'ti_av', ti_av )
2534    ENDIF
2535
[4757]2536    CALL rd_mpi_io_check_array( 'rh_av' , found = array_found )
2537    IF ( array_found )  THEN
2538       IF ( .NOT. ALLOCATED( rh_av ) )  ALLOCATE( ti_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2539       CALL rrd_mpi_io( 'rh_av', rh_av )
2540    ENDIF
2541
[4518]2542    CALL rd_mpi_io_check_array( 'u_center_av' , found = array_found )
2543    IF ( array_found )  THEN
2544       IF ( .NOT. ALLOCATED( u_center_av ) )  ALLOCATE( u_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2545       CALL rrd_mpi_io( 'u_center_av', u_center_av )
2546    ENDIF
2547
2548    CALL rd_mpi_io_check_array( 'uu_av' , found = array_found )
2549    IF ( array_found )  THEN
2550       IF ( .NOT. ALLOCATED( uu_av ) )  ALLOCATE( uu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2551       CALL rrd_mpi_io( 'uu_av', uu_av )
2552    ENDIF
2553
2554    CALL rd_mpi_io_check_array( 'uv_10m_av' , found = array_found )
2555    IF ( array_found )  THEN
2556       IF ( .NOT. ALLOCATED( uv_10m_av ) )  ALLOCATE( uv_10m_av(nysg:nyng,nxlg:nxrg) )
2557       CALL rrd_mpi_io( 'uv_10m_av', uv_10m_av )
2558    ENDIF
2559
2560    CALL rd_mpi_io_check_array( 'v_center_av' , found = array_found )
2561    IF ( array_found )  THEN
2562       IF ( .NOT. ALLOCATED( v_center_av ) )  ALLOCATE( v_center_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2563       CALL rrd_mpi_io( 'v_center_av', v_center_av )
2564    ENDIF
2565
2566    CALL rd_mpi_io_check_array( 'vv_av' , found = array_found )
2567    IF ( array_found )  THEN
2568       IF ( .NOT. ALLOCATED( vv_av ) )  ALLOCATE( vv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2569       CALL rrd_mpi_io( 'vv_av', vv_av )
2570    ENDIF
2571
2572    CALL rd_mpi_io_check_array( 'ww_av' , found = array_found )
2573    IF ( array_found )  THEN
2574       IF ( .NOT. ALLOCATED( ww_av ) )  ALLOCATE( ww_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2575       CALL rrd_mpi_io( 'ww_av', ww_av )
2576    ENDIF
2577
2578    CALL rd_mpi_io_check_array( 'wu_av' , found = array_found )
2579    IF ( array_found )  THEN
2580       IF ( .NOT. ALLOCATED( wu_av ) )  ALLOCATE( wu_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2581       CALL rrd_mpi_io( 'wu_av', wu_av )
2582    ENDIF
2583
2584    CALL rd_mpi_io_check_array( 'wv_av' , found = array_found )
2585    IF ( array_found )  THEN
2586       IF ( .NOT. ALLOCATED( wv_av ) )  ALLOCATE( wv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2587       CALL rrd_mpi_io( 'wv_av', wv_av )
2588    ENDIF
2589
2590    CALL rd_mpi_io_check_array( 'wtheta_av' , found = array_found )
2591    IF ( array_found )  THEN
2592       IF ( .NOT. ALLOCATED( wtheta_av ) )  ALLOCATE( wtheta_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2593       CALL rrd_mpi_io( 'wtheta_av', wtheta_av )
2594    ENDIF
2595
2596    CALL rd_mpi_io_check_array( 'wq_av' , found = array_found )
2597    IF ( array_found )  THEN
2598       IF ( .NOT. ALLOCATED( wq_av ) )  ALLOCATE( wq_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2599       CALL rrd_mpi_io( 'wq_av', wq_av )
2600    ENDIF
2601
[4866]2602    CALL rd_mpi_io_check_array( 'ws_av' , found = array_found )
2603    IF ( array_found )  THEN
2604       IF ( .NOT. ALLOCATED( ws_av ) )  ALLOCATE( ws_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2605       CALL rrd_mpi_io( 'ws_av', ws_av )
2606    ENDIF
2607
[4518]2608    CALL rd_mpi_io_check_array( 'wspeed_av' , found = array_found )
2609    IF ( array_found )  THEN
2610       IF ( .NOT. ALLOCATED( wspeed_av ) )  ALLOCATE( wspeed_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
2611       CALL rrd_mpi_io( 'wspeed_av', wspeed_av )
2612    ENDIF
2613
2614 END SUBROUTINE doq_rrd_local_mpi
2615
[4583]2616!--------------------------------------------------------------------------------------------------!
[4518]2617! Description:
2618! ------------
[4039]2619!> Subroutine writes local (subdomain) restart data
[4583]2620!--------------------------------------------------------------------------------------------------!
[4039]2621 SUBROUTINE doq_wrd_local
[3994]2622
[4039]2623    IMPLICIT NONE
2624
[4583]2625
[4517]2626    IF ( TRIM( restart_data_format_output ) == 'fortran_binary' )  THEN
[4331]2627
[4517]2628       IF ( ALLOCATED( pt_2m_av ) )  THEN
2629          CALL wrd_write_string( 'pt_2m_av' )
2630          WRITE ( 14 )  pt_2m_av
2631       ENDIF
[4331]2632
[4757]2633       IF ( ALLOCATED( rh_av ) )  THEN
2634          CALL wrd_write_string( 'rh_av' )
2635          WRITE ( 14 )  rh_av
2636       ENDIF
2637
[4861]2638       IF ( ALLOCATED( ta_av ) )  THEN
2639          CALL wrd_write_string( 'ta_av' )
2640          WRITE ( 14 )  ta_av
2641       ENDIF
2642
[4517]2643       IF ( ALLOCATED( ti_av ) )  THEN
2644          CALL wrd_write_string( 'ti_av' )
2645          WRITE ( 14 )  ti_av
2646       ENDIF
[4331]2647
[4518]2648       IF ( ALLOCATED( u_center_av ) )  THEN
2649          CALL wrd_write_string( 'u_center_av' )
2650          WRITE ( 14 )  u_center_av
2651       ENDIF
2652
[4517]2653       IF ( ALLOCATED( uu_av ) )  THEN
2654          CALL wrd_write_string( 'uu_av' )
2655          WRITE ( 14 )  uu_av
2656       ENDIF
[4331]2657
[4517]2658       IF ( ALLOCATED( uv_10m_av ) )  THEN
2659          CALL wrd_write_string( 'uv_10m_av' )
2660          WRITE ( 14 )  uv_10m_av
2661       ENDIF
[4331]2662
[4518]2663       IF ( ALLOCATED( v_center_av ) )  THEN
2664          CALL wrd_write_string( 'v_center_av' )
2665          WRITE ( 14 )  v_center_av
2666       ENDIF
2667
[4517]2668       IF ( ALLOCATED( vv_av ) )  THEN
2669          CALL wrd_write_string( 'vv_av' )
2670          WRITE ( 14 )  vv_av
2671       ENDIF
[4039]2672
[4517]2673       IF ( ALLOCATED( ww_av ) )  THEN
2674          CALL wrd_write_string( 'ww_av' )
2675          WRITE ( 14 )  ww_av
2676       ENDIF
[4039]2677
[4517]2678       IF ( ALLOCATED( wu_av ) )  THEN
2679          CALL wrd_write_string( 'wu_av' )
2680          WRITE ( 14 )  wu_av
2681       ENDIF
[4431]2682
[4517]2683       IF ( ALLOCATED( wv_av ) )  THEN
2684          CALL wrd_write_string( 'wv_av' )
2685          WRITE ( 14 )  wv_av
2686       ENDIF
[4431]2687
[4517]2688       IF ( ALLOCATED( wtheta_av ) )  THEN
2689          CALL wrd_write_string( 'wtheta_av' )
2690          WRITE ( 14 )  wtheta_av
2691       ENDIF
[4351]2692
[4517]2693       IF ( ALLOCATED( wq_av ) )  THEN
2694          CALL wrd_write_string( 'wq_av' )
2695          WRITE ( 14 )  wq_av
2696       ENDIF
[4431]2697
[4866]2698       IF ( ALLOCATED( ws_av ) )  THEN
2699          CALL wrd_write_string( 'ws_av' )
2700          WRITE ( 14 )  ws_av
2701       ENDIF
2702
[4517]2703       IF ( ALLOCATED( wspeed_av ) )  THEN
2704          CALL wrd_write_string( 'wspeed_av' )
2705          WRITE ( 14 )  wspeed_av
2706       ENDIF
[4431]2707
[4535]2708    ELSEIF ( restart_data_format_output(1:3) == 'mpi' )  THEN
[4517]2709
[4518]2710       IF ( ALLOCATED( pt_2m_av ) )     CALL wrd_mpi_io( 'pt_2m_av', pt_2m_av )
[4757]2711       IF ( ALLOCATED( rh_av ) )        CALL wrd_mpi_io( 'rh_av', rh_av )
[4861]2712       IF ( ALLOCATED( ta_av ) )        CALL wrd_mpi_io( 'ta_av', ta_av )
[4518]2713       IF ( ALLOCATED( ti_av ) )        CALL wrd_mpi_io( 'ti_av', ti_av )
2714       IF ( ALLOCATED( u_center_av ) )  CALL wrd_mpi_io( 'u_center_av', u_center_av )
2715       IF ( ALLOCATED( uu_av ) )        CALL wrd_mpi_io( 'uu_av', uu_av )
2716       IF ( ALLOCATED( uv_10m_av ) )    CALL wrd_mpi_io( 'uv_10m_av', uv_10m_av )
2717       IF ( ALLOCATED( vv_av ) )        CALL wrd_mpi_io( 'vv_av', vv_av )
2718       IF ( ALLOCATED( v_center_av ) )  CALL wrd_mpi_io( 'v_center_av', v_center_av )
2719       IF ( ALLOCATED( wtheta_av ) )    CALL wrd_mpi_io( 'wtheta_av', wtheta_av )
2720       IF ( ALLOCATED( wq_av ) )        CALL wrd_mpi_io( 'wq_av', wq_av )
[4866]2721       IF ( ALLOCATED( ws_av ) )        CALL wrd_mpi_io( 'ws_av', ws_av )
[4518]2722       IF ( ALLOCATED( wspeed_av ) )    CALL wrd_mpi_io( 'wspeed_av', wspeed_av )
2723       IF ( ALLOCATED( wu_av ) )        CALL wrd_mpi_io( 'wu_av', wu_av )
2724       IF ( ALLOCATED( wv_av ) )        CALL wrd_mpi_io( 'wv_av', wv_av )
2725       IF ( ALLOCATED( ww_av ) )        CALL wrd_mpi_io( 'ww_av', ww_av )
[4517]2726
[4431]2727    ENDIF
2728
[4039]2729 END SUBROUTINE doq_wrd_local
2730
[4431]2731
[3994]2732 END MODULE diagnostic_output_quantities_mod
Note: See TracBrowser for help on using the repository browser.