Index: palm/trunk/SCRIPTS/NCL/.ncl_preferences
===================================================================
--- palm/trunk/SCRIPTS/NCL/.ncl_preferences	(revision 174)
+++ palm/trunk/SCRIPTS/NCL/.ncl_preferences	(revision 175)
@@ -1,77 +1,129 @@
-___Parameter list________________________________________________________________________________
-___Please write your values/strings(wihtout " ") below the corresponding parameter line__________
-_________________________________________________________________________________________________
-___parameter name____example____________________________________default value_____________used by
-_________________________________________________________________________________________________
+___Parameter list_______________________________________________________________________________________
+___Please write your values/strings(wihtout " ") below the corresponding parameter line_________________
+________________________________________________________________________________________________________
+___parameter name__________example____________________________________default value______________used by
+________________________________________________________________________________________________________
 
-___'file_in'_________[path/name.nc; i.e. ~/.../example.nc]______default = input file______[ALL]
-input file
-___'format_out'______[ncgm; ps; eps; epsi; pdf; x11]____________default = x11_____________[ALL]
+___'file_1'________________[path/name.nc; i.e. ~/.../example.nc]______default = File in_________[ALL]
+File in
+___'format_out'____________[ncgm; ps; eps; epsi; pdf; x11]____________default = x11_____________[ALL]
 x11
-___'file_out'________[path/name; i.e. ~/.../test]_______________default = test____________[ALL]
+___'file_out'______________[path/name; i.e. ~/.../test]_______________default = test____________[ALL]
 test
-___start_time_step___[0;1;...]__________________________________default = t(0)____________[ALL]
+___start_time_step_________[0;1;...]__________________________________default = t(0)____________[ALL]
 t(0)
-___end_time_step_____[1;2;...]__________________________________default = t(end)__________[ALL]
+___end_time_step___________[1;2;...]__________________________________default = t(end)__________[ALL]
 t(end)
-___no_columns________[1;2;...]__________________________________default = 1_______________[ALL]
+___no_columns______________[1;2;...]__________________________________default = 1_______________[ALL]
 1
-___no_lines__________[1;2;...]__________________________________default = 2_______________[ALL]
+___no_lines________________[1;2;...]__________________________________default = 2_______________[ALL]
 2
-___'var'_____________[,E,dt,umax,;...]__________________________default = variables_______[ALL]
+___'var'___________________[,E,dt,umax,;...]__________________________default = variables_______[ALL]
 variables
-___combine___________[0] or [1]_________________________________default = 0_______________[PR]
+___combine_________________[0] or [1]_________________________________default = 0_______________[PR]
 0
-___number_comb_______[2] or [3]_________________________________default = 0_______________[PR]
+___number_comb_____________[2] or [3]_________________________________default = 0_______________[PR]
 0
-___'c_var'___________[,us2,vs2,ws2,;...]________________________default = c_variables_____[PR]
+___'c_var'_________________[,us2,vs2,ws2,;...]________________________default = c_variables_____[PR]
 c_variables
-___dash______________[0] or [1]_________________________________default = 0_______________[PR]
+___dash____________________[0] or [1]_________________________________default = 0_______________[PR+SP]
 0
-___black_____________[0] or [1]_________________________________default = 0_______________[PR]
+___black___________________[0] or [1]_________________________________default = 0_______________[PR+SP]
 0
-___min_z_____________[0;10;...]_________________________________default = 0_______________[PR]
+___min_z___________________[0;10;...]_________________________________default = 0_______________[PR]
 0
-___max_z_____________[100;1000;...]_____________________________default = max(z_u)________[PR] 
+___max_z___________________[100;1000;...]_____________________________default = max(z_u)________[PR] 
 max(z_u)
-___over______________[0] or [1]_________________________________default = 0_______________[PR+TS]
+___over____________________[0] or [1]_________________________________default = 0_______________[PR+TS]
 0
-___vector____________[0] or [1]_________________________________default = 0_______________[CS]
+___vector__________________[0] or [1]_________________________________default = 0_______________[CS]
 0
-___'vec1'____________[,u,;,v,;...]______________________________default = vec1____________[CS]
+___'vec1'__________________[,u,;,v,;...]______________________________default = vec1____________[CS]
 vec1
-___'vec2'____________[,u,;,v,;...]______________________________default = vec2____________[CS]
+___'vec2'__________________[,u,;,v,;...]______________________________default = vec2____________[CS]
 vec2
-__'plotvec'__________[,u,;,v,;...]______________________________default = plotvec_________[CS]
+__'plotvec'________________[,u,;,v,;...]______________________________default = plotvec_________[CS]
 plotvec
-___ref_mag___________[0.1;0.02;...]_____________________________default = 0.05____________[CS]
+___ref_mag_________________[0.1;0.02;...]_____________________________default = 0.05____________[CS]
 0.05
-___'mode'____________[Fill], [Line] or [Both]___________________default = Fill____________[CS]
+___'mode'__________________[Fill], [Line] or [Both]___________________default = Fill____________[CS]
 Fill
-___'sort'____________[time] or [layer]__________________________default = time____________[CS]
-time
-___'fill_mode'_______[AreaFill], [RasterFill] or [CellFill]_____default = AreaFill________[CS]
+___'sort'__________________[time] or [layer]__________________________default = layer___________[CS+SP]
+layer
+___'fill_mode'_____________[AreaFill], [RasterFill] or [CellFill]_____default = AreaFill________[CS]
 AreaFill
-___shape_____________[0] or [1]_________________________________default = 1_______________[CS]
+___shape___________________[0] or [1]_________________________________default = 1_______________[CS]
 1
-___xyc_______________[0] or [1]_________________________________default = 0_______________[CS]
+___xyc_____________________[0] or [1]_________________________________default = 0_______________[CS]
 0
-___xzc_______________[0] or [1]_________________________________default = 0_______________[CS]
+___xzc_____________________[0] or [1]_________________________________default = 0_______________[CS]
 0
-___yzc_______________[0] or [1]_________________________________default = 0_______________[CS]
+___yzc_____________________[0] or [1]_________________________________default = 0_______________[CS]
 0
-___xs________________[0,1,...]__________________________________default = x0______________[CS+PR]
+___xs______________________[0,1,...]__________________________________default = x0______________[CS+PR]
 x0
-___xe________________[1;2;...]__________________________________default = xdim____________[CS+PR]
+___xe______________________[1;2;...]__________________________________default = xdim____________[CS+PR]
 xdim
-___ys________________[0;1;...]__________________________________default = y0______________[CS]
+___ys______________________[0;1;...]__________________________________default = y0______________[CS]
 y0
-___ye________________[1;2;...]__________________________________default = ydim____________[CS]
+___ye______________________[1;2;...]__________________________________default = ydim____________[CS]
 ydim
-___zs________________[0;1;...]__________________________________default = z0______________[CS]
+___zs______________________[0;1;...]__________________________________default = z0______________[CS]
 z0
-___ze________________[1;2;...]__________________________________default = zdim____________[CS]
+___ze______________________[1;2;...]__________________________________default = zdim____________[CS]
 zdim
-___prof3d____________[0] or [1]_________________________________default = 0_______________[PR]
+___prof3d__________________[0] or [1]_________________________________default = 0_______________[PR]
 0
-___end of file ".ncl_preferences"_____________________________________________________________
+___logy____________________[0] or [1]_________________________________default = 0_______________[PR+SP]
+0
+___norm____________________[1;5.4;233;...]____________________________default = 1_______________[PR+SP]
+1
+___no_files________________[1;2,3;...]________________________________default = 1_______________[PR]
+1
+___'name_legend_1'_________[i.e. File 2]______________________________default = File 1__________[PR]
+File 1
+___'file_2'________________[path/name.nc; i.e. ~/.../example.2.nc]____default = File in_________[PR]
+File in
+___'name_legend_2'_________[i.e. File 2]______________________________default = File 2__________[PR]
+File 2
+___'file_3'________________[path/name.nc; i.e. ~/.../example.3.nc]____default = File in_________[PR]
+File in
+___'name_legend_3'_________[i.e. File 3]______________________________default = File 3__________[PR]
+File 3
+___'file_4'________________[path/name.nc; i.e. ~/.../example.4.nc]____default = File in_________[PR]
+File in
+___'name_legend_4'_________[i.e. File 4]______________________________default = File 4__________[PR]
+File 4
+___'file_5'________________[path/name.nc; i.e. ~/.../example.5.nc]____default = File in_________[PR]
+File in
+___'name_legend_5'_________[i.e. File 5]______________________________default = File 5__________[PR]
+File 5
+___'file_6'________________[path/name.nc; i.e. ~/.../example.6.nc]____default = File in_________[PR]
+File in
+___'name_legend_6'_________[i.e. File 6]______________________________default = File 6__________[PR]
+File 6
+___'file_7'________________[path/name.nc; i.e. ~/.../example.7.nc]____default = File in_________[PR]
+File in
+___'name_legend_7'_________[i.e. File 7]______________________________default = File 7__________[PR]
+File 7
+___'file_8'________________[path/name.nc; i.e. ~/.../example.8.nc]____default = File in_________[PR]
+File in
+___'name_legend_8'_________[i.e. File 8]______________________________default = File 8__________[PR]
+File 8
+___'file_9'________________[path/name.nc; i.e. ~/.../example.9.nc]____default = File in_________[PR]
+File in
+___'name_legend_9'_________[i.e. File 9]______________________________default = File 9__________[PR]
+File 9
+___'file_10'_______________[path/name.nc; i.e. ~/.../example.10.nc]___default = File in_________[PR]
+File in
+___'name_legend_10'________[i.e. File 10]_____________________________default = File 10_________[PR]
+File 10
+___'file_11'_______________[path/name.nc; i.e. ~/.../example.11.nc]___default = File in_________[PR]
+File in
+___'name_legend_11'________[i.e. File 11]_____________________________default = File 11_________[PR]
+File 11
+___'file_12'_______________[path/name.nc; i.e. ~/.../example.12.nc]___default = File in_________[PR]
+File in
+___'name_legend_12'________[i.e. File 12]_____________________________default = File 12_________[PR]
+File 12
+___end of file ".ncl_preferences"_______________________________________________________________________
Index: palm/trunk/SCRIPTS/NCL/cross_sections.ncl
===================================================================
--- palm/trunk/SCRIPTS/NCL/cross_sections.ncl	(revision 174)
+++ palm/trunk/SCRIPTS/NCL/cross_sections.ncl	(revision 175)
@@ -11,12 +11,17 @@
    ; ***************************************************
 
-   if (isfilepresent(".ncl_preferences")) then
-      parameter = asciiread(".ncl_preferences",77,"string")
+   if (isfilepresent("~/.ncl_preferences")) then
+      parameter = asciiread("~/.ncl_preferences",129,"string")
       delete(parameter@_FillValue)
    else
