Index: palm/trunk/SOURCE/netcdf_data_input_mod.f90
===================================================================
--- palm/trunk/SOURCE/netcdf_data_input_mod.f90	(revision 3554)
+++ palm/trunk/SOURCE/netcdf_data_input_mod.f90	(revision 3556)
@@ -20,5 +20,5 @@
 ! Current revisions:
 ! -----------------
-! 
+! Some formatting adjustment
 ! 
 ! Former revisions:
@@ -946,18 +946,30 @@
                            input_file_atts%rotation_angle, .TRUE. )
 
-       CALL get_attribute( id_mod, input_file_atts%author_char,         input_file_atts%author,         .TRUE., no_abort=.FALSE. )
-       CALL get_attribute( id_mod, input_file_atts%contact_person_char, input_file_atts%contact_person, .TRUE., no_abort=.FALSE. )
-       CALL get_attribute( id_mod, input_file_atts%institution_char,    input_file_atts%institution,    .TRUE., no_abort=.FALSE. )
-       CALL get_attribute( id_mod, input_file_atts%acronym_char,        input_file_atts%acronym,        .TRUE., no_abort=.FALSE. )
-
-       CALL get_attribute( id_mod, input_file_atts%campaign_char, input_file_atts%campaign, .TRUE., no_abort=.FALSE. )
-       CALL get_attribute( id_mod, input_file_atts%location_char, input_file_atts%location, .TRUE., no_abort=.FALSE. )
-       CALL get_attribute( id_mod, input_file_atts%site_char,     input_file_atts%site,     .TRUE., no_abort=.FALSE. )
-
-       CALL get_attribute( id_mod, input_file_atts%source_char,     input_file_atts%source,     .TRUE., no_abort=.FALSE. )
-       CALL get_attribute( id_mod, input_file_atts%references_char, input_file_atts%references, .TRUE., no_abort=.FALSE. )
-       CALL get_attribute( id_mod, input_file_atts%keywords_char,   input_file_atts%keywords,   .TRUE., no_abort=.FALSE. )
-       CALL get_attribute( id_mod, input_file_atts%licence_char,    input_file_atts%licence,    .TRUE., no_abort=.FALSE. )
-       CALL get_attribute( id_mod, input_file_atts%comment_char,    input_file_atts%comment,    .TRUE., no_abort=.FALSE. )
+       CALL get_attribute( id_mod, input_file_atts%author_char,                &
+                           input_file_atts%author, .TRUE., no_abort=.FALSE. )
+       CALL get_attribute( id_mod, input_file_atts%contact_person_char,        &
+                           input_file_atts%contact_person, .TRUE., no_abort=.FALSE. )
+       CALL get_attribute( id_mod, input_file_atts%institution_char,           &
+                           input_file_atts%institution,    .TRUE., no_abort=.FALSE. )
+       CALL get_attribute( id_mod, input_file_atts%acronym_char,               &
+                           input_file_atts%acronym,        .TRUE., no_abort=.FALSE. )
+
+       CALL get_attribute( id_mod, input_file_atts%campaign_char,              &
+                           input_file_atts%campaign, .TRUE., no_abort=.FALSE. )
+       CALL get_attribute( id_mod, input_file_atts%location_char,              &
+                           input_file_atts%location, .TRUE., no_abort=.FALSE. )
+       CALL get_attribute( id_mod, input_file_atts%site_char,                  &
+                           input_file_atts%site,     .TRUE., no_abort=.FALSE. )
+
+       CALL get_attribute( id_mod, input_file_atts%source_char,                &
+                           input_file_atts%source,     .TRUE., no_abort=.FALSE. )
+       CALL get_attribute( id_mod, input_file_atts%references_char,            &
+                           input_file_atts%references, .TRUE., no_abort=.FALSE. )
+       CALL get_attribute( id_mod, input_file_atts%keywords_char,              &
+                           input_file_atts%keywords,   .TRUE., no_abort=.FALSE. )
+       CALL get_attribute( id_mod, input_file_atts%licence_char,               &
+                           input_file_atts%licence,    .TRUE., no_abort=.FALSE. )
+       CALL get_attribute( id_mod, input_file_atts%comment_char,               &
+                           input_file_atts%comment,    .TRUE., no_abort=.FALSE. )
 !
 !--    Read coordinate reference system if available
Index: palm/trunk/SOURCE/surface_mod.f90
===================================================================
--- palm/trunk/SOURCE/surface_mod.f90	(revision 3554)
+++ palm/trunk/SOURCE/surface_mod.f90	(revision 3556)
@@ -21,5 +21,5 @@
 ! Current revisions:
 ! ------------------
