Index: palm/trunk/DOC/tec/message_identifiers
===================================================================
--- palm/trunk/DOC/tec/message_identifiers	(revision 263)
+++ palm/trunk/DOC/tec/message_identifiers	(revision 264)
@@ -566,4 +566,8 @@
                        v_max   = ... m/s   k=... j=...  i=...            
                        w_max   = ... m/s   k=... j=...  i=...
+    PA0313            illegal value for parameter particle_color: ...
+    PA0314            illegal value for parameter particle_dvrpsize: ...
+    PA0315            color_interval(2) <= color_interval(1)
+    PA0316            dvrpsize_interval(2) <= dvrpsize_interval(1)
     
 
Index: palm/trunk/SOURCE/CURRENT_MODIFICATIONS
===================================================================
--- palm/trunk/SOURCE/CURRENT_MODIFICATIONS	(revision 263)
+++ palm/trunk/SOURCE/CURRENT_MODIFICATIONS	(revision 264)
@@ -1,6 +1,8 @@
 New:
 ---
-Clipping of dvrp output implemented. (init_dvrp, data_output_dvrp,
-user_data_output_dvrp)
+Clipping of dvrp output implemented. Default colourtable for particles
+implemented, particle attributes (color, dvrp_size) can be set with new
+parameters particle_color, particle_dvrpsize, color_interval,
+dvrpsize_interval (init_dvrp, data_output_dvrp, modules, user_data_output_dvrp)
 
 Polygon reduction for topography and ground plate isosurface. Reduction level
@@ -19,4 +21,5 @@
 check_parameters, data_output_dvrp, init_dvrp, init_grid, init_3d_model, header, modules, package_parin, parin, prognostic_equations, read_var_list, user_check_parameters, user_data_output_dvrp, user_header, user_init_grid, write_var_list
 
+New: set_particle_attributes
 
 Changed:
Index: palm/trunk/SOURCE/Makefile
===================================================================
--- palm/trunk/SOURCE/Makefile	(revision 263)
+++ palm/trunk/SOURCE/Makefile	(revision 264)
@@ -4,5 +4,5 @@
 # Actual revisions:
 # -----------------
-# 
+# +set_particle_attributes
 #
 # Former revisions:
@@ -74,9 +74,10 @@
 	production_e.f90 prognostic_equations.f90 random_function.f90 \
 	random_gauss.f90 read_3d_binary.f90 read_var_list.f90 run_control.f90 \
-	singleton.f90 sor.f90 spline_x.f90 spline_y.f90 spline_z.f90 \
-	sum_up_3d_data.f90 surface_coupler.f90 swap_timelevel.f90 \
-	temperton_fft.f90 time_integration.f90 time_to_string.f90 timestep.f90 \
-	timestep_scheme_steering.f90 transpose.f90 user_3d_data_averaging.f90 \
-	user_actions.f90 user_additional_routines.f90 user_advec_particles.f90 \
+	set_particle_attributes.f90 singleton.f90 sor.f90 spline_x.f90 \
+	spline_y.f90 spline_z.f90 sum_up_3d_data.f90 surface_coupler.f90 \
+	swap_timelevel.f90 temperton_fft.f90 time_integration.f90 \
+        time_to_string.f90 timestep.f90 timestep_scheme_steering.f90 \
+        transpose.f90 user_3d_data_averaging.f90 user_actions.f90 \
+	user_additional_routines.f90 user_advec_particles.f90 \
 	user_check_data_output.f90 user_check_data_output_pr.f90 \
 	user_check_parameters.f90 user_data_output_2d.f90 \
@@ -113,7 +114,8 @@
 	poisfft_hybrid.o poismg.o prandtl_fluxes.o pres.o print_1d.o \
 	production_e.o prognostic_equations.o random_function.o random_gauss.o \