-      print(" ")
-      print("'.ncl_preferences' is not existent")
-      print(" ")
-      exit
+      if (isfilepresent("~/palm/current_version/trunk/SRIPTS/NCL/.ncl_preferences")) then
+         parameter = asciiread("~/palm/current_version/trunk/SRIPTS/NCL/.ncl_preferences",129,"string")
+         delete(parameter@_FillValue)
+      else
+         print(" ")
+         print("'.ncl_preferences' is not present in '~/palm/current_version/trunk/SRIPTS/NCL/'")
+         print(" ")
+         exit
+      end if
    end if
 
@@ -25,21 +30,42 @@
    ; ***************************************************
    
-   if ( .not. isvar("file_in") ) then     		; path+name of input file      
-      if (parameter(7) .EQ. "input file") then
-         print(" ")
-         print("Please provide input file 'file_in = ' either in prompt or parameter_list")
+   if ( .not. isvar("file_1") ) then     		  
+      if (parameter(7) .EQ. "File in") then
+         print(" ")
+         print("Please provide 1st input file 'file_1=' either in prompt or parameter_list")
          print(" ")
          exit
       else
          file_in = parameter(7)
-      end if      
-   end if
-   if ( .not. isvar("format_out") ) then   		; format of output file 
+      end if   
+   else
+      file_in = file_1   
+   end if
+   if (.not. isfilepresent(file_in)) then
+      print(" ")
+      print("Your 1st input file: '"+file_in+"' does not exist")
+      print(" ")
+      exit
+   end if 
+
+   if ( .not. isvar("format_out") ) then   		
       format_out = "x11"
       if (parameter(9) .NE. "x11") then
          format_out = parameter(9)  
+         if (format_out .NE. "x11" .AND. format_out .NE. "pdf" .AND. format_out .NE. "eps" .AND. format_out .NE. "ps" .AND. format_out .NE. "epsi" .AND. format_out .NE. "ncgm")then
+            print(" ")
+            print("Your 'format_out = "+format_out+"' is invalid and set to'x11'")
+            print(" ")
+         end if  
       end if 
-   end if
-   if ( .not. isvar("file_out") ) then    	 	; path+name of output file
+   else
+      if (format_out .NE. "x11" .AND. format_out .NE. "pdf" .AND. format_out .NE. "eps" .AND. format_out .NE. "ps" .AND. format_out .NE. "epsi" .AND. format_out .NE. "ncgm")then
+         print(" ")
+         print("Your 'format_out = "+format_out+"' is invalid and set to'x11'")
+         print(" ")
+      end if   
+   end if
+
+   if ( .not. isvar("file_out") ) then    	 	
       file_out = "test"
       if (parameter(11) .NE. "test") then
@@ -47,5 +73,5 @@
      end if     
    end if
-   if ( .not. isvar("no_columns") ) then    		; number of plots in one row
+   if ( .not. isvar("no_columns") ) then    		
       no_columns = 1
       if (parameter(17) .NE. "1") then
@@ -53,5 +79,5 @@
       end if
    end if
-   if ( .not. isvar("no_lines") ) then     		; number of plot-lines on one sheet
+   if ( .not. isvar("no_lines") ) then     		
       no_lines = 2
       if (parameter(19) .NE. "2") then
@@ -59,73 +85,145 @@
       end if
    end if 
-   if ( .not. isvar("sort") ) then   			; sort of plots
-      sort = "time"
-      if (parameter(51) .NE. "time") then
+
+   if ( .not. isvar("sort") ) then   			
+      sort = "layer"
+      if (parameter(51) .NE. "layer") then
          sort = parameter(51)  
+         if (sort .NE. "time") then
+            print(" ")
+            print("Your 'sort'= "+sort+" is invalid and set to 'layer'")
+            print(" ")
+            sort = "layer"  
+         end if
+      end if
+   else
+      if (sort .NE. "time" .OR. sort .NE. "layer")then
+         print(" ")
+         print("Your 'sort'= "+sort+" is invalid and set to 'layer'")
+         print(" ")
+         sort = "layer"    
       end if
    end if 
-   if ( .not. isvar("mode") ) then   			; pattern of contour plots
+
+   if ( .not. isvar("mode") ) then   			
       mode = "Fill"
       if (parameter(49) .NE. "Fill") then
          mode = parameter(49)  
-      end if
-   end if
-   if ( .not. isvar("fill_mode") ) then 		; pattern of filling 
+         if (mode .NE. "Line" .AND. mode .NE. "Both")then
+            print(" ")
+            print("Your 'mode'= "+mode+" is invalid and set to 'Fill'")
+            print(" ")
+            mode = "Fill"
+         end if
+      end if
+   else
+      if (mode .NE. "Line" .AND. mode .NE. "Both")then
+         print(" ")
+         print("Your 'mode'= "+mode+" is invalid and set to 'Fill'")
+         print(" ")
+         mode = "Fill"
+      end if   
+   end if
+
+   if ( .not. isvar("fill_mode") ) then 		
       fill_mode = "AreaFill"
       if (parameter(53) .NE. "AreaFill") then
-         mode = parameter(53)  
+         fill_mode = parameter(53)  
+         if (fill_mode .NE. "RasterFill" .AND. fill_mode .NE. "CellFill")then
+            print(" ")
+            print("Your 'fill_mode'= "+fill_mode+" is invalid and set to 'AreaFill'")
+            print(" ")
+            fill_mode = "AreaFill"
+         end if
       end if 
-   end if
-   if ( .not. isvar("shape") ) then   			; keeping of aspect ratio
+   else
+      if (fill_mode .NE. "RasterFill" .AND. fill_mode .NE. "CellFill")then
+         print(" ")
+         print("Your 'fill_mode'= "+fill_mode+" is invalid and set to 'AreaFill'")
+         print(" ")
+         fill_mode = "AreaFill"
+      end if
+   end if
+
+   if ( .not. isvar("shape") ) then   			
       shape = 1
       if (parameter(55) .NE. "1") then
          shape = stringtointeger(parameter(55)) 
-         if (stringtointeger(parameter(55)) .NE. 0) then
-            print(" ")
-            print("Please set 'shape' to 0 or 1")
-            print(" ")
-            exit
+         if (shape .NE. 0) then
+            print(" ")
+            print("Your 'shape'= "+shape+" is invalid and set to 1")
+            print(" ")
+            shape = 1
          end if 
-      end if  
-   end if  
-   if ( .not. isvar("var") ) then     			; set output of all variables
+      end if
+   else
+      if (shape .NE. 0) then
+         print(" ")
+         print("Your 'shape'= "+shape+" is invalid and set to 1")
+         print(" ")
+         shape = 1
+      end if      
+   end if 
+ 
+   if ( .not. isvar("var") ) then     			
       check = True
    end if
 
-   if ( .not. isvar("xyc") ) then   			; turn xy cross-section on or off
+   if ( .not. isvar("xyc") ) then   			
       xyc = 0
       if (parameter(57) .NE. "0") then
          xyc = stringtointeger(parameter(57)) 
-         if (stringtointeger(parameter(57)) .NE. 1) then
-            print(" ")
-            print("Please set 'xyc' to 0 or 1")
-            print(" ")
-            exit
+         if (xyc .NE. 1) then
+            print(" ")
+            print("Your 'xyc'= "+xyc+" is invalid and set to 0")
+            print(" ")
+            xyc = 0
          end if 
       end if
+   else
+      if (xyc .NE. 1) then
+         print(" ")
+         print("Your 'xyc'= "+xyc+" is invalid and set to 0")
+         print(" ")
+         xyc = 0
+      end if 
    end if 
-   if ( .not. isvar("xzc") ) then   			; turn xz cross-section on or off 
+   if ( .not. isvar("xzc") ) then   			
       xzc = 0
       if (parameter(59) .NE. "0") then
          xzc = stringtointeger(parameter(59)) 
-         if (stringtointeger(parameter(59)) .NE. 1) then
-            print(" ")
-            print("Please set 'xzc' to 0 or 1")
-            print(" ")
-            exit
+         if (xzc .NE. 1) then
+            print(" ")
+            print("Your 'xzc'= "+xzc+" is invalid and set to 0")
+            print(" ")
+            xzc = 0
          end if 
       end if
-   end if
-   if ( .not. isvar("yzc") ) then   			; turn yz cross-section on or off 
+   else
+      if (xzc .NE. 1) then
+         print(" ")
+         print("Your 'xzc'= "+xzc+" is invalid and set to 0")
+         print(" ")
+         xzc = 0
+      end if
+   end if
+   if ( .not. isvar("yzc") ) then   			
       yzc = 0
       if (parameter(61) .NE. "0") then
          yzc = stringtointeger(parameter(61)) 
-         if (stringtointeger(parameter(61)) .NE. 1) then
-            print(" ")
-            print("Please set 'yzc' to 0 or 1")
-            print(" ")
-            exit
+         if (yzc .NE. 1) then
+            print(" ")
+            print("Your 'yzc'= "+yzc+" is invalid and set to 0")
+            print(" ")
+            yzc = 0
          end if 
       end if
+   else
+      if (yzc .NE. 1) then
+         print(" ")
+         print("Your 'yzc'= "+yzc+" is invalid and set to 0")
+         print(" ")
+         yzc = 0
+      end if 
    end if
 
@@ -160,5 +258,5 @@
       end if
    end if
-   if ( .not. isvar("vector") ) then   			; sort of plots
+   if ( .not. isvar("vector") ) then   		
       vector = 0
       if (parameter(39) .NE. "0") then
@@ -172,5 +270,5 @@
       end if
    end if
-   if ( .not. isvar("ref_mag") ) then   			; sort of plots
+   if ( .not. isvar("ref_mag") ) then   			
       ref_mag = 0.05
       if (parameter(47) .NE. "0.05") then
@@ -272,5 +370,5 @@
    cs_res@gsnYAxisIrregular2Linear = True  
 
-   if( shape .eq. 1 ) then          			; keep the shape
+   if( shape .eq. 1 ) then          			
       cs_res@gsnShape = True
    end if
@@ -296,5 +394,5 @@
    cs_resP = True
    cs_resP@txString               = f@title
-   cs_resP@txFuncCode             = "~"         	; necessary for title
+   cs_resP@txFuncCode             = "~"         	
    cs_resP@txFontHeightF          = .02
    
@@ -304,5 +402,5 @@
       cs_res@gsnSpreadColors      = True
       cs_res@cnFillMode           = fill_mode
-      cs_res@lbOrientation        = "Vertical"  	; vertical label bar
+      cs_res@lbOrientation        = "Vertical"  	
       cs_res@cnLinesOn            = False	
       cs_res@cnLineLabelsOn       = False
@@ -313,5 +411,5 @@
       cs_res@gsnSpreadColors      = True
       cs_res@cnFillMode           = fill_mode
-      cs_res@lbOrientation        = "Vertical"  	; vertical label bar
+      cs_res@lbOrientation        = "Vertical"  
       cs_res@cnLinesOn            = True
       cs_res@cnLineLabelsOn       = True
@@ -450,6 +548,6 @@
 
    print(" ")
