Index: /palm/trunk/SOURCE/land_surface_model_mod.f90
===================================================================
--- /palm/trunk/SOURCE/land_surface_model_mod.f90	(revision 2722)
+++ /palm/trunk/SOURCE/land_surface_model_mod.f90	(revision 2723)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! Bugfix for spinups (end_time was increased twice in case of LSM + USM runs)
+! 
+! 2718 2018-01-02 08:49:38Z maronga
 ! Corrected "Former revisions" section
 ! 
@@ -1976,14 +1979,4 @@
 
 
-!        IF ( surf_t_surface_p%var_1d(m) < 270.0_wp .OR. surf_t_surface_p%var_1d(m) > 350.0_wp )  THEN
-!           PRINT*, "myid: ", myid
-!           PRINT*, "m: ", m
-!           PRINT*, "i,j,k: ", i, j, k
-!           PRINT*, "pt_p: ", surf_t_surface_p%var_1d(m)
-!           PRINT*, "f_shf: ", f_shf
-!           PRINT*, "f_qsws: ", f_qsws
-!         ENDIF
-
-
 !        pt(k+k_off,j+j_off,i+i_off) = surf_t_surface_p%var_1d(m) / exn  !is actually no air temperature
        surf%pt_surface(m)          = surf_t_surface_p%var_1d(m) / exn
@@ -4527,9 +4520,11 @@
        IF ( spinup_time > 0.0_wp )  THEN
           coupling_start_time = spinup_time
-          end_time = end_time + spinup_time
           IF ( spinup_pt_mean == 9999999.9_wp )  THEN
              spinup_pt_mean = pt_surface
           ENDIF
-          spinup = .TRUE.
+          IF ( .NOT. spinup )  THEN
+             end_time = end_time + spinup_time
+             spinup = .TRUE.
+          ENDIF
        ENDIF
 
Index: /palm/trunk/SOURCE/radiation_model_mod.f90
===================================================================
--- /palm/trunk/SOURCE/radiation_model_mod.f90	(revision 2722)
+++ /palm/trunk/SOURCE/radiation_model_mod.f90	(revision 2723)
@@ -25,4 +25,8 @@
 ! -----------------
 ! $Id$
+! Bugfix in calculation of rad_lw_out (clear-sky). First grid level was used
+! instead of the surface value
+! 
+! 2718 2018-01-02 08:49:38Z maronga
 ! Corrected "Former revisions" section
 ! 
@@ -2370,9 +2374,6 @@
 
              INTEGER(iwp) ::  i         !< index x-direction
-             INTEGER(iwp) ::  ioff      !< offset between surface element and adjacent grid point along x
              INTEGER(iwp) ::  j         !< index y-direction
-             INTEGER(iwp) ::  joff      !< offset between surface element and adjacent grid point along y
              INTEGER(iwp) ::  k         !< index z-direction
-             INTEGER(iwp) ::  koff      !< offset between surface element and adjacent grid point along z
              INTEGER(iwp) ::  m         !< running index for surface elements
 
@@ -2385,8 +2386,5 @@
 
                 k = nzut
-!
-!-- MS: Why k+1 ? 
-!-- MS: @Mohamed: emissivity belongs now to surface type with 3 different values for each 
-!--               surface element (due to tile approach).
+
                 exn1 = ( hyp(k+1) / 100000.0_wp )**0.286_wp
 
@@ -2409,10 +2407,4 @@
 !--          element.
              ELSE
-!
-!--             Determine index offset between surface element and adjacent 
-!--             atmospheric grid point (depends on surface orientation). 
-                ioff = surf%ioff
-                joff = surf%joff
-                koff = surf%koff
 
                 DO  m = 1, surf%ns
@@ -2429,4 +2421,7 @@
 !--                no weighted averaging is performed ( only one surface type per
 !--                surface element ). 
+!--                ATTENTION: when radiation interactions are switched on the 
+!--                calculated fluxes below are not actually used as they are 
+!--                overwritten in radiation_interaction.
                    IF ( ALLOCATED( surf%frac ) )  THEN
 
@@ -2441,6 +2436,5 @@
                                            )                                      &
                                            * sigma_sb                             &
-                                           * ( pt(k+koff,j+joff,i+ioff) * exn )**4
-
+                                           * ( surf%pt_surface(m) * exn )**4
 
                       surf%rad_lw_out_change_0(m) =                               &
@@ -2449,5 +2443,5 @@
                                          + surf%frac(2,m) * surf%emissivity(2,m)  &
                                          ) * 3.0_wp * sigma_sb                    &
-                                         * ( pt(k+koff,j+joff,i+ioff) * exn )** 3
+                                         * ( surf%pt_surface(m) * exn )** 3
 
                    ELSE
@@ -2457,10 +2451,9 @@
                       surf%rad_lw_out(m) = surf%emissivity(0,m)                   &
                                            * sigma_sb                             &
-                                           * ( pt(k+koff,j+joff,i+ioff) * exn )**4
-
+                                           * ( surf%pt_surface(m) * exn )**4
 
                       surf%rad_lw_out_change_0(m) = surf%emissivity(0,m)          &
                                            * 3.0_wp * sigma_sb                    &
