Index: /palm/trunk/SOURCE/Makefile
===================================================================
--- /palm/trunk/SOURCE/Makefile	(revision 1832)
+++ /palm/trunk/SOURCE/Makefile	(revision 1833)
@@ -20,5 +20,5 @@
 # Current revisions:
 # ------------------
-# 
+# spectrum renamed spectra_mod, depencies for spectra changed
 # 
 # Former revisions:
@@ -290,5 +290,5 @@
 	random_function.f90 random_gauss.f90 random_generator_parallel.f90 \
 	read_3d_binary.f90 read_var_list.f90 run_control.f90 \
-	set_slicer_attributes_dvrp.f90 singleton.f90 sor.f90 spectrum.f90 \
+	set_slicer_attributes_dvrp.f90 singleton.f90 sor.f90 spectra_mod.f90 \
 	subsidence.f90 sum_up_3d_data.f90 \
 	surface_coupler.f90 surface_layer_fluxes.f90 swap_timelevel.f90 temperton_fft.f90 \
@@ -361,5 +361,5 @@
 check_parameters.o: modules.o mod_kinds.o land_surface_model_mod.o \
 	netcdf_interface.o plant_canopy_model_mod.o pmc_interface.o radiation_model_mod.o \
-	subsidence.o
+	spectra_mod.o subsidence.o
 close_file.o: modules.o mod_kinds.o netcdf_interface.o
 compute_vpt.o: modules.o mod_kinds.o
@@ -376,5 +376,5 @@
    netcdf_interface.o mod_particle_attributes.o
 data_output_spectra.o: modules.o cpulog.o mod_kinds.o netcdf_interface.o \
-   spectrum.o
+   spectra_mod.o
 data_output_tseries.o: modules.o cpulog.o mod_kinds.o netcdf_interface.o
 data_output_2d.o: modules.o cpulog.o mod_kinds.o mod_particle_attributes.o\
@@ -399,6 +399,6 @@
 global_min_max.o: modules.o mod_kinds.o
 header.o: modules.o cpulog.o mod_kinds.o netcdf_interface.o land_surface_model_mod.o\
-plant_canopy_model_mod.o pmc_handle_communicator.o pmc_interface.o \
-   radiation_model_mod.o spectrum.o subsidence.o
+   plant_canopy_model.o pmc_handle_communicator.o pmc_interface.o \
+   radiation_model_mod.o spectra_mod.o subsidence.o 
 inflow_turbulence.o: modules.o cpulog.o mod_kinds.o
 init_1d_model.o: modules.o mod_kinds.o
@@ -414,5 +414,5 @@
 init_masks.o: modules.o mod_kinds.o netcdf_interface.o
 init_ocean.o: modules.o eqn_state_seawater.o mod_kinds.o
-init_pegrid.o: modules.o mod_kinds.o
+init_pegrid.o: modules.o mod_kinds.o spectra_mod.o
 init_pt_anomaly.o: modules.o mod_kinds.o
 init_rankine.o: modules.o mod_kinds.o
@@ -456,12 +456,22 @@
 mod_particle_attributes.o: mod_particle_attributes.f90 mod_kinds.o
 netcdf_interface.o: netcdf_interface.f90 modules.o mod_kinds.o \
-   land_surface_model_mod.o spectrum.o
+   land_surface_model_mod.o spectra_mod.o
 nudging.o: modules.o cpulog.o mod_kinds.o
+<<<<<<< .mine
+package_parin.o: modules.o mod_kinds.o \
+   plant_canopy_model.o radiation_model.o
+=======
 package_parin.o: modules.o mod_kinds.o spectrum.o
+>>>>>>> .r1832
 palm.o: modules.o cpulog.o ls_forcing.o mod_kinds.o nudging.o\
 	pmc_interface.o surface_layer_fluxes.o
 parin.o: modules.o cpulog.o land_surface_model_mod.o mod_kinds.o netcdf_interface.o \
+<<<<<<< .mine
+   pmc_interface.o progress_bar.o spectra_mod.o
+plant_canopy_model.o: modules.o mod_kinds.o
+=======
    plant_canopy_model_mod.o pmc_interface.o progress_bar.o radiation_model_mod.o 
 plant_canopy_model_mod.o: modules.o mod_kinds.o
+>>>>>>> .r1832
 pmc_interface.o: modules.o mod_kinds.o pmc_client.o pmc_general.o\
 	pmc_handle_communicator.o pmc_mpi_wrapper.o pmc_server.o
@@ -490,13 +500,21 @@
 random_gauss.o: mod_kinds.o random_function.o random_generator_parallel.o
 random_generator_parallel.o: mod_kinds.o
+<<<<<<< .mine
+read_3d_binary.o: modules.o cpulog.o mod_kinds.o land_surface_model_mod.o \
+   radiation_model.o random_function.o random_generator_parallel.o \
+   spectra_mod.o
+read_var_list.o: modules.o mod_kinds.o netcdf_interface.o plant_canopy_model.o \
+   spectra_mod.o
+=======
 read_3d_binary.o: modules.o cpulog.o mod_kinds.o land_surface_model_mod.o\
                   radiation_model_mod.o random_function.o\
                   random_generator_parallel.o
 read_var_list.o: modules.o mod_kinds.o netcdf_interface.o plant_canopy_model_mod.o
+>>>>>>> .r1832
 run_control.o: modules.o cpulog.o mod_kinds.o
 set_slicer_attributes_dvrp.o: modules.o mod_kinds.o
 singleton.o: mod_kinds.o singleton.f90
 sor.o: modules.o mod_kinds.o
-spectrum.o: spectrum.f90 modules.o mod_kinds.o cpulog.o fft_xy.o
+spectra_mod.o: spectra_mod.f90 modules.o mod_kinds.o cpulog.o fft_xy.o
 subsidence.o: modules.o mod_kinds.o
 sum_up_3d_data.o: modules.o cpulog.o mod_kinds.o land_surface_model_mod.o\
@@ -511,5 +529,5 @@
 	ls_forcing.o mod_kinds.o nudging.o pmc_interface.o production_e.o \
 	prognostic_equations.o progress_bar.o radiation_model_mod.o \
-	spectrum.o user_actions.o surface_layer_fluxes.o
+	spectra_mod.o user_actions.o surface_layer_fluxes.o
 time_to_string.o: mod_kinds.o
 timestep.o: modules.o cpulog.o mod_kinds.o
@@ -544,9 +562,15 @@
 user_parin.o: modules.o mod_kinds.o user_module.o
 user_read_restart_data.o: modules.o mod_kinds.o user_module.o
-user_spectra.o: modules.o mod_kinds.o spectrum.o user_module.o
+user_spectra.o: modules.o mod_kinds.o spectra_mod.o user_module.o
 user_statistics.o: modules.o mod_kinds.o netcdf_interface.o user_module.o
 wall_fluxes.o: modules.o mod_kinds.o
+<<<<<<< .mine
+write_3d_binary.o: modules.o cpulog.o mod_kinds.o land_surface_model_mod.o \
+	radiation_model.o random_function.o random_generator_parallel.o \
+	spectra_mod.o
+write_var_list.o: modules.o mod_kinds.o netcdf_interface.o plant_canopy_model.o\
+   spectra_mod.o =======
 write_3d_binary.o: modules.o cpulog.o mod_kinds.o land_surface_model_mod.o\
                    radiation_model_mod.o random_function.o\
                    random_generator_parallel.o
-write_var_list.o: modules.o mod_kinds.o netcdf_interface.o plant_canopy_model_mod.o 
+write_var_list.o: modules.o mod_kinds.o netcdf_interface.o plant_canopy_model_mod.o >>>>>>> .r1832
Index: /palm/trunk/SOURCE/check_parameters.f90
===================================================================
--- /palm/trunk/SOURCE/check_parameters.f90	(revision 1832)
+++ /palm/trunk/SOURCE/check_parameters.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! -----------------
-! 
+! check of spectra quantities moved to spectra_mod
 ! 
 ! Former revisions:
@@ -377,4 +377,6 @@
         ONLY: radiation, radiation_check_data_output,                          &
               radiation_check_data_output_pr, radiation_check_parameters
+    USE spectra_mod,                                                           &
+        ONLY:  calculate_spectra, spectra_check_parameters
     USE statistics
     USE subsidence_mod
@@ -1042,5 +1044,8 @@
     
 !
-!-- When land surface model is used, peform addtional checks
+!-- Additional checks for spectra
+    IF ( calculate_spectra )  CALL spectra_check_parameters
+!
+!-- When land surface model is used, perform additional checks
     IF ( land_surface )  CALL lsm_check_parameters
 
@@ -2009,5 +2014,5 @@
 !
 !-- Set the default intervals for data output, if necessary
-!-- NOTE: dt_dosp has already been set in package_parin
+!-- NOTE: dt_dosp has already been set in spectra_parin
     IF ( dt_data_output /= 9999999.9_wp )  THEN
        IF ( dt_dopr           == 9999999.9_wp )  dt_dopr           = dt_data_output
@@ -2027,6 +2032,4 @@
     IF ( skip_time_dopr    == 9999999.9_wp )                                   &
                                        skip_time_dopr    = skip_time_data_output
-    IF ( skip_time_dosp    == 9999999.9_wp )                                   &
-                                       skip_time_dosp    = skip_time_data_output
     IF ( skip_time_do2d_xy == 9999999.9_wp )                                   &
                                        skip_time_do2d_xy = skip_time_data_output
@@ -2045,6 +2048,5 @@
 
 !
-!-- Check the average intervals (first for 3d-data, then for profiles and
-!-- spectra)
+!-- Check the average intervals (first for 3d-data, then for profiles)
     IF ( averaging_interval > dt_data_output_av )  THEN
        WRITE( message_string, * )  'averaging_interval = ',                    &
@@ -2061,14 +2063,4 @@
              averaging_interval_pr, ' must be <= dt_dopr = ', dt_dopr
        CALL message( 'check_parameters', 'PA0086', 1, 2, 0, 6, 0 )
-    ENDIF
-
-    IF ( averaging_interval_sp == 9999999.9_wp )  THEN
-       averaging_interval_sp = averaging_interval
-    ENDIF
-
-    IF ( averaging_interval_sp > dt_dosp )  THEN
-       WRITE( message_string, * )  'averaging_interval_sp = ',                 &
-             averaging_interval_sp, ' must be <= dt_dosp = ', dt_dosp
-       CALL message( 'check_parameters', 'PA0087', 1, 2, 0, 6, 0 )
     ENDIF
 
Index: /palm/trunk/SOURCE/data_output_spectra.f90
===================================================================
--- /palm/trunk/SOURCE/data_output_spectra.f90	(revision 1832)
+++ /palm/trunk/SOURCE/data_output_spectra.f90	(revision 1833)
@@ -19,5 +19,6 @@
 ! Current revisions:
 ! ------------------
-!
+! spectrum renamed spectra_mod, spectra related variables moved to spectra_mod,
+! routines data_output_spectra_x/y removed
 !
 ! Former revisions:
@@ -78,6 +79,5 @@
 #if defined( __netcdf )
     USE control_parameters,                                                    &
-        ONLY:  average_count_sp, averaging_interval_sp, dosp_time_count,       &
-               message_string, run_description_header,                         &
+        ONLY:  message_string, run_description_header,                         &
                time_since_reference_point
 
@@ -94,9 +94,8 @@
     USE pegrid
 
-    USE spectrum,                                                              &
-        ONLY:  comp_spectra_level, data_output_sp, spectra_direction
-
-    USE statistics,                                                            &
-        ONLY:  spectrum_x, spectrum_y
+    USE spectra_mod,                                                           &
+        ONLY:  average_count_sp, averaging_interval_sp, comp_spectra_level,    &
+               data_output_sp, dosp_time_count, spectra_direction, spectrum_x, &
+               spectrum_y
 
 
@@ -226,7 +225,4 @@
         ONLY:  pi
 
-    USE control_parameters,                                                    &
-        ONLY:  dosp_time_count
-
     USE grid_variables,                                                        &
         ONLY:  dx, dy
@@ -243,9 +239,7 @@
                netcdf_handle_error
 
-    USE spectrum,                                                              &
-        ONLY:  n_sp_x, n_sp_y
-
-    USE statistics,                                                            &
-        ONLY:  spectrum_x, spectrum_y
+    USE spectra_mod,                                                           &
+        ONLY:  dosp_time_count, n_sp_x, n_sp_y, spectrum_x, spectrum_y
+
 
     IMPLICIT NONE
@@ -302,418 +296,2 @@
 #endif
  END SUBROUTINE output_spectra_netcdf