-! 
+! Surface restoring in restarts commented. Some formatting. 
 ! 
 ! Former revisions:
@@ -566,7 +566,7 @@
                                                     !< in the domain (required to activiate RTM)
 
-    LOGICAL ::  surf_bulk_cloud_model = .FALSE.                          !< use cloud microphysics
-    LOGICAL ::  surf_microphysics_morrison = .FALSE.                  !< use 2-moment Morrison (add. prog. eq. for nc and qc)
-    LOGICAL ::  surf_microphysics_seifert = .FALSE.                   !< use 2-moment Seifert and Beheng scheme
+    LOGICAL ::  surf_bulk_cloud_model = .FALSE.        !< use cloud microphysics
+    LOGICAL ::  surf_microphysics_morrison = .FALSE.   !< use 2-moment Morrison (add. prog. eq. for nc and qc)
+    LOGICAL ::  surf_microphysics_seifert = .FALSE.    !< use 2-moment Seifert and Beheng scheme
 
 
@@ -613,8 +613,6 @@
 !-- Public subroutines and functions
     PUBLIC get_topography_top_index, get_topography_top_index_ji, init_bc,     &
-           init_surfaces,                                                      &
-           init_surface_arrays, surface_rrd_local,                     & 
-           surface_restore_elements, surface_wrd_local,               &
-           surface_last_actions
+           init_surfaces, init_surface_arrays, surface_rrd_local,              & 
+           surface_restore_elements, surface_wrd_local, surface_last_actions
 
 
@@ -845,5 +843,5 @@
                       terrain  = BTEST( wall_flags_0(k,j-1,i), 5 )  .OR.       &
                                  topo_no_distinct
-                      building = BTEST( wall_flags_0(k,j-1,i), 6 )   .OR.       &
+                      building = BTEST( wall_flags_0(k,j-1,i), 6 )   .OR.      &
                                  topo_no_distinct
                       IF (  land_surface  .AND.  terrain )  THEN
@@ -859,7 +857,7 @@
 !--                   Unclassifified surface-grid point. Give error message.
                       ELSE 
-                         WRITE( message_string, * )                           &
-                                          'Unclassified northward-facing ' // &
-                                          'surface element at '//             &
+                         WRITE( message_string, * )                            &
+                                          'Unclassified northward-facing ' //  &
+                                          'surface element at '//              &
                                           'grid point (k,j,i) = ', k, j, i
                          CALL message( 'surface_mod', 'PA0999', 1, 2, 0, 6, 0 ) 
@@ -888,7 +886,7 @@
 !--                   Unclassifified surface-grid point. Give error message.
                       ELSE 
-                         WRITE( message_string, * )                           &
-                                          'Unclassified southward-facing ' // &
-                                          'surface element at '//             &
+                         WRITE( message_string, * )                            &
+                                          'Unclassified southward-facing ' //  &
+                                          'surface element at '//              &
                                           'grid point (k,j,i) = ', k, j, i
                          CALL message( 'surface_mod', 'PA0999', 1, 2, 0, 6, 0 ) 
@@ -917,7 +915,7 @@
 !--                   Unclassifified surface-grid point. Give error message.
                       ELSE 
-                         WRITE( message_string, * )                           &
-                                          'Unclassified eastward-facing ' //  &
-                                          'surface element at '//             &
+                         WRITE( message_string, * )                            &
+                                          'Unclassified eastward-facing ' //   &
+                                          'surface element at '//              &
                                           'grid point (k,j,i) = ', k, j, i
                          CALL message( 'surface_mod', 'PA0999', 1, 2, 0, 6, 0 ) 
@@ -946,7 +944,7 @@
 !--                   Unclassifified surface-grid point. Give error message.
                       ELSE 
-                         WRITE( message_string, * )                           &
-                                          'Unclassified westward-facing ' //  &
-                                          'surface element at '//             &
+                         WRITE( message_string, * )                            &
+                                          'Unclassified westward-facing ' //   &
+                                          'surface element at '//              &
                                           'grid point (k,j,i) = ', k, j, i
                          CALL message( 'surface_mod', 'PA0999', 1, 2, 0, 6, 0 ) 
@@ -1056,9 +1054,10 @@
 
 #if defined( __parallel )
-       CALL MPI_ALLREDUCE( num_surf_v_l, num_surf_v, 1, MPI_INTEGER, MPI_SUM, comm2d, ierr)
+       CALL MPI_ALLREDUCE( num_surf_v_l, num_surf_v, 1, MPI_INTEGER,           &
+                           MPI_SUM, comm2d, ierr)
 #else
        num_surf_v = num_surf_v_l
 #endif