-	read_3d_binary.o read_var_list.o run_control.o singleton.o sor.o \
-	spline_x.o spline_y.o spline_z.o sum_up_3d_data.o surface_coupler.o \
-	swap_timelevel.o temperton_fft.o time_integration.o time_to_string.o \
+	read_3d_binary.o read_var_list.o run_control.o \
+	set_particle_attributes.o singleton.o sor.o spline_x.o spline_y.o \
+        spline_z.o sum_up_3d_data.o surface_coupler.o swap_timelevel.o \
+        temperton_fft.o time_integration.o time_to_string.o \
         timestep.o timestep_scheme_steering.o transpose.o \
 	user_3d_data_averaging.o user_actions.o user_additional_routines.o \
@@ -253,4 +255,5 @@
 read_var_list.o: modules.o
 run_control.o: modules.o
+set_particle_attributes.o: modules.o
 singleton.o: singleton.f90
 sor.o: modules.o
Index: palm/trunk/SOURCE/advec_particles.f90
===================================================================
--- palm/trunk/SOURCE/advec_particles.f90	(revision 263)
+++ palm/trunk/SOURCE/advec_particles.f90	(revision 264)
@@ -4,4 +4,5 @@
 ! Current revisions:
 ! -----------------
+! Particle attributes are set with new routine set_particle_attributes.
 ! Output of NetCDF messages with aid of message handling routine.
 ! Output of messages replaced by message handling routine
@@ -3360,4 +3361,8 @@
 
 !
+!-- Set particle attributes
+    CALL set_particle_attributes
+
+!
 !-- Set particle attributes defined by the user
     CALL user_particle_attributes
Index: palm/trunk/SOURCE/check_parameters.f90
===================================================================
--- palm/trunk/SOURCE/check_parameters.f90	(revision 263)
+++ palm/trunk/SOURCE/check_parameters.f90	(revision 264)
@@ -4,4 +4,5 @@
 ! Actual revisions:
 ! -----------------
+! Check particle_color, particle_dvrpsize, color_interval, dvrpsize_interval
 ! Bugfix: pressure included for profile output
 ! Check pressure gradient conditions
@@ -90,4 +91,5 @@
     USE constants
     USE control_parameters
+    USE dvrp_variables
     USE grid_variables
     USE indices
@@ -2916,4 +2918,33 @@
 
 !
+!-- Check particle attributes
+    IF ( particle_color /= 'none' )  THEN
+       IF ( particle_color /= 'absuv'  .AND.  particle_color /= 'pt*'  .AND.  &
+            particle_color /= 'z' )  THEN
+          message_string = 'illegal value for parameter particle_color: ' // &
+                           TRIM( particle_color)
+          CALL message( 'check_parameters', 'PA0313', 1, 2, 0, 6, 0 )
+       ELSE
+          IF ( color_interval(2) <= color_interval(1) )  THEN
+             message_string = 'color_interval(2) <= color_interval(1)'
+             CALL message( 'check_parameters', 'PA0315', 1, 2, 0, 6, 0 )
+          ENDIF
+       ENDIF
+    ENDIF
+
+    IF ( particle_dvrpsize /= 'none' )  THEN
+       IF ( particle_dvrpsize /= 'absw' )  THEN
+          message_string = 'illegal value for parameter particle_dvrpsize:' // &
+                           ' ' // TRIM( particle_color)
+          CALL message( 'check_parameters', 'PA0314', 1, 2, 0, 6, 0 )
+       ELSE
+          IF ( dvrpsize_interval(2) <= dvrpsize_interval(1) )  THEN
+             message_string = 'dvrpsize_interval(2) <= dvrpsize_interval(1)'
+             CALL message( 'check_parameters', 'PA0316', 1, 2, 0, 6, 0 )
+          ENDIF
+       ENDIF
+    ENDIF
+
+!
 !-- Check &userpar parameters
     CALL user_check_parameters
Index: palm/trunk/SOURCE/data_output_dvrp.f90
===================================================================
--- palm/trunk/SOURCE/data_output_dvrp.f90	(revision 263)
+++ palm/trunk/SOURCE/data_output_dvrp.f90	(revision 264)
@@ -32,5 +32,6 @@
 ! Current revisions:
 ! -----------------
