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

Last change on this file since 1980 was 1980, checked in by suehring, 8 years ago

Bugfix, in order to steer user-defined output, setting flag found explicitly

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