-   print("Output from "+t_all(start_time_step)/3600+"h = "+t_all(start_time_step)+"s => index = "+start_time_step)
-   print("       till "+t_all(end_time_step)/3600+"h = "+t_all(end_time_step)+"s => index = "+end_time_step)
+   print("Output of time steps from "+t_all(start_time_step)/3600+" h = "+t_all(start_time_step)+" s => index = "+start_time_step)
+   print("                     till "+t_all(end_time_step)/3600+" h = "+t_all(end_time_step)+" s => index = "+end_time_step)
    print(" ")
   
@@ -735,5 +833,5 @@
    end if
   
-   if ( .not. isvar("zs") ) then     			; start of z-coordinate range
+   if ( .not. isvar("zs") ) then     			
       zs = 0
       if (parameter(71) .NE. "z0") then
@@ -1283,9 +1381,7 @@
    end do
 
-   if (no_var .EQ. 0) then
+   if (n .EQ. 0) then
       print(" ")
-      print("Please select a variable 'var=' or use the default value")
-      print(" ")
-      print("Your selection '"+var+"' does not exist on the input file")
+      print("The variables 'var=Â°"+var+"Â°' do not exist on your input file")
       print(" ")
       exit
Index: palm/trunk/SCRIPTS/NCL/profiles.ncl
===================================================================
--- palm/trunk/SCRIPTS/NCL/profiles.ncl	(revision 174)
+++ palm/trunk/SCRIPTS/NCL/profiles.ncl	(revision 175)
@@ -9,12 +9,17 @@
    ; ***************************************************
 
-   if (isfilepresent(".ncl_preferences")) then
-      parameter = asciiread(".ncl_preferences",77,"string")
+   if (isfilepresent("~/.ncl_preferences")) then
+      parameter = asciiread("~/.ncl_preferences",129,"string")
       delete(parameter@_FillValue)
    else
-      print(" ")
-      print("'.ncl_preferences' is not existent")
-      print(" ")
-      exit
+      if (isfilepresent("~/palm/current_version/trunk/SRIPTS/NCL/.ncl_preferences")) then
+         parameter = asciiread("~/palm/current_version/trunk/SRIPTS/NCL/.ncl_preferences",129,"string")
+         delete(parameter@_FillValue)
+      else
+         print(" ")
+         print("'.ncl_preferences' is not present in '~/palm/current_version/trunk/SRIPTS/NCL/'")
+         print(" ")
+         exit
+      end if
    end if
  
@@ -22,22 +27,275 @@
    ; set up default parameter values and strings if not assigned in prompt or parameter list
    ; ***************************************************
+  
+   if ( .not. isvar("no_files") ) then 
+      no_files = 1  
+      if (parameter(81) .NE. "1") then
+         no_files = stringtointeger(parameter(81))
+         if (no_files .LT. 1 .OR. no_files .GT. 12) then
+            print(" ")
+            print("Please assign 'no_files' between 1 and 12")  
+            print(" ")
+            exit
+         end if
+      end if
+   else
+      if (no_files .LT. 1 .OR. no_files .GT. 12) then
+         print(" ")
+         print("Please assign 'no_files' between 1 and 12")  
+         print(" ")
+         exit
+      end if
+   end if
+
+   file_in = new(no_files,string)
  
-   if ( .not. isvar("file_in") ) then     		; path+name of input file      
-      if (parameter(7) .EQ. "input file") then
-         print(" ")
-         print("Please provide input file 'file_in = ' either in prompt or parameter_list")
+   if ( .not. isvar("file_1") ) then     		     
+      if (parameter(7) .EQ. "File in") then
+         print(" ")
+         print("Please provide 1st input file 'file_1=' either in prompt or parameter_list")
          print(" ")
          exit
       else
-         file_in = parameter(7)
-      end if      
+         file_in(0) = parameter(7)
+      end if   
+   else
+      file_in(0) = file_1   
    end if
-   if ( .not. isvar("format_out") ) then   		; format of output file 
+   if (.not. isfilepresent(file_in(0))) then
+      print(" ")
+      print("Your 1st input file: '"+file_in(0)+"' does not exist")
+      print(" ")
+      exit
+   end if 
+
+   if (no_files .GT. 1) then
+      if (.not. isvar("file_2"))then
+         file_in(1) = parameter(85)
+         if (file_in(1) .EQ. "File in")
+            print(" ")
+            print("Please provide 2nd input file 'file_2=' either in prompt or parameter_list or change 'no_files' to 1")
+            print(" ")
+            exit
+         end if
+      else
+         file_in(1) = file_2
+      end if
+      if (.not. isfilepresent(file_in(1))) then
+         print(" ")
+         print("Your 2nd input file: '"+file_in(1)+"' does not exist")
+         print(" ")
+         exit
+      end if
+   end if 
+   if (no_files .GT. 2) then
+      if (.not. isvar("file_3"))then
+         file_in(2) = parameter(89)
+         if (file_in(2) .EQ. "File in")
+            print(" ")
+            print("Please provide 3rd input file 'file_3=' either in prompt or parameter_list or change 'no_files' to 2")
+            print(" ")
+            exit
+         end if
+      else
+         file_in(2) = file_3
+      end if
+      if (.not. isfilepresent(file_in(2))) then
+         print(" ")
+         print("Your 3rd input file: '"+file_in(2)+"' does not exist")
+         print(" ")
+         exit
+      end if
+   end if 
+   if (no_files .GT. 3) then
+      if (.not. isvar("file_4"))then
+         file_in(3) = parameter(93)
+         if (file_in(3) .EQ. "File in")
+            print(" ")
+            print("Please provide 4th input file 'file_4=' either in prompt or parameter_list or change 'no_files' to 3")
+            print(" ")
+            exit
+         end if
+      else
+         file_in(3) = file_4
+      end if
+      if (.not. isfilepresent(file_in(3))) then
+         print(" ")
+         print("Your 4th input file: '"+file_in(3)+"' does not exist")
+         print(" ")
+         exit
+      end if 
+   end if 
+   if (no_files .GT. 4) then
+      if (.not. isvar("file_5"))then
+         file_in(4) = parameter(97)
+         if (file_in(4) .EQ. "File in")
+            print(" ")
+            print("Please provide 5th input file 'file_5=' either in prompt or parameter_list or change 'no_files' to 4")
+            print(" ")
+            exit
+         end if
+      else
+         file_in(4) = file_5
+      end if
+      if (.not. isfilepresent(file_in(4))) then
+         print(" ")
+         print("Your 5th input file: '"+file_in(4)+"' does not exist")
+         print(" ")
+         exit
+      end if
+   end if 
+   if (no_files .GT. 5) then
+      if (.not. isvar("file_6"))then
+         file_in(5) = parameter(101)
+         if (file_in(5) .EQ. "File in")
+            print(" ")
+            print("Please provide 6th input file 'file_6=' either in prompt or parameter_list or change 'no_files' to 5")
+            print(" ")
+            exit
+         end if
+      else
+         file_in(5) = file_6
+      end if
+      if (.not. isfilepresent(file_in(5))) then
+         print(" ")
+         print("Your 6th input file: '"+file_in(5)+"' does not exist")
+         print(" ")
+         exit
+      end if
+   end if 
+   if (no_files .GT. 6) then
+      if (.not. isvar("file_7"))then
+         file_in(6) = parameter(105)
+         if (file_in(6) .EQ. "File in")
+            print(" ")
+            print("Please provide 7th input file 'file_7=' either in prompt or parameter_list or change 'no_files' to 6")
+            print(" ")
+            exit
+         end if
+      else
+         file_in(6) = file_7
+      end if
+      if (.not. isfilepresent(file_in(6))) then
+         print(" ")
+         print("Your 7th input file: '"+file_in(6)+"' does not exist")
+         print(" ")
+         exit
+      end if
+   end if 
+   if (no_files .GT. 7) then
+      if (.not. isvar("file_8"))then
+         file_in(7) = parameter(109)
+         if (file_in(7) .EQ. "File in")
+            print(" ")
+            print("Please provide 8th input file 'file_8=' either in prompt or parameter_list or change 'no_files' to 7")
+            print(" ")
+            exit
+         end if
+      else
+         file_in(7) = file_8
+      end if
+      if (.not. isfilepresent(file_in(7))) then
+         print(" ")
+         print("Your 8th input file: '"+file_in(7)+"' does not exist")
+         print(" ")
+         exit
+      end if
+   end if 
+   if (no_files .GT. 8) then
+      if (.not. isvar("file_9"))then
+         file_in(8) = parameter(113)
+         if (file_in(8) .EQ. "File in")
+            print(" ")
+            print("Please provide 9th input file 'file_9=' either in prompt or parameter_list or change 'no_files' to 8")
+            print(" ")
+            exit
+         end if
+      else
+         file_in(8) = file_9
+      end if 
+      if (.not. isfilepresent(file_in(8))) then
+         print(" ")
+         print("Your 9th input file: '"+file_in(8)+"' does not exist")
+         print(" ")
+         exit
+      end if
+   end if 
+   if (no_files .GT. 9) then
+      if (.not. isvar("file_10"))then
+         file_in(9) = parameter(117)
+         if (file_in(9) .EQ. "File in")
+            print(" ")
+            print("Please provide 10th input file 'file_10=' either in prompt or parameter_list or change 'no_files' to 9")
+            print(" ")
+            exit
+         end if
+      else
+         file_in(9) = file_10
+      end if
+      if (.not. isfilepresent(file_in(9))) then
+         print(" ")
+         print("Your 10th input file: '"+file_in(9)+"' does not exist")
+         print(" ")
+         exit
+      end if
+   end if 
+   if (no_files .GT. 10) then
+      if (.not. isvar("file_11"))then
+         file_in(10) = parameter(121)
+         if (file_in(10) .EQ. "File in")
+            print(" ")
+            print("Please provide 11th input file 'file_11=' either in prompt or parameter_list or change 'no_files' to 10")
+            print(" ")
+            exit
+         end if
+      else
+         file_in(10) = file_11
+      end if
+      if (.not. isfilepresent(file_in(10))) then
+         print(" ")
+         print("Your 11th input file: '"+file_in(10)+"' does not exist")
+         print(" ")
+         exit
+      end if
+   end if 
+   if (no_files .GT. 11) then
+      if (.not. isvar("file_12"))then
+         file_in(11) = parameter(125)
+         if (file_in(11) .EQ. "File in")
+            print(" ")
+            print("Please provide 12th input file 'file_12=' either in prompt or parameter_list or change 'no_files' to 11")
+            print(" ")
+            exit
+         end if
+      else
+         file_in(11) = file_12
+      end if
+      if (.not. isfilepresent(file_in(11))) then
+         print(" ")
+         print("Your 12th input file: '"+file_in(11)+"' does not exist")
+         print(" ")
+         exit
+      end if
+   end if 
+
+   if ( .not. isvar("format_out") ) then   		
       format_out = "x11"
       if (parameter(9) .NE. "x11") then
          format_out = parameter(9)  
+         if (format_out .NE. "x11" .OR. format_out .NE. "pdf" .OR. format_out .NE. "eps" .OR. format_out .NE. "ps" .OR. format_out .NE. "epsi" .OR. format_out .NE. "ncgm")then
+            print(" ")
+            print("Your 'format_out = "+format_out+"' is invalid and set to'x11'")
+            print(" ")
+         end if  
       end if 