-! Clipping of dvr-output implemented,
+! Clipping of dvr-output implemented, using a default colourtable for
+! particles,
 ! output of messages replaced by message handling routine.
 ! TEST: different colours for isosurfaces
@@ -168,16 +169,14 @@
 !--       number of particles and tails to be plotted; otherwise, all
 !--       particles/tails are plotted. dvrp_mask is used to mark the partikles.
-!          IF ( .NOT. use_particle_tails )  THEN
-             ALLOCATE( dvrp_mask(number_of_particles) )
-!          ELSE
-!             ALLOCATE( dvrp_mask(number_of_tails*maximum_number_of_tailpoints) )
-!          ENDIF
-          dvrp_mask = .TRUE.
+          ALLOCATE( dvrp_mask(number_of_particles) )
+
           IF ( dvrp_total_overlap )  THEN
-             dvrp_nop = number_of_particles
-             dvrp_not = number_of_tails
+             dvrp_mask = .TRUE.
+             dvrp_nop  = number_of_particles
+             dvrp_not  = number_of_tails
           ELSE
-             dvrp_nop = 0
-             dvrp_not = 0
+             dvrp_mask = .FALSE.
+             dvrp_nop  = 0
+             dvrp_not  = 0
              IF ( dvrp_overlap )  THEN
                 IF ( .NOT. use_particle_tails )  THEN
@@ -188,6 +187,5 @@
                            jp >= nys_dvrp  .AND.  jp <= nyn_dvrp )  THEN
                          dvrp_nop = dvrp_nop + 1
-                      ELSE
-                         dvrp_mask(n) = .FALSE.
+                         dvrp_mask(n) = .TRUE.
                       ENDIF
                    ENDDO
@@ -202,6 +200,5 @@
                               jp >= nys_dvrp  .AND.  jp <= nyn_dvrp )  THEN
                             dvrp_not = dvrp_not + 1
-                         ELSE
-                            dvrp_mask(n) = .FALSE.
+                            dvrp_mask(n) = .TRUE.
                          ENDIF
                       ENDIF
@@ -243,4 +240,5 @@
              i = 0
              k = 0
+
              DO  n = 1, number_of_particles
                 nn = particles(n)%tail_id
@@ -275,7 +273,9 @@
              CALL user_dvrp_coltab( 'particles', 'none' )
 
-             CALL DVRP_COLORTABLE_HLS( m-1, 0, interval_values_dvrp,     &
-                                       interval_h_dvrp, interval_l_dvrp, &
-                                       interval_s_dvrp, interval_a_dvrp )
+             CALL DVRP_COLORTABLE_HLS( m-1, 0, interval_values_dvrp_prt, &
+                                       interval_h_dvrp_prt,              &
+                                       interval_l_dvrp_prt,              &
+                                       interval_s_dvrp_prt,              &
+                                       interval_a_dvrp_prt )
 
              IF ( .NOT. use_particle_tails )  THEN
@@ -452,7 +452,20 @@
              END SELECT
 
+          ELSE
+!
+!--          No overlap of clipping domain with the current subdomain
+             DO  i = nxl_dvrp, nxr_dvrp+1
+                DO  j = nys_dvrp, nyn_dvrp+1
+                   DO  k = nzb, nz_do3d
+                      local_pf(i,j,k) = 0.0
+                   ENDDO
+                ENDDO
+             ENDDO
+
           ENDIF
 
           IF ( mode_dvrp(m)(1:10) == 'isosurface' )  THEN
+             WRITE(9,*) '#8'
+             CALL local_flush( 9 )
 !
 !--          DVRP-Calls for plotting isosurfaces:
@@ -477,4 +490,6 @@
                 CALL DVRP_MATERIAL_RGB( m-1, 1, tmp_r, tmp_g, tmp_b, tmp_t )
              ENDIF
+             WRITE(9,*) '#8.1'
+             CALL local_flush( 9 )
 
 !