-                                           * ( pt(k+koff,j+joff,i+ioff) * exn )** 3
+                                           * ( surf%pt_surface(m) * exn )** 3
 
                    ENDIF
@@ -2635,5 +2628,5 @@
                                            )                                      &
                                          * sigma_sb                               &
-                                         * ( pt(k+koff,j+joff,i+ioff) * exn )**4
+                                         * ( surf%pt_surface(m) * exn )**4
 
                       surf%rad_sw_in(m) = ( surf%rad_net(m) - surf%rad_lw_in(m)   &
@@ -2653,5 +2646,5 @@
                       surf%rad_lw_out(m) = surf%emissivity(0,m)                   &
                                          * sigma_sb                               &
-                                         * ( pt(k+koff,j+joff,i+ioff) * exn )**4
+                                         * ( surf%pt_surface(m) * exn )**4
 
                       surf%rad_sw_in(m) = ( surf%rad_net(m) - surf%rad_lw_in(m)   &
@@ -5363,5 +5356,5 @@
                                             surfinlw(i) - surfoutlw(i)
 !
-!--        westward-facding
+!--        westward-facing
            ELSEIF ( surfl(1,i) == iwest_l )  THEN
               surf_lsm_v(3)%rad_sw_in(m)  = surfinsw(i) 
Index: /palm/trunk/SOURCE/time_integration_spinup.f90
===================================================================
--- /palm/trunk/SOURCE/time_integration_spinup.f90	(revision 2722)
+++ /palm/trunk/SOURCE/time_integration_spinup.f90	(revision 2723)
@@ -25,4 +25,10 @@
 ! -----------------
 ! $Id$
+! Bugfix: array rad_sw_in no longer exists and is thus removed from RUN_CONTROL
+! output.
+! Added output of XY and 3D data during spinup.
+! Bugfix: time step in LSM and USM was set to dt_3d instead of dt_spinup
+! 
+! 2718 2018-01-02 08:49:38Z maronga
 ! Corrected "Former revisions" section
 ! 
@@ -62,11 +68,11 @@
                coupling_start_time, current_timestep_number,                   &
                data_output_during_spinup, disturbance_created, dopr_n, do_sum, &
-               dt_averaging_input_pr, dt_dopr, dt_dots, dt_run_control,        &
-               dt_spinup, humidity, intermediate_timestep_count,               &
+               dt_averaging_input_pr, dt_dopr, dt_dots, dt_do2d_xy, dt_do3d, dt_run_control,        &
+               dt_spinup, dt_3d, humidity, intermediate_timestep_count,               &
                intermediate_timestep_count_max, land_surface,                  &
-               simulated_time, simulated_time_chr,      &
-               skip_time_dopr, spinup, spinup_pt_amplitude, spinup_pt_mean,    &
-               spinup_time, timestep_count, timestep_scheme, time_dopr,        &
-               time_dopr_av, time_dots, time_run_control,                      &
+               simulated_time, simulated_time_chr,                             &
+               skip_time_dopr, skip_time_do2d_xy, skip_time_do3d, spinup, spinup_pt_amplitude, &
+               spinup_pt_mean, spinup_time, timestep_count, timestep_scheme,   &
+               time_dopr, time_dopr_av, time_dots, time_do2d_xy, time_do3d, time_run_control,           &
                time_since_reference_point, urban_surface
 
@@ -130,4 +136,5 @@
 
     REAL(wp) ::  pt_spinup   !< temporary storage of temperature
+    REAL(wp) ::  dt_save     !< temporary storage for time step
                   
     REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  pt_save   !< temporary storage of temperature
@@ -137,4 +144,7 @@
     CALL exchange_horiz( pt, nbgp )   
     pt_save = pt
+
+    dt_save = dt_3d
+    dt_3d   = dt_spinup
 
     CALL location_message( 'starting spinup-sequence', .TRUE. )
@@ -299,5 +309,5 @@
                == intermediate_timestep_count_max )  THEN
 
-               time_radiation = time_radiation + dt_spinup
+               time_radiation = time_radiation + dt_3d
 
              IF ( time_radiation >= dt_radiation .OR. force_radiation_call )   &
@@ -327,14 +337,20 @@
 !--    Increase simulation time and output times
        current_timestep_number_spinup = current_timestep_number_spinup + 1
-       simulated_time             = simulated_time   + dt_spinup
+       simulated_time             = simulated_time   + dt_3d
        simulated_time_chr         = time_to_string( simulated_time )
        time_since_reference_point = simulated_time - coupling_start_time
 
        IF ( data_output_during_spinup )  THEN
-          time_dots          = time_dots        + dt_spinup
+          IF ( simulated_time >= skip_time_do2d_xy )  THEN
+             time_do2d_xy       = time_do2d_xy     + dt_3d
+          ENDIF
+          IF ( simulated_time >= skip_time_do3d    )  THEN
+             time_do3d          = time_do3d        + dt_3d
+          ENDIF
+          time_dots          = time_dots        + dt_3d
           IF ( simulated_time >= skip_time_dopr )  THEN
-             time_dopr       = time_dopr        + dt_spinup
-          ENDIF
-          time_run_control   = time_run_control + dt_spinup
+             time_dopr       = time_dopr        + dt_3d
+          ENDIF
+          time_run_control   = time_run_control + dt_3d
 
 !
@@ -354,9 +370,9 @@
              ( dt_dopr - time_dopr ) <= averaging_interval_pr  .AND.           &
              simulated_time >= skip_time_dopr )  THEN