+   else
+      if (format_out .NE. "x11" .OR. format_out .NE. "pdf" .OR. format_out .NE. "eps" .OR. format_out .NE. "ps" .OR. format_out .NE. "epsi" .OR. format_out .NE. "ncgm")then
+         print(" ")
+         print("Your 'format_out = "+format_out+"' is invalid and set to'x11'")
+         print(" ")
+      end if  
    end if
-   if ( .not. isvar("file_out") ) then    	 	; path+name of output file
+
+   if ( .not. isvar("file_out") ) then    	 	
       file_out = "test"
       if (parameter(11) .NE. "test_ts") then
@@ -45,5 +303,5 @@
      end if     
    end if
-   if ( .not. isvar("no_columns") ) then    		; number of plots in one row
+   if ( .not. isvar("no_columns") ) then    		
       no_columns = 1
       if (parameter(17) .NE. "1") then
@@ -63,10 +321,31 @@
          if (stringtointeger(parameter(23)) .NE. 1) then
             print(" ")
-            print("Please set 'combine' to 0 or 1")
-            print(" ")
-            exit
-         end if
-      end if
+            print("Your 'combine'= "+combine+" is invalid and set to 0")
+            print(" ")
+            combine = 0
+         else
+            if (no_files .GT. 1) then
+               print(" ")
+               print("If you have more than one input file you cannot combine variables: combine is set to 0")
+               print(" ")
+               combine = 0
+            end if 
+         end if
+      end if
+   else
+      if (combine .NE. 0 .AND. combine .NE. 1)then
+         print(" ")
+         print("Your 'combine'= "+combine+" is invalid and set to 0")
+         print(" ")
+         combine = 0
+      end if
+      if (no_files .GT. 1) then
+         combine = 0
+         print(" ")
+         print("If you have more than one input file you cannot combine variables: combine is set to 0")
+         print(" ")
+      end if       
    end if
+
    if (combine .EQ. 1) then
       if( .not. isvar("c_var") ) then
@@ -81,4 +360,5 @@
       end if
    end if
+
    if ( .not. isvar("black") ) then  		 	; color of lines
       black = 0
@@ -87,11 +367,19 @@
          if (stringtointeger(parameter(31)) .NE. 1) then
             print(" ")
-            print("Please set 'black' to 0 or 1")
-            print(" ")
-            exit
-         end if
+            print("Your 'black'= "+black+" is invalid and set to 0")
+            print(" ")
+            black = 0
+         end if
+      end if
+   else
+      if (black .NE. 0 .AND. black .NE. 1)then
+         print(" ")
+         print("Your 'black'= "+black+" is invalid and set to 0")
+         print(" ")
+         black = 0
       end if
    end if
-   if ( .not. isvar("dash") ) then   		 	; pattern of lines
+
+   if ( .not. isvar("dash") ) then   		 	
       dash = 0
       if (parameter(29) .NE. "0") then
@@ -99,10 +387,18 @@
          if (stringtointeger(parameter(29)) .NE. 1) then
             print(" ")
-            print("Please set 'dash' to 0 or 1")
-            print(" ")
-            exit
+            print("Your 'dash'= "+dash+" is invalid and set to 0")
+            print(" ")
+            dash = 0
          end if  
       end if
+   else
+      if (dash .NE. 0 .AND. dash .NE. 1)then
+         print(" ")
+         print("Your 'dash'= "+dash+" is invalid and set to 0")
+         print(" ")
+         dash = 0
+      end if
    end if
+
    if ( .not. isvar("over") ) then     		        ; switches overlaying plots on
       over = 0
@@ -111,8 +407,22 @@
          if (stringtointeger(parameter(37)) .NE. 1) then
             print(" ")
-            print("Please set 'over' to 0 or 1")
-            print(" ")
-            exit
+            print("Your 'over'= "+over+" is invalid and set to 0")
+            print(" ")
+            over = 0
+         else
+            if (no_files .GT. 1) then
+               over = 0
+               print(" ")
+               print("If you have more than one input file you cannot overlay variables: over is set to 0")
+               print(" ")
+            end if
          end if   
+      end if
+   else
+      if (over .NE. 0 .AND. over .NE. 1)then
+         print(" ")
+         print("Your 'over'= "+over+" is invalid and set to 0")
+         print(" ")
+         over = 0
       end if
    end if
@@ -129,17 +439,89 @@
          end if   
       end if
+   else
+      if (prof3d .NE. 0 .AND. prof3d .NE. 1)then
+         print(" ")
+         print("Your 'prof3d'= "+prof3d+" is invalid; please set 'prof3d' to 0 or 1")
+         print(" ")
+         exit
+      end if
    end if
-  
+
+   if (.not. isvar("logy"))then
+      logy = 0
+      if (stringtointeger(parameter(77)) .NE. 0) then
+         logy = stringtointeger(parameter(77))
+         if (stringtointeger(parameter(77)) .NE. 1) then
+            print(" ")
+            print("Your 'logy'= "+logy+" is invalid and set to 0")
+            print(" ")
+            logy = 0
+         end if   
+      end if
+   else
+      if (logy .NE. 0 .AND. logy .NE. 1)then
+         print(" ")
+         print("Your 'logy'= "+logy+" is invalid and set to 0")
+         print(" ")
+         logy = 0
+      end if    
+   end if
+
+   if ( .not. isvar("norm") ) then     		
+      norm = 1.0
+      if (parameter(79) .NE. "1") then
+         norm = stringtofloat(parameter(79)) 
+         if (stringtofloat(parameter(79)) .EQ. 0) then 
+            print(" ")
+            print("You cannot normalise with 0, 'norm' is set to 1")
+            print(" ")
+            norm = 1.0
+         end if
+      end if
+   else
+      if (norm .EQ. 0) then 
+         print(" ")
+         print("You cannot normalise with 0, 'norm' is set to 1")
+         print(" ")
+         norm = 1.0
+      end if
+   end if 
+
    ; ***************************************************
    ; open input file
    ; ***************************************************
 
-   f=addfile( file_in,"r")
-
+   do nof=0,no_files-1
+
+   f=addfile( file_in(nof),"r")
+   
    vNam = getfilevarnames(f) 
+   if (nof .EQ. 0)then
+      vNam0=vNam
+   end if
+   if (nof .NE. 0)then
+      if (dim0 .NE. dim)then
+         print(" ")
+         print("Your Input files does not contain the same variables")
+         print(" ")
+         exit
+      else
+         do i=0,dim0-1
+            if (vNam0(i) .NE. vNam(i))then
+               print(" ")
+               print("Your Input files does not contain the same variables")
+               print(" ")
+               exit
+            end if
+         end do
+      end if
+   end if
+   nof=nof+1
    print(" ")
-   print("Variable on netCDF file: " + vNam)
+   print("Variable in input file "+nof+": " + vNam)
    print(" ")
+   nof=nof-1
    dim = dimsizes(vNam)
+   dim0=dim
   
    if (dim .EQ. 0) then
@@ -149,6 +531,19 @@
    end if
 
+   ;-----------------------------------------------------below steps only for first file -> nof=0
+   if (nof .EQ. 0) then
+
+      plot = new(dim,graphic)
+      plot_ = new(dim,graphic)
+      if (no_files .GT. 1) then
+         multi_plot = new((/no_files,dim/),graphic)
+         max_nof = new((/no_files,dim/),float)
+         min_nof = new((/no_files,dim/),float)
+         name    = new((/no_files,dim/),string)
+         unit_   = new((/no_files,dim/),string)
+      end if
+
    if (prof3d .EQ. 0) then
-
+ 
    do varn = dim-1,0,1
       if ( isStrSubset( vNam(varn), "time") .OR. isStrSubset( vNam(varn), "NORM")) then
@@ -156,5 +551,5 @@
          continue
       end if
-      if (vNam(varn) .EQ. "u" .OR. isStrSubset(vNam(varn), "u_"))
+      if (vNam(varn) .EQ. "u" .OR. isStrSubset(vNam(varn), "u_"))then
          z_u = f->$vNam(varn+1)$
          break
@@ -274,5 +669,5 @@
          continue
       end if 
-      if (vNam(varn) .EQ. "w" .OR. isStrSubset(vNam(varn), "w_"))
+      if (vNam(varn) .EQ. "w" .OR. isStrSubset(vNam(varn), "w_"))then
          z_w = f->$vNam(varn+1)$
          break
@@ -452,4 +847,5 @@
 
    if (zu .EQ. 0 .AND. zw .EQ. 0) then
+      co=0
       do varn=0,dim-1      
          if ( isStrSubset( vNam(varn), "time") .OR. isStrSubset( vNam(varn), "NORM")) then
@@ -466,5 +862,6 @@
             end if
          end if
-         if (check)
+         if (check)then
+            co=co+1
             z = f->$vNam(varn+1)$
             dimz = dimsizes(z)
@@ -472,5 +869,12 @@
          end if   
       end do
+      if (co .EQ. 0) then
+         print(" ")
+         print("The variables 'var=Â°"+var+"Â°' do not exist on your input file")
+         print(" ")
+         exit
+      end if
    end if
+
    if (isvar("z_u") ) then 
       dimz  = dimsizes(z_u)
@@ -484,9 +888,9 @@
 
       do varn = dim-1,0,1
-         if (vNam(varn) .EQ. "zu_3d")
+         if (vNam(varn) .EQ. "zu_3d")then
             z_u = f->zu_3d  
             dimz  = dimsizes(z_u)         
          else 
-            if (vNam(varn) .EQ. "zw_3d")
+            if (vNam(varn) .EQ. "zw_3d")then
                z_w = f->zw_3d
                dimz  = dimsizes(z_w)
@@ -501,12 +905,13 @@
    delta_t=t_all(nt-1)/nt
 
+
    ; ****************************************************	
    ; start of time step and different types of mistakes that could be done
    ; ****************************************************
-
+   
    if ( .not. isvar("start_time_step") ) then   	
       start_time_step=t_all(0)/3600
       if (parameter(13) .NE. "t(0)") then
-         if (stringtodouble(parameter(13)) .GT. t_all(nt-1)/3600)
+         if (stringtodouble(parameter(13)) .GT. t_all(nt-1)/3600)then
             print(" ")
             print("'start_time_step' = "+ parameter(13) +"h is greater than last time step = " + t_all(nt-1)+"s = "+t_all(nt-1)/3600+"h")
@@ -516,5 +921,5 @@
             exit
          end if
-         if (stringtofloat(parameter(13)) .LT. t_all(0)/3600)
+         if (stringtofloat(parameter(13)) .LT. t_all(0)/3600)then
             print(" ")
             print("'start_time_step' = "+ parameter(13) +"h is lower than first time step = " + t_all(0)+"s = "+t_all(0)/3600+"h")           
@@ -524,5 +929,5 @@
       end if
    else 
-      if (start_time_step .GT. t_all(nt-1)/3600)
+      if (start_time_step .GT. t_all(nt-1)/3600)then
          print(" ")
          print("'start_time_step' = "+ start_time_step +"h is greater than last time step = " + t_all(nt-1)+"s = "+t_all(nt-1)/3600+"h")