@@ -482,13 +497,21 @@
              CALL DVRP_DATA( m-1, local_pf, 1, nx_dvrp, ny_dvrp, nz_dvrp, &
                              cyclic_dvrp, cyclic_dvrp, cyclic_dvrp )
+             WRITE(9,*) '#8.2'
+             CALL local_flush( 9 )
 
              c_size_x = vc_size_x;  c_size_y = vc_size_y;  c_size_z = vc_size_z
              CALL DVRP_CLUSTER_SIZE( m-1, c_size_x, c_size_y, c_size_z )
+             WRITE(9,*) '#8.3'
+             CALL local_flush( 9 )
 
              c_mode   = vc_mode 
              CALL DVRP_CLUSTERING_MODE( m-1, c_mode )
+             WRITE(9,*) '#8.4'
+             CALL local_flush( 9 )
 
              gradient_normals = vc_gradient_normals
              CALL DVRP_GRADIENTNORMALS( m-1, gradient_normals )
+             WRITE(9,*) '#8.5'
+             CALL local_flush( 9 )
 !
 !--          A seperate procedure for setting vc_alpha will be in the next
@@ -496,4 +519,6 @@
              tmp_c_alpha = vc_alpha 
              CALL DVRP_THRESHOLD( -(m-1)-1, tmp_c_alpha )
+             WRITE(9,*) '#8.6'
+             CALL local_flush( 9 )
 
              IF ( dvrp_overlap )  THEN
@@ -502,11 +527,25 @@
                 tmp_th = 1.0   ! nothing is plotted because array values are 0
              ENDIF
+             WRITE(9,*) '#8.7'
+             CALL local_flush( 9 )
 
              CALL DVRP_THRESHOLD( m-1, tmp_th )
+             WRITE(9,*) '#8.8'
+             WRITE(9,*) 'nxl_dvrp=',nxl_dvrp, ' nxl=',nxl
+             WRITE(9,*) 'nxr_dvrp=',nxr_dvrp, ' nxr=',nxr
+             WRITE(9,*) 'nys_dvrp=',nys_dvrp, ' nys=',nys
+             WRITE(9,*) 'nyn_dvrp=',nyn_dvrp, ' nyn=',nyn
+             CALL local_flush( 9 )
              CALL DVRP_VISUALIZE( m-1, 21, dvrp_filecount )
+             WRITE(9,*) '#8.9'
+             CALL local_flush( 9 )
 
              CALL cpu_log( log_point_s(26), 'dvrp_isosurface', 'stop' )
+                      WRITE(9,*) '#9'
+                      CALL local_flush( 9 )
 
           ELSEIF ( mode_dvrp(m)(1:6) == 'slicer' )  THEN
+                      WRITE(9,*) '#10'
+                      CALL local_flush( 9 )
 !
 !--          DVRP-Calls for plotting slicers:
@@ -547,4 +586,6 @@
 
              CALL cpu_log( log_point_s(27), 'dvrp_slicer', 'stop' )
+                      WRITE(9,*) '#11'
+                      CALL local_flush( 9 )
 
           ENDIF
@@ -552,4 +593,6 @@
           DEALLOCATE( local_pf )
 
+                      WRITE(9,*) '#12'
+                      CALL local_flush( 9 )
        ELSEIF ( mode_dvrp(m)(1:9) == 'pathlines' ) THEN
 
@@ -588,4 +631,6 @@
 
     ENDDO
+                      WRITE(9,*) '#13'
+                      CALL local_flush( 9 )
 
     dvrp_filecount = dvrp_filecount + 1
Index: palm/trunk/SOURCE/init_dvrp.f90
===================================================================
--- palm/trunk/SOURCE/init_dvrp.f90	(revision 263)
+++ palm/trunk/SOURCE/init_dvrp.f90	(revision 264)
@@ -4,4 +4,5 @@
 ! Current revisions:
 ! -----------------
+! Definition of a colortable to be used for particles.
 ! Output names are changed: surface=groundplate, buildings=topography
 ! Output of messages replaced by message handling routine.