-       IF ( num_surf_v > 0 ) vertical_surfaces_exist = .TRUE.
+       IF ( num_surf_v > 0 )  vertical_surfaces_exist = .TRUE.
         
 
@@ -2900,5 +2899,5 @@
           ENDDO
 !
-!--       Gather start- and end indices
+!--       Recalculate start- and end indices for gathered surface type.
           start_index_h(l) = 1                                        
           DO  i = nxl, nxr
@@ -2906,5 +2905,5 @@
 
                 surf_h(l)%start_index(j,i) = start_index_h(l)
-                surf_h(l)%end_index(j,i)   = surf_h(l)%start_index(j,i) -1
+                surf_h(l)%end_index(j,i)   = surf_h(l)%start_index(j,i) - 1
 
                 DO  m = surf_def_h(l)%start_index(j,i),                        &
@@ -2928,6 +2927,8 @@
           ENDDO
        ENDDO
-
-
+!
+!--    Treat vertically orientated surface. Again, gather data from different
+!--    surfaces types but identical orientation (e.g. northward-facing) onto
+!--    one surface type which is output afterwards. 
        mm(0:3) = 1
        DO  l = 0, 3
@@ -3137,5 +3138,5 @@
           ENDDO
 !
-!--       Gather start- and end indices
+!--       Recalculate start- and end-indices for gathered surface type
           start_index_v(l) = 1                                        
           DO  i = nxl, nxr
@@ -3163,9 +3164,12 @@
 
        ENDDO
-
-
+!
+!--    Output strings for the total number of upward / downward-facing surfaces 
+!--    on subdomain.
        CALL wrd_write_string( 'ns_h_on_file' )
        WRITE ( 14 )  ns_h_on_file
-
+!
+!--    Output strings for the total number of north/south/east/westward-facing surfaces 
+!--    on subdomain.
        CALL wrd_write_string( 'ns_v_on_file' )
        WRITE ( 14 )  ns_v_on_file
@@ -3173,5 +3177,6 @@
 !
 !--    Write required restart data.
-!--    Start with horizontal surfaces (upward-, downward-facing, and model top)
+!--    Start with horizontal surfaces (upward-, downward-facing, and model top).
+!--    Always start with %start_index followed by %end_index
        DO  l = 0, 2
           WRITE( dum, '(I1)')  l
@@ -3310,5 +3315,6 @@
        ENDDO
 !
-!--    Write vertical surfaces
+!--    Write vertical surfaces.
+!--    Always start with %start_index followed by %end_index.
        DO  l = 0, 3
           WRITE( dum, '(I1)')  l
@@ -3480,24 +3486,24 @@
        INTEGER(iwp)       ::  m           !< running index for surface elements, refers to gathered array encompassing all surface types
        INTEGER(iwp)       ::  mm          !< running index for surface elements, refers to individual surface types
-       INTEGER(iwp)       ::  ii               !< running index over input files
-       INTEGER(iwp)       ::  kk               !< running index over previous input files covering current local domain
-       INTEGER(iwp)       ::  nxlc             !< index of left boundary on current subdomain
-       INTEGER(iwp)       ::  nxlf             !< index of left boundary on former subdomain 
-       INTEGER(iwp)       ::  nxl_on_file      !< index of left boundary on former local domain 
-       INTEGER(iwp)       ::  nxrc             !< index of right boundary on current subdomain
-       INTEGER(iwp)       ::  nxrf             !< index of right boundary on former subdomain
-       INTEGER(iwp)       ::  nxr_on_file      !< index of right boundary on former local domain  
-       INTEGER(iwp)       ::  nync             !< index of north boundary on current subdomain
-       INTEGER(iwp)       ::  nynf             !< index of north boundary on former subdomain
-       INTEGER(iwp)       ::  nyn_on_file      !< index of norht boundary on former local domain  
-       INTEGER(iwp)       ::  nysc             !< index of south boundary on current subdomain 
-       INTEGER(iwp)       ::  nysf             !< index of south boundary on former subdomain
-       INTEGER(iwp)       ::  nys_on_file      !< index of south boundary on former local domain  
-
-       INTEGER(iwp), SAVE  ::  l           !< index variable for surface type
-
-       LOGICAL                         ::  surf_match_def     !< flag indicating that surface element is of default type
-       LOGICAL                         ::  surf_match_lsm     !< flag indicating that surface element is of natural type
-       LOGICAL                         ::  surf_match_usm     !< flag indicating that surface element is of urban type
+       INTEGER(iwp)       ::  ii          !< running index over input files
+       INTEGER(iwp)       ::  kk          !< running index over previous input files covering current local domain
+       INTEGER(iwp)       ::  nxlc        !< index of left boundary on current subdomain
+       INTEGER(iwp)       ::  nxlf        !< index of left boundary on former subdomain 
+       INTEGER(iwp)       ::  nxl_on_file !< index of left boundary on former local domain 
+       INTEGER(iwp)       ::  nxrc        !< index of right boundary on current subdomain
+       INTEGER(iwp)       ::  nxrf        !< index of right boundary on former subdomain
+       INTEGER(iwp)       ::  nxr_on_file !< index of right boundary on former local domain  
+       INTEGER(iwp)       ::  nync        !< index of north boundary on current subdomain
+       INTEGER(iwp)       ::  nynf        !< index of north boundary on former subdomain
+       INTEGER(iwp)       ::  nyn_on_file !< index of norht boundary on former local domain  
+       INTEGER(iwp)       ::  nysc        !< index of south boundary on current subdomain 
+       INTEGER(iwp)       ::  nysf        !< index of south boundary on former subdomain
+       INTEGER(iwp)       ::  nys_on_file !< index of south boundary on former local domain  
+
+       INTEGER(iwp), SAVE ::  l           !< index variable for surface type
+
+       LOGICAL            ::  surf_match_def     !< flag indicating that surface element is of default type
+       LOGICAL            ::  surf_match_lsm     !< flag indicating that surface element is of natural type
+       LOGICAL            ::  surf_match_usm     !< flag indicating that surface element is of urban type
 
        LOGICAL, INTENT(OUT) ::  found