-
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> @todo Missing subroutine description.
-!------------------------------------------------------------------------------!
- SUBROUTINE data_output_spectra_x( m, cranz, pr, frame_written )
-
-    USE arrays_3d,                                                             &
-        ONLY:  zu, zw
-    USE constants,                                                             &
-        ONLY:  pi
-
-    USE control_parameters,                                                    &
-        ONLY:  averaging_interval_sp, run_description_header, simulated_time_chr
-
-    USE grid_variables,                                                        &
-        ONLY:  dx
-
-    USE indices,                                                               &
-        ONLY:  nx
-
-    USE kinds
-
-    USE pegrid
-
-    USE statistics,                                                            &
-        ONLY:  spectrum_x
-
-    USE spectrum,                                                              &
-        ONLY:  comp_spectra_level, header_char, lstyles, klist_x, n_sp_x,      &
-               plot_spectra_level, utext_char, ytext_char
-
-    IMPLICIT NONE
-
-    CHARACTER (LEN=30) ::  atext !<
-    
-    INTEGER(iwp)       ::  i     !<
-    INTEGER(iwp)       ::  j     !<
-    INTEGER(iwp)       ::  k     !<
-    INTEGER(iwp)       ::  m     !<
-    INTEGER(iwp)       ::  pr    !<
-    
-    LOGICAL            ::  frame_written   !<
-    
-    REAL(wp)           ::  frequency = 0.0_wp !<
-!
-!-- Variables needed for PROFIL-namelist
-    CHARACTER (LEN=80) ::  rtext                !<
-    CHARACTER (LEN=80) ::  utext                !<
-    CHARACTER (LEN=80) ::  xtext = 'k in m>->1' !<
-    CHARACTER (LEN=80) ::  ytext                !<
-
-    INTEGER(iwp)       ::  cranz       !<
-    INTEGER(iwp)       ::  labforx = 3 !<
-    INTEGER(iwp)       ::  labfory = 3 !<
-    INTEGER(iwp)       ::  legpos  = 3 !<
-    INTEGER(iwp)       ::  timodex = 1 !<
-    
-    INTEGER(iwp), DIMENSION(1:100) ::  cucol  = 1      !<
-    INTEGER(iwp), DIMENSION(1:100) ::  klist  = 999999 !<
-    INTEGER(iwp), DIMENSION(1:100) ::  lstyle = 0      !<
-    
-    LOGICAL ::  datleg = .TRUE. !<
-    LOGICAL ::  grid = .TRUE.   !<
-    LOGICAL ::  lclose = .TRUE. !<
-    LOGICAL ::  rand = .TRUE.   !<
-    LOGICAL ::  swap = .TRUE.   !<
-    LOGICAL ::  twoxa = .TRUE.  !<
-    LOGICAL ::  xlog = .TRUE.   !<
-    LOGICAL ::  ylog = .TRUE.   !<
-    
-    REAL(wp) ::  gwid = 0.1_wp    !<
-    REAL(wp) ::  rlegfak = 0.7_wp !<
-    REAL(wp) ::  uxmin            !<
-    REAL(wp) ::  uxmax            !<
-    REAL(wp) ::  uymin            !<
-    REAL(wp) ::  uymax            !<
-     
-    REAL(wp), DIMENSION(1:100) ::  lwid = 0.6_wp !<
-    REAL(wp), DIMENSION(100)   ::  uyma          !<
-    REAL(wp), DIMENSION(100)   ::  uymi          !<
-
-    NAMELIST /RAHMEN/  cranz, datleg, rtext, swap
-    NAMELIST /CROSS/   rand, cucol, grid, gwid, klist, labforx, labfory,      &
-                       legpos, lclose, lstyle, lwid, rlegfak, timodex, utext, &
-                       uxmin, uxmax, uymin, uymax, twoxa, xlog, xtext, ylog,  &
-                       ytext
-
-
-    rtext = '\0.5 ' // run_description_header
-
-!
-!-- Open parameter- and data-file
-    CALL check_open( 81 )
-    CALL check_open( 82 )
-
-!
-!-- Write file header,
-!-- write RAHMEN-parameters (pr=3: w-array is on zw, other arrys on zu,
-!-- pr serves as an index for output of strings (axis-labels) of the
-!-- different quantities u, v, w, pt and q)
-    DO  k = 1, n_sp_x
-       IF ( k < 100 )  THEN
-          IF ( pr == 3 )  THEN
-             WRITE ( 82, 100 )  '#', k, header_char( pr ),        &
-                                INT( zw(comp_spectra_level(k)) ), &
-                                simulated_time_chr
-          ELSE
-             WRITE ( 82, 100 )  '#', k, header_char( pr ),        &
-                                INT( zu(comp_spectra_level(k)) ), &
-                                simulated_time_chr
-          ENDIF
-       ELSE
-          IF ( pr == 3 )  THEN
-             WRITE ( 82, 101 )  '#', k, header_char( pr ),        &
-                                INT( zw(comp_spectra_level(k)) ), &
-                                simulated_time_chr
-          ELSE
-             WRITE ( 82, 101 )  '#', k, header_char( pr ),        &
-                                INT( zu(comp_spectra_level(k)) ), &
-                                simulated_time_chr
-          ENDIF
-       ENDIF
-    ENDDO
-
-    IF ( .NOT. frame_written )  THEN
-       WRITE ( 81, RAHMEN )
-       frame_written = .TRUE.
-    ENDIF
-
-!
-!-- Write all data and calculate uymi and uyma. They serve to calculate
-!-- the CROSS-parameters uymin and uymax
-    uymi = 999.999_wp; uyma = -999.999_wp
-    DO  i = 1, nx/2
-       frequency = 2.0_wp * pi * i / ( dx * ( nx + 1 ) )
-       WRITE ( 82, 102 )  frequency, ( frequency * spectrum_x(i,k,m), k = 1, &
-                          n_sp_x )
-       DO  k = 1, n_sp_x
-          uymi(k) = MIN( uymi(k), frequency * spectrum_x(i,k,m) )
-          uyma(k) = MAX( uyma(k), frequency * spectrum_x(i,k,m) )
-       ENDDO
-    ENDDO
-
-!
-!-- Determine CROSS-parameters
-    cucol(1:n_sp_x)  = (/ ( k, k = 1, n_sp_x ) /)
-    lstyle(1:n_sp_x) = (/ ( lstyles(k), k = 1, n_sp_x ) /)
-
-!
-!-- Calculate klist-values from the available comp_spectra_level values
-    i = 1; k = 1
-    DO WHILE ( i <= 100  .AND.  plot_spectra_level(i) /= 999999 )
-       DO WHILE ( k <= n_sp_x  .AND. &
-                  plot_spectra_level(i) >= comp_spectra_level(k) )
-          IF ( plot_spectra_level(i) == comp_spectra_level(k) )  THEN
-             klist(i) = k + klist_x
-          ELSE
-             uymi(k) =  999.999_wp
-             uyma(k) = -999.999_wp
-          ENDIF
-          k = k + 1
-       ENDDO
-       i = i + 1
-    ENDDO
-    uymi(k:n_sp_x) =  999.999_wp
-    uyma(k:n_sp_x) = -999.999_wp
-    utext = 'x'//utext_char( pr )
-    IF ( averaging_interval_sp /= 0.0_wp ) THEN
-       WRITE ( atext, 104 )  averaging_interval_sp
-       utext = TRIM(utext) // ',  ' // TRIM( atext )
-    ENDIF
-    uxmin = 0.8_wp * 2.0_wp * pi        / ( dx * ( nx + 1 ) )
-    uxmax = 1.2_wp * 2.0_wp * pi * nx/2 / ( dx * ( nx + 1 ) )
-    uymin = 0.8_wp * MIN (  999.999_wp, MINVAL ( uymi ) )
-    uymax = 1.2_wp * MAX ( -999.999_wp, MAXVAL ( uyma ) )
-    ytext = ytext_char( pr )
-
-!
-!-- Output of CROSS-parameters
-    WRITE ( 81, CROSS )
-
-!
-!-- Increase counter by the number of profiles written in the actual block
-    klist_x = klist_x + n_sp_x
-
-!
-!-- Write end-mark
-    WRITE ( 82, 103 )
-
-!
-!-- Close parameter- and data-file
-    CALL close_file( 81 )
-    CALL close_file( 82 )
-
-!
-!-- Formats
-100 FORMAT (A,I1,1X,A,1X,I4,'m ',A)
-101 FORMAT (A,I2,1X,A,1X,I4,'m ',A)
-102 FORMAT (E15.7,100(1X,E15.7))
-103 FORMAT ('NEXT')
-104 FORMAT ('time averaged over',F7.1,' s')
-
- END SUBROUTINE data_output_spectra_x
-
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> @todo Missing subroutine description.
-!------------------------------------------------------------------------------!
- SUBROUTINE data_output_spectra_y( m, cranz, pr, frame_written )
-
-    USE arrays_3d,                                                             &
-        ONLY:  zu, zw
- 
-    USE constants,                                                             &
-        ONLY:  pi
-
-    USE control_parameters,                                                    &
-        ONLY:  averaging_interval_sp, run_description_header, simulated_time_chr
-
-    USE grid_variables,                                                        &
-        ONLY:  dy
-
-    USE indices,                                                               &
-        ONLY:  ny
-
-    USE kinds
-
-    USE pegrid
-    
-    USE statistics,                                                             &
-        ONLY:  spectrum_y
-
-    USE spectrum,                                                              &
-        ONLY:  comp_spectra_level, header_char, klist_y, lstyles, n_sp_y,      &
-               plot_spectra_level, utext_char, ytext_char
-
-    IMPLICIT NONE
-
-    
-    CHARACTER (LEN=30) ::  atext !<
-    
-    INTEGER(iwp)       ::  i     !<
-    INTEGER(iwp)       ::  j     !<
-    INTEGER(iwp)       ::  k     !<
-    INTEGER(iwp)       ::  m     !<
-    INTEGER(iwp)       ::  pr    !<
-    
-    LOGICAL            :: frame_written   !<
-    
-    REAL(wp)           :: frequency = 0.0_wp !<
-
-!
-!-- Variables needed for PROFIL-namelist
-    CHARACTER (LEN=80) ::  rtext                !<
-    CHARACTER (LEN=80) ::  utext                !<
-    CHARACTER (LEN=80) ::  xtext = 'k in m>->1' !<
-    CHARACTER (LEN=80) ::  ytext                !<
-
-    INTEGER(iwp) ::  cranz       !<
-    INTEGER(iwp) ::  labforx = 3 !<
-    INTEGER(iwp) ::  labfory = 3 !<
-    INTEGER(iwp) ::  legpos  = 3 !<
-    INTEGER(iwp) ::  timodex = 1 !<
-    
-    INTEGER(iwp), DIMENSION(1:100) ::  cucol  = 1      !<
-    INTEGER(iwp), DIMENSION(1:100) ::  klist  = 999999 !<
-    INTEGER(iwp), DIMENSION(1:100) ::  lstyle = 0      !<
-    
-    LOGICAL ::  datleg = .TRUE. !<
-    LOGICAL ::  grid = .TRUE.   !<
-    LOGICAL ::  lclose = .TRUE. !<
-    LOGICAL ::  rand = .TRUE.   !<
-    LOGICAL ::  swap = .TRUE.   !<
-    LOGICAL ::  twoxa = .TRUE.  !<
-    LOGICAL ::  xlog = .TRUE.   !<
-    LOGICAL ::  ylog = .TRUE.   !<
-    
-    REAL(wp) ::  gwid = 0.1_wp     !<
-    REAL(wp) ::  rlegfak = 0.7_wp  !<
-    REAL(wp) ::  uxmin             !<
-    REAL(wp) ::  uxmax             !<
-    REAL(wp) ::  uymin             !<
-    REAL(wp) ::  uymax             !<
-    
-    REAL(wp), DIMENSION(1:100) ::  lwid = 0.6_wp !<
-    
-    REAL(wp), DIMENSION(100)   ::  uyma          !<
-    REAL(wp), DIMENSION(100)   ::  uymi          !<
-
-    NAMELIST /RAHMEN/  cranz, datleg, rtext, swap
-    NAMELIST /CROSS/   rand, cucol, grid, gwid, klist, labforx, labfory,      &
-                       legpos, lclose, lstyle, lwid, rlegfak, timodex, utext, &
-                       uxmin, uxmax, uymin, uymax, twoxa, xlog, xtext, ylog,  &
-                       ytext
-
-
-    rtext = '\0.5 ' // run_description_header
-
-!
-!-- Open parameter- and data-file
-    CALL check_open( 83 )
-    CALL check_open( 84 )
-
-!
-!-- Write file header,
-!-- write RAHMEN-parameters (pr=3: w-array is on zw, other arrys on zu,
-!-- pr serves as an index for output of strings (axis-labels) of the
-!-- different quantities u, v, w, pt and q)
-    DO  k = 1, n_sp_y
-       IF ( k < 100 )  THEN
-          IF ( pr == 3 ) THEN
-             WRITE ( 84, 100 )  '#', k, header_char( pr ),        &
-                                INT( zw(comp_spectra_level(k)) ), &
-                                simulated_time_chr
-          ELSE
-             WRITE ( 84, 100 )  '#', k, header_char( pr ),        &
-                                INT( zu(comp_spectra_level(k)) ), &
-                                simulated_time_chr
-          ENDIF
-       ELSE
-          IF ( pr == 3 )  THEN
-             WRITE ( 84, 101 )  '#', k, header_char( pr ),        &
-                                INT( zw(comp_spectra_level(k)) ), &
-                                simulated_time_chr
-          ELSE
-             WRITE ( 84, 101 )  '#', k, header_char( pr ),        &
-                                INT( zu(comp_spectra_level(k)) ), &
-                                simulated_time_chr
-          ENDIF
-       ENDIF
-    ENDDO
-
-    IF ( .NOT. frame_written )  THEN
-       WRITE ( 83, RAHMEN )
-       frame_written = .TRUE.
-    ENDIF
-
-!
-!-- Write all data and calculate uymi and uyma. They serve to calculate
-!-- the CROSS-parameters uymin and uymax
-    uymi = 999.999_wp; uyma = -999.999_wp
-    DO  j = 1, ny/2
-       frequency = 2.0_wp * pi * j / ( dy * ( ny + 1 ) )
-       WRITE ( 84, 102 ) frequency, ( frequency * spectrum_y(j,k,m), &
-                                      k = 1, n_sp_y ) 
-       DO k = 1, n_sp_y
-          uymi(k) = MIN( uymi(k), frequency * spectrum_y(j,k,m) )
-          uyma(k) = MAX( uyma(k), frequency * spectrum_y(j,k,m) )
-       ENDDO
-    ENDDO
-
-!
-!-- Determine CROSS-parameters
-    cucol(1:n_sp_y)  = (/ ( k, k = 1, n_sp_y ) /)
-    lstyle(1:n_sp_y) = (/ ( lstyles(k), k = 1, n_sp_y ) /)
-
-!
-!-- Calculate klist-values from the available comp_spectra_level values
-    j = 1; k = 1
-    DO WHILE ( j <= 100  .AND.  plot_spectra_level(j) /= 999999 )
-       DO WHILE ( k <= n_sp_y  .AND. &
-                  plot_spectra_level(j) >= comp_spectra_level(k) )
-          IF ( plot_spectra_level(j) == comp_spectra_level(k) )  THEN
-             klist(j) = k + klist_y
-          ELSE
-             uymi(k) =  999.999_wp
-             uyma(k) = -999.999_wp
-          ENDIF
-          k = k + 1
-       ENDDO
-       j = j + 1
-    ENDDO
-    uymi(k:n_sp_y) =  999.999_wp
-    uyma(k:n_sp_y) = -999.999_wp
-    utext = 'y'//utext_char( pr )
-    IF ( averaging_interval_sp /= 0.0_wp )  THEN
-       WRITE ( atext, 104 )  averaging_interval_sp
-       utext = TRIM(utext) // ',  ' // TRIM( atext )
-    ENDIF
-    uxmin = 0.8_wp * 2.0_wp * pi        / ( dy * ( ny + 1 ) )
-    uxmax = 1.2_wp * 2.0_wp * pi * ny/2 / ( dy * ( ny + 1 ) )
-    uymin = 0.8_wp * MIN (  999.999_wp, MINVAL ( uymi ) )
-    uymax = 1.2_wp * MAX ( -999.999_wp, MAXVAL ( uyma ) )
-    ytext = ytext_char( pr )
-
-!
-!-- Output CROSS-parameters
-    WRITE ( 83, CROSS )
-
-!
-!-- Increase counter by the number of profiles written in the actual block
-    klist_y = klist_y + n_sp_y
-
-!
-!-- Write end-mark
-    WRITE ( 84, 103 )  
-
-!
-!-- Close parameter- and data-file
-    CALL close_file( 83 )
-    CALL close_file( 84 )
-
-!
-!-- Formats
-100 FORMAT (A,I1,1X,A,1X,I4,'m ',A)
-101 FORMAT (A,I2,1X,A,1X,I4,'m ',A)
-102 FORMAT (E15.7,100(1X,E15.7))
-103 FORMAT ('NEXT')
-104 FORMAT ('time averaged over',F7.1,' s')
-
- END SUBROUTINE data_output_spectra_y
Index: /palm/trunk/SOURCE/header.f90
===================================================================
--- /palm/trunk/SOURCE/header.f90	(revision 1832)
+++ /palm/trunk/SOURCE/header.f90	(revision 1833)
@@ -19,5 +19,6 @@
 ! Current revisions:
 ! -----------------