-             time_dopr_av = time_dopr_av + dt_spinup
+             time_dopr_av = time_dopr_av + dt_3d
              IF ( time_dopr_av >= dt_averaging_input_pr )  THEN
                 do_sum = .TRUE.
                 time_dopr_av = MOD( time_dopr_av,                              &
-                               MAX( dt_averaging_input_pr, dt_spinup ) )
+                               MAX( dt_averaging_input_pr, dt_3d ) )
              ENDIF
           ENDIF
@@ -367,5 +383,5 @@
           IF ( time_dopr >= dt_dopr )  THEN
              IF ( dopr_n /= 0 )  CALL data_output_profiles
-             time_dopr = MOD( time_dopr, MAX( dt_dopr, dt_spinup ) )
+             time_dopr = MOD( time_dopr, MAX( dt_dopr, dt_3d ) )
              time_dopr_av = 0.0_wp    ! due to averaging (see above)
           ENDIF
@@ -375,6 +391,21 @@
           IF ( time_dots >= dt_dots )  THEN
              CALL data_output_tseries
-             time_dots = MOD( time_dots, MAX( dt_dots, dt_spinup ) )
-          ENDIF
+             time_dots = MOD( time_dots, MAX( dt_dots, dt_3d ) )
+          ENDIF
+
+!
+!--       2d-data output (cross-sections)
+          IF ( time_do2d_xy >= dt_do2d_xy )  THEN
+             CALL data_output_2d( 'xy', 0 )
+             time_do2d_xy = MOD( time_do2d_xy, MAX( dt_do2d_xy, dt_3d ) )
+          ENDIF
+
+!
+!--       3d-data output (volume data)
+          IF ( time_do3d >= dt_do3d )  THEN
+             CALL data_output_3d( 0 )
+             time_do3d = MOD( time_do3d, MAX( dt_do3d, dt_3d ) )
+          ENDIF
+
 
        ENDIF
@@ -389,5 +420,5 @@
 !           IF ( time_run_control >= dt_run_control )  THEN
 !              time_run_control = MOD( time_run_control,                         &
-!                                      MAX( dt_run_control, dt_spinup ) )
+!                                      MAX( dt_run_control, dt_3d ) )
 !           ENDIF
 !        ENDIF
@@ -408,5 +439,5 @@
 !
 !--       Write some general information about the spinup in run control file
-          WRITE ( 15, 101 )  current_timestep_number_spinup, simulated_time_chr, dt_spinup, pt_spinup, rad_sw_in(0,nysg,nxlg)
+          WRITE ( 15, 101 )  current_timestep_number_spinup, simulated_time_chr, dt_3d, pt_spinup
 !
 !--       Write buffer contents to disc immediately
@@ -423,4 +454,8 @@
     pt_p(:,:,:) = pt_save
 
+!
+!-- Reset time step
+    dt_3d = dt_save
+
     DEALLOCATE(pt_save)
 
@@ -431,7 +466,7 @@
 !-- Formats
 100 FORMAT (///'Spinup control output:'/  &
-            '----------------------------------------'// &
-            'ITER.  HH:MM:SS    DT   PT(z_MO)     SWD'/   &
-            '----------------------------------------')
+            '--------------------------------'// &
+            'ITER.  HH:MM:SS    DT   PT(z_MO)'/   &
+            '--------------------------------')
 101 FORMAT (I5,2X,A9,1X,F6.2,3X,F6.2,2X,F6.2)
 
Index: /palm/trunk/SOURCE/urban_surface_mod.f90
===================================================================
--- /palm/trunk/SOURCE/urban_surface_mod.f90	(revision 2722)
+++ /palm/trunk/SOURCE/urban_surface_mod.f90	(revision 2723)
@@ -26,4 +26,7 @@
 ! -----------------
 ! $Id$
+! Bugfix for spinups (end_time was increased twice in case of LSM + USM runs)
+! 
+! 2720 2018-01-02 16:27:15Z kanani
 ! Correction of comment
 ! 
@@ -5107,9 +5110,11 @@
        IF ( spinup_time > 0.0_wp )  THEN
           coupling_start_time = spinup_time
-          end_time = end_time + spinup_time
           IF ( spinup_pt_mean == 9999999.9_wp )  THEN
              spinup_pt_mean = pt_surface
           ENDIF
-          spinup = .TRUE.
+          IF ( .NOT. spinup )  THEN
+             end_time = end_time + spinup_time
+             spinup = .TRUE.
+          ENDIF
        ENDIF
 