@@ -3506,6 +3512,6 @@
        LOGICAL, SAVE ::  vertical_surface   !< flag indicating vertical surfaces
 
-       TYPE(surf_type), DIMENSION(0:2), SAVE ::  surf_h             !< horizontal surface type on file
-       TYPE(surf_type), DIMENSION(0:3), SAVE ::  surf_v             !< vertical surface type on file
+       TYPE(surf_type), DIMENSION(0:2), SAVE ::  surf_h !< horizontal surface type on file
+       TYPE(surf_type), DIMENSION(0:3), SAVE ::  surf_v !< vertical surface type on file
 
 
@@ -3513,5 +3519,7 @@
 
        SELECT CASE ( restart_string(1:length) )
-
+!
+!--       Read the number of horizontally orientated surface elements and 
+!--       allocate arrays
           CASE ( 'ns_h_on_file' )
              IF ( kk == 1 )  THEN
@@ -3524,5 +3532,5 @@
                 IF ( ALLOCATED( surf_h(2)%start_index ) )                      &
                    CALL deallocate_surface_attributes_h_top( surf_h(2) )       
-
+!
 !--             Allocate memory for number of surface elements on file. 
 !--             Please note, these number is not necessarily the same as 
@@ -3535,5 +3543,5 @@
                                         nys_on_file, nyn_on_file,              &
                                         nxl_on_file, nxr_on_file )
-
+!
 !--             Horizontal downward facing
                 surf_h(1)%ns = ns_h_on_file(1)
@@ -3541,5 +3549,5 @@
                                         nys_on_file, nyn_on_file,              &
                                         nxl_on_file, nxr_on_file )
-
+!
 !--             Model top
                 surf_h(2)%ns = ns_h_on_file(2)
@@ -3555,5 +3563,7 @@
 
              ENDIF   
-
+!
+!--       Read the number of vertically orientated surface elements and 
+!--       allocate arrays
           CASE ( 'ns_v_on_file' )
              IF ( kk == 1 ) THEN
@@ -3565,5 +3575,4 @@
                 ENDDO
 
-!--                Vertical surfaces
                 DO  l = 0, 3
                    surf_v(l)%ns = ns_v_on_file(l)
@@ -3574,5 +3583,6 @@
 
              ENDIF
-
+!
+!--       Read start and end indices of surface elements at each (ji)-gridpoint
           CASE ( 'surf_h(0)%start_index' )
              IF ( kk == 1 )                                                    &
@@ -3584,4 +3594,6 @@
              horizontal_surface = .TRUE.
              vertical_surface   = .FALSE.
+!
+!--       Read specific attributes
           CASE ( 'surf_h(0)%us' )         
              IF ( ALLOCATED( surf_h(0)%us )  .AND.  kk == 1 )                  &
@@ -4178,5 +4190,6 @@
        END SELECT
 !