@@ -532,5 +937,5 @@
          exit
       end if
-      if (start_time_step .LT. t_all(0)/3600)
+      if (start_time_step .LT. t_all(0)/3600)then
          print(" ")
          print("'start_time_step' = "+ start_time_step +"h is lower than first time step = " + t_all(0)+"s = "+t_all(0)/3600+"h")
@@ -540,7 +945,5 @@
    end if
    start_time_step = start_time_step*3600
-   ;if (start_time_step .GE. t_all(0) .AND. start_time_step .LT. t_all(1)+delta_t/2)then
-   ;   st=0
-   ;end if
+
    do i=0,nt-1     
       if (start_time_step .GE. t_all(i)-delta_t/2 .AND. start_time_step .LT. t_all(i)+delta_t/2)then
@@ -557,5 +960,5 @@
       end_time_step = t_all(nt-1)/3600
       if (parameter(15) .NE. "t(end)") then
-         if (stringtodouble(parameter(15)) .GT. t_all(nt-1)/3600)
+         if (stringtodouble(parameter(15)) .GT. t_all(nt-1)/3600)then
             print(" ")
             print("'end_time_step' = "+ parameter(15) +"h is greater than last time step = " + t_all(nt-1)+"s = "+t_all(nt-1)/3600+"h")
@@ -565,5 +968,5 @@
             exit
          end if
-         if (stringtodouble(parameter(15)) .LT. start_time_step/3600)
+         if (stringtodouble(parameter(15)) .LT. start_time_step/3600)then
             print(" ")
             print("'end_time_step' = "+ parameter(15) +"h is lower than 'start_time_step' = "+start_time_step/3600+"h")
@@ -576,5 +979,5 @@
       end if    
    else
-      if (end_time_step .GT. t_all(nt-1)/3600)
+      if (end_time_step .GT. t_all(nt-1)/3600)then
          print(" ")
          print("'end_time_step' = "+ end_time_step +"h is greater than last time step = " + t_all(nt-1)+"s = "+t_all(nt-1)/3600+"h")
@@ -584,5 +987,5 @@
          exit
       end if
-      if (end_time_step .LT. start_time_step/3600)
+      if (end_time_step .LT. start_time_step/3600)then
          print(" ")
          print("'end_time_step' = "+ end_time_step +"h is lower than 'start_time_step' = "+start_time_step+"h")
@@ -594,7 +997,5 @@
    end if
    end_time_step = end_time_step*3600
-   ;if (end_time_step .GE. t_all(0) .AND. end_time_step .LT. t_all(1)+delta_t/2)then
-   ;   et=1
-   ;end if
+
    do i=0,nt-1     
       if (end_time_step .GE. t_all(i)-delta_t/2 .AND. end_time_step .LT. t_all(i)+delta_t/2)then
@@ -608,4 +1009,9 @@
    delete(end_time_step)
    end_time_step=round(et,3)
+
+   print(" ")
+   print("Output of time steps from "+t_all(start_time_step)/3600+" h = "+t_all(start_time_step)+" s => index = "+start_time_step)
+   print("                     till "+t_all(end_time_step)/3600+" h = "+t_all(end_time_step)+" s => index = "+end_time_step)
+   print(" ")
 
    ; ****************************************************
@@ -615,12 +1021,12 @@
    legend_label=new(nt,double)
    do p=start_time_step,end_time_step
-      legend_label(p-start_time_step)=t_all(p)
+      if (t_all(p)/3600 .LT. 1) then
+         legend_label(p-start_time_step)=decimalPlaces(t_all(p)/3600,2,True)
+      else
+         legend_label(p-start_time_step)=decimalPlaces(t_all(p)/3600,0,True)
+      end if
    end do
    
    np = end_time_step-start_time_step+1  
-   if ( black .eq. 0 ) then
-      res=True      
-      res@xyLineColors = ispan(2,237,235/np)
-   end if
 
    ; ***************************************************
@@ -628,4 +1034,5 @@
    ; ***************************************************
 
+   res                         = True
    res@gsnDraw                 = False
    res@gsnFrame                = False
@@ -648,14 +1055,32 @@
    res@pmLegendOrthogonalPosF  = -1.0
    res@pmLegendWidthF          = 0.12
-   res@pmLegendHeightF         = 0.3
+   res@pmLegendHeightF         = 0.04*(end_time_step-start_time_step+1)
    res@lgLabelFontHeightF     = .02 
    res@txFontHeightF      = 0.02
    res@tiXAxisFontHeightF = 0.02
    res@tiYAxisFontHeightF = 0.02
-   res@tiXAxisString      = " "    
-   res@tiYAxisString      = "Height [z]"
-   
+   res@tiXAxisString      = " " 
+   if ( black .eq. 0 ) then  
+      res@xyLineColors = ispan(2,237,235/np)
+   end if
+   if (norm .EQ. 1)then
+      res@tiYAxisString      = "Height [m]"
+   else   
+      res@tiYAxisString      = "Height / "+norm+" [m]"
+   end if
+    
+   if (logy .EQ. 1) then
+      res@trYLog = True 
+   end if
+
    if ( dash .eq. 0 ) then
-      res@xyMonoDashPattern       = True  
+      res@xyMonoDashPattern       = True 
+   else
+      if (no_files .GT. 1)
+         res@xyMonoDashPattern       = True   
+         print(" ") 
+         print("If you use more than one file, patterns for different timesteps cannot be used")
+         print(" ") 
+      end if 
    end if 
 
@@ -670,7 +1095,4 @@
    ; ***************************************************
 
-   plot = new(dim,graphic)
-   plot_ = new(dim,graphic)
-
    if (combine .EQ. 1) then
       if ( .not. isvar("number_comb") ) then        	
@@ -702,5 +1124,198 @@
    wks=gsn_open_wks(format_out,file_out)
    gsn_define_colormap(wks,"rainbow+white") 
- 
+
+   ; ***************************************************
+   ; set up minimum and maximum height
+   ; *************************************************** 
+
+   if (logy .EQ. 1)then
+      if (.not. isvar("min_z"))then
+         if (isvar("z_u"))then
+            min_z=z_u(1)
+         else
+            min_z=z_w(1)
+         end if
+         if (stringtofloat(parameter(33)) .NE. z_u(1) .OR. stringtofloat(parameter(33)) .NE. z_w(1)) then
+            if (stringtofloat(parameter(33)) .EQ. 0) then
+               if (isvar("z_u"))then
+                  min_z=z_u(1)
+               else
+                  min_z=z_w(1)
+               end if  
+            else
+               if (stringtofloat(parameter(33)) .GE. max(z_u) ) then
+                  print(" ")
+                  print("Minimum of height ('min_z'="+stringtofloat(parameter(33))+") is greater than available heights (="+max(z_u)+")")
+                  print(" ")
+                  exit
+               end if
+               if (stringtofloat(parameter(33)) .LT. z_u(1)) then
+                  print(" ")
+                  print("Begin height 'min_z' at least at level k=1 (="+z_u(1)+"m) due to the logarithmic scale of the y-axis")
+                  print(" ")
+                  exit
+               end if 
+               min_z=stringtofloat(parameter(33))
+            end if
+         end if
+      else
+         if (min_z .GE. max(z_u) ) then
+            print(" ")
+            print("Minimum of height ('min_z'="+min_z+") is greater than available heights (="+max(z_u)+")")
+            print(" ")
+            exit
+         end if
+         if (min_z .LT. z_u(1) ) then
+            print(" ")
+            print("Begin height 'min_z' at least at level k=1 (="+z_u(1)+"m) due to the logarithmic scale of the y-axis")
+            print(" ")
+            exit
+         end if
+      end if
+   else
+      if (.not. isvar("min_z"))then
+         min_z=0.0     
+         if (stringtofloat(parameter(33)) .NE. 0) then
+            if (stringtofloat(parameter(33)) .GE. max(z_u) ) then
+               print(" ")
+               print("Minimum of height ('min_z'="+stringtofloat(parameter(33))+") is greater than available heights (="+max(z_u)+")")
+               print(" ")
+               exit
+            end if
+            if (stringtofloat(parameter(33)) .LT. 0 ) then
+               print(" ")
+               print("Begin minimum of height 'min_z' with 0")
+               print(" ")
+               exit
+            end if 
+            min_z=stringtofloat(parameter(33))
+         end if
+      else
+         if (min_z .GE. max(z_u) ) then
+            print(" ")
+            print("Minimum of height ('min_z'="+min_z+") is greater than available heights (="+max(z_u)+")")
+            print(" ")
+            exit
+         end if
+         if (min_z .LT. 0 ) then
+            print(" ")
+            print("Begin minimum of height 'min_z' with 0")
+            print(" ")
+            exit
+         end if
+      end if
+   end if
+
+   if (.not. isvar("max_z"))then
+      max_z=max(z_u)
+      if ((parameter(35)) .NE. "max(z_u)") then
+         if (stringtofloat(parameter(35)) .GE. max(z_u) ) then
+            print(" ")
+            print("Maximum of height ('max_z'="+parameter(35)+") is greater than available heights (="+max(z_u)+")")
+            print(" ")
+            exit
+         end if
+         if (stringtofloat(parameter(35)) .LE. 0 ) then
+            print(" ")
+            print("Maximum of height 'max_z' should be at least 1")
+            print(" ")
+            exit
+         end if 
+         max_z=stringtofloat(parameter(35))
+      end if
+   else
+      if (max_z .GE. max(z_u) ) then
+         print(" ")
+         print("Maximum of height ('max_z'="+max_z+") is greater than available heights (="+max(z_u)+")")
+         print(" ")
+         exit
+      end if
+      if (max_z .LE. 0 ) then
+         print(" ")
+         print("Maximum of height 'max_z' should be at least 1")
+         print(" ")
+         exit
+      end if
+   end if
+
+   min_z=min_z/norm
+   max_z=max_z/norm
+
+   ; ***************************************************
+   ; read data and create plots
+   ; ***************************************************
+      
+   do ti = start_time_step, end_time_step
+      if( t_all(ti) .lt. 10^36) then
+         start_time_step = ti
+         break
+      end if
+   end do  
+   
+   if (logy .EQ. 1) then
+      data   = new((/dim,(end_time_step-start_time_step)+1,dimz-1/),float) 
+      data_0 = new((/(end_time_step-start_time_step)+1,dimz-1/),float)
+      data_0 = 0.1
+      t      = new((/(end_time_step-start_time_step)+1,dimz-1/),float) 
+      t      = 0.0
+      unit   = new(dim,string)
+      if (isvar("z_u"))then
+         if (typeof(z_u) .EQ. "double")then
+            z_v    = new((/dim,dimz/),double)
+            z_     = new((/dim,dimz-1/),double)
+         else
+            if (typeof(z_u) .EQ. "float")then 
+               z_v    = new((/dim,dimz/),float)
+               z_     = new((/dim,dimz-1/),float)
+            end if
+         end if
+      else
+         if (isvar("z_w"))then
+            if (typeof(z_w) .EQ. "double")then
+               z_v    = new((/dim,dimz/),double)
+               z_     = new((/dim,dimz-1/),double)
+            else
+               if (typeof(z_w) .EQ. "float")then 
+                  z_v    = new((/dim,dimz/),float)
+                  z_     = new((/dim,dimz-1/),float)
+               end if
+            end if
+         end if
+      end if
+   else
+      data   = new((/dim,(end_time_step-start_time_step)+1,dimz/),float) 
+      data_0 = new((/(end_time_step-start_time_step)+1,dimz/),float)
+      data_0 = 0.0
+      t      = new((/(end_time_step-start_time_step)+1,dimz/),float) 
+      t      = 0.0
+      unit   = new(dim,string)
+      if (isvar("z_u"))then
+         if (typeof(z_u) .EQ. "double")then
+            z_v    = new((/dim,dimz/),double)
+            z_     = new((/dim,dimz/),double)
+         else
+            if (typeof(z_u) .EQ. "float")then 
+               z_v    = new((/dim,dimz/),float)
+               z_     = new((/dim,dimz/),float)
+            end if
+         end if
+      else
+         if (isvar("z_w"))then
+            if (typeof(z_w) .EQ. "double")then
+               z_v    = new((/dim,dimz/),double)
+               z_     = new((/dim,dimz/),double)
+            else
+               if (typeof(z_w) .EQ. "float")then 
+                  z_v    = new((/dim,dimz/),float)
+                  z_     = new((/dim,dimz/),float)
+               end if
+            end if
+         end if
+      end if
+   end if
+
+   end if
+   ;------------------------------------------------------ above steps only for first file
+
    ; ***************************************************
    ; indicate plot number 