@@ -112,8 +113,12 @@
        dvrp_overlap       = .FALSE.
        dvrp_total_overlap = .FALSE.
-       nxl_dvrp = nx_dvrp_l
-       nxr_dvrp = nxl_dvrp + 4
-       nys_dvrp = ny_dvrp_s
-       nyn_dvrp = nys_dvrp + 4
+!       nxl_dvrp = nx_dvrp_l
+!       nxr_dvrp = nxl_dvrp + 4
+!       nys_dvrp = ny_dvrp_s
+!       nyn_dvrp = nys_dvrp + 4
+       nxl_dvrp = nxl
+       nxr_dvrp = MIN( nxl+4, nxr )
+       nys_dvrp = nys
+       nyn_dvrp = MIN( nys+4, nyn )
 
     ENDIF
@@ -250,6 +255,6 @@
           tmp_3 = clip_dvrp_s * superelevation_y
           tmp_4 = 0.0
-          tmp_5 = clip_dvrp_r * superelevation_x
-          tmp_6 = clip_dvrp_n * superelevation_y
+          tmp_5 = (clip_dvrp_r+dx) * superelevation_x
+          tmp_6 = (clip_dvrp_n+dy) * superelevation_y
           tmp_7 = zu(nz_do3d) * superelevation
           CALL DVRP_BOUNDINGBOX( m-1, 1, tmp_1, tmp_2, tmp_3, tmp_4, tmp_5, &
@@ -644,4 +649,23 @@
        ENDIF
 
+       IF ( mode_dvrp(m) == 'particles' )  THEN
+!
+!--       Define a default colourtable for particles
+          DO  i = 1, 11
+             interval_values_dvrp_prt(1,i) = i - 1.0
+             interval_values_dvrp_prt(2,i) = REAL( i )
+             interval_h_dvrp_prt(:,i) = 270.0 - ( i - 1.0 ) * 9.0
+          ENDDO
+
+          DO  j = 12, 22
+             interval_values_dvrp_prt(1,i) = i - 1.0
+             interval_values_dvrp_prt(2,i) = REAL( i )
+             interval_h_dvrp_prt(:,i) = 70.0 - ( i - 12.0 ) * 9.5
+          ENDDO
+
+          dvrp_colortable_entries_prt = 22
+
+       ENDIF
+
        m = m + 1
 
Index: palm/trunk/SOURCE/modules.f90
===================================================================
--- palm/trunk/SOURCE/modules.f90	(revision 263)
+++ palm/trunk/SOURCE/modules.f90	(revision 264)
@@ -5,6 +5,7 @@
 ! Actual revisions:
 ! -----------------
-! +clip_dvrp_*, cluster_size, dvrp_overlap, dvrp_total_overlap,
-! local_dvrserver_running, n*_dvrp in dvrp_variables
+! +clip_dvrp_*, cluster_size, color_interval, dvrpsize_interval, dvrp_overlap,
+! dvrp_total_overlap, local_dvrserver_running, n*_dvrp, interval_*_dvrp_prt,
+! particle_color, particle_dvrpsize in dvrp_variables
 ! +canyon_height, canyon_width_x, canyon_width_y, canyon_wall_left, 
 ! canyon_wall_south, conserve_volume_flow_mode, dp_external, dp_level_b, 
@@ -527,5 +528,6 @@
 !------------------------------------------------------------------------------!
 
-    CHARACTER (LEN=10) ::  dvrp_output = 'rtsp'
+    CHARACTER (LEN=10) ::  dvrp_output = 'rtsp', particle_color = 'none', &
+                           particle_dvrpsize = 'none'
 
     CHARACTER (LEN=20), DIMENSION(10) ::  mode_dvrp = &
@@ -538,5 +540,6 @@
                            dvrp_username  = ' '
 
-    INTEGER ::  cluster_size = 1, dvrp_colourtable_entries = 4, islice_dvrp,   &
+    INTEGER ::  cluster_size = 1, dvrp_colortable_entries = 4,                 &
+                dvrp_colortable_entries_prt = 22, islice_dvrp,                 &
                 nx_dvrp, nxl_dvrp, nxr_dvrp, ny_dvrp, nyn_dvrp, nys_dvrp,      &
                 nz_dvrp, pathlines_fadeintime = 5, pathlines_fadeouttime = 5,  &
@@ -557,4 +560,7 @@
                 superelevation_y = 1.0, vc_alpha = 38.0
 
+    REAL, DIMENSION(2) ::  color_interval = (/ 0.0, 1.0 /), &
+                           dvrpsize_interval = (/ 0.0, 1.0 /)
+
 #if defined( __decalpha )
     REAL, DIMENSION(2,10)  ::  slicer_range_limits_dvrp = (/ &
@@ -563,15 +569,25 @@
                                 -1.0, 1.0, -1.0, 1.0 /)
 
-    REAL/4), DIMENSION(2,100) ::  interval_values_dvrp, interval_h_dvrp =      &
+    REAL(4), DIMENSION(2,100) ::  interval_values_dvrp, interval_h_dvrp =      &
                                   (/ 270.0, 225.0, 225.0, 180.0, 70.0, 25.0,   &
                                      25.0, -25.0, ( 0.0, i9 = 1, 192 ) /),     &
                                   interval_l_dvrp = 0.5, interval_s_dvrp = 1.0,&