-!--    Redistribute surface elements on its respective type. 
+!--    Redistribute surface elements on its respective type. Start with 
+!--    horizontally orientated surfaces.
        IF ( horizontal_surface  .AND.                                          &
             .NOT. INDEX( restart_string(1:length), '%start_index' ) /= 0 )     &
@@ -4187,16 +4200,26 @@
              jc = nysc
              DO  j = nysf, nynf
-
+!
+!--             Determine type of surface element, i.e. default, natural, 
+!--             urban, at current grid point. 
                 surf_match_def  = surf_def_h(l)%end_index(jc,ic) >=            &
                                   surf_def_h(l)%start_index(jc,ic)
                 surf_match_lsm  = ( surf_lsm_h%end_index(jc,ic)  >=            &
                                     surf_lsm_h%start_index(jc,ic) )            &
-                            .AND.  l == 0 
+                             .AND.  l == 0 
                 surf_match_usm  = ( surf_usm_h%end_index(jc,ic)  >=            &
                                     surf_usm_h%start_index(jc,ic) )            &
-                            .AND.  l == 0 
-                           
+                             .AND.  l == 0 
+!
+!--             Write restart data onto default-type surfaces if required.
                 IF ( surf_match_def )  THEN
+!
+!--                Set the start index for the local surface element 
                    mm = surf_def_h(l)%start_index(jc,ic)
+!
+!--                For index pair (j,i) on file loop from start to end index,
+!--                and in case the local surface element mm is smaller than
+!--                the local end index, assign the respective surface data
+!--                to this element. 
                    DO  m = surf_h(l)%start_index(j,i),                         &
                            surf_h(l)%end_index(j,i)
@@ -4207,5 +4230,9 @@
                    ENDDO
                 ENDIF
-
+!
+!--             Same for natural-type surfaces. Please note, it is implicitly
+!--             assumed that natural surface elements are below urban 
+!--             urban surface elements if there are several horizontal surfaces
+!--             at (j,i). An example would be bridges. 
                 IF ( surf_match_lsm )  THEN
                    mm = surf_lsm_h%start_index(jc,ic)
@@ -4218,5 +4245,6 @@
                    ENDDO
                 ENDIF
-
+!
+!--             Same for urban-type surfaces
                 IF ( surf_match_usm )  THEN
                    mm = surf_usm_h%start_index(jc,ic)
@@ -4241,5 +4269,7 @@
              jc = nysc
              DO  j = nysf, nynf
-
+!
+!--             Determine type of surface element, i.e. default, natural, 
+!--             urban, at current grid point. 
                 surf_match_def  = surf_def_v(l)%end_index(jc,ic) >=            &
                                   surf_def_v(l)%start_index(jc,ic)
@@ -4248,7 +4278,15 @@
                 surf_match_usm  = surf_usm_v(l)%end_index(jc,ic) >=            &
                                   surf_usm_v(l)%start_index(jc,ic)
-                                  
+!
+!--             Write restart data onto default-type surfaces if required.               
                 IF ( surf_match_def )  THEN
+!
+!--                Set the start index for the local surface element 
                    mm = surf_def_v(l)%start_index(jc,ic)
+!
+!--                For index pair (j,i) on file loop from start to end index,
+!--                and in case the local surface element mm is smaller than
+!--                the local end index, assign the respective surface data
+!--                to this element. 
                    DO  m = surf_v(l)%start_index(j,i),                         &
                            surf_v(l)%end_index(j,i)
@@ -4259,5 +4297,12 @@
                    ENDDO
                 ENDIF
-
+!
+!--             Same for natural-type surfaces. Please note, it is implicitly
+!--             assumed that natural surface elements are below urban 
+!--             urban surface elements if there are several vertical surfaces
+!--             at (j,i). An example a terrain elevations with a building on 
+!--             top. So far, initialization of urban surfaces below natural
+!--             surfaces on the same (j,i) is not possible, so that this case
+!--             cannot occur. 
                 IF ( surf_match_lsm )  THEN
                    mm = surf_lsm_v(l)%start_index(jc,ic)
@@ -4293,5 +4338,5 @@
 ! Description:
 ! ------------
-!> Restores surfacle elements back on its respective type.
+!> Restores surface elements back on its respective type.
 !------------------------------------------------------------------------------!
           SUBROUTINE restore_surface_elements( surf_target, m_target,          &
@@ -4511,5 +4556,5 @@
 ! Description:
 ! ------------
-!> Routine maps surface data read from file after restart - 1D arrays
+!> Routine maps surface data read from file after restart - 1D arrays. 
 !------------------------------------------------------------------------------!
     SUBROUTINE surface_restore_elements_1d( surf_target, surf_file,            &