@@ -712,91 +1327,4 @@
       n = 0
    end if
-
-   ; ***************************************************
-   ; set up minimum and maximum height
-   ; ***************************************************
-
-   if (.not. isvar("min_z"))
-      min_z=0
-      if (stringtointeger(parameter(33)) .NE. 0) then
-         if (stringtointeger(parameter(33)) .GE. max(z_u) ) then
-            print(" ")
-            print("Minimum of height ('min_z'="+stringtointeger(parameter(33))+") is greater than available heights (="+max(z_u)+")")
-            print(" ")
-            exit
-         end if
-         if (stringtointeger(parameter(33)) .LT. 0 ) then
-            print(" ")
-            print("Begin minimum of height 'min_z' with 0")
-            print(" ")
-            exit
-         end if 
-         min_z=stringtointeger(parameter(33))
-      end if
-   else
-      if (min_z .GE. max(z_u) ) then
-         print(" ")
-         print("Minimum of height ('min_z'="+min_z+") is greater than available heights (="+max(z_u)+")")
-         print(" ")
-         exit
-      end if
-      if (min_z .LT. 0 ) then
-         print(" ")
-         print("Begin minimum of height 'min_z' with 0")
-         print(" ")
-         exit
-      end if
-   end if
-
-   if (.not. isvar("max_z"))
-      max_z=max(z_u)
-      if ((parameter(35)) .NE. "max(z_u)") then
-         if (stringtofloat(parameter(35)) .GE. max(z_u) ) then
-            print(" ")
-            print("Maximum of height ('max_z'="+parameter(35)+") is greater than available heights (="+max(z_u)+")")
-            print(" ")
-            exit
-         end if
-         if (stringtointeger(parameter(35)) .LE. 0 ) then
-            print(" ")
-            print("Maximum of height 'max_z' should be at least 1")
-            print(" ")
-            exit
-         end if 
-         max_z=stringtointeger(parameter(35))
-      end if
-   else
-      if (max_z .GE. max(z_pr) ) then
-         print(" ")
-         print("Maximum of height ('max_z'="+max_z+") is greater than available heights (="+max(z_u)+")")
-         print(" ")
-         exit
-      end if
-      if (max_z .LE. 0 ) then
-         print(" ")
-         print("Maximum of height 'max_z' should be at least 1")
-         print(" ")
-         exit
-      end if
-   end if
-
-   ; ***************************************************
-   ; read data and create plots
-   ; ***************************************************
-      
-   do ti = start_time_step, end_time_step
-      if( t_all(ti) .lt. 10^36) then
-         start_time_step = ti
-         break
-      end if
-   end do  
- 
-   data   = new((/dim,(end_time_step-start_time_step)+1,dimz/),float) 
-   data_0 = new((/(end_time_step-start_time_step)+1,dimz/),float)
-   data_0 = 0.0
-   t      = new((/(end_time_step-start_time_step)+1,dimz/),float) 
-   t      = 0.0
-   unit   = new(dim,string)
-   z_v    = new((/dim,dimz/),double)
 
    if (over .EQ. 1) then 
@@ -976,15 +1504,28 @@
       if (combine .EQ. 1) then         
          com=isStrSubset(c_var,","+vNam(varn)+"," )      
-         if (com)      
+         if (com) then      
             if (prof3d .EQ. 0) then
                temp = f->$vNam(varn)$
-               data(varn,:,:) = temp(start_time_step:end_time_step,0:dimz-1)
+               if (logy .EQ. 1) then
+                  data(varn,:,:) = temp(start_time_step:end_time_step,1:dimz-1)
+               else
+                  data(varn,:,:) = temp(start_time_step:end_time_step,0:dimz-1)
+               end if
             else
-               do i=0,dimz-1
-                  do j=start_time_step,end_time_step
-                     temp= f->$vNam(varn)$(j,i,:,:)
-                     data(varn,j-start_time_step,i) = dim_avg_Wrap(dim_avg_Wrap(temp(:,:)))
+               if (logy .EQ. 1) then
+                  do i=1,dimz-1
+                     do j=start_time_step,end_time_step
+                        temp= f->$vNam(varn)$(j,i,:,:)
+                        data(varn,j-start_time_step,i-1) = dim_avg_Wrap(dim_avg_Wrap(temp(:,:)))
+                     end do
                   end do
-               end do
+               else
+                  do i=0,dimz-1
+                     do j=start_time_step,end_time_step
+                        temp= f->$vNam(varn)$(j,i,:,:)
+                        data(varn,j-start_time_step,i) = dim_avg_Wrap(dim_avg_Wrap(temp(:,:)))
+                     end do
+                  end do
+               end if
                print(" ")
                print("Variable for combine '"+vNam(varn)+"' is read")
@@ -1006,4 +1547,11 @@
       if(check) then
 
+         if (isStrSubset(vNam(varn),"_0" ))then
+            print(" ")
+            print("If you have Outputs of statistic regions you cannot overlay variables; 'over' is set to 0")
+            print(" ")
+            over = 0
+         end if
+        
          count_var=count_var+1
 
@@ -1011,10 +1559,19 @@
             temp = f->$vNam(varn)$
          else
-            do i=0,dimz-1
-               do j=start_time_step,end_time_step
-                  temp = f->$vNam(varn)$(j,i,:,:)
-                  data(varn,j-start_time_step,i) = dim_avg_Wrap(dim_avg_Wrap(temp(:,:)))
+             if (logy .EQ. 1) then
+               do i=1,dimz-1
+                  do j=start_time_step,end_time_step
+                     temp= f->$vNam(varn)$(j,i,:,:)
+                     data(varn,j-start_time_step,i-1) = dim_avg_Wrap(dim_avg_Wrap(temp(:,:)))
+                  end do
                end do
-            end do
+            else
+               do i=0,dimz-1
+                  do j=start_time_step,end_time_step
+                     temp= f->$vNam(varn)$(j,i,:,:)
+                     data(varn,j-start_time_step,i) = dim_avg_Wrap(dim_avg_Wrap(temp(:,:)))
+                  end do
+               end do
+            end if
             print(" ")
             print("Variable '"+vNam(varn)+"' is read")
@@ -1026,21 +1583,47 @@
     
          if (prof3d .EQ. 0) then
-            z = f->$vNam(varn+1)$ 
-            unit(varn) = temp@units
-            data(varn,:,:) = temp(start_time_step:end_time_step,0:dimz-1) 
+            if (logy .EQ. 1) then
+               z = f->$vNam(varn+1)$(1:dimz-1) 
+               unit(varn) = temp@units
+               data(varn,:,:) = temp(start_time_step:end_time_step,1:dimz-1) 
+            else
+               z = f->$vNam(varn+1)$ 
+               unit(varn) = temp@units
+               data(varn,:,:) = temp(start_time_step:end_time_step,:) 
+            end if
          else 
             do i=0,b-1            
-               if (isStrSubset( a(i),"zu_3d" ))
+               if (isStrSubset( a(i),"zu_3d" ))then
                   z_v(varn,:) = z_u
-                  z = z_v(varn,:)
-               else
-                  if (isStrSubset( a(i),"zw_3d" ))
+                  if (logy .EQ. 1) then
+                     z = z_v(varn,1:dimz-1)
+                  else
+                     z = z_v(varn,:)
+                  end if 
+               else
+                  if (isStrSubset( a(i),"zw_3d" ))then
                      z_v(varn,:) = z_w
-                     z = z_v(varn,:)
+                     if (logy .EQ. 1) then
+                        z = z_v(varn,1:dimz-1)
+                     else
+                        z = z_v(varn,:)
+                     end if
                   end if                   
                end if
             end do            
          end if
+
+         if (data(varn,0,0) .GT. 1.E35) then
+            data(varn,0,:) = data(varn,1,:)
+         end if
+         if (nof .EQ. 0) then
+            z_(n,:)=z/norm
+            z    = z_(n,:)
+         else
+            z=z/norm
+         end if
+    
          if (over .EQ. 0) then  
+            res@xyDashPattern      = nof
             res@gsnLeftString      = vNam(varn)
             res@gsnRightString     = unit(varn)
@@ -1064,5 +1647,5 @@
             else 
                res@trXMaxF            = xe  
-            end if                     
+            end if         
             plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)  
          end if       
@@ -1073,5 +1656,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_u = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_u = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1095,5 +1678,5 @@
                      res@trXMaxF            = xe  
                   end if                
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res)  
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)  
                end if
             end if
@@ -1103,5 +1686,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_v = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_v = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1125,5 +1708,5 @@
                      res@trXMaxF            = xe  
                   end if                
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if  
             end if
@@ -1133,5 +1716,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 2
-                  plot_w = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_w = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1155,5 +1738,5 @@
                      res@trXMaxF            = xe  
                   end if           
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if   
             end if
@@ -1164,5 +1747,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_pt = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_pt = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1186,5 +1769,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if  
             end if
@@ -1194,5 +1777,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_vpt = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_vpt = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1216,5 +1799,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if 
             end if
@@ -1224,5 +1807,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 2
-                  plot_lpt = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_lpt = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1246,5 +1829,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if 
             end if
@@ -1255,5 +1838,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_q = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_q = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1277,5 +1860,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -1285,5 +1868,5 @@
                if (over .EQ. 1) then 
                   res@xyDashPattern  = 1
-                  plot_qv = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_qv = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1307,5 +1890,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if
             end if