-                                  interval_a_dvrp = 0.0
+                                  interval_a_dvrp = 0.0,                       &
+                                  interval_values_dvrp_prt,                    &
+                                  interval_h_dvrp_prt,                         &
+                                  (/ 270.0, 225.0, 225.0, 180.0, 70.0, 25.0,   &
+                                     25.0, -25.0, ( 0.0, i9 = 1, 192 ) /),     &
+                                  interval_l_dvrp_prt = 0.5,                   &
+                                  interval_s_dvrp_prt = 1.0,                   &
+                                  interval_a_dvrp_prt = 0.0
 #else
     REAL, DIMENSION(2,10)  ::  slicer_range_limits_dvrp
 
-    REAL(4), DIMENSION(2,100) ::  interval_values_dvrp, interval_h_dvrp,       &
-                                  interval_l_dvrp = 0.5, interval_s_dvrp = 1.0,&
-                                  interval_a_dvrp = 0.0
+    REAL(4), DIMENSION(2,100) ::  interval_values_dvrp,                       &
+                                  interval_values_dvrp_prt, interval_h_dvrp,  &
+                                  interval_h_dvrp_prt, interval_l_dvrp = 0.5, &
+                                  interval_l_dvrp_prt, interval_s_dvrp = 1.0, &
+                                  interval_s_dvrp_prt, interval_a_dvrp = 0.0, &
+                                  interval_a_dvrp_prt
 
     DATA  slicer_range_limits_dvrp / -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, &
Index: palm/trunk/SOURCE/package_parin.f90
===================================================================
--- palm/trunk/SOURCE/package_parin.f90	(revision 263)
+++ palm/trunk/SOURCE/package_parin.f90	(revision 264)
@@ -4,5 +4,6 @@
 ! Actual revisions:
 ! -----------------
-! +clip_dvrp_*, cluster_size in dvrp_graphics_par
+! +clip_dvrp_*, cluster_size, color_interval, dvrpsize_interval, particle_color,
+! dvrpsize_color in dvrp_graphics_par
 ! Variables for dvrp-mode pathlines added
 !
@@ -46,15 +47,17 @@
 
     NAMELIST /dvrp_graphics_par/  clip_dvrp_l, clip_dvrp_n, clip_dvrp_r,       &