-!
+! spectrum renamed spectra_mod, output of spectra related quantities moved to
+! spectra_mod
 !
 ! Former revisions:
@@ -310,7 +311,6 @@
         ONLY:  radiation, radiation_header
     
-    USE spectrum,                                                              &
-        ONLY:  comp_spectra_level, data_output_sp, plot_spectra_level,         &
-               spectra_direction
+    USE spectra_mod,                                                           &
+        ONLY:  calculate_spectra, spectra_header
 
     IMPLICIT NONE
@@ -1631,22 +1631,6 @@
 
 !
-!-- Spectra output
-    IF ( dt_dosp /= 9999999.9_wp )  THEN
-       WRITE ( io, 370 )
-
-       output_format = netcdf_data_format_string
-       IF ( netcdf_deflate == 0 )  THEN
-          WRITE ( io, 344 )  output_format
-       ELSE
-          WRITE ( io, 354 )  TRIM( output_format ), netcdf_deflate
-       ENDIF
-       WRITE ( io, 371 )  dt_dosp
-       IF ( skip_time_dosp /= 0.0_wp )  WRITE ( io, 339 )  skip_time_dosp
-       WRITE ( io, 372 )  ( data_output_sp(i), i = 1,10 ),     &
-                          ( spectra_direction(i), i = 1,10 ),  &
-                          ( comp_spectra_level(i), i = 1,100 ), &
-                          ( plot_spectra_level(i), i = 1,100 ), &
-                          averaging_interval_sp, dt_averaging_input_pr
-    ENDIF
+!-- Output of spectra related quantities
+    IF ( calculate_spectra )  CALL spectra_header( io )
 
     WRITE ( io, 99 )
@@ -2089,22 +2073,4 @@
 367 FORMAT ('       Polygon reduction for topography: cluster_size = ', I1)
 #endif
-370 FORMAT ('    Spectra:')
-371 FORMAT ('       Output every ',F7.1,' s'/)
-372 FORMAT ('       Arrays:     ', 10(A5,',')/                         &
-            '       Directions: ', 10(A5,',')/                         &
-            '       height levels  k = ', 20(I3,',')/                  &
-            '                          ', 20(I3,',')/                  &
-            '                          ', 20(I3,',')/                  &
-            '                          ', 20(I3,',')/                  &
-            '                          ', 19(I3,','),I3,'.'/           &
-            '       height levels selected for standard plot:'/        &
-            '                      k = ', 20(I3,',')/                  &
-            '                          ', 20(I3,',')/                  &
-            '                          ', 20(I3,',')/                  &
-            '                          ', 20(I3,',')/                  &
-            '                          ', 19(I3,','),I3,'.'/           &
-            '       Time averaged over ', F7.1, ' s,' /                &
-            '       Profiles for the time averaging are taken every ', &
-                    F6.1,' s')
 400 FORMAT (//' Physical quantities:'/ &
               ' -------------------'/)
Index: /palm/trunk/SOURCE/init_3d_model.f90
===================================================================
--- /palm/trunk/SOURCE/init_3d_model.f90	(revision 1832)
+++ /palm/trunk/SOURCE/init_3d_model.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! ------------------
-!
+! initialization of spectra quantities moved to spectra_mod
 !
 ! Former revisions:
@@ -309,7 +309,7 @@
     USE statistics,                                                            &
         ONLY:  hom, hom_sum, mean_surface_level_height, pr_palm, rmask,        &
-               spectrum_x, spectrum_y, statistic_regions, sums, sums_divnew_l, &
-               sums_divold_l, sums_l, sums_l_l, sums_up_fraction_l,            &
-               sums_wsts_bc_l, ts_value, var_d, weight_pres, weight_substep
+               statistic_regions, sums, sums_divnew_l, sums_divold_l, sums_l,  &
+               sums_l_l, sums_up_fraction_l, sums_wsts_bc_l, ts_value,         &
+               weight_pres, weight_substep
  
     USE surface_layer_fluxes_mod,                                              &
@@ -578,14 +578,4 @@
     ENDIF
 
-    IF ( dt_dosp /= 9999999.9_wp )  THEN
-       ALLOCATE( spectrum_x( 1:nx/2, 1:10, 1:10 ),                             &
-                 spectrum_y( 1:ny/2, 1:10, 1:10 ) )
-       spectrum_x = 0.0_wp
-       spectrum_y = 0.0_wp
-
-       ALLOCATE( var_d(nzb:nzt+1) )
-       var_d = 0.0_wp
-    ENDIF
-
 !
 !-- 1D-array for large scale subsidence velocity
Index: /palm/trunk/SOURCE/init_pegrid.f90
===================================================================
--- /palm/trunk/SOURCE/init_pegrid.f90	(revision 1832)
+++ /palm/trunk/SOURCE/init_pegrid.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! ------------------
-! 
+! spectra related variables moved to spectra_mod
 !
 ! Former revisions:
@@ -141,5 +141,5 @@
     USE control_parameters,                                                    &
         ONLY:  bc_lr, bc_ns, coupling_mode, coupling_mode_remote,              &
-               coupling_topology, dt_dosp, gathered_size, grid_level,          &
+               coupling_topology, gathered_size, grid_level,                   &
                grid_level_count, host, inflow_l, inflow_n, inflow_r, inflow_s, &
                io_blocks, io_group, maximum_grid_level,                        &
@@ -163,5 +163,8 @@
       
     USE pegrid
- 
+
+    USE spectra_mod,                                                           &
+        ONLY:  dt_dosp
+
     USE transpose_indices,                                                     &
         ONLY:  nxl_y, nxl_yd, nxl_z, nxr_y, nxr_yd, nxr_z, nyn_x, nyn_z, nys_x,&
Index: /palm/trunk/SOURCE/modules.f90
===================================================================
--- /palm/trunk/SOURCE/modules.f90	(revision 1832)
+++ /palm/trunk/SOURCE/modules.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! ------------------
-!
+! spectra parameter moved to spectra module
 !
 ! Former revisions:
@@ -663,5 +663,5 @@
     CHARACTER (LEN=20), DIMENSION(0:1,100) ::  do2d = ' ', do3d = ' '
 
-    INTEGER(iwp) ::  abort_mode = 1, average_count_pr = 0, average_count_sp = 0, &
+    INTEGER(iwp) ::  abort_mode = 1, average_count_pr = 0, &
                      average_count_3d = 0, current_timestep_number = 0, &
                      coupling_topology = 0, &
@@ -669,5 +669,5 @@
                      disturbance_level_ind_t, doav_n = 0, dopr_n = 0, &
                      dopr_time_count = 0, dopts_time_count = 0, &
-                     dosp_time_count = 0, dots_time_count = 0, &
+                     dots_time_count = 0, &
                      do2d_xy_n = 0, do2d_xz_n = 0, do2d_yz_n = 0, do3d_avs_n = 0, &
                      dp_level_ind_b = 0, dvrp_filecount = 0, &
@@ -781,6 +781,5 @@
                  alpha_surface = 0.0_wp, atmos_ocean_sign = 1.0_wp, &
                  averaging_interval = 0.0_wp, averaging_interval_pr = 9999999.9_wp, &
-                 averaging_interval_sp = 9999999.9_wp, bc_pt_t_val, bc_q_t_val, &
-                 bottom_salinityflux = 0.0_wp, &
+                 bc_pt_t_val, bc_q_t_val, bottom_salinityflux = 0.0_wp, &
                  building_height = 50.0_wp, building_length_x = 50.0_wp, &
                  building_length_y = 50.0_wp, building_wall_left = 9999999.9_wp, &
@@ -799,5 +798,5 @@
                  dt_data_output_av = 9999999.9_wp, dt_disturb = 9999999.9_wp, &
                  dt_dopr = 9999999.9_wp, dt_dopr_listing = 9999999.9_wp, &
-                 dt_dopts = 9999999.9_wp, dt_dosp = 9999999.9_wp, dt_dots = 9999999.9_wp, &
+                 dt_dopts = 9999999.9_wp, dt_dots = 9999999.9_wp, &
                  dt_do2d_xy = 9999999.9_wp, dt_do2d_xz = 9999999.9_wp, &
                  dt_do2d_yz = 9999999.9_wp, dt_do3d = 9999999.9_wp, dt_dvrp = 9999999.9_wp, &
@@ -828,5 +827,5 @@
                  simulated_time = 0.0_wp, simulated_time_at_begin, sin_alpha_surface, &
                  skip_time_data_output = 0.0_wp, skip_time_data_output_av = 9999999.9_wp,&
-                 skip_time_dopr = 9999999.9_wp, skip_time_dosp = 9999999.9_wp, &
+                 skip_time_dopr = 9999999.9_wp, &
                  skip_time_do2d_xy = 9999999.9_wp, skip_time_do2d_xz = 9999999.9_wp, &
                  skip_time_do2d_yz = 9999999.9_wp, skip_time_do3d = 9999999.9_wp, &
@@ -1291,6 +1290,5 @@
     REAL(wp), DIMENSION(:), ALLOCATABLE       ::  mean_surface_level_height,      &
                                                   sums_divnew_l, sums_divold_l,   &
-                                                  var_d, weight_substep,          &
-                                                  weight_pres
+                                                  weight_substep, weight_pres
     REAL(wp), DIMENSION(:,:), ALLOCATABLE     ::  sums, sums_wsts_bc_l, ts_value, &
                                                   sums_wsus_ws_l, sums_wsvs_ws_l, &
@@ -1304,7 +1302,7 @@
                                                   sums_ls_l
                                               
-    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE   ::  hom_sum, rmask, spectrum_x,     &
-                                                  spectrum_y, sums_l, sums_l_l,   &
-                                                  sums_up_fraction_l
+    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE   ::  hom_sum, rmask, sums_l,      &
+                                                  sums_l_l, sums_up_fraction_l
+
     REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE ::  hom
 
Index: /palm/trunk/SOURCE/netcdf_interface.f90
===================================================================
--- /palm/trunk/SOURCE/netcdf_interface.f90	(revision 1832)
+++ /palm/trunk/SOURCE/netcdf_interface.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! ------------------
-!
+! spectrum renamed spectra_mod
 !
 ! Former revisions:
@@ -348,7 +348,7 @@
     USE control_parameters,                                                    &
         ONLY:  averaging_interval, averaging_interval_pr,                      &
-               averaging_interval_sp, data_output_pr,  domask,  dopr_n,        &
+               data_output_pr,  domask,  dopr_n,        &
                dopr_time_count, dopts_time_count, dots_time_count,             &
-               dosp_time_count, do2d, do2d_xz_time_count, do3d,                &
+               do2d, do2d_xz_time_count, do3d,                &
                do2d_yz_time_count, dt_data_output_av, dt_do2d_xy, dt_do2d_xz,  &
                dt_do2d_yz, dt_do3d, mask_size, do2d_xy_time_count,             &