@@ -1315,5 +1898,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 2
-                  plot_ql = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_ql = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1337,5 +1920,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -1346,5 +1929,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_e = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_e = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1368,5 +1951,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if 
             end if
@@ -1376,5 +1959,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_es = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_es = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1398,5 +1981,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if
             end if
@@ -1407,5 +1990,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_km = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_km = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1429,5 +2012,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if 
             end if
@@ -1437,5 +2020,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_kh = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_kh = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1459,5 +2042,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if
             end if
@@ -1468,5 +2051,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_wpup = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wpup = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1490,5 +2073,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -1498,5 +2081,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_wsus = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wsus = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1520,5 +2103,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if
             end if
@@ -1528,5 +2111,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 2 
-                  plot_wu = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wu = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1550,5 +2133,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -1559,5 +2142,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_wpvp = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wpvp = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1581,5 +2164,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -1589,5 +2172,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_wsvs = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wsvs = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1611,5 +2194,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if
             end if
@@ -1619,5 +2202,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 2 
-                  plot_wv = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wv = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1641,5 +2224,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -1650,5 +2233,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_wpptp = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wpptp = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1672,5 +2255,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -1680,5 +2263,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_wspts = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wspts = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1702,5 +2285,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if
             end if
@@ -1710,5 +2293,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 2 
-                  plot_wpt = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wpt = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1732,5 +2315,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -1741,5 +2324,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_wsptsBC = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wsptsBC = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1763,5 +2346,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if 
             end if             
@@ -1771,5 +2354,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_wptBC = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wptBC = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1793,5 +2376,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if
             end if
@@ -1802,5 +2385,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_wpvptp = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wpvptp = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1824,5 +2407,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -1832,5 +2415,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_wsvpts = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wsvpts = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1854,5 +2437,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if
             end if
@@ -1862,5 +2445,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 2 
-                  plot_wvpt = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wvpt = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1884,5 +2467,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -1893,5 +2476,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_wpqp = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wpqp = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1915,5 +2498,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -1923,5 +2506,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_wsqs = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wsqs = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1945,5 +2528,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if
             end if
@@ -1953,5 +2536,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 2
-                  plot_wq = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wq = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -1975,5 +2558,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -1984,5 +2567,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0 
-                  plot_wpqvp = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wpqvp = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2006,5 +2589,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if 
             end if
@@ -2014,5 +2597,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_wsqvs = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wsqvs = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2036,5 +2619,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if
             end if
@@ -2044,5 +2627,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 2
-                  plot_wqv = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wqv = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2066,5 +2649,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -2075,5 +2658,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_wpsp = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wpsp = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2097,5 +2680,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if 
             end if 
@@ -2105,5 +2688,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1 
-                  plot_wsss = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wsss = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2127,5 +2710,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if
             end if
@@ -2135,5 +2718,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 2
-                  plot_ws = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_ws = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2157,5 +2740,5 @@
                      res@trXMaxF            = xe  
                   end if
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if 
             end if
@@ -2166,5 +2749,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_wpsap = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wpsap = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2188,5 +2771,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if 
             end if
@@ -2196,5 +2779,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_wssas = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wssas = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2218,5 +2801,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if
             end if
@@ -2226,5 +2809,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 2
-                  plot_wsa = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wsa = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2248,5 +2831,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -2257,5 +2840,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_us2 = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_us2 = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2279,5 +2862,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if 
             end if
@@ -2287,5 +2870,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_vs2 = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_vs2 = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2309,5 +2892,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if
             end if
@@ -2317,5 +2900,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 2 
-                  plot_ws2 = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_ws2 = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2339,5 +2922,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if
@@ -2348,5 +2931,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 0
-                  plot_wsususodz = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_wsususodz = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2370,5 +2953,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if  
             end if
@@ -2378,5 +2961,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 1
-                  plot_wspsodz = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot_wspsodz = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2400,5 +2983,5 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_u,res)
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
                end if 
             end if
@@ -2408,5 +2991,5 @@
                if (over .EQ. 1) then
                   res@xyDashPattern  = 2 
-                  plot_wpeodz = gsn_csm_xy(wks,data(varn,:,:),z_w,res)
+                  plot_wpeodz = gsn_csm_xy(wks,data(varn,:,:),z,res)
                else
                   res@gsnLeftString      = vNam(varn)
@@ -2430,7 +3013,14 @@
                      res@trXMaxF            = xe  
                   end if 
-                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
                end if
             end if 
+         if (no_files .GT. 1) then
+            multi_plot(nof,n)=plot(n) 
+            max_nof(nof,n)=max(data(varn,:,:))
+            min_nof(nof,n)=min(data(varn,:,:))
+            name(nof,n)   =vNam(varn)
+            unit_(nof,n)  =unit(varn)
+         end if
          if (over .EQ. 0) then
             n=n+1  
@@ -2442,4 +3032,131 @@
       end if          
    end do
+   if (no_files .GT. 1) then
+      delete(vNam)
+   end if
+   end do
+
+   if (count_var .EQ. 0) then
+      print(" ")
+      print("The variables 'var=Â°"+var+"Â°' do not exist on your input file")
+      print(" ")
+      exit
+   end if
+   
+   if (no_files .GT. 1) then 
+      multi_legend=new(12,string)
+      string_len=new(12,integer)
+      multi_dash=new(no_files,string)
+      if (.not. isvar("name_legend_1")) 
+         multi_legend(0)="  "+parameter(83) 
+      else
+         multi_legend(0)="  "+name_legend_1
+      end if 
+      string_len(0)=strlen(multi_legend(0))
+      if (.not. isvar("name_legend_2")) 
+         multi_legend(1)="  "+parameter(87) 
+      else
+         multi_legend(1)="  "+name_legend_2
+      end if 
+      string_len(1)=strlen(multi_legend(1))
+      if (.not. isvar("name_legend_3")) 
+         multi_legend(2)="  "+parameter(91) 
+      else
+         multi_legend(2)="  "+name_legend_3
+      end if 
+      string_len(2)=strlen(multi_legend(2))
+      if (.not. isvar("name_legend_4")) 
+         multi_legend(3)="  "+parameter(95) 
+      else
+         multi_legend(3)="  "+name_legend_4
+      end if 
+      string_len(3)=strlen(multi_legend(3))
+      if (.not. isvar("name_legend_5")) 
+         multi_legend(4)="  "+parameter(99) 
+      else
+         multi_legend(4)="  "+name_legend_5
+      end if 
+      string_len(4)=strlen(multi_legend(4))
+      if (.not. isvar("name_legend_6")) 
+         multi_legend(5)="  "+parameter(103) 
+      else
+         multi_legend(5)="  "+name_legend_6
+      end if 
+      string_len(5)=strlen(multi_legend(5))
+      if (.not. isvar("name_legend_7")) 
+         multi_legend(6)="  "+parameter(107) 
+      else
+         multi_legend(6)="  "+name_legend_7
+      end if
+      string_len(6)=strlen(multi_legend(6)) 
+      if (.not. isvar("name_legend_8")) 
+         multi_legend(7)="  "+parameter(111) 
+      else
+         multi_legend(7)="  "+name_legend_8
+      end if 
+      string_len(7)=strlen(multi_legend(7))
+      if (.not. isvar("name_legend_9")) 
+         multi_legend(8)="  "+parameter(115) 
+      else
+         multi_legend(8)="  "+name_legend_9
+      end if
+      string_len(8)=strlen(multi_legend(8)) 
+      if (.not. isvar("name_legend_10")) 
+         multi_legend(9)="  "+parameter(119) 
+      else
+         multi_legend(9)="  "+name_legend_10
+      end if 
+      string_len(9)=strlen(multi_legend(9))
+      if (.not. isvar("name_legend_11")) 
+         multi_legend(10)="  "+parameter(123) 
+      else
+         multi_legend(10)="  "+name_legend_11
+      end if 
+      string_len(10)=strlen(multi_legend(10))
+      if (.not. isvar("name_legend_12")) 
+         multi_legend(11)="  "+parameter(127) 
+      else
+         multi_legend(11)="  "+name_legend_12
+      end if 
+      string_len(11)=strlen(multi_legend(11))
+      do ml=1,no_files
+         multi_dash(ml-1)=ml-1
+      end do
+      delete(plot)
+      plot = new(dim,graphic)
+      do pl=0,n-1
+         plot0 = new(1,graphic)
+         res@trXMinF = min(min_nof(:,pl))
+         res@trXMaxF = max(max_nof(:,pl))
+         res@gsnLeftString  = name(0,pl)
+         res@gsnRightString = unit_(0,pl) 
+   
+         plot0 = gsn_csm_xy(wks,data_0(:,:),z_(pl,:),res)
+
+         ; ***************************************************
+         ; legend for combined plot
+         ; ***************************************************
+
+         lgres                    = True
+         lgMonoDashIndex          = False
+         lgres@lgLabelFont        = "helvetica"    
+         lgres@lgLabelFontHeightF = .2            
+         lgres@vpWidthF           = max(string_len)*0.012           
+         lgres@vpHeightF          = 0.04*no_files          
+         lgres@lgDashIndexes      = multi_dash(no_files-1:0)
+         lbid = gsn_create_legend(wks,no_files,multi_legend(no_files-1:0),lgres)       
+
+         amres = True
+         amres@amParallelPosF   = max(string_len)*0.01+0.78		    
+         amres@amOrthogonalPosF = -0.0315*no_files+0.431         
+         annoid1 = gsn_add_annotation(plot0,lbid,amres) 
+
+         do plo=0,no_files-1
+            overlay(plot0,multi_plot(plo,pl))
+            plot(pl)=plot0
+         end do
+         delete(plot0)
+      end do
+   end if
 
    if (count_var .EQ. 0) then
@@ -2527,15 +3244,35 @@
          end if      
 
-         if (check)
-
-            if (prof3d .EQ. 1) then
-
-            end if
-          
+         if (check)then
+
             if (prof3d .EQ. 0) then
-               z = f->$vNam(varn+1)$ 
-            else
-               z = z_v(varn,:)
-            end if
+               if (logy .EQ. 1) then
+                  z = f->$vNam(varn+1)$(1:dimz-1) 
+               else
+                  z = f->$vNam(varn+1)$               
+               end if
+            else 
+               do i=0,b-1            
+                  if (isStrSubset( a(i),"zu_3d" ))then
+                     z_v(varn,:) = z_u
+                     if (logy .EQ. 1) then
+                        z = z_v(varn,1:dimz-1)
+                     else
+                        z = z_v(varn,:)
+                     end if 
+                  else
+                     if (isStrSubset( a(i),"zw_3d" ))then
+                        z_v(varn,:) = z_w
+                        if (logy .EQ. 1) then
+                           z = z_v(varn,1:dimz-1)
+                        else
+                           z = z_v(varn,:)
+                        end if
+                     end if                   
+                  end if
+               end do            
+            end if
+
+            z=z/norm
 
             res@gsnLeftString      = vNam(varn)
@@ -2584,10 +3321,6 @@
                   else 
                      res@trXMaxF = xe  
