source: palm/trunk/SOURCE/data_output_dvrp.f90 @ 1350

Last change on this file since 1350 was 1347, checked in by heinze, 11 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 25.0 KB
RevLine 
[1320]1 MODULE dvrp_color
2
[1036]3!--------------------------------------------------------------------------------!
4! This file is part of PALM.
5!
6! PALM is free software: you can redistribute it and/or modify it under the terms
7! of the GNU General Public License as published by the Free Software Foundation,
8! either version 3 of the License, or (at your option) any later version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
[1310]17! Copyright 1997-2014 Leibniz Universitaet Hannover
[1036]18!--------------------------------------------------------------------------------!
19!
[254]20! Current revisions:
[1]21! -----------------
[1321]22!
23!
24! Former revisions:
25! -----------------
26! $Id: data_output_dvrp.f90 1347 2014-03-27 13:23:00Z maronga $
27!
[1347]28! 1346 2014-03-27 13:18:20Z heinze
29! Bugfix: REAL constants provided with KIND-attribute especially in call of
30! intrinsic function like MAX, MIN, SIGN
31!
[1321]32! 1320 2014-03-20 08:40:49Z raasch
[1320]33! ONLY-attribute added to USE-statements,
34! kind-parameters added to all INTEGER and REAL declaration statements,
35! kinds are defined in new module kinds,
36! revision history before 2012 removed,
37! comment fields (!:) to be used for variable explanations added to
38! all variable declaration statements
[1]39!
[1319]40! 1318 2014-03-17 13:35:16Z raasch
41! module interfaces removed
42!
[1037]43! 1036 2012-10-22 13:43:42Z raasch
44! code put under GPL (PALM 3.9)
45!
[829]46! 828 2012-02-21 12:00:36Z raasch
47! particle feature color renamed class
48!
[1]49! Revision 1.1  2000/04/27 06:27:17  raasch
50! Initial revision
51!
52!
53! Description:
54! ------------
55! Plot of isosurface, particles and slicers with dvrp-software
56!------------------------------------------------------------------------------!
[1036]57
58    USE dvrp_variables
[1320]59   
60    USE kinds
[1036]61
62    IMPLICIT NONE
63
64 CONTAINS
65
66    SUBROUTINE color_dvrp( value, color )
67
[1320]68       REAL(wp), INTENT(IN)  ::  value    !:
69       REAL(wp), INTENT(OUT) ::  color(4) !:
[1036]70
[1320]71       REAL(wp)              ::  scale    !:
[1036]72
[1320]73       scale = ( value - slicer_range_limits_dvrp(1,islice_dvrp) ) /           &
74               ( slicer_range_limits_dvrp(2,islice_dvrp) -                     &
[1036]75                 slicer_range_limits_dvrp(1,islice_dvrp) )
76
[1346]77       scale = MODULO( 180.0 + 180.0 * scale, 360.0_wp )
[1036]78
[1320]79       color = (/ scale, 0.5_wp, 1.0_wp, 0.0_wp /)
[1036]80
81    END SUBROUTINE color_dvrp
82
83 END MODULE dvrp_color
84
85
86 RECURSIVE SUBROUTINE data_output_dvrp
87
[1]88#if defined( __dvrp_graphics )
89
[1320]90    USE arrays_3d,                                                             &
91        ONLY:  p, pt, q, ql, ts, u, us, v, w, zu
92       
93    USE cloud_parameters,                                                      &
94        ONLY:  l_d_cp, pt_d_t
95       
96    USE constants,                                                             &
97        ONLY:  pi
98       
99    USE control_parameters,                                                    &
100        ONLY:  cloud_droplets, cloud_physics, do2d, do3d, humidity, ibc_uv_b,  &
101               message_string, nz_do3d, passive_scalar, simulated_time,        &
102               threshold
103       
104    USE cpulog,                                                                &
105        ONLY:  log_point, log_point_s, cpu_log
106       
[1]107    USE DVRP
[1320]108   
[1]109    USE dvrp_color
[1320]110       
[1]111    USE dvrp_variables
[1320]112       
113    USE grid_variables,                                                        &
114        ONLY:  dx, dy
115       
116    USE indices,                                                               &
117        ONLY:  nxl, nxr, nyn, nys, nzb
118       
119    USE kinds
120   
121    USE particle_attributes,                                                   &
122        ONLY:  maximum_number_of_tailpoints, number_of_particles,              &
123               number_of_tails, particle_advection, particle_advection_start,  &
124               particle_tail_coordinates, particles, uniform_particles,        &
125               use_particle_tails
126       
[1]127    USE pegrid
128
129    IMPLICIT NONE
130
[1320]131    CHARACTER (LEN=2) ::  section_chr      !:
132    CHARACTER (LEN=6) ::  output_variable  !:
133   
134    INTEGER(iwp) ::  c_mode           !: 
135    INTEGER(iwp) ::  c_size_x         !:
136    INTEGER(iwp) ::  c_size_y         !:
137    INTEGER(iwp) ::  c_size_z         !:
138    INTEGER(iwp) ::  dvrp_nop         !:
139    INTEGER(iwp) ::  dvrp_not         !:
140    INTEGER(iwp) ::  gradient_normals !:
141    INTEGER(iwp) ::  i                !:
142    INTEGER(iwp) ::  ip               !:
143    INTEGER(iwp) ::  j                !:
144    INTEGER(iwp) ::  jp               !:
145    INTEGER(iwp) ::  k                !:
146    INTEGER(iwp) ::  l                !:
147    INTEGER(iwp) ::  m                !:
148    INTEGER(iwp) ::  n                !:
149    INTEGER(iwp) ::  n_isosurface     !:
150    INTEGER(iwp) ::  n_slicer         !:
151    INTEGER(iwp) ::  nn               !:
152    INTEGER(iwp) ::  section_mode     !:
153    INTEGER(iwp) ::  vn               !:
154    INTEGER(iwp), DIMENSION(:), ALLOCATABLE ::  p_c  !:
155    INTEGER(iwp), DIMENSION(:), ALLOCATABLE ::  p_t  !:
[242]156
[1320]157    LOGICAL, DIMENSION(:), ALLOCATABLE ::  dvrp_mask  !:
[242]158
[1320]159    REAL(sp) ::  slicer_position  !:
160    REAL(sp) ::  tmp_alpha        !:
161    REAL(sp) ::  tmp_alpha_w      !:
162    REAL(sp) ::  tmp_b            !:
163    REAL(sp) ::  tmp_c_alpha      !:
164    REAL(sp) ::  tmp_g            !:
165    REAL(sp) ::  tmp_norm         !:
166    REAL(sp) ::  tmp_pos          !:
167    REAL(sp) ::  tmp_r            !:
168    REAL(sp) ::  tmp_t            !:
169    REAL(sp) ::  tmp_th           !:
170    REAL(sp), DIMENSION(:),     ALLOCATABLE   ::  psize  !:
171    REAL(sp), DIMENSION(:),     ALLOCATABLE   ::  p_x    !:
172    REAL(sp), DIMENSION(:),     ALLOCATABLE   ::  p_y    !:
173    REAL(sp), DIMENSION(:),     ALLOCATABLE   ::  p_z    !:
174    REAL(sp), DIMENSION(:,:,:), ALLOCATABLE   ::  local_pf  !:
175    REAL(sp), DIMENSION(:,:,:,:), ALLOCATABLE ::  local_pfi !:
[1]176
177
178    CALL cpu_log( log_point(27), 'data_output_dvrp', 'start' )
179
180!
181!-- Loop over all output modes choosed
[284]182    m            = 1
183    n_isosurface = 0  ! isosurface counter (for threshold values and color)
184    n_slicer     = 0  ! slice plane counter (for range of values)
[1]185    DO WHILE ( mode_dvrp(m) /= ' ' )
186!
187!--    Update of the steering variables
188       IF ( .NOT. lock_steering_update )  THEN
189!
190!--       Set lock to avoid recursive calls of DVRP_STEERING_UPDATE
191          lock_steering_update = .TRUE.
[210]192!          CALL DVRP_STEERING_UPDATE( m-1, data_output_dvrp )
[1]193          lock_steering_update = .FALSE.
194       ENDIF
195
196!
197!--    Determine the variable which shall be plotted (in case of slicers or
198!--    isosurfaces)
199       IF ( mode_dvrp(m)(1:10) == 'isosurface' )  THEN
[130]200          READ ( mode_dvrp(m), '(10X,I2)' )  vn
[1]201          output_variable = do3d(0,vn)
[284]202          n_isosurface = n_isosurface + 1
[1]203       ELSEIF ( mode_dvrp(m)(1:6) == 'slicer' )  THEN
[130]204          READ ( mode_dvrp(m), '(6X,I2)' )  vn
[1]205          output_variable = do2d(0,vn)
206          l = MAX( 2, LEN_TRIM( do2d(0,vn) ) )
207          section_chr = do2d(0,vn)(l-1:l)
208          SELECT CASE ( section_chr )
209             CASE ( 'xy' )
210                section_mode = 2
211                slicer_position = zu(MIN( slicer_position_dvrp(m), nz_do3d ))
212             CASE ( 'xz' )
213                section_mode = 1
214                slicer_position = slicer_position_dvrp(m) * dy
215             CASE ( 'yz' )
216                section_mode = 0
217                slicer_position = slicer_position_dvrp(m) * dx
218          END SELECT
219       ENDIF
220
221!
222!--    Select the plot mode (in case of isosurface or slicer only if user has
223!--    defined a variable which shall be plotted; otherwise do nothing)
[86]224       IF ( mode_dvrp(m)(1:9) == 'particles'  .AND.  particle_advection  .AND. &
225            simulated_time >= particle_advection_start )  THEN
[1]226!
227!--       DVRP-Calls for plotting particles:
228          CALL cpu_log( log_point_s(28), 'dvrp_particles', 'start' )
229
230!
231!--       Definition of characteristics of particle material
[210]232!          tmp_r = 0.1;  tmp_g = 0.7;  tmp_b = 0.1;  tmp_t = 0.0
233          tmp_r = 0.0;  tmp_g = 0.0;  tmp_b = 0.0;  tmp_t = 0.0
234          CALL DVRP_MATERIAL_RGB( m-1, 1, tmp_r, tmp_g, tmp_b, tmp_t )
[1]235
236!
[242]237!--       If clipping is active and if this subdomain is clipped, find out the
238!--       number of particles and tails to be plotted; otherwise, all
[246]239!--       particles/tails are plotted. dvrp_mask is used to mark the partikles.
[264]240          ALLOCATE( dvrp_mask(number_of_particles) )
241
[242]242          IF ( dvrp_total_overlap )  THEN
[264]243             dvrp_mask = .TRUE.
244             dvrp_nop  = number_of_particles
245             dvrp_not  = number_of_tails
[242]246          ELSE
[264]247             dvrp_mask = .FALSE.
248             dvrp_nop  = 0
249             dvrp_not  = 0
[242]250             IF ( dvrp_overlap )  THEN
251                IF ( .NOT. use_particle_tails )  THEN
252                   DO  n = 1, number_of_particles
253                      ip = particles(n)%x / dx
254                      jp = particles(n)%y / dy
255                      IF ( ip >= nxl_dvrp  .AND.  ip <= nxr_dvrp  .AND.  &
256                           jp >= nys_dvrp  .AND.  jp <= nyn_dvrp )  THEN
257                         dvrp_nop = dvrp_nop + 1
[264]258                         dvrp_mask(n) = .TRUE.
[242]259                      ENDIF
260                   ENDDO
261                ELSE
262                   k = 0
263                   DO  n = 1, number_of_particles
264                      IF ( particles(n)%tail_id /= 0 )  THEN
265                         k = k + 1
266                         ip = particles(n)%x / dx
267                         jp = particles(n)%y / dy
268                         IF ( ip >= nxl_dvrp  .AND.  ip <= nxr_dvrp  .AND.  &
269                              jp >= nys_dvrp  .AND.  jp <= nyn_dvrp )  THEN
270                            dvrp_not = dvrp_not + 1
[264]271                            dvrp_mask(n) = .TRUE.
[242]272                         ENDIF
273                      ENDIF
274                   ENDDO
275                ENDIF
276             ENDIF
277          ENDIF
278
279!
[1]280!--       Move particle coordinates to one-dimensional arrays
281          IF ( .NOT. use_particle_tails )  THEN
282!
283!--          All particles are output
[242]284             ALLOCATE( psize(dvrp_nop), p_t(dvrp_nop), p_c(dvrp_nop), &
285                       p_x(dvrp_nop), p_y(dvrp_nop), p_z(dvrp_nop) )
[1]286             psize = 0.0;  p_t = 0;  p_c = 0.0;  p_x = 0.0;  p_y = 0.0
[242]287             p_z   = 0.0
288             k = 0
289             DO  n = 1, number_of_particles
290                IF ( dvrp_mask(n) )  THEN
291                   k = k + 1
292                   psize(k) = particles(n)%dvrp_psize
293                   p_x(k)   = particles(n)%x * superelevation_x
294                   p_y(k)   = particles(n)%y * superelevation_y
295                   p_z(k)   = particles(n)%z * superelevation
[828]296                   p_c(k)   = particles(n)%class
[242]297                ENDIF
298             ENDDO
[1]299          ELSE
300!
301!--          Particles have a tail
[242]302             ALLOCATE( psize(dvrp_not), p_t(dvrp_not),             &
303                       p_c(dvrp_not*maximum_number_of_tailpoints), &
304                       p_x(dvrp_not*maximum_number_of_tailpoints), &
305                       p_y(dvrp_not*maximum_number_of_tailpoints), &
306                       p_z(dvrp_not*maximum_number_of_tailpoints) )
[1]307             psize = 0.0;  p_t = 0;  p_c = 0.0;  p_x = 0.0;  p_y = 0.0
308             p_z   = 0.0;
309             i = 0
310             k = 0
[264]311
[1]312             DO  n = 1, number_of_particles
313                nn = particles(n)%tail_id
[242]314                IF ( nn /= 0  .AND.  dvrp_mask(n) )  THEN
[1]315                   k = k + 1
316                   DO  j = 1, particles(n)%tailpoints
317                      i = i + 1
318                      p_x(i) = particle_tail_coordinates(j,1,nn) * &
319                                                                superelevation_x
320                      p_y(i) = particle_tail_coordinates(j,2,nn) * &
321                                                                superelevation_y
322                      p_z(i) = particle_tail_coordinates(j,3,nn) * &
323                                                                superelevation
324                      p_c(i) = particle_tail_coordinates(j,4,nn)
325                   ENDDO
326                   psize(k) = particles(n)%dvrp_psize
327                   p_t(k)   = particles(n)%tailpoints - 1
328                ENDIF               
329             ENDDO
330          ENDIF
331
332!
333!--       Compute and plot particles in dvr-format
334          IF ( uniform_particles  .AND.  .NOT. use_particle_tails )  THEN
335!
336!--          All particles have the same color. Use simple routine to set
337!--          the particle attributes (produces less output data)
338             CALL DVRP_PARTICLES( m-1, p_x, p_y, p_z, psize )
339          ELSE
340!
341!--          Set color definitions
342             CALL user_dvrp_coltab( 'particles', 'none' )
343
[264]344             CALL DVRP_COLORTABLE_HLS( m-1, 0, interval_values_dvrp_prt, &
345                                       interval_h_dvrp_prt,              &
346                                       interval_l_dvrp_prt,              &
347                                       interval_s_dvrp_prt,              &
348                                       interval_a_dvrp_prt )
[1]349
350             IF ( .NOT. use_particle_tails )  THEN
[242]351                CALL DVRP_PARTICLES( m-1, dvrp_nop, p_x, p_y, p_z, 3, psize, &
352                                     p_c, p_t )
[1]353             ELSE
[242]354                CALL DVRP_PARTICLES( m-1, dvrp_not, p_x, p_y, p_z, 15, psize, &
355                                     p_c, p_t )
[1]356             ENDIF
357          ENDIF
358
359          CALL DVRP_VISUALIZE( m-1, 3, dvrp_filecount )
360
[242]361          DEALLOCATE( dvrp_mask, psize, p_c, p_t, p_x, p_y, p_z )
[1]362
363          CALL cpu_log( log_point_s(28), 'dvrp_particles', 'stop' )
364
365
366       ELSEIF ( ( mode_dvrp(m)(1:10) == 'isosurface'  .OR.   &
367                  mode_dvrp(m)(1:6)  == 'slicer'           ) &
368                  .AND.  output_variable /= ' ' )  THEN
369
370!
371!--       Create an intermediate array, properly dimensioned for plot-output
[246]372          ALLOCATE( local_pf(nxl_dvrp:nxr_dvrp+1,nys_dvrp:nyn_dvrp+1, &
373                             nzb:nz_do3d) )
[1]374
375!
376!--       Move original array to intermediate array
[246]377          IF ( dvrp_overlap )  THEN
[1]378
[246]379             SELECT CASE ( output_variable )
380
381                CASE ( 'u', 'u_xy', 'u_xz', 'u_yz' )
382                   DO  i = nxl_dvrp, nxr_dvrp+1
383                      DO  j = nys_dvrp, nyn_dvrp+1
384                         DO  k = nzb, nz_do3d
385                            local_pf(i,j,k) = u(k,j,i)
386                         ENDDO
[1]387                      ENDDO
388                   ENDDO
389!
[246]390!--                Replace mirrored values at lower surface by real surface
391!--                values
392                   IF ( output_variable == 'u_xz'  .OR. &
393                        output_variable == 'u_yz' )  THEN
394                      IF ( ibc_uv_b == 0 )  local_pf(:,:,nzb) = 0.0
395                   ENDIF
[1]396
397
[246]398                CASE ( 'v', 'v_xy', 'v_xz', 'v_yz' )
399                   DO  i = nxl_dvrp, nxr_dvrp+1
400                      DO  j = nys_dvrp, nyn_dvrp+1
401                         DO  k = nzb, nz_do3d
402                            local_pf(i,j,k) = v(k,j,i)
403                         ENDDO
[1]404                      ENDDO
405                   ENDDO
406!
[246]407!--                Replace mirrored values at lower surface by real surface
408!--                values
409                   IF ( output_variable == 'v_xz'  .OR. &
410                        output_variable == 'v_yz' )  THEN
411                      IF ( ibc_uv_b == 0 )  local_pf(:,:,nzb) = 0.0
412                   ENDIF
[1]413
[246]414                CASE ( 'w', 'w_xy', 'w_xz', 'w_yz' )
415                   DO  i = nxl_dvrp, nxr_dvrp+1
416                      DO  j = nys_dvrp, nyn_dvrp+1
417                         DO  k = nzb, nz_do3d
418                            local_pf(i,j,k) = w(k,j,i)
419                         ENDDO
[1]420                      ENDDO
421                   ENDDO
[106]422! Averaging for Langmuir circulation
[246]423!                   DO  k = nzb, nz_do3d
424!                      DO  j = nys_dvrp+1, nyn_dvrp
425!                         DO  i = nxl_dvrp, nxr_dvrp+1
426!                            local_pf(i,j,k) = 0.25 * local_pf(i,j-1,k) + &
427!                                              0.50 * local_pf(i,j,k)   + &
428!                                              0.25 * local_pf(i,j+1,k)
429!                         ENDDO
[106]430!                      ENDDO
431!                   ENDDO
[1]432
[246]433                CASE ( 'p', 'p_xy', 'p_xz', 'p_yz' )
434                   DO  i = nxl_dvrp, nxr_dvrp+1
435                      DO  j = nys_dvrp, nyn_dvrp+1
436                         DO  k = nzb, nz_do3d
437                            local_pf(i,j,k) = p(k,j,i)
438                         ENDDO
[1]439                      ENDDO
440                   ENDDO
441
[246]442                CASE ( 'pt', 'pt_xy', 'pt_xz', 'pt_yz' )
443                   IF ( .NOT. cloud_physics ) THEN
444                      DO  i = nxl_dvrp, nxr_dvrp+1
445                         DO  j = nys_dvrp, nyn_dvrp+1
446                            DO  k = nzb, nz_do3d
447                               local_pf(i,j,k) = pt(k,j,i)
448                            ENDDO
[1]449                         ENDDO
450                      ENDDO
[246]451                   ELSE
452                      DO  i = nxl_dvrp, nxr_dvrp+1
453                         DO  j = nys_dvrp, nyn_dvrp+1
454                            DO  k = nzb, nz_do3d
455                               local_pf(i,j,k) = pt(k,j,i) + l_d_cp * &
456                                                 pt_d_t(k) * ql(k,j,i)
457                            ENDDO
[1]458                         ENDDO
459                      ENDDO
[246]460                   ENDIF
[1]461
[246]462                CASE ( 'q', 'q_xy', 'q_xz', 'q_yz' )
463                   IF ( humidity  .OR.  passive_scalar )  THEN
464                      DO  i = nxl_dvrp, nxr_dvrp+1
465                         DO  j = nys_dvrp, nyn_dvrp+1
466                            DO  k = nzb, nz_do3d
467                               local_pf(i,j,k) = q(k,j,i)
468                            ENDDO
[1]469                         ENDDO
[246]470                      ENDDO           
[254]471                   ELSE                   
[274]472                      message_string = 'if humidity/passive_scalar = '    // & 
473                            'FALSE output of ' // TRIM( output_variable ) // &
474                            'is not provided' 
475                      CALL message( 'data_output_dvrp', 'PA0183',&
476                                                                 0, 0, 0, 6, 0 )
[1]477                   ENDIF
478             
[246]479                CASE ( 'ql', 'ql_xy', 'ql_xz', 'ql_yz' )
480                   IF ( cloud_physics  .OR.  cloud_droplets )  THEN
481                      DO  i = nxl_dvrp, nxr_dvrp+1
482                         DO  j = nys_dvrp, nyn_dvrp+1
483                            DO  k = nzb, nz_do3d
484                               local_pf(i,j,k) = ql(k,j,i)
485                            ENDDO
[1]486                         ENDDO
487                      ENDDO
[254]488                   ELSE                     
[274]489                      message_string = 'if cloud_physics = FALSE '       // & 
490                                  'output of ' // TRIM( output_variable) // &
491                                  'is not provided' 
492                      CALL message( 'data_output_dvrp', 'PA0184',&
493                                                                 0, 0, 0, 6, 0 )
[1]494                   ENDIF
495
[246]496                CASE ( 'u*_xy' )
497                   DO  i = nxl_dvrp, nxr_dvrp+1
498                      DO  j = nys_dvrp, nyn_dvrp+1
499                         local_pf(i,j,nzb+1) = us(j,i)
500                      ENDDO
[1]501                   ENDDO
[246]502                   slicer_position = zu(nzb+1)
[1]503
[246]504                CASE ( 't*_xy' )
505                   DO  i = nxl_dvrp, nxr_dvrp+1
506                      DO  j = nys_dvrp, nyn_dvrp+1
507                         local_pf(i,j,nzb+1) = ts(j,i)
508                      ENDDO
[1]509                   ENDDO
[246]510                   slicer_position = zu(nzb+1)
[1]511
512
[246]513                CASE DEFAULT
[130]514!
[246]515!--                The DEFAULT case is reached either if output_variable
516!--                contains unsupported variable or if the user has coded a
517!--                special case in the user interface. There, the subroutine
518!--                user_data_output_dvrp checks which of these two conditions
519!--                applies.
520                   CALL user_data_output_dvrp( output_variable, local_pf )
[1]521
[130]522
[246]523             END SELECT
[1]524
[264]525          ELSE
526!
527!--          No overlap of clipping domain with the current subdomain
528             DO  i = nxl_dvrp, nxr_dvrp+1
529                DO  j = nys_dvrp, nyn_dvrp+1
530                   DO  k = nzb, nz_do3d
531                      local_pf(i,j,k) = 0.0
532                   ENDDO
533                ENDDO
534             ENDDO
535
[246]536          ENDIF
[1]537
538          IF ( mode_dvrp(m)(1:10) == 'isosurface' )  THEN
[392]539
[1]540!
541!--          DVRP-Calls for plotting isosurfaces:
542             CALL cpu_log( log_point_s(26), 'dvrp_isosurface', 'start' )
543
544!
[284]545!--          Definition of isosurface color
546             tmp_r = isosurface_color(1,n_isosurface)
547             tmp_g = isosurface_color(2,n_isosurface)
548             tmp_b = isosurface_color(3,n_isosurface)
549             tmp_t = 0.0
550             CALL DVRP_MATERIAL_RGB( m-1, 1, tmp_r, tmp_g, tmp_b, tmp_t )
[1]551
552!
553!--          Compute and plot isosurface in dvr-format
554             CALL DVRP_DATA( m-1, local_pf, 1, nx_dvrp, ny_dvrp, nz_dvrp, &
555                             cyclic_dvrp, cyclic_dvrp, cyclic_dvrp )
[210]556
557             c_size_x = vc_size_x;  c_size_y = vc_size_y;  c_size_z = vc_size_z
558             CALL DVRP_CLUSTER_SIZE( m-1, c_size_x, c_size_y, c_size_z )
559
560             c_mode   = vc_mode 
561             CALL DVRP_CLUSTERING_MODE( m-1, c_mode )
562
563             gradient_normals = vc_gradient_normals
564             CALL DVRP_GRADIENTNORMALS( m-1, gradient_normals )
[392]565
[210]566!
567!--          A seperate procedure for setting vc_alpha will be in the next
568!--          version of libDVRP
569             tmp_c_alpha = vc_alpha 
570             CALL DVRP_THRESHOLD( -(m-1)-1, tmp_c_alpha )
571
[246]572             IF ( dvrp_overlap )  THEN
[284]573                tmp_th = threshold(n_isosurface)
[246]574             ELSE
575                tmp_th = 1.0   ! nothing is plotted because array values are 0
576             ENDIF
577
[210]578             CALL DVRP_THRESHOLD( m-1, tmp_th )
[392]579
[210]580             CALL DVRP_VISUALIZE( m-1, 21, dvrp_filecount )
581
[1]582             CALL cpu_log( log_point_s(26), 'dvrp_isosurface', 'stop' )
583
584          ELSEIF ( mode_dvrp(m)(1:6) == 'slicer' )  THEN
[392]585
[1]586!
587!--          DVRP-Calls for plotting slicers:
588             CALL cpu_log( log_point_s(27), 'dvrp_slicer', 'start' )
589
590!
591!--          Material and color definitions
[210]592             tmp_r = 0.0;  tmp_g = 0.0;  tmp_b = 0.0;  tmp_t = 0.0
593             CALL DVRP_MATERIAL_RGB( m-1, 1, tmp_r, tmp_g, tmp_b, tmp_t )
[1]594
[284]595             n_slicer = n_slicer + 1
596
597!
598!--           Using dolorfunction has not been properly tested
599!             islice_dvrp = n_slicer
600!             CALL DVRP_COLORFUNCTION( m-1, DVRP_CM_HLS, 25,                 &
601!                                      slicer_range_limits_dvrp(:,n_slicer), &
[1]602!                                      color_dvrp )
603
[284]604!
605!--          Set interval of values defining the colortable
606             CALL set_slicer_attributes_dvrp( n_slicer )
607
608!
609!--          Create user-defined colortable
[1]610             CALL user_dvrp_coltab( 'slicer', output_variable )
611
612             CALL DVRP_COLORTABLE_HLS( m-1, 1, interval_values_dvrp,     &
613                                       interval_h_dvrp, interval_l_dvrp, &
614                                       interval_s_dvrp, interval_a_dvrp )
615
616!
617!--          Compute and plot slicer in dvr-format
618             CALL DVRP_DATA( m-1, local_pf, 1, nx_dvrp, ny_dvrp, nz_dvrp, &
619                             cyclic_dvrp, cyclic_dvrp, cyclic_dvrp )
[262]620             tmp_pos = slicer_position
621             CALL DVRP_SLICER( m-1, section_mode, tmp_pos )
[106]622
[1]623             CALL DVRP_VISUALIZE( m-1, 2, dvrp_filecount )
624
625             CALL cpu_log( log_point_s(27), 'dvrp_slicer', 'stop' )
626
627          ENDIF
628
629          DEALLOCATE( local_pf )
630
[210]631       ELSEIF ( mode_dvrp(m)(1:9) == 'pathlines' ) THEN
632
633          ALLOCATE( local_pfi(4,nxl:nxr+1,nys:nyn+1,nzb:nz_do3d) )
634          DO  i = nxl, nxr+1
635             DO  j = nys, nyn+1
636                DO  k = nzb, nz_do3d
637                   local_pfi(1,i,j,k) = u(k,j,i)
638                   local_pfi(2,i,j,k) = v(k,j,i)
639                   local_pfi(3,i,j,k) = w(k,j,i)
640                   tmp_norm           = SQRT( u(k,j,i) * u(k,j,i) + &
641                                              v(k,j,i) * v(k,j,i) + &
642                                              w(k,j,i) * w(k,j,i) )
643                   tmp_alpha          = ACOS( 0.0 * u(k,j,i) / tmp_norm + &
644                                              0.0 * v(k,j,i) / tmp_norm - &
645                                              1.0 * w(k,j,i) / tmp_norm )
646                   tmp_alpha_w        = tmp_alpha / pi * 180.0
647                   local_pfi(4,i,j,k) = tmp_alpha_w
648                ENDDO
649             ENDDO
650          ENDDO
651
652          CALL cpu_log( log_point_s(31), 'dvrp_pathlines', 'start' )
653
654          CALL DVRP_DATA( m-1, local_pfi, 4, nx_dvrp, ny_dvrp, nz_dvrp, &
655                          cyclic_dvrp, cyclic_dvrp, cyclic_dvrp )
656          CALL DVRP_VISUALIZE( m-1, 20, dvrp_filecount )
657
658          CALL cpu_log( log_point_s(31), 'dvrp_pathlines', 'stop' )
659
660          DEALLOCATE( local_pfi )
661
[1]662       ENDIF
663
664       m = m + 1
665
666    ENDDO
667
668    dvrp_filecount = dvrp_filecount + 1
669
670    CALL cpu_log( log_point(27), 'data_output_dvrp', 'stop' )
671
672#endif
673 END SUBROUTINE data_output_dvrp
Note: See TracBrowser for help on using the repository browser.