-                                  clip_dvrp_s, cluster_size, dt_dvrp,          &
-                                  dvrp_directory, dvrp_file, dvrp_host,        &
-                                  dvrp_output, dvrp_password, dvrp_username,   &
-                                  mode_dvrp, pathlines_fadeintime,             &
-                                  pathlines_fadeouttime, pathlines_linecount,  &
-                                  pathlines_maxhistory, pathlines_wavecount,   &
-                                  pathlines_wavetime, slicer_range_limits_dvrp,&
-                                  superelevation, superelevation_x,            &
-                                  superelevation_y, threshold, vc_alpha,       &
-                                  vc_gradient_normals, vc_mode, vc_size_x,     &
-                                  vc_size_y, vc_size_z
+                                  clip_dvrp_s, cluster_size, color_interval,   &
+                                  dt_dvrp, dvrpsize_interval, dvrp_directory,  &
+                                  dvrp_file, dvrp_host, dvrp_output,           &
+                                  dvrp_password, dvrp_username, mode_dvrp,     &
+                                  particle_color, particle_dvrpsize,           &
+                                  pathlines_fadeintime, pathlines_fadeouttime, &
+                                  pathlines_linecount, pathlines_maxhistory,   &
+                                  pathlines_wavecount, pathlines_wavetime,     &
+                                  slicer_range_limits_dvrp, superelevation,    &
+                                  superelevation_x, superelevation_y,          &
+                                  threshold, vc_alpha, vc_gradient_normals,    &
+                                  vc_mode, vc_size_x, vc_size_y, vc_size_z
+
     NAMELIST /particles_par/      bc_par_b, bc_par_lr, bc_par_ns, bc_par_t,    &
                                   density_ratio, radius, dt_dopts,             &