-                  end if 
-                  if (vNam(varn) .EQ. "w") then 
-                     plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
-                  else
-                     plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_u,res) 
-                  end if
+                  end if  
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -2641,5 +3374,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_u,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -2692,5 +3425,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_u,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -2744,5 +3477,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_u,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -2795,5 +3528,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_u,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -2847,5 +3580,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -2898,5 +3631,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -2949,5 +3682,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -3000,5 +3733,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -3051,5 +3784,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -3102,5 +3835,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -3153,5 +3886,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -3204,5 +3937,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -3255,5 +3988,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -3307,5 +4040,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -3359,5 +4092,5 @@
                   end if
 
-                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z_w,res) 
+                  plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res) 
 
                   ; ***************************************************
@@ -3414,5 +4147,5 @@
          if (.not. isvar("var")) then
             if (parameter(21) .NE. "variables") then
-               var=parameter(21)
+               var=parameter(21)                
                check = isStrSubset( var,","+vNam(varn)+"," )
             end if
@@ -3421,19 +4154,42 @@
          end if      
 
-         if (check)
-
+         if (check)then
+            
             if (prof3d .EQ. 0) then
-               z = f->$vNam(varn+1)$ 
-            else
-               z = z_v(varn,:)  
-            end if 
-     
+               if (logy .EQ. 1) then
+                  z = f->$vNam(varn+1)$(1:dimz-1) 
+               else
+                  z = f->$vNam(varn+1)$               
+               end if
+            else 
+               do i=0,b-1            
+                  if (isStrSubset( a(i),"zu_3d" ))then
+                     z_v(varn,:) = z_u
+                     if (logy .EQ. 1) then
+                        z = z_v(varn,1:dimz-1)
+                     else
+                        z = z_v(varn,:)
+                     end if 
+                  else
+                     if (isStrSubset( a(i),"zw_3d" ))then
+                        z_v(varn,:) = z_w
+                        if (logy .EQ. 1) then
+                           z = z_v(varn,1:dimz-1)
+                        else
+                           z = z_v(varn,:)
+                        end if
+                     end if                   
+                  end if
+               end do            
+            end if
+
+            z=z/norm
+           
             com=isStrSubset( c_var,","+vNam(varn)+"," )
         
-            if (com) 
+            if (com)then 
                co = co+1           
                if (n_o .EQ. 1) then
-                  res@xyDashPattern  = 1 
-                                  
+                  res@xyDashPattern  = 1                                   
                else            
                   if (n_o .EQ. 2) then
@@ -3441,5 +4197,5 @@
                   else
                      res@xyDashPattern  = 0
-                     res@gsnLeftString  = c_var
+                     res@gsnLeftString  = "Combined Plot of "+c_var
                      res@gsnRightString = unit(varn)
                      if (.not. isvar("xs")) then 
@@ -3478,5 +4234,13 @@
             overlay(plot_o(0),plot_o(1))
          else
-            print("'c_var' is not conform with existing variables on the input file")
+            if (.not. isvar("var") .AND. parameter(21) .EQ. "variables") then
+               print(" ")
+               print("'c_var'(= "+c_var+") is not conform with existing variables on the input file")
+               print(" ")
+            else
+               print(" ")
+               print("'c_var'(= "+c_var+") must include two variables of the general plots ('var' = "+var+")")
+               print(" ")
+            end if
          end if 
       end if
@@ -3486,5 +4250,13 @@
             overlay(plot_o(0),plot_o(2))
          else
-            print("'c_var' is not conform with existing variables on the input file")
+            if (.not. isvar("var") .AND. parameter(21) .EQ. "variables") then
+               print(" ")
+               print("'c_var'(= "+c_var+") is not conform with existing variables on the input file")
+               print(" ")
+            else
+               print(" ")
+               print("'c_var'(= "+c_var+") must include three variables of the general plots ('var' = "+var+")")
+               print(" ")
+            end if
          end if 
       end if
@@ -3522,5 +4294,5 @@
 
    if (format_out .EQ. "eps" .OR. format_out .EQ. "epsi") then
-      gsn_panel(wks,plot_(0:n-1),(/n,1/),resP)
+      gsn_panel(wks,plot_,(/n,1/),resP)
    else    
       do i = 0,n-1, no_lines*no_columns
Index: palm/trunk/SCRIPTS/NCL/timeseries.ncl
===================================================================
--- palm/trunk/SCRIPTS/NCL/timeseries.ncl	(revision 174)
+++ palm/trunk/SCRIPTS/NCL/timeseries.ncl	(revision 175)
@@ -10,12 +10,17 @@
    ; ***************************************************
    
-   if (isfilepresent(".ncl_preferences")) then
-      parameter = asciiread(".ncl_preferences",77,"string")
+   if (isfilepresent("~/.ncl_preferences")) then
+      parameter = asciiread("~/.ncl_preferences",129,"string")
       delete(parameter@_FillValue)
    else
-      print(" ")
-      print("'.ncl_preferences' is not existent")
-      print(" ")
-      exit
+      if (isfilepresent("~/palm/current_version/trunk/SRIPTS/NCL/.ncl_preferences")) then
+         parameter = asciiread("~/palm/current_version/trunk/SRIPTS/NCL/.ncl_preferences",129,"string")
+         delete(parameter@_FillValue)
+      else
+         print(" ")
+         print("'.ncl_preferences' is not present in '~/palm/current_version/trunk/SRIPTS/NCL/'")
+         print(" ")
+         exit
+      end if
    end if
 
@@ -24,20 +29,41 @@
    ; ***************************************************
   
-   if ( .not. isvar("file_in") ) then     		; path+name of input file      
-      if (parameter(7) .EQ. "input file") then
-         print(" ")
-         print("Please provide input file 'file_in = ' either in prompt or parameter_list")
+   if ( .not. isvar("file_1") ) then     		  
+      if (parameter(7) .EQ. "File in") then
+         print(" ")
+         print("Please provide 1st input file 'file_1=' either in prompt or parameter_list")
          print(" ")
          exit
       else
          file_in = parameter(7)
-      end if      
-   end if
-   if ( .not. isvar("format_out") ) then   		; format of output file 
+      end if   
+   else
+      file_in = file_1   
+   end if
+   if (.not. isfilepresent(file_in)) then
+      print(" ")
+      print("Your 1st input file: '"+file_in+"' does not exist")
+      print(" ")
+      exit
+   end if 
+
+   if ( .not. isvar("format_out") ) then   		
       format_out = "x11"
       if (parameter(9) .NE. "x11") then
          format_out = parameter(9)  
+         if (format_out .NE. "x11" .AND. format_out .NE. "pdf" .AND. format_out .NE. "eps" .AND. format_out .NE. "ps" .AND. format_out .NE. "epsi" .AND. format_out .NE. "ncgm")then
+            print(" ")
+            print("Your 'format_out = "+format_out+"' is invalid and set to'x11'")
+            print(" ")
+         end if  
       end if 
-   end if
+   else
+      if (format_out .NE. "x11" .AND. format_out .NE. "pdf" .AND. format_out .NE. "eps" .AND. format_out .NE. "ps" .AND. format_out .NE. "epsi" .AND. format_out .NE. "ncgm")then
+         print(" ")
+         print("Your 'format_out = "+format_out+"' is invalid and set to'x11'")
+         print(" ")
+      end if   
+   end if
+
    if ( .not. isvar("file_out") ) then    	 	; path+name of output file
       file_out = "test"
@@ -46,4 +72,5 @@
       end if     
    end if
+
    if ( .not. isvar("no_columns") ) then    		; number of plots in one row
       no_columns = 1
@@ -52,4 +79,5 @@
       end if
    end if
+
    if ( .not. isvar("no_lines") ) then     		; number of plot-lines on one sheet
       no_lines = 2
@@ -58,7 +86,9 @@
       end if
    end if  
+
    if ( .not. isvar("var") ) then   			; variable name
       check = True
    end if
+
    if ( .not. isvar("over") ) then     		        ; switches overlaying plots on
       over = 0
@@ -67,11 +97,18 @@
          if (stringtointeger(parameter(37)) .NE. 1) then
             print(" ")
-            print("Please set 'over' to 0 or 1")
-            print(" ")
-            exit
+            print("Your 'over'= "+over+" is invalid and set to 0")
+            print(" ")
+            over = 0
          end if   
       end if
-   end if 
-  
+   else
+      if (over .NE. 0 .AND. over .NE. 1)then
+         print(" ")
+         print("Your 'over'= "+over+" is invalid and set to 0")
+         print(" ")
+         over = 0
+      end if
+   end if
+
    ; ***************************************************
    ; open input file
@@ -220,4 +257,9 @@
    delete(end_time_step)
    end_time_step=round(et,3)
+
+   print(" ")
+   print("Output of time steps from "+t_all(start_time_step)/3600+" h = "+t_all(start_time_step)+" s => index = "+start_time_step)
+   print("                     till "+t_all(end_time_step)/3600+" h = "+t_all(end_time_step)+" s => index = "+end_time_step)
+   print(" ")
 
    t = f->time(start_time_step:end_time_step)
@@ -477,10 +519,8 @@
    if (count_var .EQ. 0) then
       print(" ")
-      print("Please select a variable 'var=' or use the default value")
-      print(" ")
-      print("Your selection '"+var+"' does not exist on the input file")
+      print("The variables 'var=Â°"+var+"Â°' do not exist on your input file")
       print(" ")
       exit
-   end if 
+   end if
 
    do varn = dim-1,0,1 
@@ -500,4 +540,11 @@
   
       if(check) then
+
+        if (isStrSubset(vNam(varn),"_0" ))then
+            print(" ")
+            print("If you have Outputs of statistic regions you cannot overlay variables; 'over' is set to 0")
+            print(" ")
+            over = 0
+         end if
 
         if (over .EQ. 1) then
Index: palm/trunk/SOURCE/CURRENT_MODIFICATIONS
===================================================================
--- palm/trunk/SOURCE/CURRENT_MODIFICATIONS	(revision 174)
+++ palm/trunk/SOURCE/CURRENT_MODIFICATIONS	(revision 175)
@@ -77,4 +77,7 @@
 Errors:
 ------
+Bugfix: Determination of the number of vertical levels for which spectra 
+shall be output - 100 replaced by 10. (netcdf)
+
 Bugfix: Lengths of those characters that are needed for the output of the 
 characteristic levels of potential temperature, the geostrophic wind, 
Index: palm/trunk/SOURCE/netcdf.f90
===================================================================
--- palm/trunk/SOURCE/netcdf.f90	(revision 174)
+++ palm/trunk/SOURCE/netcdf.f90	(revision 175)
@@ -7,4 +7,6 @@
 ! Current revisions:
 ! ------------------
+! bug fix in the determination of the number of output heights for spectra
+!
 ! user-defined spectra
 !
@@ -3057,5 +3059,5 @@
 !--       shall be output
           ns = 1
-          DO WHILE ( comp_spectra_level(ns) /= 999999  .AND.  ns <= 100 )
+          DO WHILE ( comp_spectra_level(ns) /= 999999  .AND.  ns <= 10 )
              ns = ns + 1
           ENDDO
