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

Last change on this file since 4557 was 4535, checked in by raasch, 5 years ago

bugfix for restart data format query

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