source: palm/trunk/SOURCE/data_output_3d.f90 @ 2224

Last change on this file since 2224 was 2210, checked in by kanani, 8 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 25.4 KB
RevLine 
[1682]1!> @file data_output_3d.f90
[2000]2!------------------------------------------------------------------------------!
[1036]3! This file is part of PALM.
4!
[2000]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.
[1036]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!
[2101]17! Copyright 1997-2017 Leibniz Universitaet Hannover
[2000]18!------------------------------------------------------------------------------!
[1036]19!
[254]20! Current revisions:
[1106]21! ------------------
[1961]22!
[2210]23!
[1552]24! Former revisions:
25! -----------------
26! $Id: data_output_3d.f90 2210 2017-04-19 09:35:50Z suehring $
27!
[2210]28! 2209 2017-04-19 09:34:46Z kanani
29! Added plant canopy model output
30!
[2032]31! 2031 2016-10-21 15:11:58Z knoop
32! renamed variable rho to rho_ocean and rho_av to rho_ocean_av
33!
[2012]34! 2011 2016-09-19 17:29:57Z kanani
35! Flag urban_surface is now defined in module control_parameters,
36! changed prefix for urban surface model output to "usm_",
37! introduced control parameter varnamelength for LEN of trimvar.
38!
[2008]39! 2007 2016-08-24 15:47:17Z kanani
40! Added support for new urban surface model (temporary modifications of
41! SELECT CASE ( ) necessary, see variable trimvar)
42!
[2001]43! 2000 2016-08-20 18:09:15Z knoop
44! Forced header and separation lines into 80 columns
45!
[1981]46! 1980 2016-07-29 15:51:57Z suehring
47! Bugfix, in order to steer user-defined output, setting flag found explicitly
48! to .F.
49!
[1977]50! 1976 2016-07-27 13:28:04Z maronga
51! Output of radiation quantities is now done directly in the respective module
52!
[1973]53! 1972 2016-07-26 07:52:02Z maronga
54! Output of land surface quantities is now done directly in the respective module.
55! Unnecessary directive __parallel removed.
56!
[1961]57! 1960 2016-07-12 16:34:24Z suehring
58! Scalar surface flux added
59!
[1851]60! 1849 2016-04-08 11:33:18Z hoffmann
61! prr moved to arrays_3d
62!
[1823]63! 1822 2016-04-07 07:49:42Z hoffmann
64! prr vertical dimensions set to nzb_do to nzt_do. Unused variables deleted.
65!
[1809]66! 1808 2016-04-05 19:44:00Z raasch
67! test output removed
68!
[1784]69! 1783 2016-03-06 18:36:17Z raasch
70! name change of netcdf routines and module + related changes
71!
[1746]72! 1745 2016-02-05 13:06:51Z gronemeier
73! Bugfix: test if time axis limit exceeds moved to point after call of check_open
74!
[1692]75! 1691 2015-10-26 16:17:44Z maronga
76! Added output of radiative heating rates for RRTMG
77!
[1683]78! 1682 2015-10-07 23:56:08Z knoop
79! Code annotations made doxygen readable
80!
[1586]81! 1585 2015-04-30 07:05:52Z maronga
82! Added support for RRTMG
83!
[1552]84! 1551 2015-03-03 14:18:16Z maronga
[1551]85! Added suppport for land surface model and radiation model output. In the course
86! of this action, the limits for vertical loops have been changed (from nzb and
87! nzt+1 to nzb_do and nzt_do, respectively in order to allow soil model output).
88! Moreover, a new vertical grid zs was introduced.
[1329]89!
[1360]90! 1359 2014-04-11 17:15:14Z hoffmann
91! New particle structure integrated.
92!
[1354]93! 1353 2014-04-08 15:21:23Z heinze
94! REAL constants provided with KIND-attribute
95!
[1329]96! 1327 2014-03-21 11:00:16Z raasch
97! parts concerning avs output removed,
98! -netcdf output queries
99!
[1321]100! 1320 2014-03-20 08:40:49Z raasch
[1320]101! ONLY-attribute added to USE-statements,
102! kind-parameters added to all INTEGER and REAL declaration statements,
103! kinds are defined in new module kinds,
104! old module precision_kind is removed,
105! revision history before 2012 removed,
106! comment fields (!:) to be used for variable explanations added to
107! all variable declaration statements
[674]108!
[1319]109! 1318 2014-03-17 13:35:16Z raasch
110! barrier argument removed from cpu_log,
111! module interfaces removed
112!
[1309]113! 1308 2014-03-13 14:58:42Z fricke
114! Check, if the limit of the time dimension is exceeded for parallel output
115! To increase the performance for parallel output, the following is done:
116! - Update of time axis is only done by PE0
117!
[1245]118! 1244 2013-10-31 08:16:56Z raasch
119! Bugfix for index bounds in case of 3d-parallel output
120!
[1116]121! 1115 2013-03-26 18:16:16Z hoffmann
122! ql is calculated by calc_liquid_water_content
123!
[1107]124! 1106 2013-03-04 05:31:38Z raasch
125! array_kind renamed precision_kind
126!
[1077]127! 1076 2012-12-05 08:30:18Z hoffmann
128! Bugfix in output of ql
129!
[1054]130! 1053 2012-11-13 17:11:03Z hoffmann
131! +nr, qr, prr, qc and averaged quantities
132!
[1037]133! 1036 2012-10-22 13:43:42Z raasch
134! code put under GPL (PALM 3.9)
135!
[1035]136! 1031 2012-10-19 14:35:30Z raasch
137! netCDF4 without parallel file support implemented
138!
[1008]139! 1007 2012-09-19 14:30:36Z franke
140! Bugfix: missing calculation of ql_vp added
141!
[1]142! Revision 1.1  1997/09/03 06:29:36  raasch
143! Initial revision
144!
145!
146! Description:
147! ------------
[1682]148!> Output of the 3D-arrays in netCDF and/or AVS format.
[1]149!------------------------------------------------------------------------------!
[1682]150 SUBROUTINE data_output_3d( av )
151 
[1]152
[1320]153    USE arrays_3d,                                                             &
[2031]154        ONLY:  e, nr, p, pt, prr, q, qc, ql, ql_c, ql_v, qr, rho_ocean, s, sa, tend, &
[1960]155               u, v, vpt, w
[1320]156       
[1]157    USE averaging
[1320]158       
159    USE cloud_parameters,                                                      &
[1849]160        ONLY:  l_d_cp, pt_d_t
[1320]161       
162    USE control_parameters,                                                    &
[1822]163        ONLY:  cloud_physics, do3d, do3d_no, do3d_time_count, io_blocks,       &
164               io_group, message_string, ntdim_3d, nz_do3d, psolver,           &
[2011]165               simulated_time, time_since_reference_point, urban_surface,      &
166               varnamelength
[1320]167       
168    USE cpulog,                                                                &
169        ONLY:  log_point, cpu_log
170       
171    USE indices,                                                               &
[1822]172        ONLY:  nbgp, nx, nxl, nxlg, nxr, nxrg, ny, nyn, nyng, nys, nysg, nzb
[1320]173       
174    USE kinds
175   
[1551]176    USE land_surface_model_mod,                                                &
[1972]177        ONLY: land_surface, lsm_data_output_3d, nzb_soil, nzt_soil
[1551]178
[1783]179#if defined( __netcdf )
180    USE NETCDF
181#endif
182
183    USE netcdf_interface,                                                      &
184        ONLY:  id_set_3d, id_var_do3d, id_var_time_3d, nc_stat,                &
185               netcdf_data_format, netcdf_handle_error
[1320]186       
187    USE particle_attributes,                                                   &
[1359]188        ONLY:  grid_particles, number_of_particles, particles,                 &
189               particle_advection_start, prt_count
[1320]190       
[1]191    USE pegrid
192
[2209]193    USE plant_canopy_model_mod,                                                &
194        ONLY:  pcm_data_output_3d, plant_canopy
195       
[1585]196    USE radiation_model_mod,                                                   &
[1976]197        ONLY:  radiation, radiation_data_output_3d
[1585]198
[2007]199    USE urban_surface_mod,                                                     &
[2011]200        ONLY:  nzub, nzut, usm_data_output_3d
[1585]201
[2007]202
[1]203    IMPLICIT NONE
204
[1682]205    INTEGER(iwp) ::  av        !<
206    INTEGER(iwp) ::  i         !<
207    INTEGER(iwp) ::  if        !<
208    INTEGER(iwp) ::  j         !<
209    INTEGER(iwp) ::  k         !<
210    INTEGER(iwp) ::  n         !<
211    INTEGER(iwp) ::  nzb_do    !< vertical lower limit for data output
212    INTEGER(iwp) ::  nzt_do    !< vertical upper limit for data output
[1]213
[1682]214    LOGICAL      ::  found     !<
215    LOGICAL      ::  resorted  !<
[1]216
[1682]217    REAL(wp)     ::  mean_r    !<
218    REAL(wp)     ::  s_r2      !<
219    REAL(wp)     ::  s_r3      !<
[1]220
[1682]221    REAL(sp), DIMENSION(:,:,:), ALLOCATABLE ::  local_pf  !<
[1]222
[1682]223    REAL(wp), DIMENSION(:,:,:), POINTER ::  to_be_resorted  !<
[1]224
[2011]225    CHARACTER (LEN=varnamelength) ::  trimvar  !< TRIM of output-variable string
[2007]226
[1]227!
228!-- Return, if nothing to output
229    IF ( do3d_no(av) == 0 )  RETURN
230
231    CALL cpu_log (log_point(14),'data_output_3d','start')
232
233!
234!-- Open output file.
235!-- Also creates coordinate and fld-file for AVS.
[1031]236!-- For classic or 64bit netCDF output or output of other (old) data formats,
[493]237!-- for a run on more than one PE, each PE opens its own file and
[1]238!-- writes the data of its subdomain in binary format (regardless of the format
239!-- the user has requested). After the run, these files are combined to one
240!-- file by combine_plot_fields in the format requested by the user (netcdf
[493]241!-- and/or avs).
[1031]242!-- For netCDF4/HDF5 output, data is written in parallel into one file.
[1327]243    IF ( netcdf_data_format < 5 )  THEN
244       CALL check_open( 30 )
245       IF ( myid == 0 )  CALL check_open( 106+av*10 )
[493]246    ELSE
[1327]247       CALL check_open( 106+av*10 )
[493]248    ENDIF
[1]249
250!
[1745]251!-- For parallel netcdf output the time axis must be limited. Return, if this
252!-- limit is exceeded. This could be the case, if the simulated time exceeds
253!-- the given end time by the length of the given output interval.
254    IF ( netcdf_data_format > 4 )  THEN
255       IF ( do3d_time_count(av) + 1 > ntdim_3d(av) )  THEN
256          WRITE ( message_string, * ) 'Output of 3d data is not given at t=',  &
257                                      simulated_time, '&because the maximum ', & 
258                                      'number of output time levels is ',      &
259                                      'exceeded.'
260          CALL message( 'data_output_3d', 'PA0387', 0, 1, 0, 6, 0 )
261          CALL cpu_log( log_point(14), 'data_output_3d', 'stop' )
262          RETURN
263       ENDIF
264    ENDIF
265
266!
[1031]267!-- Update the netCDF time axis
[1308]268!-- In case of parallel output, this is only done by PE0 to increase the
269!-- performance.
[1]270#if defined( __netcdf )
[1308]271    do3d_time_count(av) = do3d_time_count(av) + 1
272    IF ( myid == 0 )  THEN
[1327]273       nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_time_3d(av),           &
274                               (/ time_since_reference_point /),            &
275                               start = (/ do3d_time_count(av) /),           &
276                               count = (/ 1 /) )
[1783]277       CALL netcdf_handle_error( 'data_output_3d', 376 )
[1]278    ENDIF
279#endif
280
281!
282!-- Loop over all variables to be written.
283    if = 1
284
285    DO  WHILE ( do3d(av,if)(1:1) /= ' ' )
[2007]286
[1]287!
[2007]288!--    Temporary solution to account for data output within the new urban
289!--    surface model (urban_surface_mod.f90), see also SELECT CASE ( trimvar ).
[1]290!--    Store the array chosen on the temporary array.
[2007]291       trimvar = TRIM( do3d(av,if) )
[2011]292       IF ( urban_surface  .AND.  trimvar(1:4) == 'usm_' )  THEN
[2007]293          trimvar = 'usm_output'
294          resorted = .TRUE.
295          nzb_do   = nzub
296          nzt_do   = nzut
297       ELSE
298          resorted = .FALSE.
299          nzb_do   = nzb
300          nzt_do   = nz_do3d
301       ENDIF
[1551]302!
[1980]303!--    Set flag to steer output of radiation, land-surface, or user-defined
304!--    quantities
305       found = .FALSE.
306!
[1551]307!--    Allocate a temporary array with the desired output dimensions.
308       ALLOCATE( local_pf(nxlg:nxrg,nysg:nyng,nzb_do:nzt_do) )
309
[2007]310       SELECT CASE ( trimvar )
[1]311
312          CASE ( 'e' )
313             IF ( av == 0 )  THEN
314                to_be_resorted => e
315             ELSE
316                to_be_resorted => e_av
317             ENDIF
318
[771]319          CASE ( 'lpt' )
320             IF ( av == 0 )  THEN
321                to_be_resorted => pt
322             ELSE
323                to_be_resorted => lpt_av
324             ENDIF
325
[1053]326          CASE ( 'nr' )
327             IF ( av == 0 )  THEN
328                to_be_resorted => nr
329             ELSE
330                to_be_resorted => nr_av
331             ENDIF
332
[1]333          CASE ( 'p' )
334             IF ( av == 0 )  THEN
[727]335                IF ( psolver /= 'sor' )  CALL exchange_horiz( p, nbgp )
[1]336                to_be_resorted => p
337             ELSE
[727]338                IF ( psolver /= 'sor' )  CALL exchange_horiz( p_av, nbgp )
[1]339                to_be_resorted => p_av
340             ENDIF
341
342          CASE ( 'pc' )  ! particle concentration (requires ghostpoint exchange)
343             IF ( av == 0 )  THEN
[1359]344                IF ( simulated_time >= particle_advection_start )  THEN
345                   tend = prt_count
346                   CALL exchange_horiz( tend, nbgp )
347                ELSE
348                   tend = 0.0_wp
349                ENDIF
[667]350                DO  i = nxlg, nxrg
351                   DO  j = nysg, nyng
[1551]352                      DO  k = nzb_do, nzt_do
[1]353                         local_pf(i,j,k) = tend(k,j,i)
354                      ENDDO
355                   ENDDO
356                ENDDO
357                resorted = .TRUE.
358             ELSE
[667]359                CALL exchange_horiz( pc_av, nbgp )
[1]360                to_be_resorted => pc_av
361             ENDIF
362
[1359]363          CASE ( 'pr' )  ! mean particle radius (effective radius)
[1]364             IF ( av == 0 )  THEN
[1359]365                IF ( simulated_time >= particle_advection_start )  THEN
366                   DO  i = nxl, nxr
367                      DO  j = nys, nyn
[1551]368                         DO  k = nzb_do, nzt_do
[1359]369                            number_of_particles = prt_count(k,j,i)
370                            IF (number_of_particles <= 0)  CYCLE
371                            particles => grid_particles(k,j,i)%particles(1:number_of_particles)
372                            s_r2 = 0.0_wp
373                            s_r3 = 0.0_wp
374                            DO  n = 1, number_of_particles
375                               IF ( particles(n)%particle_mask )  THEN
376                                  s_r2 = s_r2 + particles(n)%radius**2 * &
377                                         particles(n)%weight_factor
378                                  s_r3 = s_r3 + particles(n)%radius**3 * &
379                                         particles(n)%weight_factor
380                               ENDIF
381                            ENDDO
382                            IF ( s_r2 > 0.0_wp )  THEN
383                               mean_r = s_r3 / s_r2
384                            ELSE
385                               mean_r = 0.0_wp
386                            ENDIF
387                            tend(k,j,i) = mean_r
[1]388                         ENDDO
389                      ENDDO
390                   ENDDO
[1359]391                   CALL exchange_horiz( tend, nbgp )
392                ELSE
393                   tend = 0.0_wp
394                ENDIF
[667]395                DO  i = nxlg, nxrg
396                   DO  j = nysg, nyng
[1551]397                      DO  k = nzb_do, nzt_do
[1]398                         local_pf(i,j,k) = tend(k,j,i)
399                      ENDDO
400                   ENDDO
401                ENDDO
402                resorted = .TRUE.
403             ELSE
[667]404                CALL exchange_horiz( pr_av, nbgp )
[1]405                to_be_resorted => pr_av
406             ENDIF
407
[1053]408          CASE ( 'prr' )
409             IF ( av == 0 )  THEN
410                CALL exchange_horiz( prr, nbgp )
411                DO  i = nxlg, nxrg
412                   DO  j = nysg, nyng
[1822]413                      DO  k = nzb_do, nzt_do
[1053]414                         local_pf(i,j,k) = prr(k,j,i)
415                      ENDDO
416                   ENDDO
417                ENDDO
418             ELSE
419                CALL exchange_horiz( prr_av, nbgp )
420                DO  i = nxlg, nxrg
421                   DO  j = nysg, nyng
[1822]422                      DO  k = nzb_do, nzt_do
[1053]423                         local_pf(i,j,k) = prr_av(k,j,i)
424                      ENDDO
425                   ENDDO
426                ENDDO
427             ENDIF
428             resorted = .TRUE.
429
[1]430          CASE ( 'pt' )
431             IF ( av == 0 )  THEN
432                IF ( .NOT. cloud_physics ) THEN
433                   to_be_resorted => pt
434                ELSE
[667]435                   DO  i = nxlg, nxrg
436                      DO  j = nysg, nyng
[1551]437                         DO  k = nzb_do, nzt_do
[1320]438                            local_pf(i,j,k) = pt(k,j,i) + l_d_cp *             &
439                                                          pt_d_t(k) *          &
[1]440                                                          ql(k,j,i)
441                         ENDDO
442                      ENDDO
443                   ENDDO
444                   resorted = .TRUE.
445                ENDIF
446             ELSE
447                to_be_resorted => pt_av
448             ENDIF
449
450          CASE ( 'q' )
451             IF ( av == 0 )  THEN
452                to_be_resorted => q
453             ELSE
454                to_be_resorted => q_av
455             ENDIF
[691]456
[1053]457          CASE ( 'qc' )
[1]458             IF ( av == 0 )  THEN
[1115]459                to_be_resorted => qc
[1]460             ELSE
[1115]461                to_be_resorted => qc_av
[1]462             ENDIF
463
[1053]464          CASE ( 'ql' )
465             IF ( av == 0 )  THEN
[1115]466                to_be_resorted => ql
[1053]467             ELSE
[1115]468                to_be_resorted => ql_av
[1053]469             ENDIF
470
[1]471          CASE ( 'ql_c' )
472             IF ( av == 0 )  THEN
473                to_be_resorted => ql_c
474             ELSE
475                to_be_resorted => ql_c_av
476             ENDIF
477
478          CASE ( 'ql_v' )
479             IF ( av == 0 )  THEN
480                to_be_resorted => ql_v
481             ELSE
482                to_be_resorted => ql_v_av
483             ENDIF
484
485          CASE ( 'ql_vp' )
486             IF ( av == 0 )  THEN
[1359]487                IF ( simulated_time >= particle_advection_start )  THEN
488                   DO  i = nxl, nxr
489                      DO  j = nys, nyn
[1551]490                         DO  k = nzb_do, nzt_do
[1359]491                            number_of_particles = prt_count(k,j,i)
492                            IF (number_of_particles <= 0)  CYCLE
493                            particles => grid_particles(k,j,i)%particles(1:number_of_particles)
494                            DO  n = 1, number_of_particles
495                               IF ( particles(n)%particle_mask )  THEN
496                                  tend(k,j,i) =  tend(k,j,i) +                 &
497                                                 particles(n)%weight_factor /  &
498                                                 prt_count(k,j,i)
499                               ENDIF
500                            ENDDO
[1007]501                         ENDDO
502                      ENDDO
503                   ENDDO
[1359]504                   CALL exchange_horiz( tend, nbgp )
505                ELSE
506                   tend = 0.0_wp
507                ENDIF
[1007]508                DO  i = nxlg, nxrg
509                   DO  j = nysg, nyng
[1551]510                      DO  k = nzb_do, nzt_do
[1007]511                         local_pf(i,j,k) = tend(k,j,i)
512                      ENDDO
513                   ENDDO
514                ENDDO
515                resorted = .TRUE.
[1]516             ELSE
[1007]517                CALL exchange_horiz( ql_vp_av, nbgp )
[1]518                to_be_resorted => ql_vp_av
519             ENDIF
520
[1053]521          CASE ( 'qr' )
522             IF ( av == 0 )  THEN
523                to_be_resorted => qr
524             ELSE
525                to_be_resorted => qr_av
526             ENDIF
527
[1]528          CASE ( 'qv' )
529             IF ( av == 0 )  THEN
[667]530                DO  i = nxlg, nxrg
531                   DO  j = nysg, nyng
[1551]532                      DO  k = nzb_do, nzt_do
[1]533                         local_pf(i,j,k) = q(k,j,i) - ql(k,j,i)
534                      ENDDO
535                   ENDDO
536                ENDDO
537                resorted = .TRUE.
538             ELSE
539                to_be_resorted => qv_av
540             ENDIF
541
[2031]542          CASE ( 'rho_ocean' )
[96]543             IF ( av == 0 )  THEN
[2031]544                to_be_resorted => rho_ocean
[96]545             ELSE
[2031]546                to_be_resorted => rho_ocean_av
[96]547             ENDIF
[691]548
[1]549          CASE ( 's' )
550             IF ( av == 0 )  THEN
[1960]551                to_be_resorted => s
[1]552             ELSE
[355]553                to_be_resorted => s_av
[1]554             ENDIF
[691]555
[96]556          CASE ( 'sa' )
557             IF ( av == 0 )  THEN
558                to_be_resorted => sa
559             ELSE
560                to_be_resorted => sa_av
561             ENDIF
[691]562
[1]563          CASE ( 'u' )
564             IF ( av == 0 )  THEN
565                to_be_resorted => u
566             ELSE
567                to_be_resorted => u_av
568             ENDIF
569
570          CASE ( 'v' )
571             IF ( av == 0 )  THEN
572                to_be_resorted => v
573             ELSE
574                to_be_resorted => v_av
575             ENDIF
576
577          CASE ( 'vpt' )
578             IF ( av == 0 )  THEN
579                to_be_resorted => vpt
580             ELSE
581                to_be_resorted => vpt_av
582             ENDIF
583
584          CASE ( 'w' )
585             IF ( av == 0 )  THEN
586                to_be_resorted => w
587             ELSE
588                to_be_resorted => w_av
589             ENDIF
[2007]590!             
591!--       Block of urban surface model outputs   
592          CASE ( 'usm_output' )
593             CALL usm_data_output_3d( av, do3d(av,if), found, local_pf,     &
594                                         nzb_do, nzt_do )
[1]595
596          CASE DEFAULT
[1972]597
[1]598!
[1972]599!--          Land surface quantity
600             IF ( land_surface )  THEN
601!
602!--             For soil model quantities, it is required to re-allocate local_pf
603                nzb_do = nzb_soil
604                nzt_do = nzt_soil
605
606                DEALLOCATE ( local_pf )
607                ALLOCATE( local_pf(nxlg:nxrg,nysg:nyng,nzb_do:nzt_do) )
608
609                CALL lsm_data_output_3d( av, do3d(av,if), found, local_pf )
610                resorted = .TRUE.
[1976]611
612!
613!--             If no soil model variable was found, re-allocate local_pf
614                IF ( .NOT. found )  THEN
615                   nzb_do = nzb
616                   nzt_do = nz_do3d
617
618                   DEALLOCATE ( local_pf )
619                   ALLOCATE( local_pf(nxlg:nxrg,nysg:nyng,nzb_do:nzt_do) )                 
620                ENDIF
621
[1972]622             ENDIF
623
624!
[1976]625!--          Radiation quantity
626             IF ( .NOT. found  .AND.  radiation )  THEN
627                CALL radiation_data_output_3d( av, do3d(av,if), found,         &
628                                               local_pf )
629                resorted = .TRUE.
630             ENDIF
631
632!
[2209]633!--          Plant canopy model output
634             IF ( .NOT. found  .AND.  plant_canopy )  THEN
635                CALL pcm_data_output_3d( av, do3d(av,if), found,         &
636                                               local_pf )
637                resorted = .TRUE.
638             ENDIF
639
640!
[1]641!--          User defined quantity
[1972]642             IF ( .NOT. found )  THEN
643                CALL user_data_output_3d( av, do3d(av,if), found, local_pf,    &
644                                          nzb_do, nzt_do )
645                resorted = .TRUE.
646             ENDIF
[1]647
[254]648             IF ( .NOT. found )  THEN
[1320]649                message_string =  'no output available for: ' //               &
[274]650                                  TRIM( do3d(av,if) )
[254]651                CALL message( 'data_output_3d', 'PA0182', 0, 0, 0, 6, 0 )
[1]652             ENDIF
653
654       END SELECT
655
656!
657!--    Resort the array to be output, if not done above
658       IF ( .NOT. resorted )  THEN
[667]659          DO  i = nxlg, nxrg
660             DO  j = nysg, nyng
[1551]661                DO  k = nzb_do, nzt_do
[1]662                   local_pf(i,j,k) = to_be_resorted(k,j,i)
663                ENDDO
664             ENDDO
665          ENDDO
666       ENDIF
667
668!
[1327]669!--    Output of the 3D-array
670#if defined( __parallel )
671       IF ( netcdf_data_format < 5 )  THEN
[1]672!
[1327]673!--       Non-parallel netCDF output. Data is output in parallel in
674!--       FORTRAN binary format here, and later collected into one file by
675!--       combine_plot_fields
676          IF ( myid == 0 )  THEN
677             WRITE ( 30 )  time_since_reference_point,                   &
678                           do3d_time_count(av), av
[1]679          ENDIF
[1327]680          DO  i = 0, io_blocks-1
681             IF ( i == io_group )  THEN
[1551]682                WRITE ( 30 )  nxlg, nxrg, nysg, nyng, nzb_do, nzt_do
683                WRITE ( 30 )  local_pf(:,:,nzb_do:nzt_do)
[1327]684             ENDIF
[1972]685
[1327]686             CALL MPI_BARRIER( comm2d, ierr )
[1972]687
[1327]688          ENDDO
[559]689
[1327]690       ELSE
[646]691#if defined( __netcdf )
[493]692!
[1327]693!--       Parallel output in netCDF4/HDF5 format.
694!--       Do not output redundant ghost point data except for the
695!--       boundaries of the total domain.
696          IF ( nxr == nx  .AND.  nyn /= ny )  THEN
697             nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,if),  &
[1551]698                               local_pf(nxl:nxr+1,nys:nyn,nzb_do:nzt_do),    &
699                start = (/ nxl+1, nys+1, nzb_do+1, do3d_time_count(av) /),  &
700                count = (/ nxr-nxl+2, nyn-nys+1, nzt_do-nzb_do+1, 1 /) )
[1327]701          ELSEIF ( nxr /= nx  .AND.  nyn == ny )  THEN
702             nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,if),  &
[1551]703                               local_pf(nxl:nxr,nys:nyn+1,nzb_do:nzt_do),    &
704                start = (/ nxl+1, nys+1, nzb_do+1, do3d_time_count(av) /),  &
705                count = (/ nxr-nxl+1, nyn-nys+2, nzt_do-nzb_do+1, 1 /) )
[1327]706          ELSEIF ( nxr == nx  .AND.  nyn == ny )  THEN
707             nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,if),  &
[1551]708                             local_pf(nxl:nxr+1,nys:nyn+1,nzb_do:nzt_do  ),  &
709                start = (/ nxl+1, nys+1, nzb_do+1, do3d_time_count(av) /),  &
710                count = (/ nxr-nxl+2, nyn-nys+2, nzt_do-nzb_do+1, 1 /) )
[1327]711          ELSE
712             nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,if),  &
[1551]713                                 local_pf(nxl:nxr,nys:nyn,nzb_do:nzt_do),    &
714                start = (/ nxl+1, nys+1, nzb_do+1, do3d_time_count(av) /),  &
715                count = (/ nxr-nxl+1, nyn-nys+1, nzt_do-nzb_do+1, 1 /) )
[1327]716          ENDIF
[1783]717          CALL netcdf_handle_error( 'data_output_3d', 386 )
[646]718#endif
[1327]719       ENDIF
[1]720#else
721#if defined( __netcdf )
[1327]722       nc_stat = NF90_PUT_VAR( id_set_3d(av), id_var_do3d(av,if),        &
[1551]723                         local_pf(nxl:nxr+1,nys:nyn+1,nzb_do:nzt_do),        &
[1327]724                         start = (/ 1, 1, 1, do3d_time_count(av) /),     &
[1551]725                         count = (/ nx+2, ny+2, nzt_do-nzb_do+1, 1 /) )
[1783]726       CALL netcdf_handle_error( 'data_output_3d', 446 )
[1]727#endif
728#endif
729
730       if = if + 1
731
732!
[1551]733!--    Deallocate temporary array
734       DEALLOCATE ( local_pf )
[1]735
[1551]736    ENDDO
[1]737
[1318]738    CALL cpu_log( log_point(14), 'data_output_3d', 'stop' )
[1]739
740!
741!-- Formats.
[1320]7423300 FORMAT ('variable ',I4,'  file=',A,'  filetype=unformatted  skip=',I12/   &
[1]743             'label = ',A,A)
744
745 END SUBROUTINE data_output_3d
Note: See TracBrowser for help on using the repository browser.