@@ -380,6 +380,6 @@
         ONLY:  crmax, cross_profiles, dopr_index, profile_columns, profile_rows
 
-    USE spectrum,                                                              &
-        ONLY:  comp_spectra_level, data_output_sp, spectra_direction
+    USE spectra_mod,                                                           &
+        ONLY:  averaging_interval_sp, comp_spectra_level, data_output_sp, dosp_time_count, spectra_direction
 
     USE statistics,                                                            &
Index: /palm/trunk/SOURCE/package_parin.f90
===================================================================
--- /palm/trunk/SOURCE/package_parin.f90	(revision 1832)
+++ /palm/trunk/SOURCE/package_parin.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! -----------------
-!
+! reading of spectra_par moved to spectra_mod
 ! 
 ! Former revisions:
@@ -126,6 +126,5 @@
 
     USE control_parameters,                                                    &
-        ONLY:  averaging_interval_sp, dt_data_output, dt_dopts, dt_dosp,       &
-               dt_dvrp, particle_maximum_age, skip_time_dosp,                  &
+        ONLY:  dt_data_output, dt_dopts, dt_dvrp, particle_maximum_age,        &
                threshold
 
@@ -159,8 +158,4 @@
                seed_follows_topography, use_sgs_for_particles,                 &
                vertical_particle_advection, write_particle_statistics
-
-    USE spectrum,                                                              &
-        ONLY:  comp_spectra_level, data_output_sp, plot_spectra_level,         &
-               spectra_direction
 
     IMPLICIT NONE
@@ -206,8 +201,4 @@
                                   write_particle_statistics
 