Index: palm/trunk/SOURCE/set_particle_attributes.f90
===================================================================
--- palm/trunk/SOURCE/set_particle_attributes.f90	(revision 264)
+++ palm/trunk/SOURCE/set_particle_attributes.f90	(revision 264)
@@ -0,0 +1,161 @@
+ SUBROUTINE set_particle_attributes
+
+!------------------------------------------------------------------------------!
+! Actual revisions:
+! -----------------
+! 
+!
+! Former revisions:
+! -----------------
+! $Id$
+!
+! Initial version (2009/03/21)
+!
+! Description:
+! ------------
+! This routine sets certain particle attributes depending on the values that
+! other PALM variables have at the current particle position.
+!------------------------------------------------------------------------------!
+
+    USE arrays_3d
+    USE control_parameters
+    USE cpulog
+    USE dvrp_variables
+    USE grid_variables
+    USE indices
+    USE interfaces
+    USE particle_attributes
+    USE pegrid
+    USE statistics
+
+    IMPLICIT NONE
+
+    INTEGER ::  i, j, k, n
+    REAL    ::  aa, bb, cc, dd, gg, pt_int, pt_int_l, pt_int_u, w_int, &
+                w_int_l, w_int_u, x, y
+
+
+    CALL cpu_log( log_point_s(49), 'set_particle_attributes', 'start' )
+
+!
+!-- Set particle color
+    IF ( particle_color == 'absuv' )  THEN
+
+    ELSEIF ( particle_color == 'pt*' )  THEN
+!
+!--    Set particle color depending on the resolved scale temperature
+!--    fluctuation.
+!--    First, calculate the horizontal average of the potential temperature
+!--    (This is also done in flow_statistics, but flow_statistics is called
+!--    after this routine.)
+       sums_l(:,4,0) = 0.0
+       DO  i = nxl, nxr
+          DO  j =  nys, nyn
+             DO  k = nzb, nzt+1
+                sums_l(k,4,0) = sums_l(k,4,0) + pt(k,j,i)
+             ENDDO
+          ENDDO
+       ENDDO
+
+#if defined( __parallel )
+
+       CALL MPI_ALLREDUCE( sums_l(nzb,4,0), sums(nzb,4), nzt+2-nzb, &
+                           MPI_REAL, MPI_SUM, comm2d, ierr )
+
+#else
+
+       sums(:,4) = sums_l(:,4,0)
+
+#endif
+
+       sums(:,4) = sums(:,4) / ngp_2dh(0)
+
+       DO  n = 1, number_of_particles
+!
+!--       Interpolate temperature to the current particle position
+          i = particles(n)%x * ddx
+          j = particles(n)%y * ddy
+          k = ( particles(n)%z + 0.5 * dz * atmos_ocean_sign ) / dz  &
+              + offset_ocean_nzt                     ! only exact if equidistant
+
+          x  = particles(n)%x - i * dx
+          y  = particles(n)%y - j * dy
+          aa = x**2          + y**2
+          bb = ( dx - x )**2 + y**2
+          cc = x**2          + ( dy - y )**2
+          dd = ( dx - x )**2 + ( dy - y )**2
+          gg = aa + bb + cc + dd
+
+          pt_int_l = ( ( gg - aa ) * pt(k,j,i)   + ( gg - bb ) * pt(k,j,i+1)   &
+                     + ( gg - cc ) * pt(k,j+1,i) + ( gg - dd ) * pt(k,j+1,i+1) &
+                     ) / ( 3.0 * gg ) - sums(k,4)
+
+          pt_int_u = ( ( gg-aa ) * pt(k+1,j,i)   + ( gg-bb ) * pt(k+1,j,i+1)   &
+                     + ( gg-cc ) * pt(k+1,j+1,i) + ( gg-dd ) * pt(k+1,j+1,i+1) &
+                     ) / ( 3.0 * gg ) - sums(k,4)
+
+          pt_int = pt_int_l + ( particles(n)%z - zu(k) ) / dz * &
+                              ( pt_int_u - pt_int_l )
+
+!
+!--       Limit values by the given interval and normalize to interval [0,1]
+          pt_int = MIN( pt_int, color_interval(2) )
+          pt_int = MAX( pt_int, color_interval(1) )
+
+          pt_int = ( pt_int - color_interval(1) ) / &
+                   ( color_interval(2) - color_interval(1) )
+
+!
+!--       Number of available colors is defined in init_dvrp
+          particles(n)%color = 1 + pt_int * ( dvrp_colortable_entries_prt - 1 )
+
+    ENDDO
+
+    ELSEIF ( particle_color == 'z' )  THEN
+
+    ENDIF
+
+!
+!-- Set particle size for dvrp graphics
+    IF ( particle_dvrpsize == 'absw' )  THEN
+
+       DO  n = 1, number_of_particles
+!
+!--       Interpolate w-component to the current particle position
+          i = particles(n)%x * ddx
+          j = particles(n)%y * ddy
+          k = particles(n)%z / dz
+
+          x  = particles(n)%x - i * dx
+          y  = particles(n)%y - j * dy
+          aa = x**2          + y**2
+          bb = ( dx - x )**2 + y**2
+          cc = x**2          + ( dy - y )**2
+          dd = ( dx - x )**2 + ( dy - y )**2
+          gg = aa + bb + cc + dd
+
+          w_int_l = ( ( gg - aa ) * w(k,j,i)   + ( gg - bb ) * w(k,j,i+1)   &
+                    + ( gg - cc ) * w(k,j+1,i) + ( gg - dd ) * w(k,j+1,i+1) &
+                    ) / ( 3.0 * gg )
+
+          IF ( k+1 == nzt+1 )  THEN
+             w_int = w_int_l
+          ELSE
+             w_int_u = ( ( gg-aa ) * w(k+1,j,i)   + ( gg-bb ) * w(k+1,j,i+1)   &
+                       + ( gg-cc ) * w(k+1,j+1,i) + ( gg-dd ) * w(k+1,j+1,i+1) &
+                       ) / ( 3.0 * gg )
+             w_int = w_int_l + ( particles(n)%z - zw(k) ) / dz * &
+                               ( w_int_u - w_int_l )
+          ENDIF
+
+          particles(n)%dvrp_psize = ( 0.2 + MIN( 0.6, ABS( w_int ) / 6.0 ) ) * &
+                                    dx * 1.4
+
+       ENDDO
+
+    ENDIF
+
+    CALL cpu_log( log_point_s(49), 'set_particle_attributes', 'stop' )
+
+
+ END SUBROUTINE set_particle_attributes
