source: palm/trunk/SOURCE/data_output_spectra.f90 @ 254

Last change on this file since 254 was 254, checked in by heinze, 13 years ago

Output of messages replaced by message handling routine.

  • Property svn:keywords set to Id
File size: 16.5 KB
Line 
1 SUBROUTINE data_output_spectra
2
3!------------------------------------------------------------------------------!
4! Current revisions:
5! -----------------
6! Output of messages replaced by message handling routine.
7!
8!
9! Former revisions:
10! -----------------
11! $Id: data_output_spectra.f90 254 2009-03-05 15:33:42Z heinze $
12!
13! 189 2008-08-13 17:09:26Z letzel
14! allow 100 spectra levels instead of 10 for consistency with
15! define_netcdf_header, +user-defined spectra
16!
17! February 2007
18! RCS Log replace by Id keyword, revision history cleaned up
19!
20! Revision 1.7  2006/04/11 14:56:38  raasch
21! pl_spectra renamed data_output_sp
22!
23! Revision 1.1  2001/01/05 15:14:20  raasch
24! Initial revision
25!
26!
27! Description:
28! ------------
29! Writing spectra data on file, using a special format which allows
30! plotting of these data with PROFIL-graphic-software
31!------------------------------------------------------------------------------!
32#if defined( __spectra )
33
34    USE arrays_3d
35    USE control_parameters
36    USE cpulog
37    USE interfaces
38    USE netcdf_control
39    USE pegrid
40    USE spectrum
41    USE statistics
42
43
44    IMPLICIT NONE
45
46    INTEGER :: m, pr, cranz_x, cranz_y
47    LOGICAL :: frame_x, frame_y
48
49    CALL cpu_log( log_point(31), 'data_output_spectra', 'start' )
50
51!
52!-- Output is only performed on PE0
53    IF ( myid == 0 )  THEN
54
55!
56!--    Open file for spectra output in NetCDF format
57       IF ( netcdf_output )  CALL check_open( 107 )
58
59!
60!--    Increment the counter for number of output times
61       dosp_time_count = dosp_time_count + 1
62
63#if defined( __netcdf )
64!
65!--    Update the spectra time axis
66       nc_stat = NF90_PUT_VAR( id_set_sp, id_var_time_sp, (/ simulated_time /),&
67                               start = (/ dosp_time_count /), count = (/ 1 /) )
68       IF (nc_stat /= NF90_NOERR)  CALL handle_netcdf_error( 47 )
69#endif
70
71       IF ( profil_output )  THEN
72!
73!--       Compute RAHMEN-Parameter CRANZ for x- and y-spectra separately
74          cranz_x = 0; cranz_y = 0; frame_x = .FALSE.; frame_y = .FALSE.
75
76          m = 1
77          DO WHILE ( data_output_sp(m) /= ' ' .AND. m <= 10 )
78
79             IF ( INDEX( spectra_direction(m), 'x' ) /= 0 )  THEN
80                cranz_x = cranz_x + 1
81             ENDIF
82             IF ( INDEX( spectra_direction(m), 'y' ) /= 0 )  THEN
83                cranz_y = cranz_y + 1
84             ENDIF
85
86             m = m + 1
87
88          ENDDO
89
90       ENDIF
91
92!
93!--    If necessary, calculate time average and reset average counter
94       IF ( average_count_sp == 0 )  THEN
95           message_string = 'no spectra data available'
96           CALL message( 'data_output_spectra', 'PA0186', 0, 0, 0, 6, 0 )
97       ENDIF
98       IF ( average_count_sp /= 1 )  THEN
99          spectrum_x = spectrum_x / REAL( average_count_sp )
100          spectrum_y = spectrum_y / REAL( average_count_sp )
101          average_count_sp = 0
102       ENDIF
103
104!
105!--    Loop over all spectra defined by the user
106       m = 1
107       DO WHILE ( data_output_sp(m) /= ' '  .AND.  m <= 10 )
108
109          SELECT CASE ( TRIM( data_output_sp(m) ) )
110
111             CASE ( 'u' )
112                pr = 1
113
114             CASE ( 'v' )
115                pr = 2
116
117             CASE ( 'w' )
118                pr = 3
119
120             CASE ( 'pt' )
121                pr = 4
122
123             CASE ( 'q' )
124                pr = 5
125
126             CASE DEFAULT
127!
128!--             The DEFAULT case is reached either if the parameter
129!--             data_output_sp(m) contains a wrong character string or if the
130!--             user has coded a special case in the user interface. There, the
131!--             subroutine user_spectra checks which of these two conditions
132!--             applies.
133                CALL user_spectra( 'data_output', m, pr )
134
135          END SELECT
136
137!
138!--       Output of spectra in NetCDF format
139          IF ( netcdf_output )  THEN
140!
141!--          Output of x-spectra
142             IF ( INDEX( spectra_direction(m), 'x' ) /= 0 ) THEN
143                CALL output_spectra_netcdf( m, 'x' )
144             ENDIF
145!
146!--          Output of y-spectra
147             IF ( INDEX( spectra_direction(m), 'y' ) /= 0 ) THEN
148                CALL output_spectra_netcdf( m, 'y' )
149             ENDIF
150          ENDIF
151
152!
153!--       Output of spectra in profil format
154          IF ( profil_output )  THEN
155!
156!--          Output of x-spectra
157             IF ( INDEX( spectra_direction(m), 'x' ) /= 0 ) THEN
158                CALL data_output_spectra_x( m, cranz_x, pr, frame_x )
159             ENDIF
160
161!
162!--          Output of y-spectra
163             IF ( INDEX( spectra_direction(m), 'y' ) /= 0 ) THEN
164                CALL data_output_spectra_y( m, cranz_y, pr, frame_y )
165             ENDIF
166          ENDIF
167
168!
169!--       Increase counter for next spectrum
170          m = m + 1
171
172       ENDDO
173
174!
175!--    Reset spectra values
176       spectrum_x = 0.0; spectrum_y = 0.0
177
178    ENDIF
179
180    CALL cpu_log( log_point(31), 'data_output_spectra', 'stop' )
181
182#if defined( __parallel )
183!    CALL MPI_BARRIER( comm2d, ierr )  ! really necessary
184#endif
185
186#endif
187 END SUBROUTINE data_output_spectra
188
189
190 SUBROUTINE output_spectra_netcdf( nsp, direction )
191#if defined( __netcdf )
192
193    USE constants
194    USE control_parameters
195    USE grid_variables
196    USE indices
197    USE netcdf_control
198    USE spectrum
199    USE statistics
200
201    IMPLICIT NONE
202
203    CHARACTER (LEN=1), INTENT(IN) ::  direction
204
205    INTEGER, INTENT(IN) ::  nsp
206
207    INTEGER ::  i, k
208
209    REAL ::  frequency
210
211    REAL, DIMENSION(nx/2) ::  netcdf_data_x
212    REAL, DIMENSION(ny/2) ::  netcdf_data_y
213
214
215    IF ( direction == 'x' )  THEN
216
217       DO  k = 1, n_sp_x
218
219          DO  i = 1, nx/2
220             frequency = 2.0 * pi * i / ( dx * ( nx + 1 ) )
221             netcdf_data_x(i) = frequency * spectrum_x(i,k,nsp)
222          ENDDO
223
224          nc_stat = NF90_PUT_VAR( id_set_sp, id_var_dospx(nsp), netcdf_data_x, &
225                                  start = (/ 1, k, dosp_time_count /), &
226                                  count = (/ nx/2, 1, 1 /) )
227          IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 348 )
228
229       ENDDO
230
231    ENDIF
232
233    IF ( direction == 'y' )  THEN
234
235       DO  k = 1, n_sp_y
236
237          DO  i = 1, ny/2
238             frequency = 2.0 * pi * i / ( dy * ( ny + 1 ) )
239             netcdf_data_y(i) = frequency * spectrum_y(i,k,nsp)
240          ENDDO
241
242          nc_stat = NF90_PUT_VAR( id_set_sp, id_var_dospy(nsp), netcdf_data_y, &
243                                  start = (/ 1, k, dosp_time_count /), &
244                                  count = (/ ny/2, 1, 1 /) )
245          IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 349 )
246
247       ENDDO
248
249    ENDIF
250
251#endif
252 END SUBROUTINE output_spectra_netcdf
253
254
255#if defined( __spectra )
256 SUBROUTINE data_output_spectra_x( m, cranz, pr, frame_written )
257
258    USE arrays_3d
259    USE constants
260    USE control_parameters
261    USE grid_variables
262    USE indices
263    USE pegrid
264    USE singleton
265    USE spectrum
266    USE statistics
267    USE transpose_indices
268
269    IMPLICIT NONE
270
271    CHARACTER (LEN=30) ::  atext
272    INTEGER            ::  i, j, k, m, pr
273    LOGICAL            ::  frame_written
274    REAL               ::  frequency = 0.0
275
276!
277!-- Variables needed for PROFIL-namelist
278    INTEGER                  :: cranz, labforx = 3, labfory = 3, legpos = 3, &
279                                timodex = 1
280    INTEGER, DIMENSION(1:100):: cucol = 1, klist = 999999, lstyle = 0
281    LOGICAL                  :: datleg = .TRUE., grid = .TRUE., &
282                                lclose = .TRUE., rand = .TRUE., &
283                                swap = .TRUE., twoxa = .TRUE.,  &
284                                xlog = .TRUE., ylog = .TRUE.
285    CHARACTER (LEN=80)       :: rtext, utext, xtext = 'k in m>->1', ytext
286    REAL                     :: gwid = 0.1, rlegfak = 0.7, uxmin, uxmax, &
287                                uymin, uymax
288    REAL, DIMENSION(1:100)   :: lwid = 0.6
289    REAL, DIMENSION(100)     :: uyma, uymi
290
291    NAMELIST /RAHMEN/  cranz, datleg, rtext, swap
292    NAMELIST /CROSS/   rand, cucol, grid, gwid, klist, labforx, labfory,      &
293                       legpos, lclose, lstyle, lwid, rlegfak, timodex, utext, &
294                       uxmin, uxmax, uymin, uymax, twoxa, xlog, xtext, ylog,  &
295                       ytext
296
297
298    rtext = '\0.5 ' // run_description_header
299
300!
301!-- Open parameter- and data-file
302    CALL check_open( 81 )
303    CALL check_open( 82 )
304
305!
306!-- Write file header,
307!-- write RAHMEN-parameters (pr=3: w-array is on zw, other arrys on zu,
308!-- pr serves as an index for output of strings (axis-labels) of the
309!-- different quantities u, v, w, pt and q)
310    DO  k = 1, n_sp_x
311       IF ( k < 100 )  THEN
312          IF ( pr == 3 )  THEN
313             WRITE ( 82, 100 )  '#', k, header_char( pr ),        &
314                                INT( zw(comp_spectra_level(k)) ), &
315                                simulated_time_chr
316          ELSE
317             WRITE ( 82, 100 )  '#', k, header_char( pr ),        &
318                                INT( zu(comp_spectra_level(k)) ), &
319                                simulated_time_chr
320          ENDIF
321       ELSE
322          IF ( pr == 3 )  THEN
323             WRITE ( 82, 101 )  '#', k, header_char( pr ),        &
324                                INT( zw(comp_spectra_level(k)) ), &
325                                simulated_time_chr
326          ELSE
327             WRITE ( 82, 101 )  '#', k, header_char( pr ),        &
328                                INT( zu(comp_spectra_level(k)) ), &
329                                simulated_time_chr
330          ENDIF
331       ENDIF
332    ENDDO
333
334    IF ( .NOT. frame_written )  THEN
335       WRITE ( 81, RAHMEN )
336       frame_written = .TRUE.
337    ENDIF
338
339!
340!-- Write all data and calculate uymi and uyma. They serve to calculate
341!-- the CROSS-parameters uymin and uymax
342    uymi = 999.999; uyma = -999.999
343    DO  i = 1, nx/2
344       frequency = 2.0 * pi * i / ( dx * ( nx + 1 ) )
345       WRITE ( 82, 102 )  frequency, ( frequency * spectrum_x(i,k,m), k = 1, &
346                          n_sp_x )
347       DO  k = 1, n_sp_x
348          uymi(k) = MIN( uymi(k), frequency * spectrum_x(i,k,m) )
349          uyma(k) = MAX( uyma(k), frequency * spectrum_x(i,k,m) )
350       ENDDO
351    ENDDO
352
353!
354!-- Determine CROSS-parameters
355    cucol(1:n_sp_x)  = (/ ( k, k = 1, n_sp_x ) /)
356    lstyle(1:n_sp_x) = (/ ( lstyles(k), k = 1, n_sp_x ) /)
357
358!
359!-- Calculate klist-values from the available comp_spectra_level values
360    i = 1; k = 1
361    DO WHILE ( i <= 100  .AND.  plot_spectra_level(i) /= 999999 )
362       DO WHILE ( k <= n_sp_x  .AND. &
363                  plot_spectra_level(i) >= comp_spectra_level(k) )
364          IF ( plot_spectra_level(i) == comp_spectra_level(k) )  THEN
365             klist(i) = k + klist_x
366          ELSE
367             uymi(k) =  999.999
368             uyma(k) = -999.999
369          ENDIF
370          k = k + 1
371       ENDDO
372       i = i + 1
373    ENDDO
374    uymi(k:n_sp_x) =  999.999
375    uyma(k:n_sp_x) = -999.999
376    utext = 'x'//utext_char( pr )
377    IF ( averaging_interval_sp /= 0.0 ) THEN
378       WRITE ( atext, 104 )  averaging_interval_sp
379       utext = TRIM(utext) // ',  ' // TRIM( atext )
380    ENDIF
381    uxmin = 0.8 * 2.0 * pi        / ( dx * ( nx + 1 ) )
382    uxmax = 1.2 * 2.0 * pi * nx/2 / ( dx * ( nx + 1 ) )
383    uymin = 0.8 * MIN (  999.999, MINVAL ( uymi ) )
384    uymax = 1.2 * MAX ( -999.999, MAXVAL ( uyma ) )
385    ytext = ytext_char( pr )
386
387!
388!-- Output of CROSS-parameters
389    WRITE ( 81, CROSS )
390
391!
392!-- Increase counter by the number of profiles written in the actual block
393    klist_x = klist_x + n_sp_x
394
395!
396!-- Write end-mark
397    WRITE ( 82, 103 )
398
399!
400!-- Close parameter- and data-file
401    CALL close_file( 81 )
402    CALL close_file( 82 )
403
404!
405!-- Formats
406100 FORMAT (A,I1,1X,A,1X,I4,'m ',A)
407101 FORMAT (A,I2,1X,A,1X,I4,'m ',A)
408102 FORMAT (E15.7,100(1X,E15.7))
409103 FORMAT ('NEXT')
410104 FORMAT ('time averaged over',F7.1,' s')
411
412 END SUBROUTINE data_output_spectra_x
413
414
415 SUBROUTINE data_output_spectra_y( m, cranz, pr, frame_written )
416
417    USE arrays_3d
418    USE constants
419    USE control_parameters
420    USE grid_variables
421    USE indices
422    USE pegrid
423    USE singleton
424    USE spectrum
425    USE statistics
426    USE transpose_indices
427
428    IMPLICIT NONE
429
430    CHARACTER (LEN=30) ::  atext
431    INTEGER            :: i, j, k, m, pr
432    LOGICAL            :: frame_written
433    REAL               :: frequency = 0.0
434
435!
436!-- Variables needed for PROFIL-namelist
437    INTEGER                  :: cranz, labforx = 3, labfory = 3, legpos = 3, &
438                                timodex = 1
439    INTEGER, DIMENSION(1:100):: cucol = 1, klist = 999999, lstyle = 0
440    LOGICAL                  :: datleg = .TRUE., grid = .TRUE., &
441                                lclose = .TRUE., rand = .TRUE., &
442                                swap = .TRUE., twoxa = .TRUE.,  &
443                                xlog = .TRUE., ylog = .TRUE.
444    CHARACTER (LEN=80)       :: rtext, utext, xtext = 'k in m>->1', ytext
445    REAL                     :: gwid = 0.1, rlegfak = 0.7, uxmin, uxmax, &
446                                uymin, uymax
447    REAL, DIMENSION(1:100)   :: lwid = 0.6
448    REAL, DIMENSION(100)     :: uyma, uymi
449
450    NAMELIST /RAHMEN/  cranz, datleg, rtext, swap
451    NAMELIST /CROSS/   rand, cucol, grid, gwid, klist, labforx, labfory,      &
452                       legpos, lclose, lstyle, lwid, rlegfak, timodex, utext, &
453                       uxmin, uxmax, uymin, uymax, twoxa, xlog, xtext, ylog,  &
454                       ytext
455
456
457    rtext = '\0.5 ' // run_description_header
458
459!
460!-- Open parameter- and data-file
461    CALL check_open( 83 )
462    CALL check_open( 84 )
463
464!
465!-- Write file header,
466!-- write RAHMEN-parameters (pr=3: w-array is on zw, other arrys on zu,
467!-- pr serves as an index for output of strings (axis-labels) of the
468!-- different quantities u, v, w, pt and q)
469    DO  k = 1, n_sp_y
470       IF ( k < 100 )  THEN
471          IF ( pr == 3 ) THEN
472             WRITE ( 84, 100 )  '#', k, header_char( pr ),        &
473                                INT( zw(comp_spectra_level(k)) ), &
474                                simulated_time_chr
475          ELSE
476             WRITE ( 84, 100 )  '#', k, header_char( pr ),        &
477                                INT( zu(comp_spectra_level(k)) ), &
478                                simulated_time_chr
479          ENDIF
480       ELSE
481          IF ( pr == 3 )  THEN
482             WRITE ( 84, 101 )  '#', k, header_char( pr ),        &
483                                INT( zw(comp_spectra_level(k)) ), &
484                                simulated_time_chr
485          ELSE
486             WRITE ( 84, 101 )  '#', k, header_char( pr ),        &
487                                INT( zu(comp_spectra_level(k)) ), &
488                                simulated_time_chr
489          ENDIF
490       ENDIF
491    ENDDO
492
493    IF ( .NOT. frame_written )  THEN
494       WRITE ( 83, RAHMEN )
495       frame_written = .TRUE.
496    ENDIF
497
498!
499!-- Write all data and calculate uymi and uyma. They serve to calculate
500!-- the CROSS-parameters uymin and uymax
501    uymi = 999.999; uyma = -999.999
502    DO  j = 1, ny/2
503       frequency = 2.0 * pi * j / ( dy * ( ny + 1 ) )
504       WRITE ( 84, 102 ) frequency, ( frequency * spectrum_y(j,k,m), &
505                                      k = 1, n_sp_y ) 
506       DO k = 1, n_sp_y
507          uymi(k) = MIN( uymi(k), frequency * spectrum_y(j,k,m) )
508          uyma(k) = MAX( uyma(k), frequency * spectrum_y(j,k,m) )
509       ENDDO
510    ENDDO
511
512!
513!-- Determine CROSS-parameters
514    cucol(1:n_sp_y)  = (/ ( k, k = 1, n_sp_y ) /)
515    lstyle(1:n_sp_y) = (/ ( lstyles(k), k = 1, n_sp_y ) /)
516
517!
518!-- Calculate klist-values from the available comp_spectra_level values
519    j = 1; k = 1
520    DO WHILE ( j <= 100  .AND.  plot_spectra_level(j) /= 999999 )
521       DO WHILE ( k <= n_sp_y  .AND. &
522                  plot_spectra_level(j) >= comp_spectra_level(k) )
523          IF ( plot_spectra_level(j) == comp_spectra_level(k) )  THEN
524             klist(j) = k + klist_y
525          ELSE
526             uymi(k) =  999.999
527             uyma(k) = -999.999
528          ENDIF
529          k = k + 1
530       ENDDO
531       j = j + 1
532    ENDDO
533    uymi(k:n_sp_y) =  999.999
534    uyma(k:n_sp_y) = -999.999
535    utext = 'y'//utext_char( pr )
536    IF ( averaging_interval_sp /= 0.0 )  THEN
537       WRITE ( atext, 104 )  averaging_interval_sp
538       utext = TRIM(utext) // ',  ' // TRIM( atext )
539    ENDIF
540    uxmin = 0.8 * 2.0 * pi        / ( dy * ( ny + 1 ) )
541    uxmax = 1.2 * 2.0 * pi * ny/2 / ( dy * ( ny + 1 ) )
542    uymin = 0.8 * MIN (  999.999, MINVAL ( uymi ) )
543    uymax = 1.2 * MAX ( -999.999, MAXVAL ( uyma ) )
544    ytext = ytext_char( pr )
545
546!
547!-- Output CROSS-parameters
548    WRITE ( 83, CROSS )
549
550!
551!-- Increase counter by the number of profiles written in the actual block
552    klist_y = klist_y + n_sp_y
553
554!
555!-- Write end-mark
556    WRITE ( 84, 103 ) 
557
558!
559!-- Close parameter- and data-file
560    CALL close_file( 83 )
561    CALL close_file( 84 )
562
563!
564!-- Formats
565100 FORMAT (A,I1,1X,A,1X,I4,'m ',A)
566101 FORMAT (A,I2,1X,A,1X,I4,'m ',A)
567102 FORMAT (E15.7,100(1X,E15.7))
568103 FORMAT ('NEXT')
569104 FORMAT ('time averaged over',F7.1,' s')
570
571 END SUBROUTINE data_output_spectra_y
572#endif
Note: See TracBrowser for help on using the repository browser.