-    NAMELIST /spectra_par/        averaging_interval_sp, comp_spectra_level,   &
-                                  data_output_sp, dt_dosp, plot_spectra_level, &
-                                  skip_time_dosp, spectra_direction
-
 !
 !-- Position the namelist-file at the beginning (it was already opened in
@@ -251,24 +242,3 @@
  30 CONTINUE
 
-    REWIND ( 11 )
-    line = ' '
-    DO   WHILE ( INDEX( line, '&spectra_par' ) == 0 )
-       READ ( 11, '(A)', END=40 )  line
-    ENDDO
-    BACKSPACE ( 11 )
-
-!
-!-- Read user-defined namelist
-    READ ( 11, spectra_par )
-
-!
-!-- Default setting of dt_dosp here (instead of check_parameters), because its
-!-- current value is needed in init_pegrid
-    IF ( dt_dosp == 9999999.9_wp )  dt_dosp = dt_data_output
-
- 40 CONTINUE
-
-
-
-
  END SUBROUTINE package_parin
Index: /palm/trunk/SOURCE/palm.f90
===================================================================
--- /palm/trunk/SOURCE/palm.f90	(revision 1832)
+++ /palm/trunk/SOURCE/palm.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! -----------------
-!
+! required user interface version changed
 !
 ! Former revisions:
@@ -192,5 +192,5 @@
 
     version = 'PALM 4.0'
-    user_interface_required_revision = 'r1783'
+    user_interface_required_revision = 'r1819'
 
 #if defined( __parallel )
Index: /palm/trunk/SOURCE/parin.f90
===================================================================
--- /palm/trunk/SOURCE/parin.f90	(revision 1832)
+++ /palm/trunk/SOURCE/parin.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! -----------------
-!
+! call of spectra_parin
 ! 
 ! Former revisions:
@@ -245,7 +245,10 @@
         ONLY :  batch_job
 
-     USE radiation_model_mod,                                                  &
-         ONLY: radiation_parin  
-        
+    USE radiation_model_mod,                                                   &
+        ONLY: radiation_parin  
+
+    USE spectra_mod,                                                           &
+        ONLY :  spectra_parin
+
     USE statistics,                                                            &
         ONLY:  hom, hom_sum, pr_palm, region, statistic_regions
@@ -473,9 +476,12 @@
           ENDIF
 
-          
 !
 !--       Check if land surface model is used and read &lsm_par if required
           CALL lsm_parin
- 
+
+!
+!--       Check if spectra shall be calculated and read spectra_par if required
+          CALL spectra_parin
+
 !
 !--       Check if radiation model is used and read &radiation_par if required
Index: /palm/trunk/SOURCE/pmc_client.f90
===================================================================
--- /palm/trunk/SOURCE/pmc_client.f90	(revision 1832)
+++ /palm/trunk/SOURCE/pmc_client.f90	(revision 1833)
@@ -20,5 +20,6 @@
 ! Current revisions:
 ! ------------------
-!
+! gfortran requires pointer attributes for some array declarations,
+! long line wrapped
 !
 ! Former revisions:
@@ -354,5 +355,5 @@
        IMPLICIT none
 
-       REAL(wp), INTENT(IN) ,DIMENSION(:,:) ::  array
+       REAL(wp), INTENT(IN) ,DIMENSION(:,:), POINTER ::  array
 
        INTEGER                              :: NrDims
@@ -387,5 +388,5 @@
        IMPLICIT none
 
-       REAL(wp),INTENT(IN),DIMENSION(:,:,:) ::  array
+       REAL(wp),INTENT(IN),DIMENSION(:,:,:), POINTER ::  array
 
        INTEGER                              ::  NrDims
@@ -551,5 +552,6 @@
               ar%SendBuf = c_loc(base_array_cs(ar%SendIndex))
               if(ar%SendIndex+ar%SendSize > bufsize) then
-                 write(0,'(a,i4,4i7,1x,a)') 'Client Buffer too small ',i,ar%SendIndex,ar%SendSize,ar%SendIndex+ar%SendSize,bufsize,trim(ar%name)
+                 write(0,'(a,i4,4i7,1x,a)') 'Client Buffer too small ',i,      &
+                     ar%SendIndex,ar%SendSize,ar%SendIndex+ar%SendSize,bufsize,trim(ar%name)
                  CALL MPI_Abort (MPI_COMM_WORLD, istat, ierr)
               end if
Index: /palm/trunk/SOURCE/pmc_server.f90
===================================================================
--- /palm/trunk/SOURCE/pmc_server.f90	(revision 1832)
+++ /palm/trunk/SOURCE/pmc_server.f90	(revision 1833)
@@ -20,5 +20,6 @@
 ! Current revisions:
 ! ------------------
-!
+! gfortran requires pointer attributes for some array declarations,
+! long line wrapped
 !
 ! Former revisions:
@@ -277,7 +278,8 @@
         IMPLICIT none
 
-        INTEGER,INTENT(IN)                         :: ClientId
-        REAL(wp), INTENT(IN), DIMENSION(:,:)           ::  array
-        REAL(wp), INTENT(IN), DIMENSION(:,:), OPTIONAL ::  array_2
+        INTEGER,INTENT(IN) ::  ClientId
+
+        REAL(wp), INTENT(IN), DIMENSION(:,:), POINTER           ::  array
+        REAL(wp), INTENT(IN), DIMENSION(:,:), POINTER, OPTIONAL ::  array_2
 
         INTEGER                           :: NrDims
@@ -307,7 +309,8 @@
         IMPLICIT none
 
-        INTEGER,INTENT(IN)                         :: ClientId
-        REAL(wp), INTENT(IN), DIMENSION(:,:,:)           ::  array
-        REAL(wp), INTENT(IN), DIMENSION(:,:,:), OPTIONAL ::  array_2
+        INTEGER,INTENT(IN) ::  ClientId
+
+        REAL(wp), INTENT(IN), DIMENSION(:,:,:), POINTER           ::  array
+        REAL(wp), INTENT(IN), DIMENSION(:,:,:), POINTER, OPTIONAL ::  array_2
         INTEGER,INTENT(IN)                         :: nz_cl
         INTEGER,INTENT(IN)                         :: nz
@@ -433,5 +436,6 @@
             ar%SendBuf = c_loc(base_array_sc(ar%SendIndex))
             if(ar%SendIndex+ar%SendSize > bufsize) then
-               write(0,'(a,i4,4i7,1x,a)') 'Server Buffer too small ',i,ar%SendIndex,ar%SendSize,ar%SendIndex+ar%SendSize,bufsize,trim(ar%name)
+               write(0,'(a,i4,4i7,1x,a)') 'Server Buffer too small ',i,        &
+                  ar%SendIndex,ar%SendSize,ar%SendIndex+ar%SendSize,bufsize,trim(ar%name)
                CALL MPI_Abort (MPI_COMM_WORLD, istat, ierr)
             end if
Index: /palm/trunk/SOURCE/read_3d_binary.f90
===================================================================
--- /palm/trunk/SOURCE/read_3d_binary.f90	(revision 1832)
+++ /palm/trunk/SOURCE/read_3d_binary.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! -----------------
-! 
+! statistics module replaced by spectra module
 !
 ! Former revisions:
@@ -146,5 +146,5 @@
         ONLY:  id_random_array, seq_random_array
 
-    USE statistics,                                                            &
+    USE spectra_mod,                                                           &
         ONLY:  spectrum_x, spectrum_y
 
Index: /palm/trunk/SOURCE/read_var_list.f90
===================================================================
--- /palm/trunk/SOURCE/read_var_list.f90	(revision 1832)
+++ /palm/trunk/SOURCE/read_var_list.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! ------------------
-!
+! spectra_mod added
 !
 ! Former revisions:
@@ -193,5 +193,8 @@
 
     USE radiation_model_mod,                                                   &
-        ONLY: time_radiation
+        ONLY:  time_radiation
+
+    USE spectra_mod,                                                           &
+        ONLY:  average_count_sp
 
     USE statistics,                                                            &
Index: /palm/trunk/SOURCE/spectra_mod.f90
===================================================================
--- /palm/trunk/SOURCE/spectra_mod.f90	(revision 1833)
+++ /palm/trunk/SOURCE/spectra_mod.f90	(revision 1833)
@@ -0,0 +1,852 @@
+!> @file spectra_mod.f90
+!--------------------------------------------------------------------------------!
+! This file is part of PALM.
+!
+! PALM is free software: you can redistribute it and/or modify it under the terms
+! of the GNU General Public License as published by the Free Software Foundation,
+! either version 3 of the License, or (at your option) any later version.
+!
+! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
+! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
+! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
+!
+! You should have received a copy of the GNU General Public License along with
+! PALM. If not, see <http://www.gnu.org/licenses/>.
+!
+! Copyright 1997-2016 Leibniz Universitaet Hannover
+!--------------------------------------------------------------------------------!
+!
+! Current revisions:
+! -----------------
+! file renamed, reading the spectra_par NAMELIST moved from package_parin to
+! here
+!
+! Former revisions:
+! -----------------
+! $Id$
+!
+! 1815 2016-04-06 13:49:59Z raasch
+! bugfix: preprocessor directives included for the non-parallel case
+!
+! 1808 2016-04-05 19:44:00Z raasch
+! MPI module used by default on all machines
+!
+! 1786 2016-03-08 05:49:27Z raasch
+! routine is modularized, filename renamed from calc_spectra to spectrum,
+! privious data module spectrum moved from modules.f90 to here,
+! cpp-direktives for spectra removed, immediate return if no spectra levels are
+! given
+!
+! 1682 2015-10-07 23:56:08Z knoop
+! Code annotations made doxygen readable 
+!
+! 1575 2015-03-27 09:56:27Z raasch
+! adjustments for psolver-queries
+!
+! 1511 2014-12-16 15:54:16Z suehring
+! Bugfix concerning spectra normalization
+!
+! 1431 2014-07-15 14:47:17Z suehring
+! Wavenumber-integrated spectra coincide with respective variance. 
+! 
+! 1342 2014-03-26 17:04:47Z kanani
+! REAL constants defined as wp-kinds
+!
+! 1324 2014-03-21 09:13:16Z suehring
+! Bugfix: nzb_x, nzb_yd, nyn_x, nyn_x, nzt_x, nzt_yd belong to transpose_indices 
+!
+! 1320 2014-03-20 08:40:49Z raasch
+! ONLY-attribute added to USE-statements,
+! kind-parameters added to all INTEGER and REAL declaration statements,
+! kinds are defined in new module kinds,
+! revision history before 2012 removed,
+! comment fields (!:) to be used for variable explanations added to
+! all variable declaration statements
+!
+! 1318 2014-03-17 13:35:16Z raasch
+! module interfaces removed
+!
+! 1216 2013-08-26 09:31:42Z raasch
+! resorting of array moved to separate routine resort_for_zx,
+! one argument removed from the transpose_..d routines
+!
+! 1120 2013-04-05 15:11:35Z raasch
+! bugfix: calls of fft_x|y replaced by fft_x|y_1d
+!
+! 1036 2012-10-22 13:43:42Z raasch
+! code put under GPL (PALM 3.9)
+!
+! 1003 2012-09-14 14:35:53Z raasch
+! adjustment of array tend for cases with unequal subdomain sizes removed
+!
+! Revision 1.1  2001/01/05 15:08:07  raasch
+! Initial revision
+!
+!
+! Description:
+! ------------
+!> Calculate horizontal spectra along x and y. 
+!> ATTENTION: 1d-decomposition along y still needs improvement, because in that
+!>            case the gridpoint number along z still depends on the PE number
+!>            because transpose_xz has to be used (and possibly also
+!>            transpose_zyd needs modification).
+!------------------------------------------------------------------------------!
+ MODULE spectra_mod
+
+    USE kinds
+
+    PRIVATE
+
+    CHARACTER (LEN=2),  DIMENSION(10) ::  spectra_direction = 'x'
+    CHARACTER (LEN=10), DIMENSION(10) ::  data_output_sp  = ' '
+
+    INTEGER(iwp) ::  average_count_sp = 0
+    INTEGER(iwp) ::  dosp_time_count = 0
+    INTEGER(iwp) ::  n_sp_x = 0, n_sp_y = 0
+
+    INTEGER(iwp) ::  comp_spectra_level(100) = 999999
+
+    LOGICAL ::  calculate_spectra   = .FALSE.  !< internal switch that spectra are calculated
+    LOGICAL ::  spectra_initialized = .FALSE.  !< internal switch that spectra related quantities are initialized
+
+    REAL(wp) ::  averaging_interval_sp = 9999999.9_wp  !< averaging interval for spectra output
+    REAL(wp) ::  dt_dosp = 9999999.9_wp                !< time interval for spectra output
+    REAL(wp) ::  skip_time_dosp = 9999999.9_wp         !< no output of spectra data before this interval has passed
+
+    REAL(wp), DIMENSION(:), ALLOCATABLE ::  var_d
+
+    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  spectrum_x, spectrum_y
+
+    SAVE
+
+    INTERFACE calc_spectra
+       MODULE PROCEDURE calc_spectra
+    END INTERFACE calc_spectra
+
+    INTERFACE preprocess_spectra
+       MODULE PROCEDURE preprocess_spectra
+    END INTERFACE preprocess_spectra
+
+    INTERFACE calc_spectra_x
+       MODULE PROCEDURE calc_spectra_x
+    END INTERFACE calc_spectra_x
+
+    INTERFACE calc_spectra_y
+       MODULE PROCEDURE calc_spectra_y
+    END INTERFACE calc_spectra_y
+
+    INTERFACE spectra_check_parameters
+       MODULE PROCEDURE spectra_check_parameters
+    END INTERFACE spectra_check_parameters
+
+    INTERFACE spectra_header
+       MODULE PROCEDURE spectra_header
+    END INTERFACE spectra_header
+
+    INTERFACE spectra_init
+       MODULE PROCEDURE spectra_init
+    END INTERFACE spectra_init
+
+    INTERFACE spectra_parin
+       MODULE PROCEDURE spectra_parin
+    END INTERFACE spectra_parin
+
+    PUBLIC average_count_sp, averaging_interval_sp, calc_spectra,              &
+           calculate_spectra, comp_spectra_level, data_output_sp,              &
+           dosp_time_count, dt_dosp, n_sp_x, n_sp_y, plot_spectra_level,       &
+           skip_time_dosp, spectra_check_parameters, spectra_direction,        &
+           spectra_header, spectra_init, spectra_parin, spectrum_x,            &
+           spectrum_y, var_d
+
+
+ CONTAINS
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Parin for &spectra_par for calculating spectra
+!------------------------------------------------------------------------------!
+    SUBROUTINE spectra_parin
+
+       USE control_parameters,                                                 &
+           ONLY:  dt_data_output
+
+       IMPLICIT NONE
+
+       CHARACTER (LEN=80) ::  line  !< dummy string that contains the current  &
+                                    !< line of the parameter file
+
+       NAMELIST /spectra_par/  averaging_interval_sp, comp_spectra_level,      &
+                               data_output_sp, dt_dosp, skip_time_dosp,        &
+                               spectra_direction
+
+
+!
+!--    Position the namelist-file at the beginning (it was already opened in
+!--    parin), search for the namelist-group of the package and position the
+!--    file at this line.
+       line = ' '
+
+!
+!--    Try to find the spectra package
+       REWIND ( 11 )
+       line = ' '
+       DO   WHILE ( INDEX( line, '&spectra_par' ) == 0 )
+          READ ( 11, '(A)', END=10 )  line
+       ENDDO
+       BACKSPACE ( 11 )
+
+!
+!--    Read namelist
+       READ ( 11, spectra_par )
+
+!
+!--    Default setting of dt_dosp here (instead of check_parameters), because
+!--    its current value is needed in init_pegrid
+       IF ( dt_dosp == 9999999.9_wp )  dt_dosp = dt_data_output
+
+!
+!--    Set general switch that spectra shall be calculated
+       calculate_spectra = .TRUE.
+
+ 10    CONTINUE
+
+    END SUBROUTINE spectra_parin
+
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Initialization of spectra related variables
+!------------------------------------------------------------------------------!
+    SUBROUTINE spectra_init
+
+       USE indices,                                                            &
+           ONLY:  nx, ny, nzb, nzt
+
+       IMPLICIT NONE
+
+       IF ( spectra_initialized )  RETURN
+
+       IF ( dt_dosp /= 9999999.9_wp )  THEN
+          ALLOCATE( spectrum_x( 1:nx/2, 1:10, 1:10 ),                          &
+                    spectrum_y( 1:ny/2, 1:10, 1:10 ) )
+          spectrum_x = 0.0_wp
+          spectrum_y = 0.0_wp
+
+          ALLOCATE( var_d(nzb:nzt+1) )
+          var_d = 0.0_wp
+       ENDIF
+
+       spectra_initialized = .TRUE.
+
+    END SUBROUTINE spectra_init
+
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Check spectra related quantities
+!------------------------------------------------------------------------------!
+    SUBROUTINE spectra_check_parameters
+
+       USE control_parameters,                                                 &
+           ONLY:  averaging_interval, message_string, skip_time_data_output
+
+       IMPLICIT NONE
+
+!
+!--    Check the average interval
+       IF ( averaging_interval_sp == 9999999.9_wp )  THEN
+          averaging_interval_sp = averaging_interval
+       ENDIF
+
+       IF ( averaging_interval_sp > dt_dosp )  THEN
+          WRITE( message_string, * )  'averaging_interval_sp = ',              &
+                averaging_interval_sp, ' must be <= dt_dosp = ', dt_dosp
+          CALL message( 'spectra_check_parameters', 'PA0087', 1, 2, 0, 6, 0 )
+       ENDIF
+
+!
+!--    Set the default skip time interval for data output, if necessary
+       IF ( skip_time_dosp == 9999999.9_wp )                                   &
+                                          skip_time_dosp = skip_time_data_output
+
+    END SUBROUTINE spectra_check_parameters
+
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Header output for spectra
+!>
+!> @todo Output of netcdf data format and compression level
+!------------------------------------------------------------------------------!
+    SUBROUTINE spectra_header ( io )
+
+       USE control_parameters,                                                 &
+           ONLY:  dt_averaging_input_pr
+
+!       USE netcdf_interface,                                                  &
+!           ONLY:  netcdf_data_format_string, netcdf_deflate
+
+       IMPLICIT NONE
+
+       CHARACTER (LEN=40) ::  output_format       !< internal string
+
+       INTEGER(iwp) ::  i                         !< internal counter
+       INTEGER(iwp), INTENT(IN) ::  io            !< Unit of the output file
+
+!
+!--    Spectra output
+       IF ( dt_dosp /= 9999999.9_wp )  THEN
+          WRITE ( io, 1 )
+
+!          output_format = netcdf_data_format_string
+!          IF ( netcdf_deflate == 0 )  THEN
+!             WRITE ( io, 2 )  output_format
+!          ELSE
+!             WRITE ( io, 3 )  TRIM( output_format ), netcdf_deflate
+!          ENDIF
+          WRITE ( io, 2 )  'see profiles or other quantities'
+          WRITE ( io, 4 )  dt_dosp
+          IF ( skip_time_dosp /= 0.0_wp )  WRITE ( io, 5 )  skip_time_dosp
+          WRITE ( io, 6 )  ( data_output_sp(i), i = 1,10 ),     &
+                           ( spectra_direction(i), i = 1,10 ),  &
+                           ( comp_spectra_level(i), i = 1,100 ), &
+                           averaging_interval_sp, dt_averaging_input_pr
+       ENDIF
+
+     1 FORMAT ('    Spectra:')
+     2 FORMAT ('       Output format: ',A/)
+     3 FORMAT ('       Output format: ',A, '   compressed with level: ',I1/)
+     4 FORMAT ('       Output every ',F7.1,' s'/)
+     5 FORMAT ('       No output during initial ',F8.2,' s')
+     6 FORMAT ('       Arrays:     ', 10(A5,',')/                         &
+               '       Directions: ', 10(A5,',')/                         &
+               '       height levels  k = ', 20(I3,',')/                  &
+               '                          ', 20(I3,',')/                  &
+               '                          ', 20(I3,',')/                  &
+               '                          ', 20(I3,',')/                  &
+               '                          ', 19(I3,','),I3,'.'/           &
+               '       Time averaged over ', F7.1, ' s,' /                &
+               '       Profiles for the time averaging are taken every ', &
+                    F6.1,' s')
+
+    END SUBROUTINE spectra_header
+
+
+
+    SUBROUTINE calc_spectra
+
+       USE arrays_3d,                                                          &
+           ONLY:  d, tend
+
+       USE control_parameters,                                                 &
+           ONLY:  bc_lr_cyc, bc_ns_cyc, message_string, psolver
+
+       USE cpulog,                                                             &
+           ONLY:  cpu_log, log_point
+
+       USE fft_xy,                                                             &
+           ONLY:  fft_init
+
+       USE indices,                                                            &
+           ONLY:  nxl, nxr, nyn, nys, nzb, nzt
+
+       USE kinds
+
+       USE pegrid,                                                             &
+           ONLY:  myid, pdims
+
+       IMPLICIT NONE
+
+       INTEGER(iwp) ::  m  !<
+       INTEGER(iwp) ::  pr !<
+
+
+!
+!--    Check if user gave any levels for spectra to be calculated
+       IF ( comp_spectra_level(1) == 999999 )  RETURN
+
+       CALL cpu_log( log_point(30), 'calc_spectra', 'start' )
+
+!
+!--    Initialize spectra related quantities
+       CALL spectra_init
+
+!
+!--    Initialize ffts
+       CALL fft_init
+
+!
+!--    Reallocate array d in required size
+       IF ( psolver(1:9) == 'multigrid' )  THEN
+          DEALLOCATE( d )
+          ALLOCATE( d(nzb+1:nzt,nys:nyn,nxl:nxr) )
+       ENDIF
+
+       m = 1
+       DO WHILE ( data_output_sp(m) /= ' '  .AND.  m <= 10 )
+!
+!--       Transposition from z --> x  ( y --> x in case of a 1d-decomposition
+!--       along x)
+          IF ( INDEX( spectra_direction(m), 'x' ) /= 0 )  THEN
+
+!
+!--          Calculation of spectra works for cyclic boundary conditions only
+             IF ( .NOT. bc_lr_cyc )  THEN
+
+                message_string = 'non-cyclic lateral boundaries along x do'//  &
+                                 ' not & allow calculation of spectra along x'
+                CALL message( 'calc_spectra', 'PA0160', 1, 2, 0, 6, 0 )
+             ENDIF
+
+             CALL preprocess_spectra( m, pr )
+
+#if defined( __parallel )
+             IF ( pdims(2) /= 1 )  THEN
+                CALL resort_for_zx( d, tend )
+                CALL transpose_zx( tend, d )
+             ELSE
+                CALL transpose_yxd( d, d )
+             ENDIF
+             CALL calc_spectra_x( d, pr, m )
+#else
+             message_string = 'sorry, calculation of spectra in non paral' //  &
+                              'lel mode& is still not realized'
+             CALL message( 'calc_spectra', 'PA0161', 1, 2, 0, 6, 0 )
+#endif
+
+          ENDIF
+
+!
+!--       Transposition from z --> y (d is rearranged only in case of a
+!--       1d-decomposition along x)
+          IF ( INDEX( spectra_direction(m), 'y' ) /= 0 )  THEN
+
+!
+!--          Calculation of spectra works for cyclic boundary conditions only
+             IF ( .NOT. bc_ns_cyc )  THEN
+                IF ( myid == 0 )  THEN
+                   message_string = 'non-cyclic lateral boundaries along y' // &
+                                    ' do not & allow calculation of spectr' // &
+                                    'a along y'
+                   CALL message( 'calc_spectra', 'PA0162', 1, 2, 0, 6, 0 )
+                ENDIF
+                CALL local_stop
+             ENDIF
+
+             CALL preprocess_spectra( m, pr )
+
+#if defined( __parallel )
+             CALL transpose_zyd( d, d )
+             CALL calc_spectra_y( d, pr, m )
+#else
+             message_string = 'sorry, calculation of spectra in non paral' //  &
+                              'lel mode& is still not realized'
+             CALL message( 'calc_spectra', 'PA0161', 1, 2, 0, 6, 0 )
+#endif
+
+          ENDIF
+
+!
+!--       Increase counter for next spectrum
+          m = m + 1
+          
+       ENDDO
+
+!
+!--    Increase counter for averaging process in routine plot_spectra
+       average_count_sp = average_count_sp + 1
+
+       CALL cpu_log( log_point(30), 'calc_spectra', 'stop' )
+
+    END SUBROUTINE calc_spectra
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> @todo Missing subroutine description.
+!------------------------------------------------------------------------------!
+    SUBROUTINE preprocess_spectra( m, pr )
+
+       USE arrays_3d,                                                          &
+           ONLY:  d, pt, q, u, v, w
+
+       USE indices,                                                            &
+           ONLY:  ngp_2dh, nxl, nxr, nyn, nys, nzb, nzt
+
+       USE kinds
+
+#if defined( __parallel )
+#if defined( __mpifh )
+       INCLUDE "mpif.h"
+#else
+       USE MPI
+#endif
+#endif
+       USE pegrid,                                                             &
+           ONLY:  collective_wait, comm2d, ierr
+
+       USE statistics,                                                         &
+           ONLY:  hom
+
+
+       IMPLICIT NONE
+
+       INTEGER(iwp) :: i  !<
+       INTEGER(iwp) :: j  !<
+       INTEGER(iwp) :: k  !<
+       INTEGER(iwp) :: m  !<
+       INTEGER(iwp) :: pr !<
+
+       REAL(wp), DIMENSION(nzb:nzt+1) :: var_d_l
+
+       SELECT CASE ( TRIM( data_output_sp(m) ) )
+          
+       CASE ( 'u' )
+          pr = 1
+          d(nzb+1:nzt,nys:nyn,nxl:nxr) = u(nzb+1:nzt,nys:nyn,nxl:nxr)
+       
+       CASE ( 'v' )
+          pr = 2
+          d(nzb+1:nzt,nys:nyn,nxl:nxr) = v(nzb+1:nzt,nys:nyn,nxl:nxr)
+       
+       CASE ( 'w' )
+          pr = 3
+          d(nzb+1:nzt,nys:nyn,nxl:nxr) = w(nzb+1:nzt,nys:nyn,nxl:nxr)
+       
+       CASE ( 'pt' )
+          pr = 4
+          d(nzb+1:nzt,nys:nyn,nxl:nxr) = pt(nzb+1:nzt,nys:nyn,nxl:nxr)
+       
+       CASE ( 'q' )
+          pr = 41
+          d(nzb+1:nzt,nys:nyn,nxl:nxr) = q(nzb+1:nzt,nys:nyn,nxl:nxr)
+       
+       CASE DEFAULT
+!
+!--       The DEFAULT case is reached either if the parameter data_output_sp(m)
+!--       contains a wrong character string or if the user has coded a special
+!--       case in the user interface. There, the subroutine user_spectra
+!--       checks which of these two conditions applies.
+          CALL user_spectra( 'preprocess', m, pr )
+          
+       END SELECT
+
+!
+!--    Subtract horizontal mean from the array, for which spectra have to be
+!--    calculated
+       var_d_l(:) = 0.0_wp
+       DO  i = nxl, nxr
+          DO  j = nys, nyn
+             DO  k = nzb+1, nzt
+                d(k,j,i)   = d(k,j,i) - hom(k,1,pr,0)
+                var_d_l(k) = var_d_l(k) + d(k,j,i) * d(k,j,i)
+             ENDDO
+          ENDDO
+       ENDDO
+!
+!--    Compute total variance from local variances
+       var_d(:) = 0.0_wp
+#if defined( __parallel ) 
+       IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
+       CALL MPI_ALLREDUCE( var_d_l(0), var_d(0), nzt+1-nzb, MPI_REAL, MPI_SUM, &
+                           comm2d, ierr )
+#else
+       var_d(:) = var_d_l(:)
+#endif
+       var_d(:) = var_d(:) / ngp_2dh(0)
+
+    END SUBROUTINE preprocess_spectra
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> @todo Missing subroutine description.
+!------------------------------------------------------------------------------!
+    SUBROUTINE calc_spectra_x( ddd, pr, m )
+
+       USE control_parameters,                                                 &
+           ONLY:  fft_method
+
+       USE fft_xy,                                                             &
+           ONLY:  fft_x_1d
+
+       USE grid_variables,                                                     &
+           ONLY:  dx
+
+       USE indices,                                                            &
+           ONLY:  nx, ny
+
+       USE kinds
+
+#if defined( __parallel )
+#if defined( __mpifh )
+       INCLUDE "mpif.h"
+#else
+       USE MPI
+#endif
+#endif
+       USE pegrid,                                                             &
+           ONLY:  comm2d, ierr, myid
+
+       USE transpose_indices,                                                  &
+           ONLY:  nyn_x, nys_x, nzb_x, nzt_x
+
+
+       IMPLICIT NONE
+
+       INTEGER(iwp) ::  i         !<
+       INTEGER(iwp) ::  ishape(1) !<
+       INTEGER(iwp) ::  j         !<
+       INTEGER(iwp) ::  k         !<
+       INTEGER(iwp) ::  m         !<
+       INTEGER(iwp) ::  n         !<
+       INTEGER(iwp) ::  pr        !<
+
+       REAL(wp) ::  exponent     !<
+       REAL(wp) ::  sum_spec_dum !< wavenumber-integrated spectrum
+    
+       REAL(wp), DIMENSION(0:nx) ::  work !<
+    
+       REAL(wp), DIMENSION(0:nx/2) ::  sums_spectra_l !<
+    
+       REAL(wp), DIMENSION(0:nx/2,100) ::  sums_spectra !<
+    
+       REAL(wp), DIMENSION(0:nx,nys_x:nyn_x,nzb_x:nzt_x) ::  ddd !<
+
+!
+!--    Exponent for geometric average
+       exponent = 1.0_wp / ( ny + 1.0_wp )
+
+!
+!--    Loop over all levels defined by the user
+       n = 1
+       DO WHILE ( comp_spectra_level(n) /= 999999  .AND.  n <= 100 )
+
+          k = comp_spectra_level(n)
+
+!
+!--       Calculate FFT only if the corresponding level is situated on this PE
+          IF ( k >= nzb_x  .AND.  k <= nzt_x )  THEN
+          
+             DO  j = nys_x, nyn_x
+
+                work = ddd(0:nx,j,k)
+                CALL fft_x_1d( work, 'forward' )
+
+                ddd(0,j,k) = dx * work(0)**2
+                DO  i = 1, nx/2
+                   ddd(i,j,k) = dx * ( work(i)**2 + work(nx+1-i)**2 )
+                ENDDO
+
+             ENDDO
+
+!
+!--          Local sum and geometric average of these spectra
+!--          (WARNING: no global sum should be performed, because floating
+!--          point overflow may occur)
+             DO  i = 0, nx/2
+
+                sums_spectra_l(i) = 1.0_wp
+                DO  j = nys_x, nyn_x
+                   sums_spectra_l(i) = sums_spectra_l(i) * ddd(i,j,k)**exponent
+                ENDDO
+
+             ENDDO
+          
+          ELSE
+
+             sums_spectra_l = 1.0_wp
+
+          ENDIF
+
+!
+!--       Global sum of spectra on PE0 (from where they are written on file)
+          sums_spectra(:,n) = 0.0_wp
+#if defined( __parallel )   
+          CALL MPI_BARRIER( comm2d, ierr )  ! Necessary?
+          CALL MPI_REDUCE( sums_spectra_l(0), sums_spectra(0,n), nx/2+1,       &
+                           MPI_REAL, MPI_PROD, 0, comm2d, ierr )
+#else
+          sums_spectra(:,n) = sums_spectra_l
+#endif
+!
+!--       Normalize spectra by variance
+          sum_spec_dum = SUM( sums_spectra(:,n) )
+          IF ( sum_spec_dum /= 0.0_wp )  THEN
+             sums_spectra(:,n) = sums_spectra(:,n) * var_d(k) / sum_spec_dum
+          ENDIF
+          n = n + 1
+
+       ENDDO
+       n = n - 1
+
+       IF ( myid == 0 )  THEN
+!
+!--       Sum of spectra for later averaging (see routine data_output_spectra)
+          DO  i = 1, nx/2
+             DO k = 1, n
+                spectrum_x(i,k,m) = spectrum_x(i,k,m) + sums_spectra(i,k)
+             ENDDO
+          ENDDO
+
+       ENDIF
+!
+!--    n_sp_x is needed by data_output_spectra_x
+       n_sp_x = n
+
+    END SUBROUTINE calc_spectra_x
+
+
+!------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> @todo Missing subroutine description.
+!------------------------------------------------------------------------------!
+    SUBROUTINE calc_spectra_y( ddd, pr, m )
+
+       USE control_parameters,                                                 &
+           ONLY:  fft_method
+
+       USE fft_xy,                                                             &
+           ONLY:  fft_y_1d
+
+       USE grid_variables,                                                     &
+           ONLY:  dy
+
+       USE indices,                                                            &
+           ONLY:  nx, ny
+
+       USE kinds
+
+#if defined( __parallel )
+#if defined( __mpifh )
+       INCLUDE "mpif.h"
+#else
+       USE MPI
+#endif
+#endif
+       USE pegrid,                                                             &
+           ONLY:  comm2d, ierr, myid
+
+       USE transpose_indices,                                                  &
+           ONLY:  nxl_yd, nxr_yd, nzb_yd, nzt_yd
+
+
+       IMPLICIT NONE
+
+       INTEGER(iwp) ::  i         !<
+       INTEGER(iwp) ::  j         !<
+       INTEGER(iwp) ::  jshape(1) !<
+       INTEGER(iwp) ::  k         !<
+       INTEGER(iwp) ::  m         !<
+       INTEGER(iwp) ::  n         !<
+       INTEGER(iwp) ::  pr        !<
+
+       REAL(wp) ::  exponent !<
+       REAL(wp) ::  sum_spec_dum !< wavenumber-integrated spectrum
+    
+       REAL(wp), DIMENSION(0:ny) ::  work !<
+    
+       REAL(wp), DIMENSION(0:ny/2) ::  sums_spectra_l !<
+    
+       REAL(wp), DIMENSION(0:ny/2,100) ::  sums_spectra !<
+    
+       REAL(wp), DIMENSION(0:ny,nxl_yd:nxr_yd,nzb_yd:nzt_yd) :: ddd !<
+
+
+!
+!--    Exponent for geometric average
+       exponent = 1.0_wp / ( nx + 1.0_wp )
+
+!
+!--    Loop over all levels defined by the user
+       n = 1
+       DO WHILE ( comp_spectra_level(n) /= 999999  .AND.  n <= 100 )
+
+          k = comp_spectra_level(n)
+
+!
+!--       Calculate FFT only if the corresponding level is situated on this PE
+          IF ( k >= nzb_yd  .AND.  k <= nzt_yd )  THEN
+          
+             DO  i = nxl_yd, nxr_yd
+
+                work = ddd(0:ny,i,k)
+                CALL fft_y_1d( work, 'forward' )
+
+                ddd(0,i,k) = dy * work(0)**2
+                DO  j = 1, ny/2
+                   ddd(j,i,k) = dy * ( work(j)**2 + work(ny+1-j)**2 )
+                ENDDO
+
+             ENDDO
+
+!
+!--          Local sum and geometric average of these spectra
+!--          (WARNING: no global sum should be performed, because floating
+!--          point overflow may occur)
+             DO  j = 0, ny/2
+
+                sums_spectra_l(j) = 1.0_wp
+                DO  i = nxl_yd, nxr_yd
+                   sums_spectra_l(j) = sums_spectra_l(j) * ddd(j,i,k)**exponent
+                ENDDO
+
+             ENDDO
+          
+          ELSE
+
+             sums_spectra_l = 1.0_wp
+
+          ENDIF
+
+!
+!--       Global sum of spectra on PE0 (from where they are written on file)
+          sums_spectra(:,n) = 0.0_wp
+#if defined( __parallel )   
+          CALL MPI_BARRIER( comm2d, ierr )  ! Necessary?
+          CALL MPI_REDUCE( sums_spectra_l(0), sums_spectra(0,n), ny/2+1,       &
+                           MPI_REAL, MPI_PROD, 0, comm2d, ierr )
+#else
+          sums_spectra(:,n) = sums_spectra_l
+#endif
+!
+!--       Normalize spectra by variance
+          sum_spec_dum = SUM( sums_spectra(:,n) )
+          IF ( SUM(sums_spectra(:,n)) /= 0.0_wp )  THEN
+             sums_spectra(:,n) = sums_spectra(:,n) *                           &
+                                 var_d(k) / SUM(sums_spectra(:,n))
+          ENDIF
+          n = n + 1
+
+       ENDDO
+       n = n - 1
+
+
+       IF ( myid == 0 )  THEN
+!
+!--       Sum of spectra for later averaging (see routine data_output_spectra)
+          DO  j = 1, ny/2
+             DO k = 1, n
+                spectrum_y(j,k,m) = spectrum_y(j,k,m) + sums_spectra(j,k)
+             ENDDO
+          ENDDO
+
+       ENDIF
+
+!
+!--    n_sp_y is needed by data_output_spectra_y
+       n_sp_y = n
+
+    END SUBROUTINE calc_spectra_y
+
+ END MODULE spectra_mod
Index: lm/trunk/SOURCE/spectrum.f90
===================================================================
--- /palm/trunk/SOURCE/spectrum.f90	(revision 1832)
+++ 	(revision )
@@ -1,671 +1,0 @@
-!> @file spectrum.f90
-!--------------------------------------------------------------------------------!
-! This file is part of PALM.
-!
-! PALM is free software: you can redistribute it and/or modify it under the terms
-! of the GNU General Public License as published by the Free Software Foundation,
-! either version 3 of the License, or (at your option) any later version.
-!
-! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
-! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
-! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
-!
-! You should have received a copy of the GNU General Public License along with
-! PALM. If not, see <http://www.gnu.org/licenses/>.
-!
-! Copyright 1997-2016 Leibniz Universitaet Hannover
-!--------------------------------------------------------------------------------!
-!
-! Current revisions:
-! -----------------
-!
-!
-! Former revisions:
-! -----------------
-! $Id$
-!
-! 1815 2016-04-06 13:49:59Z raasch
-! bugfix: preprocessor directives included for the non-parallel case
-!
-! 1808 2016-04-05 19:44:00Z raasch
-! MPI module used by default on all machines
-!
-! 1786 2016-03-08 05:49:27Z raasch
-! routine is modularized, filename renamed from calc_spectra to spectrum,
-! privious data module spectrum moved from modules.f90 to here,
-! cpp-direktives for spectra removed, immediate return if no spectra levels are
-! given
-!
-! 1682 2015-10-07 23:56:08Z knoop
-! Code annotations made doxygen readable 
-!
-! 1575 2015-03-27 09:56:27Z raasch
-! adjustments for psolver-queries
-!
-! 1511 2014-12-16 15:54:16Z suehring
-! Bugfix concerning spectra normalization
-!
-! 1431 2014-07-15 14:47:17Z suehring
-! Wavenumber-integrated spectra coincide with respective variance. 
-! 
-! 1342 2014-03-26 17:04:47Z kanani
-! REAL constants defined as wp-kinds
-!
-! 1324 2014-03-21 09:13:16Z suehring
-! Bugfix: nzb_x, nzb_yd, nyn_x, nyn_x, nzt_x, nzt_yd belong to transpose_indices 
-!
-! 1320 2014-03-20 08:40:49Z raasch
-! ONLY-attribute added to USE-statements,
-! kind-parameters added to all INTEGER and REAL declaration statements,
-! kinds are defined in new module kinds,
-! revision history before 2012 removed,
-! comment fields (!:) to be used for variable explanations added to
-! all variable declaration statements
-!
-! 1318 2014-03-17 13:35:16Z raasch
-! module interfaces removed
-!
-! 1216 2013-08-26 09:31:42Z raasch
-! resorting of array moved to separate routine resort_for_zx,
-! one argument removed from the transpose_..d routines
-!
-! 1120 2013-04-05 15:11:35Z raasch
-! bugfix: calls of fft_x|y replaced by fft_x|y_1d
-!
-! 1036 2012-10-22 13:43:42Z raasch
-! code put under GPL (PALM 3.9)
-!
-! 1003 2012-09-14 14:35:53Z raasch
-! adjustment of array tend for cases with unequal subdomain sizes removed
-!
-! Revision 1.1  2001/01/05 15:08:07  raasch
-! Initial revision
-!
-!
-! Description:
-! ------------
-!> Calculate horizontal spectra along x and y. 
-!> ATTENTION: 1d-decomposition along y still needs improvement, because in that
-!>            case the gridpoint number along z still depends on the PE number
-!>            because transpose_xz has to be used (and possibly also
-!>            transpose_zyd needs modification).
-!------------------------------------------------------------------------------!
- MODULE spectrum
-
-    USE kinds
-
-    PRIVATE
-
-    CHARACTER (LEN=6),  DIMENSION(1:5) ::  header_char = (/ 'PS(u) ', 'PS(v) ',&
-                                           'PS(w) ', 'PS(pt)', 'PS(q) ' /)
-    CHARACTER (LEN=2),  DIMENSION(10)  ::  spectra_direction = 'x'
-    CHARACTER (LEN=10), DIMENSION(10)  ::  data_output_sp  = ' '
-    CHARACTER (LEN=25), DIMENSION(1:5) ::  utext_char =                    &
-                                           (/ '-power spectrum of u     ', &
-                                              '-power spectrum of v     ', &
-                                              '-power spectrum of w     ', &
-                                              '-power spectrum of ^1185 ', &
-                                              '-power spectrum of q     ' /)
-    CHARACTER (LEN=39), DIMENSION(1:5) ::  ytext_char =                        &
-                                 (/ 'k ^2236 ^2566^2569<u(k) in m>2s>->2    ', &
-                                    'k ^2236 ^2566^2569<v(k) in m>2s>->2    ', &
-                                    'k ^2236 ^2566^2569<w(k) in m>2s>->2    ', &
-                                    'k ^2236 ^2566^2569<^1185(k) in m>2s>->2', &
-                                    'k ^2236 ^2566^2569<q(k) in m>2s>->2    ' /)
-
-    INTEGER(iwp) ::  klist_x = 0, klist_y = 0, n_sp_x = 0, n_sp_y = 0
-
-    INTEGER(iwp) ::  comp_spectra_level(100) = 999999,                   &
-                     lstyles(100) = (/ 0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
-                                       0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
-                                       0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
-                                       0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
-                                       0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
-                                       0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
-                                       0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
-                                       0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
-                                       0, 7, 3, 10, 1, 4, 9, 2, 6, 8,    &
-                                       0, 7, 3, 10, 1, 4, 9, 2, 6, 8 /), &
-                     plot_spectra_level(100) = 999999
-
-    REAL(wp)    ::  time_to_start_sp = 0.0_wp
-
-    PUBLIC  comp_spectra_level, data_output_sp, header_char, klist_x, klist_y, &
-            lstyles, n_sp_x, n_sp_y, plot_spectra_level, spectra_direction,    &
-            utext_char, ytext_char
-
-    SAVE
-
-    INTERFACE calc_spectra
-       MODULE PROCEDURE calc_spectra
-    END INTERFACE calc_spectra
-
-    INTERFACE preprocess_spectra
-       MODULE PROCEDURE preprocess_spectra
-    END INTERFACE preprocess_spectra
-
-    INTERFACE calc_spectra_x
-       MODULE PROCEDURE calc_spectra_x
-    END INTERFACE calc_spectra_x
-
-    INTERFACE calc_spectra_y
-       MODULE PROCEDURE calc_spectra_y
-    END INTERFACE calc_spectra_y
-
-    PUBLIC calc_spectra
-
-
- CONTAINS
-
-    SUBROUTINE calc_spectra
-
-       USE arrays_3d,                                                          &
-           ONLY:  d, tend
-
-       USE control_parameters,                                                 &
-           ONLY:  average_count_sp, bc_lr_cyc, bc_ns_cyc, message_string, psolver
-
-       USE cpulog,                                                             &
-           ONLY:  cpu_log, log_point
-
-       USE fft_xy,                                                             &
-           ONLY:  fft_init
-
-       USE indices,                                                            &
-           ONLY:  nxl, nxr, nyn, nys, nzb, nzt
-
-       USE kinds
-
-       USE pegrid,                                                             &
-           ONLY:  myid, pdims
-
-       IMPLICIT NONE
-
-       INTEGER(iwp) ::  m  !<
-       INTEGER(iwp) ::  pr !<
-
-
-!
-!--    Check if user gave any levels for spectra to be calculated
-       IF ( comp_spectra_level(1) == 999999 )  RETURN
-
-       CALL cpu_log( log_point(30), 'calc_spectra', 'start' )
-
-!
-!--    Initialize ffts
-       CALL fft_init
-
-!
-!--    Reallocate array d in required size
-       IF ( psolver(1:9) == 'multigrid' )  THEN
-          DEALLOCATE( d )
-          ALLOCATE( d(nzb+1:nzt,nys:nyn,nxl:nxr) )
-       ENDIF
-
-       m = 1
-       DO WHILE ( data_output_sp(m) /= ' '  .AND.  m <= 10 )
-!
-!--       Transposition from z --> x  ( y --> x in case of a 1d-decomposition
-!--       along x)
-          IF ( INDEX( spectra_direction(m), 'x' ) /= 0 )  THEN
-
-!
-!--          Calculation of spectra works for cyclic boundary conditions only
-             IF ( .NOT. bc_lr_cyc )  THEN
-
-                message_string = 'non-cyclic lateral boundaries along x do'//  &
-                                 ' not & allow calculation of spectra along x'
-                CALL message( 'calc_spectra', 'PA0160', 1, 2, 0, 6, 0 )
-             ENDIF
-
-             CALL preprocess_spectra( m, pr )
-
-#if defined( __parallel )
-             IF ( pdims(2) /= 1 )  THEN
-                CALL resort_for_zx( d, tend )
-                CALL transpose_zx( tend, d )
-             ELSE
-                CALL transpose_yxd( d, d )
-             ENDIF
-             CALL calc_spectra_x( d, pr, m )
-#else
-             message_string = 'sorry, calculation of spectra in non paral' //  &
-                              'lel mode& is still not realized'
-             CALL message( 'calc_spectra', 'PA0161', 1, 2, 0, 6, 0 )
-#endif
-
-          ENDIF
-
-!
-!--       Transposition from z --> y (d is rearranged only in case of a
-!--       1d-decomposition along x)
-          IF ( INDEX( spectra_direction(m), 'y' ) /= 0 )  THEN
-
-!
-!--          Calculation of spectra works for cyclic boundary conditions only
-             IF ( .NOT. bc_ns_cyc )  THEN
-                IF ( myid == 0 )  THEN
-                   message_string = 'non-cyclic lateral boundaries along y' // &
-                                    ' do not & allow calculation of spectr' // &
-                                    'a along y'
-                   CALL message( 'calc_spectra', 'PA0162', 1, 2, 0, 6, 0 )
-                ENDIF
-                CALL local_stop
-             ENDIF
-
-             CALL preprocess_spectra( m, pr )
-
-#if defined( __parallel )
-             CALL transpose_zyd( d, d )
-             CALL calc_spectra_y( d, pr, m )
-#else
-             message_string = 'sorry, calculation of spectra in non paral' //  &
-                              'lel mode& is still not realized'
-             CALL message( 'calc_spectra', 'PA0161', 1, 2, 0, 6, 0 )
-#endif
-
-          ENDIF
-
-!
-!--       Increase counter for next spectrum
-          m = m + 1
-          
-       ENDDO
-
-!
-!--    Increase counter for averaging process in routine plot_spectra
-       average_count_sp = average_count_sp + 1
-
-       CALL cpu_log( log_point(30), 'calc_spectra', 'stop' )
-
-    END SUBROUTINE calc_spectra
-
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> @todo Missing subroutine description.
-!------------------------------------------------------------------------------!
-    SUBROUTINE preprocess_spectra( m, pr )
-
-       USE arrays_3d,                                                          &
-           ONLY:  d, pt, q, u, v, w
-
-       USE indices,                                                            &
-           ONLY:  ngp_2dh, nxl, nxr, nyn, nys, nzb, nzt
-
-       USE kinds
-
-#if defined( __parallel )
-#if defined( __mpifh )
-       INCLUDE "mpif.h"
-#else
-       USE MPI
-#endif
-#endif
-       USE pegrid,                                                             &
-           ONLY:  collective_wait, comm2d, ierr
-
-       USE statistics,                                                         &
-           ONLY:  hom, var_d
-
-
-       IMPLICIT NONE
-
-       INTEGER(iwp) :: i  !<
-       INTEGER(iwp) :: j  !<
-       INTEGER(iwp) :: k  !<
-       INTEGER(iwp) :: m  !<
-       INTEGER(iwp) :: pr !<
-
-       REAL(wp), DIMENSION(nzb:nzt+1) :: var_d_l
-
-       SELECT CASE ( TRIM( data_output_sp(m) ) )
-          
-       CASE ( 'u' )
-          pr = 1
-          d(nzb+1:nzt,nys:nyn,nxl:nxr) = u(nzb+1:nzt,nys:nyn,nxl:nxr)
-       
-       CASE ( 'v' )
-          pr = 2
-          d(nzb+1:nzt,nys:nyn,nxl:nxr) = v(nzb+1:nzt,nys:nyn,nxl:nxr)
-       
-       CASE ( 'w' )
-          pr = 3
-          d(nzb+1:nzt,nys:nyn,nxl:nxr) = w(nzb+1:nzt,nys:nyn,nxl:nxr)
-       
-       CASE ( 'pt' )
-          pr = 4
-          d(nzb+1:nzt,nys:nyn,nxl:nxr) = pt(nzb+1:nzt,nys:nyn,nxl:nxr)
-       
-       CASE ( 'q' )
-          pr = 41
-          d(nzb+1:nzt,nys:nyn,nxl:nxr) = q(nzb+1:nzt,nys:nyn,nxl:nxr)
-       
-       CASE DEFAULT
-!
-!--       The DEFAULT case is reached either if the parameter data_output_sp(m)
-!--       contains a wrong character string or if the user has coded a special
-!--       case in the user interface. There, the subroutine user_spectra
-!--       checks which of these two conditions applies.
-          CALL user_spectra( 'preprocess', m, pr )
-          
-       END SELECT
-
-!
-!--    Subtract horizontal mean from the array, for which spectra have to be
-!--    calculated
-       var_d_l(:) = 0.0_wp
-       DO  i = nxl, nxr
-          DO  j = nys, nyn
-             DO  k = nzb+1, nzt
-                d(k,j,i)   = d(k,j,i) - hom(k,1,pr,0)
-                var_d_l(k) = var_d_l(k) + d(k,j,i) * d(k,j,i)
-             ENDDO
-          ENDDO
-       ENDDO
-!
-!--    Compute total variance from local variances
-       var_d(:) = 0.0_wp
-#if defined( __parallel ) 
-       IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
-       CALL MPI_ALLREDUCE( var_d_l(0), var_d(0), nzt+1-nzb, MPI_REAL, MPI_SUM, &
-                           comm2d, ierr )
-#else
-       var_d(:) = var_d_l(:)
-#endif
-       var_d(:) = var_d(:) / ngp_2dh(0)
-
-    END SUBROUTINE preprocess_spectra
-
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> @todo Missing subroutine description.
-!------------------------------------------------------------------------------!
-    SUBROUTINE calc_spectra_x( ddd, pr, m )
-
-       USE control_parameters,                                                 &
-           ONLY:  fft_method
-
-       USE fft_xy,                                                             &
-           ONLY:  fft_x_1d
-
-       USE grid_variables,                                                     &
-           ONLY:  dx
-
-       USE indices,                                                            &
-           ONLY:  nx, ny
-
-       USE kinds
-
-#if defined( __parallel )
-#if defined( __mpifh )
-       INCLUDE "mpif.h"
-#else
-       USE MPI
-#endif
-#endif
-       USE pegrid,                                                             &
-           ONLY:  comm2d, ierr, myid
-
-       USE statistics,                                                         &
-           ONLY:  spectrum_x, var_d
-
-       USE transpose_indices,                                                  &
-           ONLY:  nyn_x, nys_x, nzb_x, nzt_x
-
-
-       IMPLICIT NONE
-
-       INTEGER(iwp) ::  i         !<
-       INTEGER(iwp) ::  ishape(1) !<
-       INTEGER(iwp) ::  j         !<
-       INTEGER(iwp) ::  k         !<
-       INTEGER(iwp) ::  m         !<
-       INTEGER(iwp) ::  n         !<
-       INTEGER(iwp) ::  pr        !<
-
-       REAL(wp) ::  exponent     !<
-       REAL(wp) ::  sum_spec_dum !< wavenumber-integrated spectrum
-    
-       REAL(wp), DIMENSION(0:nx) ::  work !<
-    
-       REAL(wp), DIMENSION(0:nx/2) ::  sums_spectra_l !<
-    
-       REAL(wp), DIMENSION(0:nx/2,100) ::  sums_spectra !<
-    
-       REAL(wp), DIMENSION(0:nx,nys_x:nyn_x,nzb_x:nzt_x) ::  ddd !<
-
-!
-!--    Exponent for geometric average
-       exponent = 1.0_wp / ( ny + 1.0_wp )
-
-!
-!--    Loop over all levels defined by the user
-       n = 1
-       DO WHILE ( comp_spectra_level(n) /= 999999  .AND.  n <= 100 )
-
-          k = comp_spectra_level(n)
-
-!
-!--       Calculate FFT only if the corresponding level is situated on this PE
-          IF ( k >= nzb_x  .AND.  k <= nzt_x )  THEN
-          
-             DO  j = nys_x, nyn_x
-
-                work = ddd(0:nx,j,k)
-                CALL fft_x_1d( work, 'forward' )
-
-                ddd(0,j,k) = dx * work(0)**2
-                DO  i = 1, nx/2
-                   ddd(i,j,k) = dx * ( work(i)**2 + work(nx+1-i)**2 )
-                ENDDO
-
-             ENDDO
-
-!
-!--          Local sum and geometric average of these spectra
-!--          (WARNING: no global sum should be performed, because floating
-!--          point overflow may occur)
-             DO  i = 0, nx/2
-
-                sums_spectra_l(i) = 1.0_wp
-                DO  j = nys_x, nyn_x
-                   sums_spectra_l(i) = sums_spectra_l(i) * ddd(i,j,k)**exponent
-                ENDDO
-
-             ENDDO
-          
-          ELSE
-
-             sums_spectra_l = 1.0_wp
-
-          ENDIF
-
-!
-!--       Global sum of spectra on PE0 (from where they are written on file)
-          sums_spectra(:,n) = 0.0_wp
-#if defined( __parallel )   
-          CALL MPI_BARRIER( comm2d, ierr )  ! Necessary?
-          CALL MPI_REDUCE( sums_spectra_l(0), sums_spectra(0,n), nx/2+1,       &
-                           MPI_REAL, MPI_PROD, 0, comm2d, ierr )
-#else
-          sums_spectra(:,n) = sums_spectra_l
-#endif
-!
-!--       Normalize spectra by variance
-          sum_spec_dum = SUM( sums_spectra(:,n) )
-          IF ( sum_spec_dum /= 0.0_wp )  THEN
-             sums_spectra(:,n) = sums_spectra(:,n) * var_d(k) / sum_spec_dum
-          ENDIF
-          n = n + 1
-
-       ENDDO
-       n = n - 1
-
-       IF ( myid == 0 )  THEN
-!
-!--       Sum of spectra for later averaging (see routine data_output_spectra)
-          DO  i = 1, nx/2
-             DO k = 1, n
-                spectrum_x(i,k,m) = spectrum_x(i,k,m) + sums_spectra(i,k)
-             ENDDO
-          ENDDO
-
-       ENDIF
-!
-!--    n_sp_x is needed by data_output_spectra_x
-       n_sp_x = n
-
-    END SUBROUTINE calc_spectra_x
-
-
-!------------------------------------------------------------------------------!
-! Description:
-! ------------
-!> @todo Missing subroutine description.
-!------------------------------------------------------------------------------!
-    SUBROUTINE calc_spectra_y( ddd, pr, m )
-
-       USE control_parameters,                                                 &
-           ONLY:  fft_method
-
-       USE fft_xy,                                                             &
-           ONLY:  fft_y_1d
-
-       USE grid_variables,                                                     &
-           ONLY:  dy
-
-       USE indices,                                                            &
-           ONLY:  nx, ny
-
-       USE kinds
-
-#if defined( __parallel )
-#if defined( __mpifh )
-       INCLUDE "mpif.h"
-#else
-       USE MPI
-#endif
-#endif
-       USE pegrid,                                                             &
-           ONLY:  comm2d, ierr, myid
-
-       USE statistics,                                                         &
-           ONLY:  spectrum_y, var_d
-
-       USE transpose_indices,                                                  &
-           ONLY:  nxl_yd, nxr_yd, nzb_yd, nzt_yd
-
-
-       IMPLICIT NONE
-
-       INTEGER(iwp) ::  i         !<
-       INTEGER(iwp) ::  j         !<
-       INTEGER(iwp) ::  jshape(1) !<
-       INTEGER(iwp) ::  k         !<
-       INTEGER(iwp) ::  m         !<
-       INTEGER(iwp) ::  n         !<
-       INTEGER(iwp) ::  pr        !<
-
-       REAL(wp) ::  exponent !<
-       REAL(wp) ::  sum_spec_dum !< wavenumber-integrated spectrum
-    
-       REAL(wp), DIMENSION(0:ny) ::  work !<
-    
-       REAL(wp), DIMENSION(0:ny/2) ::  sums_spectra_l !<
-    
-       REAL(wp), DIMENSION(0:ny/2,100) ::  sums_spectra !<
-    
-       REAL(wp), DIMENSION(0:ny,nxl_yd:nxr_yd,nzb_yd:nzt_yd) :: ddd !<
-
-
-!
-!--    Exponent for geometric average
-       exponent = 1.0_wp / ( nx + 1.0_wp )
-
-!
-!--    Loop over all levels defined by the user
-       n = 1
-       DO WHILE ( comp_spectra_level(n) /= 999999  .AND.  n <= 100 )
-
-          k = comp_spectra_level(n)
-
-!
-!--       Calculate FFT only if the corresponding level is situated on this PE
-          IF ( k >= nzb_yd  .AND.  k <= nzt_yd )  THEN
-          
-             DO  i = nxl_yd, nxr_yd
-
-                work = ddd(0:ny,i,k)
-                CALL fft_y_1d( work, 'forward' )
-
-                ddd(0,i,k) = dy * work(0)**2
-                DO  j = 1, ny/2
-                   ddd(j,i,k) = dy * ( work(j)**2 + work(ny+1-j)**2 )
-                ENDDO
-
-             ENDDO
-
-!
-!--          Local sum and geometric average of these spectra
-!--          (WARNING: no global sum should be performed, because floating
-!--          point overflow may occur)
-             DO  j = 0, ny/2
-
-                sums_spectra_l(j) = 1.0_wp
-                DO  i = nxl_yd, nxr_yd
-                   sums_spectra_l(j) = sums_spectra_l(j) * ddd(j,i,k)**exponent
-                ENDDO
-
-             ENDDO
-          
-          ELSE
-
-             sums_spectra_l = 1.0_wp
-
-          ENDIF
-
-!
-!--       Global sum of spectra on PE0 (from where they are written on file)
-          sums_spectra(:,n) = 0.0_wp
-#if defined( __parallel )   
-          CALL MPI_BARRIER( comm2d, ierr )  ! Necessary?
-          CALL MPI_REDUCE( sums_spectra_l(0), sums_spectra(0,n), ny/2+1,       &
-                           MPI_REAL, MPI_PROD, 0, comm2d, ierr )
-#else
-          sums_spectra(:,n) = sums_spectra_l
-#endif
-!
-!--       Normalize spectra by variance
-          sum_spec_dum = SUM( sums_spectra(:,n) )
-          IF ( SUM(sums_spectra(:,n)) /= 0.0_wp )  THEN
-             sums_spectra(:,n) = sums_spectra(:,n) *                           &
-                                 var_d(k) / SUM(sums_spectra(:,n))
-          ENDIF
-          n = n + 1
-
-       ENDDO
-       n = n - 1
-
-
-       IF ( myid == 0 )  THEN
-!
-!--       Sum of spectra for later averaging (see routine data_output_spectra)
-          DO  j = 1, ny/2
-             DO k = 1, n
-                spectrum_y(j,k,m) = spectrum_y(j,k,m) + sums_spectra(j,k)
-             ENDDO
-          ENDDO
-
-       ENDIF
-
-!
-!--    n_sp_y is needed by data_output_spectra_y
-       n_sp_y = n
-
-    END SUBROUTINE calc_spectra_y
-
- END MODULE spectrum
Index: /palm/trunk/SOURCE/time_integration.f90
===================================================================
--- /palm/trunk/SOURCE/time_integration.f90	(revision 1832)
+++ /palm/trunk/SOURCE/time_integration.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! ------------------
-!
+! spectrum renamed spectra_mod, spectra related variables moved to spectra_mod
 !
 ! Former revisions:
@@ -207,6 +207,6 @@
     USE control_parameters,                                                    &
         ONLY:  advected_distance_x, advected_distance_y, average_count_3d,     &
-               average_count_sp, averaging_interval, averaging_interval_pr,    &
-               averaging_interval_sp, bc_lr_cyc, bc_ns_cyc, bc_pt_t_val,       &
+               averaging_interval, averaging_interval_pr,                      &
+               bc_lr_cyc, bc_ns_cyc, bc_pt_t_val,                              &
                bc_q_t_val, call_psolver_at_all_substeps, cloud_droplets,       &
                cloud_physics, constant_flux_layer, constant_heatflux,          &
@@ -217,5 +217,5 @@
                dt_coupling, dt_data_output_av, dt_disturb, dt_do2d_xy,         &
                dt_do2d_xz, dt_do2d_yz, dt_do3d, dt_domask,dt_dopts, dt_dopr,   &
-               dt_dopr_listing, dt_dosp, dt_dots, dt_dvrp, dt_run_control,     &
+               dt_dopr_listing, dt_dots, dt_dvrp, dt_run_control,              &
                end_time, first_call_lpm, galilei_transformation, humidity,     &
                intermediate_timestep_count,                                    &
@@ -229,5 +229,5 @@
                skip_time_do2d_xy, skip_time_do2d_xz, skip_time_do2d_yz,        &
                skip_time_do3d, skip_time_domask, skip_time_dopr,               &
-               skip_time_dosp, skip_time_data_output_av, sloping_surface,      &
+               skip_time_data_output_av, sloping_surface,                      &
                stop_dt, terminate_coupled, terminate_run, timestep_scheme,     &
                time_coupling, time_do2d_xy, time_do2d_xz, time_do2d_yz,        &
@@ -292,6 +292,7 @@
               skip_time_do_radiation, time_radiation
 
-    USE spectrum,                                                              &
-        ONLY: calc_spectra
+    USE spectra_mod,                                                           &
+        ONLY: average_count_sp, averaging_interval_sp, calc_spectra, dt_dosp,  &
+              skip_time_dosp
 
     USE statistics,                                                            &
Index: /palm/trunk/SOURCE/user_parin.f90
===================================================================
--- /palm/trunk/SOURCE/user_parin.f90	(revision 1832)
+++ /palm/trunk/SOURCE/user_parin.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! -----------------
-! 
+! required interface revision changed
 !
 ! Former revisions:
@@ -92,5 +92,5 @@
 !-- current revision does not match with previous revisions (e.g. if routines
 !-- have been added/deleted or if parameter lists in subroutines have been changed).
-    user_interface_current_revision = 'r1783'
+    user_interface_current_revision = 'r1819'
 
 !
Index: /palm/trunk/SOURCE/user_spectra.f90
===================================================================
--- /palm/trunk/SOURCE/user_spectra.f90	(revision 1832)
+++ /palm/trunk/SOURCE/user_spectra.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! -----------------
-! 
+! spectrum renamed spectra_mod
 ! 
 ! Former revisions:
@@ -59,5 +59,5 @@
     USE kinds
     
-    USE spectrum
+    USE spectra_mod
     
     USE statistics
Index: /palm/trunk/SOURCE/write_3d_binary.f90
===================================================================
--- /palm/trunk/SOURCE/write_3d_binary.f90	(revision 1832)
+++ /palm/trunk/SOURCE/write_3d_binary.f90	(revision 1833)
@@ -19,5 +19,5 @@
 ! Current revisions:
 ! -----------------
-!
+! statistics module replaced by spectra module
 ! 
 ! Former revisions:
@@ -134,5 +134,5 @@
         ONLY:  id_random_array, seq_random_array
         
-    USE statistics,                                                            &
+    USE spectra_mod,                                                           &
         ONLY:  spectrum_x, spectrum_y
 
Index: /palm/trunk/SOURCE/write_var_list.f90
===================================================================
--- /palm/trunk/SOURCE/write_var_list.f90	(revision 1832)
+++ /palm/trunk/SOURCE/write_var_list.f90	(revision 1833)
@@ -19,4 +19,5 @@
 ! Current revisions:
 ! -----------------
+! spectra_mod added
 !
 !
@@ -172,5 +173,8 @@
 
     USE radiation_model_mod,                                                   &
-        ONLY: time_radiation
+        ONLY:  time_radiation
+
+    USE spectra_mod,                                                           &
+        ONLY:  average_count_sp
 
     USE statistics,                                                            &
