Index: palm/trunk/SOURCE/virtual_measurement_mod.f90
===================================================================
--- palm/trunk/SOURCE/virtual_measurement_mod.f90	(revision 4705)
+++ palm/trunk/SOURCE/virtual_measurement_mod.f90	(revision 4706)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Revise setting of measurement height above ground for trajectory measurements
+!
+! 4695 2020-09-24 11:30:03Z suehring
 ! Introduce additional namelist parameters to further customize sampling in the horizontal and
 ! vertical surroundings of the original observation coordinates
@@ -325,6 +328,6 @@
        REAL(wp) ::  origin_y_obs       !< origin of the observation in UTM coordiates in y-direction
 
-       REAL(wp), DIMENSION(:), ALLOCATABLE   ::  depth         !< measurement depth in soil
-       REAL(wp), DIMENSION(:), ALLOCATABLE   ::  zar           !< measurement height above ground level
+       REAL(wp), DIMENSION(:), ALLOCATABLE   ::  depth  !< measurement depth in soil
+       REAL(wp), DIMENSION(:), ALLOCATABLE   ::  zar    !< measurement height above ground level
 
        REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  measured_vars       !< measured variables
@@ -339,4 +342,5 @@
     ! This need to be generalized
     CHARACTER(LEN=10) ::  char_fill = '_FillValue'                 !< attribute name for fill value
+    CHARACTER(LEN=1)  ::  char_height = "height"                   !< attribute name indicating height above surface (for trajectories only)
     CHARACTER(LEN=9)  ::  char_long = 'long_name'                  !< attribute name for long_name
     CHARACTER(LEN=18) ::  char_mv = "measured_variables"           !< variable name for the array with the measured variable names
@@ -822,7 +826,10 @@
 
     NAMELIST /virtual_measurement_parameters/  dt_virtual_measurement,                             &
+                                               off_pr,                                             &
+                                               off_pr_z,                                           &
+                                               off_tr,                                             &
+                                               off_tr_z,                                           &
                                                off_ts,                                             &
-                                               off_pr,                                             &
-                                               off_tr,                                             &
+                                               off_ts_z,                                           &
                                                use_virtual_measurement,                            &
                                                vm_time_start
@@ -902,10 +909,12 @@
     LOGICAL  ::  on_pe  !< flag indicating that the respective measurement coordinate is on subdomain
 
-    REAL(wp) ::  fill_eutm !< _FillValue for coordinate array E_UTM
-    REAL(wp) ::  fill_nutm !< _FillValue for coordinate array N_UTM
-    REAL(wp) ::  fill_zar  !< _FillValue for height coordinate
+    REAL(wp) ::  fill_height !< _FillValue for height coordinate (for trajectories)
+    REAL(wp) ::  fill_eutm   !< _FillValue for coordinate array E_UTM
+    REAL(wp) ::  fill_nutm   !< _FillValue for coordinate array N_UTM
+    REAL(wp) ::  fill_zar    !< _FillValue for zar coordinate
 
     REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  e_utm      !< easting UTM coordinate, temporary variable
     REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  e_utm_tmp  !< EUTM coordinate before rotation
+    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  height     !< observation height above ground (for trajectories)
     REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  n_utm      !< northing UTM coordinate, temporary variable
     REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  n_utm_tmp  !< NUTM coordinate before rotation
@@ -1052,8 +1061,10 @@
           ALLOCATE( station_h(1:vmea(l)%n_tr_st,1:dim_ntime)   )
           ALLOCATE( zar(1:vmea(l)%n_tr_st,1:dim_ntime)         )
+          IF ( vmea(l)%trajectory )  ALLOCATE( height(1:vmea(l)%n_tr_st,1:dim_ntime) )
           e_utm     = 0.0_wp
           n_utm     = 0.0_wp
           station_h = 0.0_wp
           zar       = 0.0_wp
+          IF ( vmea(l)%trajectory )  height = 0.0_wp
 
           ALLOCATE( e_utm_tmp(1:vmea(l)%n_tr_st,1:dim_ntime) )
@@ -1062,7 +1073,9 @@
 !--       Read UTM and height coordinates for all trajectories and times. Note, in case
 !--       these obtain any missing values, replace them with default _FillValues.
