source: palm/trunk/UTIL/combine_plot_fields.f90

Last change on this file was 4843, checked in by raasch, 19 months ago

local namelist parameter added to switch off the module although the respective module namelist appears in the namelist file, further copyright updates

  • Property svn:keywords set to Id
File size: 30.2 KB
Line 
1 PROGRAM combine_plot_fields
2
3!------------------------------------------------------------------------------!
4! This file is part of the PALM model system.
5!
6! PALM is free software: you can redistribute it and/or modify it under the
7! terms of the GNU General Public License as published by the Free Software
8! Foundation, either version 3 of the License, or (at your option) any later
9! version.
10!
11! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
12! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
13! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
14!
15! You should have received a copy of the GNU General Public License along with
16! PALM. If not, see <http://www.gnu.org/licenses/>.
17!
18! Copyright 1997-2021  Leibniz Universitaet Hannover
19!------------------------------------------------------------------------------!
20!
21! Current revisions:
22! -----------------
23!
24!
25! Former revisions:
26! -----------------
27! $Id: combine_plot_fields.f90 4843 2021-01-15 15:22:11Z banzhafs $
28! bugfix for processing data from atmosphere-ocean coupled runs (only atmosphere data were
29! processed, but twice)
30!
31! 3241 2018-09-12 15:02:00Z raasch
32! unused variables removed
33!
34! 2718 2018-01-02 08:49:38Z maronga
35! Corrected "Former revisions" section
36!
37! 2696 2017-12-14 17:12:51Z kanani
38! Change in file header (GPL part)
39!
40! 2669 2017-12-06 16:03:27Z raasch
41! data of 3d-nest runs are completely processed now
42!
43! 2523 2017-10-05 14:42:47Z kanani
44! Increased LEN for CHARACTER variable var_name, equal to the value in PALM
45!
46! 2512 2017-10-04 08:26:59Z raasch
47! PALM output does not contain ghost layer data any more
48! avs- and iso2d-related parts removed, handling of compressed data removed
49!
50! 2365 2017-08-21 14:59:59Z kanani
51! Vertical grid nesting implemented (SadiqHuq)
52!
53! 1809 2016-04-05 20:13:28Z raasch
54!
55! 1808 2016-04-05 19:44:00Z raasch
56! cpu measurements are done with standard FORTRAN routine on every machine
57!
58! 1551 2015-03-03 14:18:16Z maronga
59! Adjustments for data output of soil model quantities
60!
61! 1468 2014-09-24 14:06:57Z maronga
62! Adapted for use on up to 6-digit processor cores (not tested)
63!
64! 1394 2014-05-06 10:17:31Z keck
65! KIND-attribute added to all INTEGER and REAL declaration statements
66!
67! 1046 2012-11-09 14:38:45Z maronga
68! code put under GPL (PALM 3.9)
69!
70! 691 2011-03-04 08:45:30Z maronga
71! Bugfix: check for precursor ocean runs, removed typo
72!
73! 493 2010-03-01 08:30:24Z raasch
74! Exit in case of already complete NetCDF data (due to parallel output in PALM)
75! cpu measurements included
76!
77! 210 2008-11-06 08:54:02Z raasch
78! Size of pf3d adjusted to the required output size (1 gridpoint less, along
79! all three dimensions), because output of a subset of the data
80! (pf3d(nxa:nxe...) in the NF90_PUT_VAR statement caused segmentation fault
81! with the INTEL compiler.
82! Subdomain data are read into temporary arrays pf_tmp/pf3d_tmp in order to
83! avoid INTEL compiler warnings about (automatic) creation of temporary arrays
84! Bugfix: three misplaced #endif directives
85!
86! 114 2007-10-10 00:03:15Z raasch
87! Bugfix: model_string needed a default value
88!
89! Aug 07    Loop for processing of output by coupled runs, id_string does not
90!           contain modus any longer
91!
92! 18/01/06  Output of time-averaged data
93!
94! 25/05/05  Errors removed
95!
96! 26/04/05  Output in NetCDF format, iso2d and avs output only if parameter
97!           file exists
98!
99! 31/10/01  All comments and messages translated into English
100!
101! 23/02/99  Keine Bearbeitung komprimierter 3D-Daten
102! Ursprungsversion vom 28/07/97
103!
104!
105! Description:
106! ------------
107! This routine combines data of the PALM-subdomains into one file. In PALM
108! every processor element opens its own file and writes 2D- or 3D-binary-data
109! into it (different files are opened for xy-, xz-, yz-cross-sections and
110! 3D-data). For plotting or analyzing these PE-data have to be collected and
111! to be put into single files, which is done by this routine.
112! Output format is NetCDF. Additionally, a data are output in a binary format
113! readable by ISO2D-software (cross-sections) and by AVS (3D-data).
114!------------------------------------------------------------------------------!
115
116#if defined( __netcdf )
117    USE netcdf
118#endif
119
120    IMPLICIT NONE
121
122!
123!-- Local variables
124    CHARACTER (LEN=2)    ::  modus
125    CHARACTER (LEN=4)    ::  model_string
126    CHARACTER (LEN=6)    ::  id_string
127    CHARACTER (LEN=30)   ::  dimname, var_name
128    CHARACTER (LEN=40)   ::  filename
129
130    CHARACTER (LEN=2000), DIMENSION(0:1) ::  var_list
131
132    INTEGER, PARAMETER ::  iwp = 4, spk = SELECTED_REAL_KIND( 6 ), wp = 8
133
134    INTEGER(iwp) ::  av, danz, i, id, j, k, model, models, nc_stat,            &
135                     nxa, nxe, nya, nye,               &
136                     nza, nze, pos, time_step, xa, xe, xxa, xxe, ya, ya_do,    &
137                     ya_tot, ye, ye_do, ye_tot, yya, yye, za, ze
138
139    INTEGER(8)                        ::  count, count_rate
140
141    INTEGER(iwp), DIMENSION(0:1)      ::  current_level, current_var,          &
142                                          fanz, id_set, id_var_time, num_var
143    INTEGER(iwp), DIMENSION(4)        ::  id_dims_loc
144    INTEGER(iwp), DIMENSION(0:1,4)    ::  id_dims
145    INTEGER(iwp), DIMENSION(0:1,1000) ::  id_var, levels
146
147    LOGICAL  ::  found, found_ao, nest3d, netcdf_output, netcdf_parallel, netcdf_0,      &
148                 netcdf_1, vnest
149
150    REAL(wp) ::  cpu_start_time, cpu_end_time, simulated_time
151    REAL(wp),  DIMENSION(:,:), ALLOCATABLE   ::  pf, pf_tmp
152    REAL(spk), DIMENSION(:,:,:), ALLOCATABLE ::  pf3d, pf3d_tmp
153
154
155
156    PRINT*, ''
157    PRINT*, ''
158    PRINT*, '*** combine_plot_fields ***'
159
160!
161!-- Find out if a coupled or nested run has been carried out
162    INQUIRE( FILE='COUPLING_PORT_OPENED', EXIST=found_ao )
163    INQUIRE( FILE='VNESTING_PORT_OPENED', EXIST=vnest )
164    INQUIRE( FILE='3DNESTING', EXIST=nest3d )
165    IF ( found_ao )  THEN
166       models = 2
167       PRINT*, '    coupled run'
168    ELSEIF ( vnest )  THEN
169       models = 2
170       PRINT*, '    Vertically nested grid coupling'
171    ELSEIF ( nest3d )  THEN
172       OPEN( 90, FILE='3DNESTING', FORM='FORMATTED' )
173       READ ( 90, '(I2)' )  models
174       CLOSE ( 90 )
175       PRINT*, '    3d-nest run'
176       PRINT*, '    number of nest domains = ', models
177    ELSE
178       models = 1
179       PRINT*, '    uncoupled run'
180    ENDIF
181
182!
183!-- Find out if a precursor ocean run has been carried out
184    INQUIRE( FILE='PRECURSOR_OCEAN', EXIST=found )
185    IF ( found )  THEN
186       model_string = '_O'
187       PRINT*, '    precursor ocean run'
188    ELSE
189       model_string = ''
190    ENDIF
191
192!
193!-- Do everything for each model
194    DO model = 1, models
195!
196!--    Set the model string used to identify the filenames
197       IF ( found_ao  .OR.  vnest )  THEN
198          PRINT*, ''
199          IF ( model == 2 )  THEN
200             IF ( vnest )  THEN
201                model_string = '_N'
202                PRINT*, '    now combining FINE data'
203                PRINT*, '    ========================'
204             ELSE
205                model_string = '_O'
206                PRINT*, '    now combining ocean data'
207                PRINT*, '    ========================'
208             ENDIF
209          ELSE
210             IF ( vnest )  THEN
211                PRINT*, '    now combining COARSE data'
212                PRINT*, '    ============================='
213             ELSE
214                PRINT*, '    now combining atmosphere data'
215                PRINT*, '    ============================='
216             ENDIF
217          ENDIF
218       ELSEIF ( nest3d )  THEN
219          PRINT*, ''
220          PRINT*, '--> processing nest id = ', model
221          IF ( model == 1 )  THEN
222             model_string = ''
223          ELSE
224             WRITE( model_string, '(A2,I2.2)' )  '_N', model
225          ENDIF
226       ENDIF
227!
228!--    2D-arrays for ISO2D
229!--    Main loop for the three different cross-sections, starting with
230!--    xy-section
231       modus = 'XY'
232       PRINT*, ''
233       DO  WHILE ( modus == 'XY'  .OR.  modus == 'XZ'  .OR.  modus == 'YZ' )
234!
235!--       Take current time
236          CALL SYSTEM_CLOCK( count, count_rate )
237          cpu_start_time = REAL( count ) / REAL( count_rate )
238
239          netcdf_parallel = .FALSE.
240!
241!--       Check, if file from PE0 exists. If it does not exist, PALM did not
242!--       create any output for this cross-section.
243          danz = 0
244          WRITE (id_string,'(I6.6)')  danz
245          INQUIRE ( &
246               FILE='PLOT2D_'//modus//TRIM( model_string )//'_'//id_string, &
247               EXIST=found )
248!
249!--       Find out the number of files (equal to the number of PEs which
250!--       have been used in PALM) and open them
251          DO  WHILE ( found )
252
253             OPEN ( danz+110, &
254                  FILE='PLOT2D_'//modus//TRIM( model_string )//'_'//id_string, &
255                  FORM='UNFORMATTED' )
256             danz = danz + 1
257             WRITE (id_string,'(I6.6)')  danz
258             INQUIRE ( &
259                  FILE='PLOT2D_'//modus//TRIM( model_string )//'_'//id_string, &
260                  EXIST=found )
261
262          ENDDO
263
264!
265!--       Inquire whether a NetCDF file exists
266          INQUIRE( FILE='DATA_2D_'//modus//'_NETCDF'//TRIM( model_string ), &
267               EXIST=netcdf_0 )
268
269!
270!--       Inquire whether a NetCDF file for time-averaged data exists
271          INQUIRE( FILE='DATA_2D_'//modus//'_AV_NETCDF'//TRIM( model_string ),&
272               EXIST=netcdf_1 )
273
274          IF ( netcdf_0  .OR.  netcdf_1 )  THEN
275             netcdf_output = .TRUE.
276!
277!--          Inquire whether the NetCDF file is already complete (parallel
278!--          output)
279             INQUIRE( FILE='NO_COMBINE_PLOT_FIELDS_'//modus, &
280                      EXIST=netcdf_parallel )
281             IF ( netcdf_parallel )  THEN
282                netcdf_parallel = .TRUE.
283             ELSE
284                netcdf_parallel = .FALSE.
285             ENDIF
286          ELSE
287             netcdf_output = .FALSE.
288          ENDIF
289
290!
291!--       Info-output
292          PRINT*, ''
293#if defined( __netcdf )
294          IF ( netcdf_output )  THEN
295             IF ( netcdf_parallel )  THEN
296             PRINT*, '    NetCDF ' // modus // '-data are in one file ', &
297                          '(NetCDF4-format) - merging not neccessary'
298             ELSE
299                PRINT*, '    NetCDF output enabled'
300             ENDIF
301          ENDIF
302#else
303          IF ( netcdf_output )  THEN
304             PRINT*, '--- Sorry, no NetCDF support on this host'
305             netcdf_output = .FALSE.
306          ENDIF
307#endif
308          IF ( .NOT. netcdf_parallel )  THEN
309             IF ( danz /= 0 )  THEN
310                PRINT*, '    ',modus,'-section:  ', danz, ' file(s) found'
311             ELSE
312                PRINT*, '    no ', modus, '-section data available'
313             ENDIF
314          ENDIF
315
316          IF ( netcdf_output  .AND. .NOT. netcdf_parallel  .AND.  danz /= 0 ) &
317          THEN
318#if defined( __netcdf )
319             DO  av = 0, 1
320
321                IF ( av == 0  .AND.  .NOT.  netcdf_0 )  CYCLE
322                IF ( av == 1  .AND.  .NOT.  netcdf_1 )  CYCLE
323
324!
325!--             Open NetCDF dataset
326                IF ( av == 0 )  THEN
327                   filename = 'DATA_2D_'//modus//'_NETCDF' &
328                        //TRIM( model_string )
329                ELSE
330                   filename = 'DATA_2D_'//modus//'_AV_NETCDF' &
331                        //TRIM( model_string )
332                ENDIF
333                nc_stat = NF90_OPEN( filename, NF90_WRITE, id_set(av) )
334                IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 1 )
335
336!
337!--             Get the list of variables (order of variables corresponds with
338!--             the order of data on the binary file)
339                var_list(av) = ' '    ! GET_ATT does not assign trailing blanks
340                nc_stat = NF90_GET_ATT( id_set(av), NF90_GLOBAL, 'VAR_LIST', &
341                     var_list(av) )
342                IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 2 )
343
344!
345!--             Inquire id of the time coordinate variable
346                nc_stat = NF90_INQ_VARID( id_set(av), 'time', id_var_time(av) )
347                IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 3 )
348
349!
350!--             Count number of variables; there is one more semicolon in the
351!--             string than variable names
352                num_var(av) = -1
353                DO  i = 1, LEN( var_list(av) )
354                   IF ( var_list(av)(i:i) == ';' )  num_var(av) = num_var(av) +1
355                ENDDO
356
357!
358!--             Extract the variable names from the list and inquire their
359!--             NetCDF IDs
360                pos = INDEX( var_list(av), ';' )
361!
362!--             Loop over all variables
363                DO  i = 1, num_var(av)
364
365!
366!--                Extract variable name from list
367                   var_list(av) = var_list(av)(pos+1:)
368                   pos = INDEX( var_list(av), ';' )
369                   var_name = var_list(av)(1:pos-1)
370
371!
372!--                Get variable ID from name
373                   nc_stat = NF90_INQ_VARID( id_set(av), TRIM( var_name ), &
374                        id_var(av,i) )
375                   IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 4 )
376
377!
378!--                Get number of x/y/z levels for that variable
379                   nc_stat = NF90_INQUIRE_VARIABLE( id_set(av), id_var(av,i), &
380                        dimids = id_dims_loc )
381                   IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 5 )
382                   id_dims(av,:) = id_dims_loc
383
384!
385!--                Inquire dimension ID
386                   DO  j = 1, 4
387                      nc_stat = NF90_INQUIRE_DIMENSION( id_set(av), &
388                           id_dims(av,j), dimname, levels(av,i) )
389                      IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 6 )
390
391                      IF ( modus == 'XY' .AND. INDEX(dimname, 'z') /= 0 )  EXIT
392                      IF ( modus == 'XZ' .AND. INDEX(dimname, 'y') /= 0 )  EXIT
393                      IF ( modus == 'YZ' .AND. INDEX(dimname, 'x') /= 0 )  EXIT
394                   ENDDO
395
396                ENDDO
397
398             ENDDO   ! av = 0, 1
399
400#endif
401          ENDIF
402
403!
404!--       Read the arrays, as long as the end of the file is reached
405          IF ( .NOT. netcdf_parallel )  THEN
406
407             fanz          =         0
408             current_level =         1
409             current_var   = 999999999
410
411             DO  WHILE ( danz /= 0 )
412
413!
414!--             Loop over all files (reading data of the subdomains)
415                DO  id = 0, danz-1
416!
417!--                File from PE0 contains special information at the beginning,
418!--                concerning the lower and upper indices of the total-domain
419!--                used in PALM (nxa, nxe, nya, nye).
420!--                Allocate necessary arrays, open the output file and write
421!--                the coordinate informations needed by ISO2D.
422                   IF ( id == 0  .AND.  fanz(0) == 0  .AND.  fanz(1) == 0 ) THEN
423
424                      READ ( id+110 )  nxa, nxe, nya, nye
425                      ALLOCATE ( pf(nxa:nxe,nya:nye) )
426!
427!--                   Set actual domain bounds to total domain
428                      ya_do = nya
429                      ye_do = nye
430
431                   ENDIF
432!
433!--                Read output time
434                   IF ( netcdf_output  .AND.  id == 0 )  THEN
435                      READ ( id+110, END=998 )  simulated_time, time_step, av
436                   ENDIF
437!
438!--                Read subdomain indices
439                   READ ( id+110, END=998 )  xa, xe, ya, ye, ya_tot, ye_tot
440
441!
442!--                IF the PE made no output (in case that no part of the
443!--                cross-section is situated on this PE), indices have the
444!--                value -1
445                   IF ( .NOT. ( xa == -1  .AND.  xe == -1  .AND. &
446                                ya == -1  .AND.  ye == -1 ) )  THEN
447
448
449!
450!--                   Read the subdomain grid-point values
451                      ALLOCATE( pf_tmp(xa:xe,ya:ye) )
452                      READ ( id+110 )  pf_tmp
453!
454!--                   Calculate indices on atmospheric grid (e.g. for soil model
455!--                   quantities)
456                      IF ( ya /= ya_tot .OR. ye /= ye_tot )  THEN
457                         ye_do = ye - ya
458                         ya_do = ya
459                         pf(xa:xe,ya_do:ye_do) = pf_tmp
460                      ELSE
461                         ye_do = nye
462                         ya_do = nya
463                         pf(xa:xe,ya:ye) = pf_tmp
464                      ENDIF
465
466                      DEALLOCATE( pf_tmp )
467                   ENDIF
468                   IF ( id == 0 )  fanz(av) = fanz(av) + 1
469
470                ENDDO
471!
472!--             Write data in NetCDF format
473                IF ( netcdf_output )  THEN
474#if defined( __netcdf )
475!
476!--                Check if a new time step has begun; if yes write data to
477!--                time axis
478                   IF ( current_var(av) > num_var(av) )  THEN
479                      current_var(av) = 1
480                      nc_stat = NF90_PUT_VAR( id_set(av), id_var_time(av), &
481                                              (/ simulated_time /),        &
482                                              start = (/ time_step /),     &
483                                              count = (/ 1 /) )
484                      IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 7 )
485                   ENDIF
486
487!
488!--                Now write the data; this is mode dependent
489                   SELECT CASE ( modus )
490
491                      CASE ( 'XY' )
492                         nc_stat = NF90_PUT_VAR( id_set(av),                   &
493                                           id_var(av,current_var(av)),         &
494                                           pf(nxa:nxe,nya:nye),                &
495                             start = (/ 1, 1, current_level(av), time_step /), &
496                                      count = (/ nxe-nxa+1, nye-nya+1, 1, 1 /) )
497                         IF ( nc_stat /= NF90_NOERR )  THEN
498                            CALL handle_netcdf_error( 8 )
499                         ENDIF
500                 
501                      CASE ( 'XZ' )
502                         nc_stat = NF90_PUT_VAR( id_set(av),                   &
503                                           id_var(av,current_var(av)),         &
504                                           pf(nxa:nxe,ya_do:ye_do),            &
505                             start = (/ 1, current_level(av), 1, time_step /), &
506                                      count = (/ nxe-nxa+1, 1, ye_do-ya_do+1, 1 /) )
507                         IF ( nc_stat /= NF90_NOERR )  THEN
508                            CALL handle_netcdf_error( 9 )
509                         ENDIF
510
511                      CASE ( 'YZ' )
512                         nc_stat = NF90_PUT_VAR( id_set(av),                   &
513                                           id_var(av,current_var(av)),         &
514                                           pf(nxa:nxe,ya_do:ye_do),            &
515                             start = (/ current_level(av), 1, 1, time_step /), &
516                                      count = (/ 1, nxe-nxa+1, ye_do-ya_do+1, 1 /) )
517                         IF ( nc_stat /= NF90_NOERR )  THEN
518                            CALL handle_netcdf_error( 10 )
519                         ENDIF
520
521                   END SELECT
522
523!
524!--                Data is written, check if max level is reached
525                   current_level(av) = current_level(av) + 1
526                   IF ( current_level(av) > levels(av,current_var(av)) )  THEN
527                      current_level(av) = 1
528                      current_var(av)   = current_var(av) + 1
529                   ENDIF
530
531#endif
532                ENDIF
533
534             ENDDO
535
536          ENDIF
537
538998       IF ( danz /= 0  .AND.  .NOT. netcdf_parallel )  THEN
539!
540!--          Print the number of the arrays processed
541             WRITE (*,'(16X,I4,A)')  fanz(0)+fanz(1), ' array(s) processed'
542             IF ( fanz(1) /= 0 )  THEN
543                WRITE (*,'(16X,I4,A)')  fanz(1), ' array(s) are time-averaged'
544             ENDIF
545
546!
547!--          Close all files and deallocate arrays
548             DO  id = 0, danz-1
549                CLOSE ( id+110 )
550             ENDDO
551             CLOSE ( 2 )
552             DEALLOCATE ( pf )
553
554!
555!--          Close the NetCDF file
556             IF ( netcdf_output )  THEN
557#if defined( __netcdf )
558                IF ( netcdf_0 )  THEN
559                   nc_stat = NF90_CLOSE( id_set(0) )
560                   IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 11 )
561                ENDIF
562                IF ( netcdf_1 )  THEN
563                   nc_stat = NF90_CLOSE( id_set(1) )
564                   IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 12 )
565                ENDIF
566#endif
567             ENDIF
568          ENDIF
569
570!
571!--       Output required cpu time
572          IF ( danz /= 0  .AND.  .NOT. netcdf_parallel )  THEN
573             CALL SYSTEM_CLOCK( count, count_rate )
574             cpu_end_time = REAL( count ) / REAL( count_rate )
575             WRITE (*,'(5X,A,F9.3,A)')  'Required cpu-time: ', &
576                                        cpu_end_time-cpu_start_time, ' sec'
577          ENDIF
578
579!
580!--       Choose the next cross-section
581          SELECT CASE ( modus )
582             CASE ( 'XY' )
583                modus = 'XZ'
584             CASE ( 'XZ' )
585                modus = 'YZ'
586             CASE ( 'YZ' )
587                modus = 'no'
588          END SELECT
589
590       ENDDO
591
592
593!
594!--    Combine the 3D-arrays
595       netcdf_parallel = .FALSE.
596
597!
598!--    Info-output
599       PRINT*, ' '
600
601!
602!--    Take current time
603       CALL SYSTEM_CLOCK( count, count_rate )
604       cpu_start_time = REAL( count ) / REAL( count_rate )
605
606!
607!--    Inquire whether a NetCDF file exists
608       INQUIRE( FILE='DATA_3D_NETCDF'//TRIM( model_string ), EXIST=netcdf_0 )
609
610!
611!--    Inquire whether a NetCDF file for time-averaged data exists
612       INQUIRE( FILE='DATA_3D_AV_NETCDF'//TRIM( model_string ), EXIST=netcdf_1 )
613
614       IF ( netcdf_0  .OR.  netcdf_1 )  THEN
615          netcdf_output = .TRUE.
616!
617!--       Inquire whether the NetCDF file is already complete (parallel output)
618          INQUIRE( FILE='NO_COMBINE_PLOT_FIELDS_3D', EXIST=netcdf_parallel )
619          IF ( netcdf_parallel )  THEN
620             netcdf_parallel = .TRUE.
621          ELSE
622             netcdf_parallel = .FALSE.
623          ENDIF
624       ELSE
625          netcdf_output = .FALSE.
626       ENDIF
627
628!
629!--    Check, if file from PE0 exists; not neccessary in case of parallel
630!--    PALM output
631       IF ( .NOT. netcdf_parallel )  THEN
632          danz = 0
633          WRITE (id_string,'(I6.6)')  danz
634          INQUIRE ( &
635               FILE='PLOT3D_DATA'//TRIM( model_string )//'_'//TRIM( id_string ),  &
636               EXIST=found )
637       ELSE
638          found = .FALSE.
639       ENDIF
640
641!
642!--    Find out the number of files and open them
643       DO  WHILE ( found )
644
645          OPEN ( danz+110, &
646               FILE='PLOT3D_DATA'//TRIM( model_string )//'_'//TRIM(id_string), &
647               FORM='UNFORMATTED')
648          danz = danz + 1
649          WRITE (id_string,'(I6.6)')  danz
650          INQUIRE ( &
651               FILE='PLOT3D_DATA'//TRIM( model_string )//'_'//TRIM(id_string), &
652               EXIST=found )
653
654       ENDDO
655
656#if defined( __netcdf )
657       IF ( netcdf_output )  THEN
658          IF ( netcdf_parallel )  THEN
659             PRINT*, '    NetCDF data are in one file (NetCDF4-format)', &
660                          ' - merging not neccessary'
661          ELSE
662             PRINT*, '    NetCDF output enabled'
663          ENDIF
664       ENDIF
665#else
666       IF ( netcdf_output )  THEN
667          PRINT*, '--- Sorry, no NetCDF support on this host'
668          netcdf_output = .FALSE.
669       ENDIF
670#endif
671       IF ( .NOT. netcdf_parallel )  THEN
672          IF ( danz /= 0 )  THEN
673             PRINT*, '    3D-data:     ', danz, ' file(s) found'
674          ELSE
675             PRINT*, '    no 3D-data file available'
676          ENDIF
677       ENDIF
678
679       IF ( netcdf_output  .AND. .NOT. netcdf_parallel  .AND.  danz /= 0 )  THEN
680#if defined( __netcdf )
681          DO  av = 0, 1
682
683             IF ( av == 0  .AND.  .NOT.  netcdf_0 )  CYCLE
684             IF ( av == 1  .AND.  .NOT.  netcdf_1 )  CYCLE
685
686!
687!--          Open NetCDF dataset
688             IF ( av == 0 )  THEN
689                filename = 'DATA_3D_NETCDF'//TRIM( model_string )
690             ELSE
691                filename = 'DATA_3D_AV_NETCDF'//TRIM( model_string )
692             ENDIF
693             nc_stat = NF90_OPEN( filename, NF90_WRITE, id_set(av) )
694             IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 13 )
695
696
697!
698!--          Get the list of variables (order of variables corresponds with the
699!--          order of data on the binary file)
700             var_list(av) = ' '    ! GET_ATT does not assign trailing blanks
701             nc_stat = NF90_GET_ATT( id_set(av), NF90_GLOBAL, 'VAR_LIST', &
702                  var_list(av) )
703             IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 14 )
704
705!
706!--          Inquire id of the time coordinate variable
707             nc_stat = NF90_INQ_VARID( id_set(av), 'time', id_var_time(av) )
708             IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 15 )
709
710!
711!--          Count number of variables; there is one more semicolon in the
712!--          string than variable names
713             num_var(av) = -1
714             DO  i = 1, LEN( var_list(av) )
715                IF ( var_list(av)(i:i) == ';' )  num_var(av) = num_var(av) + 1
716             ENDDO
717
718!
719!--          Extract the variable names from the list and inquire their NetCDF
720!--          IDs
721             pos = INDEX( var_list(av), ';' )
722!
723!--          Loop over all variables
724             DO  i = 1, num_var(av)
725
726!
727!--             Extract variable name from list
728                var_list(av) = var_list(av)(pos+1:)
729                pos = INDEX( var_list(av), ';' )
730                var_name = var_list(av)(1:pos-1)
731
732!
733!--             Get variable ID from name
734                nc_stat = NF90_INQ_VARID( id_set(av), TRIM( var_name ), &
735                     id_var(av,i) )
736                IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 16 )
737
738             ENDDO
739
740          ENDDO    ! av=0,1
741
742#endif
743       ENDIF
744
745!
746!--    Read arrays, until the end of the file is reached
747       IF ( .NOT. netcdf_parallel )  THEN
748
749          current_var = 999999999
750          fanz = 0
751          DO  WHILE ( danz /= 0 )
752
753!
754!--          Loop over all files
755             DO  id = 0, danz-1
756!
757!--             File from PE0 contains at the beginning the index bounds
758!--             of PALM's total domain.
759!--             Allocate the array for storing the total domain data
760                IF ( id == 0  .AND.  fanz(0) == 0  .AND.  fanz(1) == 0 )  THEN
761                   READ ( id+110 )  nxa, nxe, nya, nye, nza, nze
762                   ALLOCATE ( pf3d(nxa:nxe,nya:nye,nza:nze) )
763                ENDIF
764
765!
766!--             Read output time
767                IF ( netcdf_output  .AND.  id == 0 )  THEN
768                   IF ( netcdf_1 )  THEN
769                      READ ( id+110, END=999 )  simulated_time, time_step, av
770                   ELSE
771!
772!--                   For compatibility with earlier PALM versions
773                      READ ( id+110, END=999 )  simulated_time, time_step
774                      av = 0
775                   ENDIF
776                ENDIF
777
778!
779!--             Read subdomain indices and grid point values
780                READ ( id+110, END=999 )  xa, xe, ya, ye, za, ze
781                ALLOCATE( pf3d_tmp(xa:xe,ya:ye,za:ze) )
782                READ ( id+110 )  pf3d_tmp
783
784                xxa = MAX( nxa, xa )
785                xxe = MIN( nxe, xe )
786                yya = MAX( nya, ya )
787                yye = MIN( nye, ye )
788                DO  k = za, ze
789                   DO  j = yya, yye
790                      DO  i = xxa, xxe
791                         pf3d(i,j,k) = pf3d_tmp(i,j,k)
792                      ENDDO
793                   ENDDO
794                ENDDO
795
796                DEALLOCATE( pf3d_tmp )
797                IF ( id == 0 )  fanz(av) = fanz(av) + 1
798
799             ENDDO
800
801!
802!--          Write data of the total domain in NetCDF format
803             IF ( netcdf_output )  THEN
804#if defined( __netcdf )
805!
806!--             Check if a new time step has begun; if yes write data to time
807!--             axis
808                IF ( current_var(av) > num_var(av) )  THEN
809                   current_var(av) = 1
810                   nc_stat = NF90_PUT_VAR( id_set(av), id_var_time(av), &
811                                      (/ simulated_time /),&
812                                      start = (/ time_step /), count = (/ 1 /) )
813                   IF ( nc_stat /= NF90_NOERR ) CALL handle_netcdf_error( 17 )
814                ENDIF
815
816!
817!--             Now write the data
818                nc_stat = NF90_PUT_VAR( id_set(av), id_var(av,current_var(av)),&
819                                        pf3d(nxa:nxe,nya:nye,za:ze), start = (/ 1, 1, 1, time_step /),&
820                              count = (/ nxe-nxa+1, nye-nya+1, ze-za+1, 1 /) )
821                IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 18 )
822
823                current_var(av) = current_var(av) + 1
824
825#endif
826             ENDIF
827
828          ENDDO
829
830       ENDIF
831
832999    IF ( danz /= 0  .AND.  .NOT. netcdf_parallel )  THEN
833!
834!--       Print the number of arrays processed
835          WRITE (*,'(16X,I4,A)')  fanz(0)+fanz(1), ' array(s) processed'
836          IF ( fanz(1) /= 0 )  THEN
837             WRITE (*,'(16X,I4,A)')  fanz(1), ' array(s) are time-averaged'
838          ENDIF
839!
840!--       Close all files and deallocate array
841          DO  id = 0, danz-1
842             CLOSE ( id+110 )
843          ENDDO
844          CLOSE ( 2 )
845          DEALLOCATE ( pf3d )
846!
847!--       Close the NetCDF file
848          IF ( netcdf_output )  THEN
849#if defined( __netcdf )
850             IF ( netcdf_0 )  THEN
851                nc_stat = NF90_CLOSE( id_set(0) )
852                IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 19 )
853             ENDIF
854             IF ( netcdf_1 )  THEN
855                nc_stat = NF90_CLOSE( id_set(1) )
856                IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 20 )
857             ENDIF
858#endif
859          ENDIF
860
861!
862!--       Output required cpu time
863          CALL SYSTEM_CLOCK( count, count_rate )
864          cpu_end_time = REAL( count ) / REAL( count_rate )
865          WRITE (*,'(5X,A,F9.3,A)')  'Required cpu-time: ', &
866                                     cpu_end_time-cpu_start_time, ' sec'
867
868       ENDIF
869
870    ENDDO  ! models
871
872
873 CONTAINS
874
875
876    SUBROUTINE handle_netcdf_error( errno )
877!
878!--    Prints out a text message corresponding to the current NetCDF status
879
880       IMPLICIT NONE
881
882       INTEGER, INTENT(IN) ::  errno
883
884#if defined( __netcdf )
885       IF ( nc_stat /= NF90_NOERR )  THEN
886          PRINT*, '+++ combine_plot_fields  netcdf: ', av, errno, &
887                  TRIM( nf90_strerror( nc_stat ) )
888       ENDIF
889#endif
890
891    END SUBROUTINE handle_netcdf_error
892
893
894 END PROGRAM combine_plot_fields
895
896
897
Note: See TracBrowser for help on using the repository browser.