-          CALL inquire_fill_value( pids_id, char_eutm // TRIM( dum ), nofill, fill_eutm )
-          CALL inquire_fill_value( pids_id, char_nutm // TRIM( dum ), nofill, fill_nutm )
-          CALL inquire_fill_value( pids_id, char_zar // TRIM( dum ), nofill, fill_zar )
+          CALL inquire_fill_value( pids_id, char_eutm   // TRIM( dum ), nofill, fill_eutm    )
+          CALL inquire_fill_value( pids_id, char_nutm   // TRIM( dum ), nofill, fill_nutm    )
+          CALL inquire_fill_value( pids_id, char_zar    // TRIM( dum ), nofill, fill_zar     )
+          IF ( vmea(l)%trajectory )                                                                &
+             CALL inquire_fill_value( pids_id, char_height // TRIM( dum ), nofill, fill_height  )
 !
 !--       Further line is just to avoid compiler warnings. nofill might be used in future.
@@ -1073,15 +1086,17 @@
 !--       UC2-standard.
           IF ( vmea(l)%trajectory )  THEN
-             CALL get_variable( pids_id, char_eutm // TRIM( dum ), e_utm, 0, dim_ntime-1, 0,       &
+             CALL get_variable( pids_id, char_eutm   // TRIM( dum ), e_utm,  0, dim_ntime-1, 0,    &
                                 vmea(l)%n_tr_st-1 )
-             CALL get_variable( pids_id, char_nutm // TRIM( dum ), n_utm, 0, dim_ntime-1, 0,       &
+             CALL get_variable( pids_id, char_nutm   // TRIM( dum ), n_utm,  0, dim_ntime-1, 0,    &
                                 vmea(l)%n_tr_st-1 )
-             CALL get_variable( pids_id, char_zar // TRIM( dum ), zar, 0, dim_ntime-1, 0,          &
+             CALL get_variable( pids_id, char_zar    // TRIM( dum ), zar,    0, dim_ntime-1, 0,    &
+                                vmea(l)%n_tr_st-1 )
+             CALL get_variable( pids_id, char_height // TRIM( dum ), height, 0, dim_ntime-1, 0,    &
                                 vmea(l)%n_tr_st-1 )
           ELSE
-             CALL get_variable( pids_id, char_eutm // TRIM( dum ), e_utm(:,1) )
-             CALL get_variable( pids_id, char_nutm // TRIM( dum ), n_utm(:,1) )
+             CALL get_variable( pids_id, char_eutm      // TRIM( dum ), e_utm(:,1)     )
+             CALL get_variable( pids_id, char_nutm      // TRIM( dum ), n_utm(:,1)     )
              CALL get_variable( pids_id, char_station_h // TRIM( dum ), station_h(:,1) )
-             CALL get_variable( pids_id, char_zar // TRIM( dum ), zar(:,1) )
+             CALL get_variable( pids_id, char_zar       // TRIM( dum ), zar(:,1)       )
           ENDIF
 
@@ -1089,7 +1104,15 @@
           n_utm = MERGE( n_utm, vmea(l)%fillout, n_utm /= fill_nutm )
           zar   = MERGE( zar,   vmea(l)%fillout, zar   /= fill_zar  )
-!
-!--       Compute observation height above ground.
-          zar  = zar - station_h
+          IF ( vmea(l)%trajectory )                                                                &
+             height = MERGE( height, vmea(l)%fillout, height /= fill_height )
+!
+!--       Compute observation height above ground. Note, for trajectory measurements the height
+!--       above the surface is actually stored in 'height'
+          IF ( vmea(l)%trajectory )  THEN
+             zar      = height
+             fill_zar = fill_height
+          ELSE
+             zar = zar - station_h
+          ENDIF
 !
 !--       Based on UTM coordinates, check if the measurement station or parts of the trajectory are
@@ -1288,4 +1311,5 @@
           IF ( ALLOCATED( zar  )      )  DEALLOCATE( vmea(l)%dim_t )
           IF ( ALLOCATED( zar  )      )  DEALLOCATE( zar  )
+          IF ( ALLOCATED( height )    )  DEALLOCATE( height )
           IF ( ALLOCATED( station_h ) )  DEALLOCATE( station_h )
 
