source: palm/trunk/SCRIPTS/NCL/profiles.ncl @ 958

Last change on this file since 958 was 958, checked in by hoffmann, 12 years ago

minor changes in palmplot pr

File size: 140.6 KB
Line 
1load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
2load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
3load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
4
5;***************************************************
6; load .ncl.config or .ncl.config.default
7;***************************************************
8
9function which_script()
10local script
11begin
12   script="profiles"
13   return(script)
14end
15   
16if (isfilepresent("$PALM_BIN/../../.ncl.config")) then
17   loadscript("$PALM_BIN/../../.ncl.config")
18else
19  if (isfilepresent("$PALM_BIN/NCL/.ncl.config.default")) then
20     loadscript( "$PALM_BIN/NCL/.ncl.config.default")
21  else
22      palm_bin_path = getenv("PALM_BIN")
23      print(" ")
24      print("Neither the personal configuration file '.ncl.config' exists in")
25      print("~/palm/current_version")
26      print("nor the default configuration file '.ncl.config.default' "+\
27            "exists in")
28      print(palm_bin_path + "/NCL")
29      print(" ")
30      exit
31   end if
32end if   
33
34   
35begin
36
37   ;***************************************************
38   ; Retrieving the NCL version used
39   ;***************************************************
40   
41   ncl_version_ch = systemfunc("ncl -V")
42   ncl_version    = stringtofloat(ncl_version_ch)
43
44   ;***************************************************
45   ; Retrieving the double quote character
46   ;***************************************************
47   
48   dq=str_get_dq()
49
50   ;***************************************************
51   ; set up default parameter values and strings
52   ;***************************************************
53
54   if (no_files .LT. 1 .OR. no_files .GT. 6) then
55      print(" ")
56      print("Assign 'no_files' between 1 and 6") 
57      print(" ")
58      exit
59   end if
60   
61   file_in   = new(no_files,string)
62   file_in_1 = new(no_files,logical)
63   start_f_t = new(no_files,integer)
64   end_f_t   = new(no_files,integer)
65   
66   if (file_1 .EQ. "File in") then
67      print(" ")
68      print("Declare 1st input file 'file_1=' in '.ncl.config' or prompt")
69      print(" ")
70      exit
71   else
72      file_in(0) = file_1   
73   end if     
74   file_in_1(0) = False
75   if (isStrSubset(file_in(0), ".nc"))then
76      start_f = -2
77      end_f = -2
78      file_in_1(0) = True     
79   end if 
80   if (start_f .EQ. -1)then
81      print(" ")
82      print("'start_f' must be one of the cyclic numbers (at least 0) "+\
83            "of your input file(s)")
84      print(" ") 
85      exit
86   end if
87   if (end_f .EQ. -1)then
88      print(" ")
89      print("'end_f' must be one of the cyclic numbers (at least 0) of "+\
90            "your input file(s)")
91      print(" ") 
92      exit
93   end if           
94   start_f_t(0) = start_f
95   end_f_t(0) = end_f   
96 
97   if (no_files .GT. 1) then
98      if (file_2 .EQ. "File in") then
99         print(" ")
100         print("Declare 2nd input file 'file_2=' in '.ncl.config' or prompt")
101         print(" ")
102         exit
103      else
104         file_in(1) = file_2   
105      end if     
106      file_in_1(1) = False
107      if (isStrSubset(file_in(1), ".nc"))then
108         start_f_2 = -2
109         end_f_2 = -2
110         file_in_1(1) = True     
111      end if
112      if (start_f_2 .EQ. -1)then
113         print(" ")
114         print("'start_f_2' must be one of the cyclic numbers (at least 0) "+\
115              "of your input file(s)")
116         print(" ") 
117         exit
118      end if
119      if (end_f_2 .EQ. -1)then
120         print(" ")
121         print("'end_f_2' must be one of the cyclic numbers (at least 0) "+\
122               "of your input file(s)")
123         print(" ") 
124         exit
125      end if     
126      start_f_t(1) = start_f_2
127      end_f_t(1) = end_f_2   
128   end if   
129   
130   if (no_files .GT. 2) then
131      if (file_3 .EQ. "File in") then
132         print(" ")
133         print("Declare 3rd input file 'file_3=' in '.ncl.config' or prompt")
134         print(" ")
135         exit
136      else
137         file_in(2) = file_3   
138      end if     
139      file_in_1(2) = False
140      if (isStrSubset(file_in(2), ".nc"))then
141         start_f_3 = -2
142         end_f_3 = -2
143         file_in_1(2) = True     
144      end if
145      if (start_f_3 .EQ. -1)then
146         print(" ")
147         print("'start_f_3' must be one of the cyclic numbers (at least 0) "+\
148               "of your input file(s)")
149         print(" ") 
150         exit
151      end if
152      if (end_f_3 .EQ. -1)then
153         print(" ")
154         print("'end_f_3' must be one of the cyclic numbers (at least 0) "+\
155               "of your input file(s)")
156         print(" ") 
157         exit
158      end if     
159      start_f_t(2) = start_f_3
160      end_f_t(2) = end_f_3 
161   end if
162   
163   if (no_files .GT. 3) then
164      if (file_4 .EQ. "File in") then
165         print(" ")
166         print("Declare 4th input file 'file_4=' in '.ncl.config' or prompt")
167         print(" ")
168         exit
169      else
170         file_in(3) = file_4   
171      end if     
172      file_in_1(3) = False
173      if (isStrSubset(file_in(3), ".nc"))then
174         start_f_4 = -2
175         end_f_4 = -2
176         file_in_1(3) = True     
177      end if
178      if (start_f_4 .EQ. -1)then
179         print(" ")
180         print("'start_f_4' must be one of the cyclic numbers (at least 0) "+\
181               "of your input file(s)")
182         print(" ") 
183         exit
184      end if
185      if (end_f_4 .EQ. -1)then
186         print(" ")
187         print("'end_f_4' must be one of the cyclic numbers (at least 0) "+\
188               "of your input file(s)")
189         print(" ") 
190         exit
191      end if     
192      start_f_t(3) = start_f_4
193      end_f_t(3) = end_f_4
194   end if
195   
196   if (no_files .GT. 4) then
197      if (file_5 .EQ. "File in") then
198         print(" ")
199         print("Declare 5th input file 'file_5=' in '.ncl.config' or prompt")
200         print(" ")
201         exit
202      else
203         file_in(4) = file_5   
204      end if     
205      file_in_1(4) = False
206      if (isStrSubset(file_in(4), ".nc"))then
207         start_f_5 = -2
208         end_f_5 = -2
209         file_in_1(4) = True     
210      end if
211      if (start_f_5 .EQ. -1)then
212         print(" ")
213         print("'start_f_5' must be one of the cyclic numbers (at least 0) "+\
214               "of your input file(s)")
215         print(" ") 
216         exit
217      end if
218      if (end_f_5 .EQ. -1)then
219         print(" ")
220         print("'end_f_5' must be one of the cyclic numbers (at least 0) "+\
221               "of your input file(s)")
222         print(" ") 
223         exit
224      end if     
225      start_f_t(4) = start_f_5
226      end_f_t(4) = end_f_5   
227   end if
228   
229   if (no_files .GT. 5) then
230      if (file_6 .EQ. "File in") then
231         print(" ")
232         print("Declare 6th input file 'file_6=' in '.ncl.config' or prompt")
233         print(" ")
234         exit
235      else
236         file_in(5) = file_6   
237      end if     
238      file_in_1(5) = False
239      if (isStrSubset(file_in(5), ".nc"))then
240         start_f_6 = -2
241         end_f_6 = -2
242         file_in_1(5) = True     
243      end if
244      if (start_f_6 .EQ. -1)then
245         print(" ")
246         print("'start_f_6' must be one of the cyclic numbers (at least 0) "+\
247               "of your input file(s)")
248         print(" ") 
249         exit
250      end if
251      if (end_f_6 .EQ. -1)then
252         print(" ")
253         print("'end_f_6' must be one of the cyclic numbers (at least 0) "+\
254               "of your input file(s)")
255         print(" ") 
256         exit
257      end if     
258      start_f_t(5) = start_f_6
259      end_f_t(5) = end_f_6   
260   end if
261
262   
263   if (format_out .NE. "x11" .AND. format_out .NE. "pdf" .AND.  \
264      format_out .NE. "eps" .AND. format_out .NE. "ps" .AND.    \
265      format_out .NE. "epsi" .AND. format_out .NE. "ncgm" .AND. \
266      format_out .NE. "png")then
267      print(" ")
268      print("'format_out = "+format_out+"' is invalid and set to'x11'")
269      print(" ")
270      format_out="x11"
271   end if 
272
273   if (ncl_version .LE. 5.1 .AND. format_out .EQ. "png") then
274      print(" ")
275      print("Output of png files not available")
276      print("png output is avaiable with NCL version 5.2.0 and higher ")
277      print("NCL version used: " + ncl_version_ch)
278      print(" ")
279      exit
280   end if
281
282   if (over .NE. 0 .AND. over .NE. 1) then
283      print(" ")
284      print("'over'= "+over+" is invalid and set to 0")
285      print(" ")
286      over = 0
287   end if
288   if (no_files .GT. 1) then
289      over = 0
290      print(" ")
291      print("If you have more than one input file - you cannot overlay "+\
292            "variables: over is set to 0")
293      print(" ")
294   end if
295
296   ;************************************************************************
297   ; check if cross_profiles is defined in file_1
298   ;************************************************************************
299
300    if (isStrSubset(file_in(0), ".nc")) then
301      firstfile=addfile(file_1,"r")
302    else
303      firstfile=addfile(file_in(0) + "." + start_f + ".nc","r")
304    end if
305
306   cross_check=isatt(firstfile,"cross_profiles")
307
308   if (var .EQ. "all") then
309      VAR_LIST = firstfile@VAR_LIST
310      no_plots = str_fields_count(VAR_LIST, ";")
311   else
312      no_plots = str_fields_count(var, ",")
313   end if
314
315   if ((.not. cross_check) .AND. (var .EQ. ",cross_profiles,")) then
316       print(" ")
317       print("cross_profiles is not defined in the NetCDF header."+\
318             "Using var='all' instead.")
319       print(" ")
320       var="all"
321   end if
322
323   wks_gsn_log = True
324
325   ;************************************************************************
326   ; check if cross_profiles can be used
327   ;************************************************************************
328
329   if ((var .EQ. ",cross_profiles,") .AND. cross_check .AND. (over .EQ. 0) .AND. (no_files .LT. 2)) then
330      c_var_log = 1
331   else
332      c_var_log = 0
333   end if
334
335   if ((no_rows .EQ. -1) .AND. isatt(firstfile,"no_rows")) then
336      no_rows = firstfile@no_rows
337   else if ((no_rows .EQ. -1) .AND. .not. isatt(firstfile,"no_rows")) then
338       print(" ")
339       print("no_rows is not defined in the NetCDF header."+\
340             "Using no_rows=3 instead.")
341       print(" ")
342      no_rows = 3
343   end if
344   end if
345
346   if ((no_columns .EQ. -1) .AND. isatt(firstfile,"no_columns")) then
347      no_columns = firstfile@no_columns
348   else if ((no_columns .EQ. -1) .AND. .not. isatt(firstfile,"no_columns")) then
349       print(" ")
350       print("no_columns is not defined in the NetCDF header."+\
351             "Using no_columns=2 instead.")
352       print(" ")
353      no_columns = 2
354   end if
355   end if
356
357   if (c_var_log .EQ. 1) then
358      c_var_str = firstfile@cross_profiles
359
360      c_var_all = str_split(c_var_str,";")
361
362      c_var_all = ","+c_var_all+","
363
364      ;************************************************************************
365      ; check if an additional combined plot is demanded by the user
366      ;************************************************************************
367
368      if ((c_var .EQ. "c_variables") .AND. (number_comb .EQ. -1) .AND. \
369             (combine .EQ. 0)) then
370         max_com_i = max(dimsizes(c_var_all)) - 1
371      else
372         max_com_i = max(dimsizes(c_var_all))
373         c_var_all_temp = c_var_all
374         delete(c_var_all)
375         c_var_all = new(max_com_i+1,string)
376         c_var_all(0:max_com_i-1) = c_var_all_temp
377         delete(c_var_all_temp)
378         c_var_all(max_com_i) = c_var
379      end if
380      number_comb_all = str_fields_count(c_var_all,",")
381
382      ;************************************************************************
383      ; currently, more than 3 plots cannot be drawn together
384      ;***********************************************************************
385
386      do com_i = 0, max_com_i
387         if(number_comb_all(com_i) .GT. 3) then
388            print(" ")
389            print("Currently, it is not possible to plot more than three"+\
390                  "profiles together. Hence, they are split to a maximum"+\
391                  "of three profiles per coordinate plane.")
392            print(" ")
393            c_var_all_temp = c_var_all
394            number_comb_all_temp = number_comb_all
395            delete(c_var_all)
396            delete(number_comb_all)
397            max_com_i = max_com_i + 1
398            c_var_all = new(max_com_i+1,string)
399            number_comb_all = new(max_com_i+1, integer)
400            if (com_i-1 .GT. 0) then
401              c_var_all(0:com_i-1) = c_var_all_temp(0:com_i-1)
402              number_comb_all(0:com_i-1) = number_comb_all_temp(0:com_i-1)
403            end if
404            c_var_long = str_split(c_var_all_temp(com_i),",")
405            c_var_all(com_i) = "," + c_var_long(0) + "," + c_var_long(1) + "," + c_var_long(2) + ","
406            number_comb_all(com_i) = 3
407            c_var_all(com_i + 1) = ","
408            do com_j = 3, max(dimsizes(c_var_long)) - 1
409               c_var_all(com_i + 1) = c_var_all(com_i + 1) + c_var_long(com_j) + ","
410            end do
411            number_comb_all(com_i+1) = max(dimsizes(c_var_long)) - 3
412            c_var_all(com_i+2:max_com_i) = c_var_all_temp(com_i+1:max_com_i-1)
413            number_comb_all(com_i+2:max_com_i) = number_comb_all_temp(com_i+1:max_com_i-1)
414            delete(c_var_all_temp)
415            delete(number_comb_all_temp)
416            delete(c_var_long)
417         end if
418      end do
419
420      c_var_plot = new(max_com_i+1,graphic)
421      c_var_end_time_step = end_time_step
422      c_var_start_time_step = start_time_step
423      number_comb_1_log = 0
424      number_comb_23_log = 0
425   else
426      max_com_i = 0
427   end if
428
429   do com_i = 0,max_com_i
430
431      ;************************************************************************
432      ; changes basic properties of the plots according to cross_profiles
433      ;************************************************************************
434
435      if (c_var_log .EQ. 1) then
436         if (number_comb_all(com_i) .EQ. 1) then
437            number_combine= -1
438            combine = 0
439            c_var = "c_variables"
440            var = c_var_all(com_i)
441            delete(end_time_step)
442            delete(start_time_step)
443            end_time_step = c_var_end_time_step
444            start_time_step= c_var_start_time_step
445            if (number_comb_1_log .EQ. 1) then
446                delete(var_char)
447            end if
448            number_comb_1_log=1
449         else
450            var = "all"
451            number_comb = number_comb_all(com_i)
452            combine = 1
453            c_var = c_var_all(com_i)
454            delete(end_time_step)
455            delete(start_time_step)
456            end_time_step = c_var_end_time_step
457            start_time_step= c_var_start_time_step
458            if ((number_comb_1_log .EQ. 1) .OR. (number_comb_23_log .EQ. 1)) then
459                delete(label)
460            end if
461            if (number_comb_23_log .EQ. 1) then
462                delete(plot_o)
463;               delete(label)
464                delete(color_o)
465                delete(mini)
466                delete(maxi)
467                delete(plot)
468            end if
469            number_comb_23_log=1
470         end if
471      end if
472
473      if (combine .NE. 0 .AND. combine .NE. 1)then
474          print(" ")
475          print("Your 'combine'= "+combine+" is invalid and set to 0")
476          print(" ")
477          combine = 0
478      end if
479      if (no_files .GT. 1) then
480          combine = 0
481          print(" ")
482          print("If you have more than one input file you cannot combine "+\
483                "variables: combine is set to 0")
484          print(" ")
485      end if   
486      if (combine .EQ. 1 .AND. number_comb .EQ. -1) then
487          print(" ")
488          print("Set 'number_comb' to 2 or 3 or combine to 0")
489          print(" ")
490          exit
491      end if
492      if (combine .EQ. 1)then
493          if (number_comb .EQ. -1) then
494            print(" ")
495            print("Set 'number_comb' to 2 or 3 or combine to 0")
496            print(" ")
497            exit
498          end if
499          if (number_comb .NE. 2 .AND. number_comb .NE. 3) then
500            print(" ")
501            print("Set 'number_comb' to 2 or 3 or combine to 0")
502            print(" ")
503            exit
504          end if
505      end if 
506      if (combine .EQ. 1 .AND. c_var .EQ. "c_variables") then
507          print(" ")
508          print("Select variables for overlaying ('c_var') or set combine to 0")
509          print(" ")
510          exit
511      end if
512
513      if (log_z .NE. 0 .AND. log_z .NE. 1)then
514          print(" ")
515          print("'log_z'= "+log_z+" is invalid and set to 0")
516          print(" ")
517          log_z = 0
518      end if   
519     
520      if (norm_z .EQ. 0) then
521          print(" ")
522          print("Normalising with 0 is not allowed, 'norm_z' is set to 1.0")
523          print(" ")
524          norm_z = 1.0
525      end if
526
527      ;***************************************************
528      ; open input file
529      ;***************************************************
530
531      do nof=0,no_files-1
532
533if (com_i .EQ. 0)
534
535      files=new(end_f_t(nof)-start_f_t(nof)+1,string)
536      if (file_in_1(nof))then
537          if (isfilepresent(file_in(nof)))then
538            files(0)=file_in(nof)
539          else
540            print(" ")
541            print("1st input file: '"+file_in(nof)+"' does not exist")
542            print(" ")
543            exit
544          end if
545      else
546          if (start_f_t(nof) .EQ. 0)then
547            if (isfilepresent(file_in(nof)+".nc"))then
548                files(0)=file_in(nof)+".nc"
549                do i=1,end_f_t(nof)
550                  if (isfilepresent(file_in(nof)+"."+i+".nc"))then   
551                      files(i)=file_in(nof)+"."+i+".nc"
552                  else
553                      print(" ")
554                      print("Input file: '"+file_in(nof)+"."+i+".nc' does not "+\
555                            "exist")
556                      print(" ")
557                      exit 
558                  end if           
559                end do         
560            else
561                print(" ")
562                print("Input file: '"+file_in(nof)+".nc' does not exist")
563                print(" ")
564                exit
565            end if
566          else
567            do i=start_f_t(nof),end_f_t(nof)
568                if (isfilepresent(file_in(nof)+"."+i+".nc"))then   
569                  files(i-start_f_t(nof))=file_in(nof)+"."+i+".nc"
570                else
571                  print(" ")
572                  print("Input file: '"+file_in(nof)+"."+i+".nc' does not exist")
573                  print(" ")
574                  exit 
575                end if
576            end do
577          end if
578      end if
579
580      f=addfiles(files,"r")
581      f_att=addfile(files(0),"r")
582      ListSetType(f,"cat")
583     
584      vNam  = getfilevarnames(f_att)
585      vType = getfilevartypes(f_att,vNam)
586
587      if ((all(vType .eq. "double"))) then ;distinction if data is double or float
588          check_vType = True
589      else
590          check_vType = False
591      end if
592     
593      if (nof .EQ. 0)then
594          vNam0=vNam
595      end if
596      if (nof .NE. 0)then
597          if (dim0 .NE. dim)then
598            print(" ")
599            print("There are 'no_files'="+no_files+" input files but they do "+\
600                  "not contain the same variables")
601            print(" ")
602            exit
603          else
604            do i=0,dim0-1
605                if (vNam0(i) .NE. vNam(i))then
606                  print(" ")
607                  print("There are 'no_files'="+no_files+" input files but "+\
608                        "they do not contain the same variables")
609                  print(" ")
610                  exit
611                end if
612            end do
613          end if
614      end if
615      nof=nof+1
616      if (com_i .EQ. 0) then
617          print(" ")
618          print("Variables in input file "+nof+":")
619          print("- "+ vNam)
620          print(" ")
621      end if
622      nof=nof-1
623      dim = dimsizes(vNam)
624      dim0=dim
625     
626      if (dim .EQ. 0) then
627          print(" ")
628          print("There is no data on file")
629          print(" ")
630      end if
631     
632      prof3d = 0
633      do varn = dim-1,0,1
634          if ( isStrSubset( vNam(varn), "zu_3d") .OR. \
635              isStrSubset( vNam(varn), "zw_3d")) then
636            prof3d = 1
637            break
638          end if
639      end do
640     
641end if
642
643      if (var .NE. "all") then
644          ;rearrange the order of the variables in vNam so that the variables
645          ;specified by "var" are in the top of vNam
646
647          vNam_static = new((/dim/),string)
648          vNam_temp   = new((/dim/),string)
649
650          var_char    = stringtocharacter(var)
651          no_char     = dimsizes(var_char)
652          comma       = 0
653         
654          do j=0,no_char-1
655            if(var_char(j) .eq. ",")
656                comma = comma + 1
657            end if
658          end do
659
660          if(comma .le. 1)
661              print(" ")
662              print("The variables 'var="+var+"'" )
663              print("do not exist on your input file;")
664              print("be sure to have one comma before and after each variable")
665              print(" ")
666              exit
667          end if
668
669          if(comma .gt. dim)
670              print(" ")
671              print("The variables 'var="+var+"'" )
672              print("do not exist on your input file;")
673              print("be sure to have one comma before and after each variable")
674              print(" ")
675              exit
676          end if
677
678          indices = new((/comma/),integer)
679          comma   = 0
680           
681          do j=0,no_char-1
682            if(var_char(j) .eq. ",")
683              indices(comma) = j
684              comma          = comma + 1
685            end if
686          end do
687
688          do j=0,comma-2
689            vNam_temp(j) = charactertostring(\
690                                      var_char(indices(j)+1:indices(j+1)-1))
691          end do
692
693          do j=0,comma-2
694            count_check = 0
695            do varn=0,dim-1
696              if(vNam_temp(j) .ne. vNam(varn))
697                count_check=count_check+1
698              end if
699            end do   
700
701            if(count_check .eq. dim)
702                print(" ")
703                print("The variables 'var="+var+"'" )
704                print("do not exist on your input file;")
705                print("be sure to have one comma before and after each variable")
706                print(" ")
707                exit
708            end if
709
710            vNam_static(j) = vNam_temp(comma - 2 - j)
711          end do
712
713          vNam_temp   = vNam_static
714          vNam_static = ""
715
716          if (prof3d .EQ. 0) then
717            i = 0
718            do j=0,dim-1,2
719                vNam_static(j)   = vNam_temp(i)
720                vNam_static(j+1) = "z" + vNam_static(j)
721                i = i+1
722                if(i .eq. comma-1)
723                    break
724                end if
725            end do
726
727            delete(vNam_temp)
728            count_variable = 2*(comma-1)
729          else
730            i = 0
731              do j=0,dim-1,1
732                vNam_static(j) = vNam_temp(i)
733                i = i+1
734                if(i .eq. comma-1)
735                    break
736                end if
737              end do
738
739              delete(vNam_temp)
740              count_variable = comma-1
741          end if
742
743          counter  = 0
744          counter2 = 0
745
746          do j=0,dim-1
747            counter = 0
748            do i = 0, count_variable - 1       
749                if(vNam_static(i) .ne. vNam(j))
750                  counter = counter + 1
751                end if
752            end do
753         
754            if (counter .eq. count_variable)
755                vNam_static(count_variable + counter2) = vNam(j)
756                counter2 = counter2 + 1
757            end if
758          end do
759       
760          vNam = vNam_static
761          delete(vNam_static)
762        end if
763
764      ;---------below steps only for first file -> nof=0
765      if (nof .EQ. 0) then
766
767          plot = new(dim,graphic)
768          plot_ = new(dim,graphic)
769          if (no_files .GT. 1) then
770            multi_plot = new((/no_files,dim/),graphic)
771            if (check_vType) then
772              max_nof = new((/no_files,dim/),double)
773              min_nof = new((/no_files,dim/),double)
774            else
775              max_nof = new((/no_files,dim/),float)
776              min_nof = new((/no_files,dim/),float)
777            end if
778            name    = new((/no_files,dim/),string)
779            unit_   = new((/no_files,dim/),string)
780          end if
781
782      if (prof3d .EQ. 0) then
783   
784      do varn=0,dim-1
785          if ( isStrSubset( vNam(varn), "time") .OR. \
786              isStrSubset( vNam(varn), "NORM")) then
787            continue
788          end if
789          if (vNam(varn) .EQ. "u" .OR. isStrSubset(vNam(varn), "u_"))then
790            z_u = f_att->$vNam(varn+1)$
791            break
792          else
793            if (vNam(varn) .EQ. "v" .OR. isStrSubset(vNam(varn), "v_"))then 
794                z_u = f_att->$vNam(varn+1)$
795                break
796            else
797                if(vNam(varn) .EQ. "pt" .OR. isStrSubset(vNam(varn), "pt_"))then
798                  z_u = f_att->$vNam(varn+1)$
799                  break
800                else
801                  if(vNam(varn) .EQ. "vpt" .OR. isStrSubset(vNam(varn), "vpt_"))then
802                      z_u = f_att->$vNam(varn+1)$
803                      break
804                  else
805                      if(vNam(varn) .EQ. "lpt" .OR. isStrSubset(vNam(varn), "lpt_"))then
806                        z_u = f_att->$vNam(varn+1)$
807                        break
808                      else   
809                        if(vNam(varn) .EQ. "q" .OR. isStrSubset(vNam(varn), "q_") )then
810                            z_u = f_att->$vNam(varn+1)$
811                            break
812                        else
813                            if(vNam(varn) .EQ. "qv" .OR. isStrSubset(vNam(varn), "qv_"))then
814                              z_u = f_att->$vNam(varn+1)$
815                              break
816                            else
817                              if(vNam(varn) .EQ. "ql" .OR. isStrSubset(vNam(varn), "ql_"))then
818                                  z_u = f_att->$vNam(varn+1)$
819                                  break
820                              else
821                                  if(vNam(varn) .EQ. "rho" .OR. isStrSubset(vNam(varn), "rho_"))then
822                                    z_u = f_att->$vNam(varn+1)$
823                                    break
824                                  else
825                                    if(vNam(varn) .EQ. "s" .OR. isStrSubset(vNam(varn), "s_") )then
826                                        z_u = f_att->$vNam(varn+1)$
827                                        break
828                                    else
829                                        if(vNam(varn) .EQ. "sa" .OR. isStrSubset(vNam(varn), "sa_"))then
830                                          z_u = f_att->$vNam(varn+1)$
831                                          break
832                                        else
833                                          if(vNam(varn) .EQ. "e" .OR. isStrSubset(vNam(varn), "e_"))then
834                                           
835                                              z_u = f_att->$vNam(varn+1)$
836                                              break
837                                          else
838                                              if(vNam(varn) .EQ. "es" .OR. \
839                                              isStrSubset(vNam(varn), "es_") .OR. vNam(varn) .EQ. "e*" .OR. isStrSubset(vNam(varn), "e*_"))then
840                                                z_u = f_att->$vNam(varn+1)$
841                                                break
842                                              else
843                                                if(vNam(varn) .EQ. "km" .OR. isStrSubset(vNam(varn), "km_"))then
844                                                    z_u = f_att->$vNam(varn+1)$
845                                                    break
846                                                else
847                                                    if(vNam(varn) .EQ. "kh" .OR. isStrSubset(vNam(varn), "kh_"))then
848                                                      z_u = f_att->$vNam(varn+1)$
849                                                      break
850                                                    else
851                                                      if(vNam(varn) .EQ. "l" .OR. isStrSubset(vNam(varn), "l_"))then
852                                                          z_u = f_att->$vNam(varn+1)$
853                                                          break
854                                                      else
855                                                          if(vNam(varn) .EQ. "us2" .OR. isStrSubset(vNam(varn), "us2_")\
856                                                            .OR. vNam(varn) .EQ. "u*2" .OR. isStrSubset(vNam(varn), "u*2_"))then
857                                                            z_u = f_att->$vNam(varn+1)$
858                                                            break
859                                                          else
860                                                            if(vNam(varn) .EQ. "vs2" .OR. isStrSubset(vNam(varn), "vs2_")\
861                                                                .OR. vNam(varn) .EQ. "v*2" .OR. isStrSubset(vNam(varn), "v*2_"))then
862                                                                z_u = f_att->$vNam(varn+1)$
863                                                                break
864                                                            else
865                                                                if(vNam(varn) .EQ. "pts2" .OR. isStrSubset(vNam(varn), "pts2_")\
866                                                                  .OR. vNam(varn) .EQ. "pt*2" .OR. isStrSubset(vNam(varn), "pt*2_"))then
867                                                                  z_u = f_att->$vNam(varn+1)$
868                                                                  break
869                                                                else
870                                                                  if(vNam(varn) .EQ. "wsususodz" .OR. isStrSubset(vNam(varn), "wsususodz_")\
871                                                                      .OR. vNam(varn) .EQ. "w*u*u*:dz" .OR. isStrSubset(vNam(varn), "w*u*u*:dz_"))then
872                                                                      z_u = f_att->$vNam(varn+1)$
873                                                                      break
874                                                                  else
875                                                                      if(vNam(varn) .EQ. "wspsodz" .OR. isStrSubset(vNam(varn), "wspsodz_")\
876                                                                        .OR. vNam(varn) .EQ. "w*p*:dz" .OR. isStrSubset(vNam(varn), "w*p*:dz_"))then
877                                                                        z_u = f_att->$vNam(varn+1)$
878                                                                        break
879                                                                      else
880                                                                        if(vNam(varn) .EQ. "wpeodz" .OR. isStrSubset(vNam(varn), "wpeodz_")\
881                                                                            .OR. vNam(varn) .EQ. "w"+dq+"e:dz" .OR. isStrSubset(vNam(varn), "w"+dq+"e:dz_"))then
882                                                                            z_u = f_att->$vNam(varn+1)$
883                                                                            break   
884                                                                        else
885                                                                            if(vNam(varn) .EQ. "qs2" .OR. isStrSubset(vNam(varn), "qs2_")\
886                                                                              .OR. vNam(varn) .EQ. "q*2" .OR. isStrSubset(vNam(varn), "q*2_"))then
887                                                                              z_u = f_att->$vNam(varn+1)$
888                                                                              break 
889                                                                          end if               
890                                                                        end if
891                                                                      end if
892                                                                  end if
893                                                                end if
894                                                            end if
895                                                          end if
896                                                      end if
897                                                    end if
898                                                end if
899                                              end if
900                                          end if
901                                        end if
902                                    end if
903                                  end if
904                              end if
905                            end if
906                        end if
907                      end if
908                  end if
909                end if
910            end if
911          end if
912      end do
913
914      do varn=0,dim-1
915          if ( isStrSubset( vNam(varn), "time") .OR. isStrSubset( vNam(varn), "NORM")) then
916            continue
917          end if
918          if (vNam(varn) .EQ. "w" .OR. isStrSubset(vNam(varn), "w_"))then
919            z_w = f_att->$vNam(varn+1)$
920            break
921          else
922            if (vNam(varn) .EQ. "wpup" .OR. isStrSubset(vNam(varn), "wpup_")\
923                .OR. vNam(varn) .EQ. "w"+dq+"u"+dq .OR. isStrSubset(vNam(varn), "w"+dq+"u"+dq+"_"))then
924                z_w = f_att->$vNam(varn+1)$
925                break
926            else
927                if(vNam(varn) .EQ. "wsus" .OR. isStrSubset(vNam(varn), "wsus_")\
928                  .OR. vNam(varn) .EQ. "w*u*" .OR. isStrSubset(vNam(varn), "w*u*_"))then
929                  z_w = f_att->$vNam(varn+1)$
930                  break
931                else
932                  if(vNam(varn) .EQ. "wu" .OR. isStrSubset(vNam(varn), "wu_"))then
933                      z_w = f_att->$vNam(varn+1)$
934                      break
935                  else
936                      if(vNam(varn) .EQ. "wpvp" .OR. isStrSubset(vNam(varn), "wpvp_")\
937                        .OR. vNam(varn) .EQ. "w"+dq+"v"+dq .OR. isStrSubset(vNam(varn), "w"+dq+"v"+dq+"_"))then
938                        z_w = f_att->$vNam(varn+1)$
939                        break
940                      else   
941                        if(vNam(varn) .EQ. "wsvs" .OR. isStrSubset(vNam(varn), "wsvs_")\
942                            .OR. vNam(varn) .EQ. "w*v*" .OR. isStrSubset(vNam(varn), "w*v*_"))then
943                            z_w = f_att->$vNam(varn+1)$
944                            break
945                        else
946                            if(vNam(varn) .EQ. "wv" .OR. isStrSubset(vNam(varn), "wv_"))then
947                              z_w = f_att->$vNam(varn+1)$
948                              break
949                            else
950                              if(vNam(varn) .EQ. "wptpp" .OR. isStrSubset(vNam(varn), "wptpp_")\
951                                  .OR. vNam(varn) .EQ. "w"+dq+"pt"+dq .OR. isStrSubset(vNam(varn), "w"+dq+"pt"+dq+"_"))then
952                                  z_w = f_att->$vNam(varn+1)$
953                                  break
954                              else
955                                  if(vNam(varn) .EQ. "wspts" .OR. isStrSubset(vNam(varn), "wspts_")\
956                                    .OR. vNam(varn) .EQ. "w*pt*" .OR. isStrSubset(vNam(varn), "w*pt*_"))then
957                                    z_w = f_att->$vNam(varn+1)$
958                                    break
959                                  else
960                                    if(vNam(varn) .EQ. "wpt" .OR. isStrSubset(vNam(varn), "wpt_"))then
961                                        z_w = f_att->$vNam(varn+1)$
962                                        break
963                                    else
964                                        if(vNam(varn) .EQ. "wsptsBC" .OR. isStrSubset(vNam(varn), "wsptsBC_")\
965                                          .OR. vNam(varn) .EQ. "w*pt*BC" .OR. isStrSubset(vNam(varn), "w*pt*BC_"))then
966                                          z_w = f_att->$vNam(varn+1)$
967                                          break
968                                        else
969                                          if(vNam(varn) .EQ. "wptBC" .OR. isStrSubset(vNam(varn), "wptBC_"))then
970                                              z_w = f_att->$vNam(varn+1)$
971                                              break
972                                          else
973                                              if(vNam(varn) .EQ. "wpvptp" .OR. isStrSubset(vNam(varn), "wpvptp_")\
974                                                  .OR. vNam(varn) .EQ. "w"+dq+"vpt"+dq .OR. isStrSubset(vNam(varn), "w"+dq+"vpt"+dq+"_"))then
975                                                z_w = f_att->$vNam(varn+1)$
976                                                break
977                                              else
978                                                if(vNam(varn) .EQ. "wsvpts" .OR. isStrSubset(vNam(varn), "wsvpts_")\
979                                                  .OR. vNam(varn) .EQ. "w*vpt*" .OR. isStrSubset(vNam(varn), "w*vpt*_"))then
980                                                    z_w = f_att->$vNam(varn+1)$
981                                                    break
982                                                else
983                                                    if(vNam(varn) .EQ. "wvpt" .OR. isStrSubset(vNam(varn), "wvpt_"))then
984                                                      z_w = f_att->$vNam(varn+1)$
985                                                      break
986                                                    else
987                                                      if(vNam(varn) .EQ. "wpqp" .OR. isStrSubset(vNam(varn), "wpqp_")\
988                                                          .OR. vNam(varn) .EQ. "w"+dq+"q"+dq .OR. isStrSubset(vNam(varn), "w"+dq+"q"+dq+"_"))then
989                                                          z_w = f_att->$vNam(varn+1)$
990                                                          break
991                                                      else
992                                                          if(vNam(varn) .EQ. "wsqs" .OR. isStrSubset(vNam(varn), "wsqs_")\
993                                                            .OR. vNam(varn) .EQ. "w*q*" .OR. isStrSubset(vNam(varn), "w*q*_"))then
994                                                            z_w = f_att->$vNam(varn+1)$
995                                                            break
996                                                          else
997                                                            if(vNam(varn) .EQ. "wq" .OR. isStrSubset(vNam(varn), "wq_"))then
998                                                                z_w = f_att->$vNam(varn+1)$
999                                                                break
1000                                                            else
1001                                                                if(vNam(varn) .EQ. "wpqvp" .OR. isStrSubset(vNam(varn), "wpqvp_")\
1002                                                                  .OR. vNam(varn) .EQ. "w"+dq+"qv"+dq .OR. isStrSubset(vNam(varn), "w"+dq+"qv"+dq+"_"))then
1003                                                                  z_w = f_att->$vNam(varn+1)$
1004                                                                  break
1005                                                                else
1006                                                                  if(vNam(varn) .EQ. "wsqvs" .OR. isStrSubset(vNam(varn), "wsqvs_")\
1007                                                                      .OR. vNam(varn) .EQ. "w*qv*" .OR. isStrSubset(vNam(varn), "w*qv*_"))then
1008                                                                      z_w = f_att->$vNam(varn+1)$
1009                                                                      break
1010                                                                  else
1011                                                                      if(vNam(varn) .EQ. "wqv" .OR. isStrSubset(vNam(varn), "wqv_"))then
1012                                                                        z_w = f_att->$vNam(varn+1)$
1013                                                                        break
1014                                                                      else
1015                                                                        if(vNam(varn) .EQ. "wpsp" .OR. isStrSubset(vNam(varn), "wpsp_")\
1016                                                                            .OR. vNam(varn) .EQ. "w"+dq+"s"+dq .OR. isStrSubset(vNam(varn), "w"+dq+"s"+dq+"_"))then
1017                                                                            z_w = f_att->$vNam(varn+1)$
1018                                                                            break
1019                                                                        else
1020                                                                            if(vNam(varn) .EQ. "wsss" .OR. isStrSubset(vNam(varn), "wsss_")\
1021                                                                                .OR. vNam(varn) .EQ. "w*s*" .OR. isStrSubset(vNam(varn), "w*s*_"))then
1022                                                                              z_w = f_att->$vNam(varn+1)$
1023                                                                              break
1024                                                                            else
1025                                                                              if(vNam(varn) .EQ. "ws" .OR. isStrSubset(vNam(varn), "ws_"))then
1026                                                                                  z_w = f_att->$vNam(varn+1)$
1027                                                                                  break
1028                                                                              else
1029                                                                                  if(vNam(varn) .EQ. "wpsap" .OR. isStrSubset(vNam(varn), "wpsap_")\
1030                                                                                      .OR. vNam(varn) .EQ. "w"+dq+"sa"+dq .OR. isStrSubset(vNam(varn), "w"+dq+"sa"+dq+"_"))then
1031                                                                                    z_w = f_att->$vNam(varn+1)$
1032                                                                                    break
1033                                                                                  else
1034                                                                                    if(vNam(varn) .EQ. "wssas" .OR. isStrSubset(vNam(varn), "wssas_")\
1035                                                                                        .OR. vNam(varn) .EQ. "w*sa*" .OR. isStrSubset(vNam(varn), "w*sa*_"))then
1036                                                                                        z_w = f_att->$vNam(varn+1)$
1037                                                                                        break
1038                                                                                    else
1039                                                                                        if(vNam(varn) .EQ. "wsa" .OR. isStrSubset(vNam(varn), "wsa_"))then
1040                                                                                          z_w = f_att->$vNam(varn+1)$
1041                                                                                          break
1042                                                                                        else
1043                                                                                          if(vNam(varn) .EQ. "wses" .OR. isStrSubset(vNam(varn), "wses_")\
1044                                                                                              .OR. vNam(varn) .EQ. "w*e*" .OR. isStrSubset(vNam(varn), "w*e*_"))then
1045                                                                                              z_w = f_att->$vNam(varn+1)$
1046                                                                                              break
1047                                                                                          else
1048                                                                                              if(vNam(varn) .EQ. "ws2" .OR. isStrSubset(vNam(varn), "ws2_")\
1049                                                                                                .OR. vNam(varn) .EQ. "w*2" .OR. isStrSubset(vNam(varn), "w*2_"))then
1050                                                                                                z_w = f_att->$vNam(varn+1)$
1051                                                                                                break
1052                                                                                              else
1053                                                                                                if(vNam(varn) .EQ. "ws3" .OR. isStrSubset(vNam(varn), "ws3_")\
1054                                                                                                    .OR. vNam(varn) .EQ. "w*3" .OR. isStrSubset(vNam(varn), "w*3_"))then
1055                                                                                                    z_w = f_att->$vNam(varn+1)$
1056                                                                                                    break
1057                                                                                                else
1058                                                                                                    if(vNam(varn) .EQ. "Sw" .OR. isStrSubset(vNam(varn), "Sw_"))then
1059                                                                                                      z_w = f_att->$vNam(varn+1)$
1060                                                                                                      break
1061                                                                                                    else
1062                                                                                                      if(vNam(varn) .EQ. "ws2pts".OR. isStrSubset(vNam(varn), "ws2pts_")\
1063                                                                                                          .OR. vNam(varn) .EQ. "w*2pt*" .OR. isStrSubset(vNam(varn), "w*2pt*_"))then
1064                                                                                                          z_w = f_att->$vNam(varn+1)$
1065                                                                                                          break
1066                                                                                                      else
1067                                                                                                          if(vNam(varn) .EQ. "wspts2" .OR. isStrSubset(vNam(varn), "wspts2_")\
1068                                                                                                            .OR. vNam(varn) .EQ. "w*pt*2" .OR. isStrSubset(vNam(varn), "w*pt*2_"))then
1069                                                                                                            z_w = f_att->$vNam(varn+1)$
1070                                                                                                            break                                           
1071                                                                                                          end if
1072                                                                                                      end if
1073                                                                                                    end if
1074                                                                                                end if
1075                                                                                              end if
1076                                                                                          end if
1077                                                                                        end if
1078                                                                                    end if
1079                                                                                  end if
1080                                                                              end if
1081                                                                            end if   
1082                                                                        end if
1083                                                                      end if
1084                                                                  end if
1085                                                                end if
1086                                                            end if
1087                                                          end if
1088                                                      end if
1089                                                    end if
1090                                                end if
1091                                              end if
1092                                          end if
1093                                        end if
1094                                    end if
1095                                  end if
1096                              end if
1097                            end if
1098                        end if
1099                      end if
1100                  end if
1101                end if
1102            end if
1103          end if
1104      end do
1105     
1106      if (.not. isvar("z_u") .AND. .not. isvar ("z_w")) then
1107          co=0
1108          do varn=0,dim-1     
1109            if ( isStrSubset( vNam(varn), "time") .OR. \
1110                  isStrSubset( vNam(varn), "NORM")) then
1111                check = False
1112            else
1113                check = True
1114                if (var .NE. "all") then
1115                  check = isStrSubset( var,","+vNam(varn)+"," )
1116                end if         
1117            end if
1118            if (check)then
1119                co=co+1
1120                z = f_att->$vNam(varn+1)$
1121                if (getvardims(z) .EQ. "zu")then
1122                  z_u = z
1123                else
1124                  if (getvardims(z) .EQ. "zw")then
1125                      z_w = z
1126                  end if
1127                end if
1128                dimz = dimsizes(z)
1129                break
1130            end if   
1131          end do
1132          if (co .EQ. 0) then
1133            print(" ")
1134            print("The variables 'var="+var+"'" )
1135            print("do not exist on your input file;")
1136            print("be sure to have one comma before and after each variable")
1137            print(" ")
1138            exit           
1139          end if
1140      end if
1141
1142      if (isvar("z_u") ) then
1143          dimz  = dimsizes(z_u)
1144      else
1145          if (isvar("z_w"))then
1146            dimz  = dimsizes(z_w)
1147          end if
1148      end if
1149
1150      else
1151
1152          do varn = dim-1,0,1
1153            if (vNam(varn) .EQ. "zu_3d")then
1154                z_u = f_att->zu_3d 
1155                dimz  = dimsizes(z_u)         
1156            else
1157                if (vNam(varn) .EQ. "zw_3d")then
1158                  z_w = f_att->zw_3d
1159                  dimz  = dimsizes(z_w)
1160                end if
1161            end if
1162            if (vNam(varn) .EQ. "x")then
1163                x = f_att->x
1164                dimx=dimsizes(x)
1165            else
1166                if (vNam(varn) .EQ. "xu")then
1167                  x = f_att->xu
1168                  dimx=dimsizes(x)
1169                end if
1170            end if
1171            if (vNam(varn) .EQ. "y")then
1172                y = f_att->y
1173                dimy=dimsizes(y)
1174            else
1175                if (vNam(varn) .EQ. "yv")then
1176                  y = f_att->yv
1177                  dimy=dimsizes(y)
1178                end if
1179            end if
1180          end do
1181
1182      end if
1183
1184      if(.not. isvar("z_u") .AND. .not. isvar("z_w"))then
1185          print(" ")
1186          print("Program aborts - there are no z-variables available")
1187          print("Be sure if 'plot_3d' is set correctly")
1188          print(" ")
1189          exit
1190      end if
1191     
1192
1193      t_all = f[:]->time
1194      nt    = dimsizes(t_all)
1195
1196     
1197      if (nt .EQ. 1)then
1198          delta_t=t_all(nt-1)/nt
1199      else
1200          delta_t=(t_all(nt-1)-t_all(0))/(nt-1)
1201      end if
1202
1203      ; ****************************************************   
1204      ; start of time step and different types of mistakes that could be done
1205      ; ****************************************************
1206     
1207      if (start_time_step .EQ. -1.) then   
1208          delete(start_time_step)       
1209          start_time_step=t_all(0)/3600   
1210      else
1211          if (start_time_step .GT. t_all(nt-1)/3600)then
1212            print(" ")
1213            print("'start_time_step' = "+ start_time_step +"h is greater "+\
1214                  "than last time step = " \
1215                  + t_all(nt-1)+"s = "+t_all(nt-1)/3600+"h")
1216            print(" ")
1217            print("Select another 'start_time_step'")
1218            print(" ")
1219            exit
1220          end if
1221          if (start_time_step .LT. t_all(0)/3600)then
1222            print(" ")
1223            print("'start_time_step' = "+ start_time_step +"h is lower "+\
1224                  "than first time step = " + t_all(0)+"s = "+t_all(0)/3600+"h")
1225            print(" ")
1226            exit
1227          end if
1228      end if
1229
1230      do i=0,nt-1   
1231          if (start_time_step .GE. (t_all(i)-delta_t/2)/3600 .AND.\
1232              start_time_step .LT. (t_all(i)+delta_t/2)/3600)then
1233            st=i
1234            break
1235          else
1236            st=0
1237          end if
1238      end do
1239     
1240      ; ****************************************************
1241      ; end of time step and different types of mistakes that could be done
1242      ; ****************************************************
1243
1244      if (end_time_step .EQ. -1.) then         
1245          delete(end_time_step)
1246          end_time_step = t_all(nt-1)/3600 
1247      else
1248          if (end_time_step .GT. t_all(nt-1)/3600)then
1249            print(" ")
1250            print("'end_time_step' = "+ end_time_step +"h is greater "+\
1251                  "than last time step = " +\
1252                    t_all(nt-1)+"s = "+t_all(nt-1)/3600+"h")
1253            print(" ")
1254            print("Select another 'end_time_step'") 
1255            print(" ")
1256            exit
1257          end if
1258          if (end_time_step .LT. start_time_step/3600)then
1259            print(" ")
1260            print("'end_time_step' = "+ end_time_step +"h is lower "+\
1261                  "than 'start_time_step' = "+start_time_step+"h")
1262            print(" ")
1263            print("Select another 'start_time_step' or 'end_time_step'")
1264            print(" ")
1265            exit
1266          end if
1267      end if
1268
1269      do i=0,nt-1     
1270          if (end_time_step .GE. (t_all(i)-delta_t/2)/3600 .AND.\
1271              end_time_step .LT. (t_all(i)+delta_t/2)/3600)then
1272            et=i
1273            break
1274          else
1275            et=0
1276          end if
1277      end do
1278     
1279      delete(start_time_step)
1280      start_time_step=round(st,3)
1281      delete(end_time_step)
1282      end_time_step=round(et,3)
1283
1284      no_time = end_time_step-start_time_step+1
1285     
1286      ; ****************************************************
1287      ; time_stride and different types of mistakes that could be done
1288      ; ****************************************************
1289   
1290      if (time_stride .LT. 1) then
1291          print(" ")
1292          print("'time_stride' has to be positive and is set to 1")
1293          print(" ")
1294          time_stride = 1
1295      end if
1296
1297      if (time_stride .GT. no_time) then
1298          print(" ")
1299          print("'time_stride' is greater than number of available "+\
1300              "time steps,")
1301          print("=> 'time_stride' is set to 1")
1302          time_stride = 1
1303      end if
1304
1305      ti_in = ispan(start_time_step,end_time_step,time_stride) ;ti_in contents
1306                                                                ;the time indices
1307                                                                ;to plot
1308      np    = dimsizes(ti_in) 
1309
1310      if (com_i .EQ. max_com_i) then
1311          print(" ")
1312          print("Output of time steps from "+t_all(start_time_step)/3600+\
1313            " h = "+t_all(start_time_step)+" s => index = "+start_time_step)
1314          print("                     till "+t_all(ti_in(np-1))/3600+" h = "\
1315            +t_all(ti_in(np-1))+" s => index = "+end_time_step)
1316          print("                     with temporal stride = "+time_stride)
1317          print(" ")
1318      end if
1319
1320      ; ****************************************************
1321      ; set up legend and colors
1322      ; ****************************************************
1323     
1324      legend_label=new(np,string)
1325      do p=0, np-1
1326          legend_label(p)=sprintf("%6.2f", t_all(ti_in(p))/3600)
1327      end do
1328
1329      ; ***************************************************
1330      ; set up recourses
1331      ; ***************************************************
1332
1333      res                         = True
1334      res@gsnDraw                 = False
1335      res@gsnFrame                = False
1336      res@txFont                  = "helvetica"
1337      res@tiMainFont              = "helvetica"
1338      res@tiXAxisFont             = "helvetica"
1339      res@tiYAxisFont             = "helvetica"
1340      res@tmXBLabelFont           = "helvetica"
1341      res@tmYLLabelFont           = "helvetica"
1342      res@lgLabelFont             = "helvetica"
1343      res@tmLabelAutoStride       = True
1344
1345
1346      if (legend .EQ. 1)then
1347          res@pmLegendDisplayMode  = "Always"
1348      end if
1349
1350      res@pmLegendSide            = "Right"
1351      res@xyExplicitLegendLabels  = legend_label
1352      res@pmLegendWidthF          = 0.12
1353      res@pmLegendHeightF         = 0.05*np
1354      res@pmLegendParallelPosF    = 1.0
1355      res@pmLegendOrthogonalPosF  = -0.2625;-1.0
1356      res@lgLabelFontHeightF      = font_size_legend
1357      res@lgTitleString           = "Time (h)"
1358      res@lgTitleFontHeightF      = font_size   
1359      res@txFontHeightF           = font_size
1360      res@tiXAxisFontHeightF      = font_size
1361      res@tiYAxisFontHeightF      = font_size
1362      res@tmXBLabelFontHeightF    = font_size
1363      res@tmYLLabelFontHeightF    = font_size
1364      res@tiXAxisString           = " "
1365      res@lgJustification         = "TopLeft"
1366
1367      if ( black .eq. 0 ) then 
1368          res@xyLineColors = -(ispan(-237,-2,235/np))
1369      end if
1370      if (norm_z .EQ. 1)then
1371          res@tiYAxisString = "Height (m)"
1372      else   
1373          res@tiYAxisString = "Height / "+norm_z+" (m)"
1374      end if
1375       
1376      if (log_z .EQ. 1) then
1377          res@trYLog = True
1378      end if
1379
1380      if (dash .EQ. 0 ) then
1381          res@xyMonoDashPattern = True
1382      else
1383          res@xyMonoDashPattern = False
1384          if (no_files .GT. 1)
1385            res@xyMonoDashPattern = True 
1386            print(" ")
1387            print("If you use more than one file, patterns for different "+\
1388                  "timesteps cannot be used")
1389            print(" ")
1390          end if       
1391      end if
1392
1393      res@tmXBMinorPerMajor = 4
1394      res@tmYLMinorPerMajor = 4
1395
1396      resP                            = True
1397      resP@gsnMaximize                = True
1398      resP@gsnPanelXWhiteSpacePercent = 3.0
1399      resP@gsnPanelYWhiteSpacePercent = 2.0
1400      resP@txFont                     = "helvetica"
1401      resP@txString                   = f_att@title
1402      resP@txFuncCode                 = "~"
1403      resP@txFontHeightF              = 0.015
1404
1405      ; ***************************************************
1406      ; set up graphics for plot
1407      ; ***************************************************
1408
1409      if (combine .EQ. 1) then
1410          plot_o = new(number_comb,graphic) 
1411          label=new(number_comb,string)
1412          color_o=new(number_comb,integer)
1413
1414          if (check_vType) then
1415            mini=new(number_comb,double)
1416            maxi=new(number_comb,double)
1417          else
1418            mini=new(number_comb,float)
1419            maxi=new(number_comb,float)
1420          end if
1421      end if
1422   
1423      if ( format_out .EQ. "pdf" .OR. format_out .EQ. "ps" ) then
1424          format_out@wkPaperSize = "A4"
1425      end if
1426      if ( format_out .EQ. "png" ) then
1427          format_out@wkWidth  = 1000
1428          format_out@wkHeight = 1000
1429      end if
1430
1431      if (wks_gsn_log) then
1432          wks=gsn_open_wks(format_out,file_out)
1433          gsn_define_colormap(wks,"rainbow+white")
1434          wks_gsn_log = False
1435      end if
1436
1437      ; ***************************************************
1438      ; set up minimum and maximum height
1439      ; ***************************************************
1440
1441      if (log_z .EQ. 1)then
1442          if (min_z .EQ. -1)then
1443            if (isvar("z_u"))then
1444                min_z=z_u(1)
1445            else
1446                min_z=z_w(1)
1447            end if       
1448          else
1449            if (isvar("z_w"))then
1450                if (min_z .GE. max(z_w) ) then
1451                  print(" ")
1452                  print("Minimum of height ('min_z'="+min_z+") is greater "+\
1453                        "than available heights (="+max(z_w)+")")
1454                  print(" ")
1455                  exit
1456                end if     
1457            else
1458                if (min_z .GE. max(z_u) ) then
1459                  print(" ")
1460                  print("Minimum of height ('min_z'="+min_z+") is greater "+\
1461                        "than available heights (="+max(z_u)+")")
1462                  print(" ")
1463                  exit
1464                end if
1465            end if
1466            if (isvar("z_u"))then   
1467                if (min_z .LT. z_u(1) ) then
1468                  print(" ")
1469                  print("Begin height 'min_z' at least at level k=1 (="+\
1470                        z_u(1)+"m) due to the logarithmic scale of the y-axis")
1471                  print(" ")
1472                  exit
1473                end if
1474            else
1475                if (min_z .LT. z_w(1) ) then
1476                  print(" ")
1477                  print("Begin height 'min_z' at least at level k=1 (="+\
1478                        z_w(1)+"m) due to the logarithmic scale of the y-axis")
1479                  print(" ")
1480                  exit
1481                end if   
1482            end if
1483          end if
1484      else
1485          if (isvar("z_u"))then
1486            if (min_z .EQ. -1)then
1487                min_z=z_u(0)
1488            end if
1489          else
1490            if (min_z .EQ. -1)then
1491                min_z=z_w(0)
1492            end if   
1493          end if
1494     
1495          if (isvar("z_w"))then
1496            if (min_z .GE. max(z_w) ) then
1497                print(" ")
1498                print("Minimum of height ('min_z'="+min_z+") is greater than "+\
1499                      "available heights (="+max(z_w)+")")
1500                print(" ")
1501                exit
1502            end if
1503          else
1504            if (min_z .EQ. -1)then
1505                min_z=z_u(0)
1506            end if
1507            if (min_z .GE. max(z_u) ) then
1508                print(" ")
1509                print("Minimum of height ('min_z'="+min_z+") is greater than "+\
1510                      "available heights (="+max(z_u)+")")
1511                print(" ")
1512                exit
1513            end if
1514          end if
1515      end if
1516     
1517      if (isvar("z_w"))then
1518          if (max_z .EQ. -1)then
1519            max_z=max(z_w)
1520          end if
1521      else
1522          if (max_z .EQ. -1)then
1523            max_z=max(z_u)
1524          end if   
1525      end if
1526     
1527      if (isvar("z_w"))then
1528          if (max_z .GT. max(z_w) ) then
1529            print(" ")
1530            print("Maximum of height ('max_z'="+max_z+") is greater than "+\
1531                  "available heights (="+max(z_w)+")")
1532            print(" ")
1533            exit
1534          end if
1535      end if
1536
1537      min_z=min_z/norm_z
1538      max_z=max_z/norm_z
1539
1540      ; ***************************************************
1541      ; read data and create plots
1542      ; ***************************************************
1543         
1544      do ti =0,np-1
1545          if( t_all(ti_in(ti)) .lt. 10^36) then
1546            start_time_step = ti
1547            break
1548          end if
1549      end do 
1550     
1551      if (log_z .EQ. 1) then     
1552          if (check_vType) then
1553            data   = new((/dim,np,dimz-1/),double)
1554            data_0 = new((/np,dimz-1/),double)
1555          else
1556            data   = new((/dim,np,dimz-1/),float)
1557            data_0 = new((/np,dimz-1/),float)
1558          end if
1559          data@_FillValue=9.96921e+36
1560          data_0 = 0.1
1561          t      = new((/np,dimz-1/),float)
1562          t      = 0.0
1563          unit   = new(dim,string)
1564          if (isvar("z_u"))then
1565            if (typeof(z_u) .EQ. "double")then
1566                z_v    = new((/dim,dimz/),double)
1567                z_     = new((/dim,dimz-1/),double)
1568            else
1569                if (typeof(z_u) .EQ. "float")then
1570                  z_v    = new((/dim,dimz/),float)
1571                  z_     = new((/dim,dimz-1/),float)
1572                end if
1573            end if
1574          else
1575            if (isvar("z_w"))then
1576                if (typeof(z_w) .EQ. "double")then
1577                  z_v    = new((/dim,dimz/),double)
1578                  z_     = new((/dim,dimz-1/),double)
1579                else
1580                  if (typeof(z_w) .EQ. "float")then
1581                      z_v    = new((/dim,dimz/),float)
1582                      z_     = new((/dim,dimz-1/),float)
1583                  end if
1584                end if
1585            end if
1586          end if
1587      else
1588          if (check_vType) then
1589            data   = new((/dim,np,dimz/),double)
1590            data_0 = new((/np,dimz/),double)
1591          else
1592            data   = new((/dim,np,dimz/),float)
1593            data_0 = new((/np,dimz/),float)
1594          end if
1595          data@_FillValue=9.96921e+36     
1596          data_0 = 0.0
1597          t      = new((/np,dimz/),float)
1598          t      = 0.0
1599          unit   = new(dim,string)
1600          if (isvar("z_u"))then
1601            if (typeof(z_u) .EQ. "double")then
1602                z_v    = new((/dim,dimz/),double)
1603                z_     = new((/dim,dimz/),double)
1604            else
1605                if (typeof(z_u) .EQ. "float")then
1606                  z_v    = new((/dim,dimz/),float)
1607                  z_     = new((/dim,dimz/),float)
1608                end if
1609            end if
1610          else
1611            if (isvar("z_w"))then
1612                if (typeof(z_w) .EQ. "double")then
1613                  z_v    = new((/dim,dimz/),double)
1614                  z_     = new((/dim,dimz/),double)
1615                else
1616                  if (typeof(z_w) .EQ. "float")then
1617                      z_v    = new((/dim,dimz/),float)
1618                      z_     = new((/dim,dimz/),float)
1619                  end if
1620                end if
1621            end if
1622          end if
1623      end if
1624
1625      end if
1626      ;-------above steps only for first file
1627
1628      ; ***************************************************
1629      ; indicate plot number
1630      ; ***************************************************
1631     
1632      if (combine .EQ. 1) then
1633          n = 1
1634      else
1635          n = 0
1636      end if
1637
1638      if (over .EQ. 1) then
1639          plot_u         = gsn_csm_xy(wks,t,data_0(:,:),res)
1640          plot_v         = gsn_csm_xy(wks,t,data_0(:,:),res)
1641          plot_w         = gsn_csm_xy(wks,t,data_0(:,:),res)
1642          plot_pt        = gsn_csm_xy(wks,t,data_0(:,:),res)     
1643          plot_vpt       = gsn_csm_xy(wks,t,data_0(:,:),res)
1644          plot_lpt       = gsn_csm_xy(wks,t,data_0(:,:),res)
1645          plot_q         = gsn_csm_xy(wks,t,data_0(:,:),res)
1646          plot_qv        = gsn_csm_xy(wks,t,data_0(:,:),res)
1647          plot_ql        = gsn_csm_xy(wks,t,data_0(:,:),res)
1648          plot_rho       = gsn_csm_xy(wks,t,data_0(:,:),res)
1649          plot_s         = gsn_csm_xy(wks,t,data_0(:,:),res)
1650          plot_sa        = gsn_csm_xy(wks,t,data_0(:,:),res)
1651          plot_e         = gsn_csm_xy(wks,t,data_0(:,:),res)
1652          plot_es        = gsn_csm_xy(wks,t,data_0(:,:),res)
1653          plot_km        = gsn_csm_xy(wks,t,data_0(:,:),res)
1654          plot_kh        = gsn_csm_xy(wks,t,data_0(:,:),res)
1655          plot_l         = gsn_csm_xy(wks,t,data_0(:,:),res)     
1656          plot_wpup      = gsn_csm_xy(wks,t,data_0(:,:),res)
1657          plot_wsus      = gsn_csm_xy(wks,t,data_0(:,:),res)
1658          plot_wu        = gsn_csm_xy(wks,t,data_0(:,:),res)
1659          plot_wpvp      = gsn_csm_xy(wks,t,data_0(:,:),res)
1660          plot_wsvs      = gsn_csm_xy(wks,t,data_0(:,:),res)
1661          plot_wv        = gsn_csm_xy(wks,t,data_0(:,:),res)
1662          plot_wpptp     = gsn_csm_xy(wks,t,data_0(:,:),res)
1663          plot_wspts     = gsn_csm_xy(wks,t,data_0(:,:),res)
1664          plot_wpt       = gsn_csm_xy(wks,t,data_0(:,:),res)
1665          plot_wsptsBC   = gsn_csm_xy(wks,t,data_0(:,:),res)
1666          plot_wptBC     = gsn_csm_xy(wks,t,data_0(:,:),res)
1667          plot_wpvptp    = gsn_csm_xy(wks,t,data_0(:,:),res)
1668          plot_wsvpts    = gsn_csm_xy(wks,t,data_0(:,:),res)
1669          plot_wvpt      = gsn_csm_xy(wks,t,data_0(:,:),res)
1670          plot_wpqp      = gsn_csm_xy(wks,t,data_0(:,:),res)
1671          plot_wsqs      = gsn_csm_xy(wks,t,data_0(:,:),res)
1672          plot_wq        = gsn_csm_xy(wks,t,data_0(:,:),res)
1673          plot_wpqvp     = gsn_csm_xy(wks,t,data_0(:,:),res)
1674          plot_wsqvs     = gsn_csm_xy(wks,t,data_0(:,:),res)
1675          plot_wqv       = gsn_csm_xy(wks,t,data_0(:,:),res)
1676          plot_wpsp      = gsn_csm_xy(wks,t,data_0(:,:),res)
1677          plot_wsss      = gsn_csm_xy(wks,t,data_0(:,:),res)
1678          plot_ws        = gsn_csm_xy(wks,t,data_0(:,:),res)
1679          plot_wpsap     = gsn_csm_xy(wks,t,data_0(:,:),res)
1680          plot_wssas     = gsn_csm_xy(wks,t,data_0(:,:),res)
1681          plot_wsa       = gsn_csm_xy(wks,t,data_0(:,:),res)
1682          plot_wses      = gsn_csm_xy(wks,t,data_0(:,:),res)
1683          plot_us2       = gsn_csm_xy(wks,t,data_0(:,:),res)
1684          plot_vs2       = gsn_csm_xy(wks,t,data_0(:,:),res)
1685          plot_ws2       = gsn_csm_xy(wks,t,data_0(:,:),res)
1686          plot_pts2      = gsn_csm_xy(wks,t,data_0(:,:),res)
1687          plot_ws3       = gsn_csm_xy(wks,t,data_0(:,:),res)
1688          plot_Sw        = gsn_csm_xy(wks,t,data_0(:,:),res)
1689          plot_ws2pts    = gsn_csm_xy(wks,t,data_0(:,:),res)
1690          plot_wspts2    = gsn_csm_xy(wks,t,data_0(:,:),res)
1691          plot_wsususodz = gsn_csm_xy(wks,t,data_0(:,:),res)
1692          plot_wspsodz   = gsn_csm_xy(wks,t,data_0(:,:),res)
1693          plot_wpeodz    = gsn_csm_xy(wks,t,data_0(:,:),res)
1694     
1695          if (check_vType) then
1696            miniu         =  100000.d
1697            maxiu         = -100000.d
1698            miniv         =  100000.d
1699            maxiv         = -100000.d
1700            miniw         =  100000.d
1701            maxiw         = -100000.d
1702            minipt        =  100000.d
1703            maxipt        = -100000.d
1704            minivpt       =  100000.d
1705            maxivpt       = -100000.d
1706            minilpt       =  100000.d
1707            maxilpt       = -100000.d
1708            miniq         =  100000.d
1709            maxiq         = -100000.d
1710            miniqv        =  100000.d
1711            maxiqv        = -100000.d
1712            miniql        =  100000.d
1713            maxiql        = -100000.d
1714            minie         =  100000.d
1715            maxie         = -100000.d
1716            minies        =  100000.d
1717            maxies        = -100000.d
1718            minikm        =  100000.d
1719            maxikm        = -100000.d
1720            minikh        =  100000.d
1721            maxikh        = -100000.d
1722            miniwpup      =  100000.d
1723            maxiwpup      = -100000.d
1724            miniwsus      =  100000.d
1725            maxiwsus      = -100000.d
1726            miniwu        =  100000.d
1727            maxiwu        = -100000.d
1728            miniwpvp      =  100000.d
1729            maxiwpvp      = -100000.d
1730            miniwsvs      =  100000.d
1731            maxiwsvs      = -100000.d
1732            miniwv        =  100000.d
1733            maxiwv        = -100000.d
1734            miniwpptp     =  100000.d
1735            maxiwpptp     = -100000.d
1736            miniwspts     =  100000.d
1737            maxiwspts     = -100000.d
1738            miniwpt       =  100000.d
1739            maxiwpt       = -100000.d
1740            miniwsptsBC   =  100000.d
1741            maxiwsptsBC   = -100000.d
1742            miniwptBC     =  100000.d
1743            maxiwptBC     = -100000.d
1744            miniwpvptp    =  100000.d
1745            maxiwpvptp    = -100000.d
1746            miniwsvpts    =  100000.d
1747            maxiwsvpts    = -100000.d
1748            miniwvpt      =  100000.d
1749            maxiwvpt      = -100000.d
1750            miniwpqp      =  100000.d
1751            maxiwpqp      = -100000.d
1752            miniwsqs      =  100000.d
1753            maxiwsqs      = -100000.d
1754            miniwq        =  100000.d
1755            maxiwq        = -100000.d
1756            miniwpqvp     =  100000.d
1757            maxiwpqvp     = -100000.d
1758            miniwsqvs     =  100000.d
1759            maxiwsqvs     = -100000.d
1760            miniwqv       =  100000.d
1761            maxiwqv       = -100000.d
1762            miniwpsp      =  100000.d
1763            maxiwpsp      = -100000.d
1764            miniwsss      =  100000.d
1765            maxiwsss      = -100000.d
1766            miniws        =  100000.d
1767            maxiws        = -100000.d
1768            miniwpsap     =  100000.d
1769            maxiwpsap     = -100000.d
1770            miniwssas     =  100000.d
1771            maxiwssas     = -100000.d
1772            miniwsa       =  100000.d
1773            maxiwsa       = -100000.d
1774            minius2       =  100000.d
1775            maxius2       = -100000.d
1776            minivs2       =  100000.d
1777            maxivs2       = -100000.d
1778            miniws2       =  100000.d
1779            maxiws2       = -100000.d
1780            miniwsususodz =  100000.d
1781            maxiwsususodz = -100000.d
1782            miniwspsodz   =  100000.d
1783            maxiwspsodz   = -100000.d
1784            miniwpeodz    =  100000.d
1785            maxiwpeodz    = -100000.d
1786          else
1787            miniu         =  100000.
1788            maxiu         = -100000.
1789            miniv         =  100000.
1790            maxiv         = -100000.
1791            miniw         =  100000.
1792            maxiw         = -100000.
1793            minipt        =  100000.
1794            maxipt        = -100000.
1795            minivpt       =  100000.
1796            maxivpt       = -100000.
1797            minilpt       =  100000.
1798            maxilpt       = -100000.
1799            miniq         =  100000.
1800            maxiq         = -100000.
1801            miniqv        =  100000.
1802            maxiqv        = -100000.
1803            miniql        =  100000.
1804            maxiql        = -100000.
1805            minie         =  100000.
1806            maxie         = -100000.
1807            minies        =  100000.
1808            maxies        = -100000.
1809            minikm        =  100000.
1810            maxikm        = -100000.
1811            minikh        =  100000.
1812            maxikh        = -100000.
1813            miniwpup      =  100000.
1814            maxiwpup      = -100000.
1815            miniwsus      =  100000.
1816            maxiwsus      = -100000.
1817            miniwu        =  100000.
1818            maxiwu        = -100000.
1819            miniwpvp      =  100000.
1820            maxiwpvp      = -100000.
1821            miniwsvs      =  100000.
1822            maxiwsvs      = -100000.
1823            miniwv        =  100000.
1824            maxiwv        = -100000.
1825            miniwpptp     =  100000.
1826            maxiwpptp     = -100000.
1827            miniwspts     =  100000.
1828            maxiwspts     = -100000.
1829            miniwpt       =  100000.
1830            maxiwpt       = -100000.
1831            miniwsptsBC   =  100000.
1832            maxiwsptsBC   = -100000.
1833            miniwptBC     =  100000.
1834            maxiwptBC     = -100000.
1835            miniwpvptp    =  100000.
1836            maxiwpvptp    = -100000.
1837            miniwsvpts    =  100000.
1838            maxiwsvpts    = -100000.
1839            miniwvpt      =  100000.
1840            maxiwvpt      = -100000.
1841            miniwpqp      =  100000.
1842            maxiwpqp      = -100000.
1843            miniwsqs      =  100000.
1844            maxiwsqs      = -100000.
1845            miniwq        =  100000.
1846            maxiwq        = -100000.
1847            miniwpqvp     =  100000.
1848            maxiwpqvp     = -100000.
1849            miniwsqvs     =  100000.
1850            maxiwsqvs     = -100000.
1851            miniwqv       =  100000.
1852            maxiwqv       = -100000.
1853            miniwpsp      =  100000.
1854            maxiwpsp      = -100000.
1855            miniwsss      =  100000.
1856            maxiwsss      = -100000.
1857            miniws        =  100000.
1858            maxiws        = -100000.
1859            miniwpsap     =  100000.
1860            maxiwpsap     = -100000.
1861            miniwssas     =  100000.
1862            maxiwssas     = -100000.
1863            miniwsa       =  100000.
1864            maxiwsa       = -100000.
1865            minius2       =  100000.
1866            maxius2       = -100000.
1867            minivs2       =  100000.
1868            maxivs2       = -100000.
1869            miniws2       =  100000.
1870            maxiws2       = -100000.
1871            miniwsususodz =  100000.
1872            maxiwsususodz = -100000.
1873            miniwspsodz   =  100000.
1874            maxiwspsodz   = -100000.
1875            miniwpeodz    =  100000.
1876            maxiwpeodz    = -100000.
1877          end if
1878
1879      end if
1880
1881      if (prof3d .EQ. 1)then
1882
1883      if (end_x .EQ. -1) then
1884          end_x=dimx-2
1885      end if
1886      if (end_y .EQ. -1)then
1887          end_y=dimy-2
1888      end if
1889      if (start_x .LT. 0)then
1890          print(" ")
1891          print("'start_x' is lower than 0 and set to 0")
1892          print(" ")
1893          start_x=0
1894      end if
1895      if (start_x .GT. dimx-1)then
1896          print(" ")
1897          print("'start_x' is greater than available x-range and set to "+\
1898                "maximum of x-range (excluding ghostpoint)")
1899          print(" ")
1900          start_x=dimx-2
1901      end if
1902      if (end_x .EQ. dimx-1)then
1903          print(" ")
1904          print("'end_x' = "+end_x+" and includes the ghostpoint")
1905          print(" ")
1906      end if
1907      if (end_x .GT. dimx-1)then
1908          print(" ")
1909          print("'end_x' = "+end_x+" is greater than available x-range and set "+\
1910                "to maximum of x-range (excluding ghostpoint)")
1911          print(" ")
1912          end_x=dimx-2
1913      end if
1914      if (end_x .LT. start_x)then
1915          print(" ")
1916          print("'end_x' = "+end_x+" is lower than 'start_x' = "+start_x+\
1917                " and set to maximum of x-range (excluding ghostpoint)")
1918          print(" ")
1919          end_x=dimx-2
1920      end if
1921      if (start_y .LT. 0)then
1922          print(" ")
1923          print("'start_y' is lower than 0 and set to 0")
1924          print(" ")
1925          start_y=0
1926      end if
1927      if (start_y .GT. dimy-1)then
1928          print(" ")
1929          print("'start_y' is greater than available y-range and set to "+\
1930                "maximum of y-range (excluding ghostpoint)")
1931          print(" ")
1932          start_x=dimy-2
1933      end if
1934      if (end_y .EQ. dimy-1)then
1935          print(" ")
1936          print("'end_y' = "+end_y+" and includes the ghostpoint")
1937          print(" ")
1938      end if
1939      if (end_y .GT. dimy-1)then
1940          print(" ")
1941          print("'end_y' = "+end_y+" is greater than available y-range and "+\
1942                "set to maximum of y-range (excluding ghostpoint)")
1943          print(" ")
1944          end_x=dimy-2
1945      end if
1946      if (end_y .LT. start_y)then
1947          print(" ")
1948          print("'end_y' = "+end_y+" is lower than 'start_y' = "+start_y+\
1949                " and set to maximum of y-range (excluding ghostpoint)")
1950          print(" ")
1951          end_y=dimy-2
1952      end if
1953     
1954      end if
1955     
1956      n_o=0
1957      count_var=0
1958
1959      res@xyDashPattern = 1*nof
1960
1961      over_remind = False
1962      if ( over .eq. 1)then
1963        over_remind = True
1964      end if
1965
1966      do varn = 0,dim-1
1967
1968          check = True
1969         
1970          if (prof3d .EQ. 0) then
1971            if ( isStrSubset( vNam(varn), "time") .OR. \
1972                  isStrSubset( vNam(varn), "NORM")) then
1973                check = False
1974            end if
1975          else
1976            if ( isStrSubset( vNam(varn), "time") .OR.  \
1977                  isStrSubset( vNam(varn), "zusi") .OR.  \
1978                  isStrSubset( vNam(varn), "zwwi") .OR.  \
1979                  isStrSubset( vNam(varn), "x") .OR.     \
1980                  isStrSubset( vNam(varn), "xu") .OR.    \
1981                  isStrSubset( vNam(varn), "y") .OR.     \
1982                  isStrSubset( vNam(varn), "yv") .OR.    \
1983                  isStrSubset( vNam(varn), "zu_3d") .OR. \
1984                  isStrSubset( vNam(varn), "zw_3d")) then
1985                check = False
1986            end if
1987          end if
1988
1989          if (var .NE. "all") then
1990            check = isStrSubset( var,","+vNam(varn)+"," )
1991          end if
1992
1993          if (combine .EQ. 1) then         
1994            com=isStrSubset(c_var,","+vNam(varn)+"," )     
1995            if (com) then     
1996                if (prof3d .EQ. 0) then             
1997                  temp     = f[:]->$vNam(varn)$       
1998                  temp_att = f_att->$vNam(varn)$
1999                  if (log_z .EQ. 1) then
2000                      do j=0,np-1
2001                        data(varn,j,:) = temp(ti_in(j),1:dimz-1)
2002                      end do
2003                  else
2004                      do j=0,np-1
2005                        data(varn,j,:) = temp(ti_in(j),0:dimz-1)
2006                      end do
2007                  end if
2008                else
2009                  if (log_z .EQ. 1) then
2010                      do i=1,dimz-1
2011                        do j=0,np-1
2012                            temp = f[:]->$vNam(varn)$
2013                            temp_att = f_att->$vNam(varn)$
2014                            data_temp = temp(ti_in(j),i,\
2015                                            start_y:end_y,start_x:end_x)
2016                            data(varn,j,i-1) = dim_avg_Wrap(\
2017                                                      dim_avg_Wrap(data_temp))
2018                        end do
2019                      end do
2020                  else
2021                      do i=0,dimz-1
2022                        do j=0,np-1
2023                            temp = f[:]->$vNam(varn)$
2024                            temp_att = f_att->$vNam(varn)$
2025                            data_temp = temp(ti_in(j),i,\
2026                                                start_y:end_y,start_x:end_x)
2027                            data(varn,j,i) = dim_avg_Wrap(dim_avg_Wrap(data_temp))
2028                        end do
2029                      end do
2030                  end if
2031                  print(" ")
2032                  print("Variable for combine '"+vNam(varn)+"' is read")
2033                  print(" ")
2034                end if                 
2035                unit(varn) = temp_att@units
2036                if (n_o .GT. number_comb-1) then
2037                  print(" ")
2038                  print("Set 'number_comb' to the number of overlaying "+\
2039                        "variables ('c_var' = "+c_var+")")
2040                  print(" ")
2041                  exit
2042                end if
2043                if (abs(min(data(varn,:,:))) .GT. 10)then
2044                  min_value = abs(0.01*min(data(varn,:,:)))
2045                  max_value = abs(0.01*max(data(varn,:,:)))
2046                else
2047                  if (abs(min(data(varn,:,:))) .LT. 0.01 .AND. \
2048                      abs(max(data(varn,:,:))) .GT. 0.01)then
2049                      min_value = abs(0.1*max(data(varn,:,:)))
2050                      max_value = abs(0.1*max(data(varn,:,:)))
2051                  else
2052                      if (abs(max(data(varn,:,:))) .LT. 0.01 .AND. \
2053                          abs(min(data(varn,:,:))) .GT. 0.01)then
2054                        min_value = abs(0.1*min(data(varn,:,:)))
2055                        max_value = abs(0.1*min(data(varn,:,:)))
2056                      else
2057                        min_value = abs(0.1*min(data(varn,:,:)))
2058                        max_value = abs(0.1*max(data(varn,:,:)))
2059                      end if
2060                  end if
2061                end if
2062                if (min(data(varn,:,:)) .EQ. 0 .AND. \
2063                    max(data(varn,:,:)) .EQ. 0)then
2064                  min_value = 0.1
2065                  max_value = 0.1
2066                end if
2067                mini(n_o)=min(data(varn,:,:))
2068                maxi(n_o)=max(data(varn,:,:))
2069                n_o=n_o+1
2070            end if
2071          end if
2072
2073          if(check) then
2074           
2075            count_var=count_var+1
2076
2077            if (prof3d .EQ. 0) then       
2078                temp = f[:]->$vNam(varn)$
2079                temp_att = f_att->$vNam(varn)$
2080            else
2081                if (log_z .EQ. 1) then
2082                  do i=1,dimz-1
2083                      do j=0,np-1
2084                        temp= f[:]->$vNam(varn)$
2085                        temp_att = f_att->$vNam(varn)$
2086                        data_temp = temp(ti_in(j),i,start_y:end_y,start_x:end_x)
2087                        data(varn,j,i-1) = dim_avg_Wrap(dim_avg_Wrap(data_temp))
2088                        delete(data_temp)
2089                      end do
2090                  end do
2091                else
2092                  do i=0,dimz-1
2093                      do j=0,np-1
2094                        temp= f[:]->$vNam(varn)$
2095                        temp_att = f_att->$vNam(varn)$
2096                        data_temp = temp(ti_in(j),i,start_y:end_y,start_x:end_x)
2097                        data_temp!0 = "t"
2098                        data_temp!1 = "z"
2099                        data(varn,j,i) = dim_avg_Wrap(dim_avg_Wrap(data_temp))
2100                        delete(data_temp)
2101                      end do
2102                  end do
2103                end if
2104                print(" ")
2105                print("Variable '"+vNam(varn)+"' is read")
2106                print(" ")
2107                unit(varn) = temp_att@units
2108                a=getvardims(temp_att)
2109                b=dimsizes(a)
2110            end if
2111
2112            if (prof3d .EQ. 0) then
2113                if (log_z .EQ. 1) then
2114                  z = f_att->$vNam(varn+1)$(1:dimz-1)
2115                  unit(varn) = temp_att@units
2116                  do j=0,np-1
2117                      data(varn,j,:) = temp(ti_in(j),1:dimz-1)
2118                  end do
2119                else
2120                  z = f_att->$vNam(varn+1)$
2121                  unit(varn) = temp_att@units
2122                  do j=0,np-1
2123                      data(varn,j,:) = temp(ti_in(j),:)
2124                  end do
2125                end if
2126            else
2127                do i=0,b-1           
2128                  if (isStrSubset( a(i),"zu_3d" ))then
2129                      z_v(varn,:) = z_u
2130                      if (log_z .EQ. 1) then
2131                        z = z_v(varn,1:dimz-1)
2132                      else
2133                        z = z_v(varn,:)
2134                      end if
2135                  else
2136                      if (isStrSubset( a(i),"zw_3d" ))then
2137                        z_v(varn,:) = z_w
2138                        if (log_z .EQ. 1) then
2139                            z = z_v(varn,1:dimz-1)
2140                        else
2141                            z = z_v(varn,:)
2142                        end if
2143                      end if                   
2144                  end if
2145                end do           
2146            end if
2147
2148            if (isStrSubset(data@long_name," SR " ) ) then
2149              over = 0
2150            end if
2151           
2152            if (nof .EQ. 0) then
2153                z_(n,:)=z/norm_z
2154                z    = z_(n,:)
2155            else
2156                z=z/norm_z
2157            end if
2158
2159            delta_z = z(2) - z(1)
2160
2161            max_z_int=doubletoint(max_z/delta_z)
2162            min_z_int=doubletoint(min_z/delta_z)
2163
2164            if(max_z_int .eq. min_z_int)
2165                print(" ")
2166                print("Please increase 'max_z' or decrease 'min_z' so that "+\
2167                      "there are")
2168                print("at least two layers for the z-axis to plot")
2169                print(" ")
2170                exit
2171            end if
2172
2173            if(min_z_int .gt. max_z_int)
2174                print(" ")
2175                print("'min_z' is greater than 'max_z',")
2176                print("please change this")
2177                print(" ")
2178                exit
2179            end if
2180
2181            if (max_z_int .ge. dimz-1)
2182                max_z_int = dimz-1
2183                if (log_z .EQ. 1) then
2184                max_z_int = max_z_int - 1
2185                end if
2186            end if
2187
2188            if (min_z_int .lt. 0)
2189                min_z_int = 0
2190                if (log_z .EQ. 1) then
2191                  min_z_int = min_z_int + 1
2192                end if
2193            end if         
2194
2195            ;data can contain missing values in case of output of t=0h (inital profiles)
2196            ;where no output is possible
2197            n_not_ismissing = num(.not.ismissing(data(varn,:,:)))
2198
2199            if (n_not_ismissing .GT. 0 ) then   
2200   
2201              if (abs(min(data(varn,:,min_z_int:max_z_int))) .GT. 10)then
2202                  min_value = abs(0.001*min(data(varn,:,min_z_int:max_z_int)))
2203                  max_value = abs(0.001*max(data(varn,:,min_z_int:max_z_int)))
2204              else
2205                  if (abs(min(data(varn,:,min_z_int:max_z_int))) .LT. 0.01 .AND. \
2206                      abs(max(data(varn,:,min_z_int:max_z_int))) .GT. 0.01)then
2207                    min_value = abs(0.1*max(data(varn,:,min_z_int:max_z_int)))
2208                    max_value = abs(0.1*max(data(varn,:,min_z_int:max_z_int)))
2209                  else
2210                    if (abs(max(data(varn,:,min_z_int:max_z_int))) .LT. 0.01 .AND.\
2211                        abs(min(data(varn,:,min_z_int:max_z_int))) .GT. 0.01)then
2212                        min_value = abs(0.1*min(data(varn,:,min_z_int:max_z_int)))
2213                        max_value = abs(0.1*min(data(varn,:,min_z_int:max_z_int)))
2214                    else
2215                        min_value = abs(0.1*min(data(varn,:,min_z_int:max_z_int)))
2216                        max_value = abs(0.1*max(data(varn,:,min_z_int:max_z_int)))
2217                    end if
2218                  end if
2219              end if
2220           
2221              if (min(data(varn,:,min_z_int:max_z_int)) .EQ. 0 .AND. \
2222                  max(data(varn,:,min_z_int:max_z_int)) .EQ. 0)then
2223                  min_value = 0.1
2224                  max_value = 0.1
2225              end if
2226
2227            else     
2228              print(" ")
2229              print(vNam(varn) + " contains only missing values")
2230              print(" ")
2231            end if   
2232
2233            if (over .EQ. 0) then 
2234                res@gsnLeftString      = " "
2235                res@tiXAxisString      = "("+unit(varn)+")"
2236                res@gsnRightString     = " "
2237                res@trYMinF            = min_z
2238                res@trYMaxF            = max_z
2239                if (xs .EQ. -1) then
2240                  res@trXMinF            = min(data(varn,:,min_z_int:max_z_int))-\
2241                                                                          min_value
2242                else
2243                  res@trXMinF            = xs     
2244                end if
2245                if (xe .EQ. -1) then
2246                  res@trXMaxF            = max(data(varn,:,min_z_int:max_z_int))+\
2247                                                                          max_value
2248                else
2249                  res@trXMaxF            = xe 
2250                end if         
2251
2252                if (c_var_log .EQ. 1) then
2253                  if ((mod(com_i, no_rows * no_columns) .EQ. 0) .AND. (no_files .LT. 2))then
2254                    res@pmLegendDisplayMode  = "Always"
2255                  else
2256                    res@pmLegendDisplayMode  = "Never"
2257                  end if
2258                else
2259                  if ((mod(no_plots-1-n+combine, no_rows * no_columns) .EQ. 0) .AND. (no_files .LT. 2))then
2260                    res@pmLegendDisplayMode  = "Always"
2261                  else
2262                    res@pmLegendDisplayMode  = "Never"
2263                  end if
2264                end if
2265
2266                plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
2267
2268            end if
2269     
2270       
2271            if (vNam(varn) .EQ. "u" ) then
2272                  miniu=min(data(varn,:,min_z_int:max_z_int))-min_value
2273                  maxiu=max(data(varn,:,min_z_int:max_z_int))+max_value
2274                  if (over .EQ. 1) then
2275                      res@xyDashPattern  = 0
2276                      plot_u = gsn_csm_xy(wks,data(varn,:,:),z,res)
2277                  else
2278                      res@gsnLeftString      = " "
2279                      res@tiXAxisString      = "("+unit(varn)+")"
2280                      res@gsnRightString     = " "                 
2281                      if (xs .EQ. -1) then 
2282                        res@trXMinF            = miniu
2283                      else
2284                        res@trXMinF            = xs
2285                      end if
2286                      if (xe .EQ. -1) then       
2287                        res@trXMaxF            = maxiu
2288                      else
2289                        res@trXMaxF            = xe 
2290                      end if               
2291                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
2292                  end if
2293                end if
2294                if (vNam(varn) .EQ. "v") then
2295                  miniv=min(data(varn,:,min_z_int:max_z_int))-min_value
2296                  maxiv=max(data(varn,:,min_z_int:max_z_int))+max_value
2297                  if (over .EQ. 1) then
2298                      res@xyMonoDashPattern = True
2299                      res@xyDashPattern  = 1
2300                      plot_v = gsn_csm_xy(wks,data(varn,:,:),z,res)
2301                  else
2302                      res@gsnLeftString      = " "
2303                      res@tiXAxisString      = "("+unit(varn)+")"
2304                      res@gsnRightString     = " "                 
2305                      if (xs .EQ. -1) then
2306                        res@trXMinF            = miniv
2307                      else
2308                        res@trXMinF            = xs
2309                      end if
2310                      if (xe .EQ. -1) then
2311                        res@trXMaxF            = maxiv
2312                      else
2313                        res@trXMaxF            = xe 
2314                      end if               
2315                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2316                  end if 
2317                end if
2318                if (vNam(varn) .EQ. "w") then
2319                  miniw=min(data(varn,:,min_z_int:max_z_int))-min_value
2320                  maxiw=max(data(varn,:,min_z_int:max_z_int))+max_value
2321                  if (over .EQ. 1) then
2322                      res@xyDashPattern = 2
2323                      plot_w = gsn_csm_xy(wks,data(varn,:,:),z,res)
2324                  else
2325                      res@gsnLeftString      = " "
2326                      res@tiXAxisString      = "("+unit(varn)+")"
2327                      res@gsnRightString     = " "
2328                      if (xs .EQ. -1) then
2329                        res@trXMinF            = miniw
2330                      else
2331                        res@trXMinF            = xs
2332                      end if
2333                      if (xe .EQ. -1) then
2334                        res@trXMaxF            = maxiw
2335                      else
2336                        res@trXMaxF            = xe 
2337                      end if           
2338                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2339                  end if   
2340                end if
2341
2342                if (vNam(varn) .EQ. "pt") then
2343                  minipt=min(data(varn,:,min_z_int:max_z_int))-min_value
2344                  maxipt=max(data(varn,:,min_z_int:max_z_int))+max_value
2345                  if (over .EQ. 1) then
2346                      res@xyDashPattern  = 0
2347                      plot_pt = gsn_csm_xy(wks,data(varn,:,:),z,res)
2348                  else
2349                      res@gsnLeftString      = " "
2350                      res@tiXAxisString      = "("+unit(varn)+")"
2351                      res@gsnRightString     = " "
2352                      if (xs .EQ. -1) then
2353                        res@trXMinF            = minipt
2354                      else
2355                        res@trXMinF            = xs
2356                      end if
2357                      if (xe .EQ. -1) then       
2358                        res@trXMaxF            = maxipt
2359                      else
2360                        res@trXMaxF            = xe 
2361                      end if
2362                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2363                  end if 
2364                end if
2365                if (vNam(varn) .EQ. "vpt") then
2366                  minivpt=min(data(varn,:,min_z_int:max_z_int))-min_value
2367                  maxivpt=max(data(varn,:,min_z_int:max_z_int))+max_value
2368                  if (over .EQ. 1) then
2369                      res@xyDashPattern  = 1
2370                      plot_vpt = gsn_csm_xy(wks,data(varn,:,:),z,res)
2371                  else
2372                      res@gsnLeftString      = " "
2373                      res@tiXAxisString      = "("+unit(varn)+")"
2374                      res@gsnRightString     = " "
2375                      if (xs .EQ. -1) then
2376                        res@trXMinF            = minivpt
2377                      else
2378                        res@trXMinF            = xs
2379                      end if
2380                      if (xe .EQ. -1) then       
2381                        res@trXMaxF            = maxivpt
2382                      else
2383                        res@trXMaxF            = xe 
2384                      end if
2385                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2386                  end if
2387                end if
2388                if (vNam(varn) .EQ. "lpt") then
2389                  minilpt=min(data(varn,:,min_z_int:max_z_int))-min_value
2390                  maxilpt=max(data(varn,:,min_z_int:max_z_int))+max_value
2391                  if (over .EQ. 1) then
2392                      res@xyDashPattern  = 2
2393                      plot_lpt = gsn_csm_xy(wks,data(varn,:,:),z,res)
2394                  else
2395                      res@gsnLeftString      = " "
2396                      res@tiXAxisString      = "("+unit(varn)+")"
2397                      res@gsnRightString     = " "
2398                      if (xs .EQ. -1) then
2399                        res@trXMinF            = minilpt
2400                      else
2401                        res@trXMinF            = xs
2402                      end if
2403                      if (xe .EQ. -1) then       
2404                        res@trXMaxF            = maxilpt
2405                      else
2406                        res@trXMaxF            = xe 
2407                      end if
2408                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2409                  end if
2410                end if
2411
2412                if (vNam(varn) .EQ. "q") then
2413                  miniq=min(data(varn,:,min_z_int:max_z_int))-min_value
2414                  maxiq=max(data(varn,:,min_z_int:max_z_int))+max_value
2415                  if (over .EQ. 1) then
2416                      res@xyDashPattern  = 0
2417                      plot_q = gsn_csm_xy(wks,data(varn,:,:),z,res)
2418                  else
2419                      res@gsnLeftString      = " "
2420                      res@tiXAxisString      = "("+unit(varn)+")"
2421                      res@gsnRightString     = " "
2422                      if (xs .EQ. -1) then
2423                        res@trXMinF            = miniq
2424                      else
2425                        res@trXMinF            = xs
2426                      end if
2427                      if (xe .EQ. -1) then       
2428                        res@trXMaxF            = maxiq
2429                      else
2430                        res@trXMaxF            = xe 
2431                      end if
2432                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2433                  end if
2434                end if
2435                if (vNam(varn) .EQ. "qv") then
2436                  miniqv=min(data(varn,:,min_z_int:max_z_int))-min_value
2437                  maxiqv=max(data(varn,:,min_z_int:max_z_int))+max_value
2438                  if (over .EQ. 1) then
2439                      res@xyDashPattern  = 1
2440                      plot_qv = gsn_csm_xy(wks,data(varn,:,:),z,res)
2441                  else
2442                      res@gsnLeftString      = " "
2443                      res@tiXAxisString      = "("+unit(varn)+")"
2444                      res@gsnRightString     = " "
2445                      if (xs .EQ. -1) then
2446                        res@trXMinF            = miniqv
2447                      else
2448                        res@trXMinF            = xs
2449                      end if
2450                      if (xe .EQ. -1) then         
2451                        res@trXMaxF            = maxiqv
2452                      else
2453                        res@trXMaxF            = xe 
2454                      end if
2455                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2456                  end if
2457                end if
2458                if (vNam(varn) .EQ. "ql") then
2459                  miniql=min(data(varn,:,min_z_int:max_z_int))-min_value
2460                  maxiql=max(data(varn,:,min_z_int:max_z_int))+max_value
2461                  if (over .EQ. 1) then
2462                      res@xyDashPattern  = 2
2463                      plot_ql = gsn_csm_xy(wks,data(varn,:,:),z,res)
2464                  else
2465                      res@gsnLeftString      = " "
2466                      res@tiXAxisString      = "("+unit(varn)+")"
2467                      res@gsnRightString     = " "
2468                      if (xs .EQ. -1) then
2469                        res@trXMinF            = miniql
2470                      else
2471                        res@trXMinF            = xs
2472                      end if
2473                      if (xe .EQ. -1) then       
2474                        res@trXMaxF            = maxiql
2475                      else
2476                        res@trXMaxF            = xe 
2477                      end if
2478                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2479                  end if
2480                end if
2481
2482                if (vNam(varn) .EQ. "e") then
2483                  minie=min(data(varn,:,min_z_int:max_z_int))-min_value
2484                  maxie=max(data(varn,:,min_z_int:max_z_int))+max_value
2485                  if (over .EQ. 1) then
2486                      res@xyDashPattern  = 0
2487                      plot_e = gsn_csm_xy(wks,data(varn,:,:),z,res)
2488                  else
2489                      res@gsnLeftString      = " "
2490                      res@tiXAxisString      = "("+unit(varn)+")"
2491                      res@gsnRightString     = " "
2492                      if (xs .EQ. -1) then
2493                        res@trXMinF            = minie
2494                      else
2495                        res@trXMinF            = xs
2496                      end if
2497                      if (xe .EQ. -1) then       
2498                        res@trXMaxF            = maxie
2499                      else
2500                        res@trXMaxF            = xe 
2501                      end if
2502                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2503                  end if
2504                end if
2505                if (vNam(varn) .EQ. "es" .OR. vNam(varn) .EQ. "e*") then
2506                  minies=min(data(varn,:,min_z_int:max_z_int))-min_value
2507                  maxies=max(data(varn,:,min_z_int:max_z_int))+max_value
2508                  if (over .EQ. 1) then
2509                      res@xyDashPattern  = 1
2510                      plot_es = gsn_csm_xy(wks,data(varn,:,:),z,res)
2511                  else
2512                      res@gsnLeftString      = " "
2513                      res@tiXAxisString      = "("+unit(varn)+")"
2514                      res@gsnRightString     = " "
2515                      if (xs .EQ. -1) then
2516                        res@trXMinF            = minies
2517                      else
2518                        res@trXMinF            = xs
2519                      end if
2520                      if (xe .EQ. -1) then       
2521                        res@trXMaxF            = maxies
2522                      else
2523                        res@trXMaxF            = xe 
2524                      end if
2525                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2526                  end if
2527                end if
2528
2529                if (vNam(varn) .EQ. "km") then
2530                  minikm=min(data(varn,:,min_z_int:max_z_int))-min_value
2531                  maxikm=max(data(varn,:,min_z_int:max_z_int))+max_value
2532                  if (over .EQ. 1) then
2533                      res@xyDashPattern  = 0
2534                      plot_km = gsn_csm_xy(wks,data(varn,:,:),z,res)
2535                  else
2536                      res@gsnLeftString      = " "
2537                      res@tiXAxisString      = "("+unit(varn)+")"
2538                      res@gsnRightString     = " "
2539                      if (xs .EQ. -1) then
2540                        res@trXMinF            = minikm
2541                      else
2542                        res@trXMinF            = xs
2543                      end if
2544                      if (xe .EQ. -1) then     
2545                        res@trXMaxF            = maxikm
2546                      else
2547                        res@trXMaxF            = xe 
2548                      end if
2549                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2550                  end if
2551                end if
2552                if (vNam(varn) .EQ. "kh") then
2553                  minikh=min(data(varn,:,min_z_int:max_z_int))-min_value
2554                  maxikh=max(data(varn,:,min_z_int:max_z_int))+max_value
2555                  if (over .EQ. 1) then
2556                      res@xyDashPattern  = 1
2557                      plot_kh = gsn_csm_xy(wks,data(varn,:,:),z,res)
2558                  else
2559                      res@gsnLeftString      = " "
2560                      res@tiXAxisString      = "("+unit(varn)+")"
2561                      res@gsnRightString     = " "
2562                      if (xs .EQ. -1) then
2563                        res@trXMinF            = minikh
2564                      else
2565                        res@trXMinF            = xs
2566                      end if
2567                      if (xe .EQ. -1) then     
2568                        res@trXMaxF            = maxikh
2569                      else
2570                        res@trXMaxF            = xe 
2571                      end if
2572                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2573                  end if
2574                end if
2575
2576                if (vNam(varn) .EQ. "wpup" .OR. vNam(varn) .EQ. "w"+dq+"u"+dq) then
2577                  miniwpup=min(data(varn,:,min_z_int:max_z_int))-min_value
2578                  maxiwpup=max(data(varn,:,min_z_int:max_z_int))+max_value
2579                  if (over .EQ. 1) then
2580                      res@xyDashPattern  = 0
2581                      plot_wpup = gsn_csm_xy(wks,data(varn,:,:),z,res)
2582                  else
2583                      res@gsnLeftString      = " "
2584                      res@tiXAxisString      = "("+unit(varn)+")"
2585                      res@gsnRightString     = " "
2586                      if (xs .EQ. -1) then
2587                        res@trXMinF            = miniwpup
2588                      else
2589                        res@trXMinF            = xs
2590                      end if
2591                      if (xe .EQ. -1) then
2592                        res@trXMaxF            = maxiwpup
2593                      else
2594                        res@trXMaxF            = xe 
2595                      end if
2596                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2597                  end if
2598                end if
2599                if (vNam(varn) .EQ. "wsus" .OR. vNam(varn) .EQ. "w*u*") then
2600                  miniwsus=min(data(varn,:,min_z_int:max_z_int))-min_value
2601                  maxiwsus=max(data(varn,:,min_z_int:max_z_int))+max_value
2602                  if (over .EQ. 1) then
2603                      res@xyDashPattern  = 1
2604                      plot_wsus = gsn_csm_xy(wks,data(varn,:,:),z,res)
2605                  else
2606                      res@gsnLeftString      = " "
2607                      res@tiXAxisString      = "("+unit(varn)+")"
2608                      res@gsnRightString     = " "
2609                      if (xs .EQ. -1) then
2610                        res@trXMinF            = miniwsus
2611                      else
2612                        res@trXMinF            = xs
2613                      end if
2614                      if (xe .EQ. -1) then     
2615                        res@trXMaxF            = maxiwsus
2616                      else
2617                        res@trXMaxF            = xe 
2618                      end if
2619                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2620                  end if
2621                end if
2622                if (vNam(varn) .EQ. "wu") then
2623                  miniwu=min(data(varn,:,min_z_int:max_z_int))-min_value
2624                  maxiwu=max(data(varn,:,min_z_int:max_z_int))+max_value
2625                  if (over .EQ. 1) then
2626                      res@xyDashPattern  = 2
2627                      plot_wu = gsn_csm_xy(wks,data(varn,:,:),z,res)
2628                  else
2629                      res@gsnLeftString      = " "
2630                      res@tiXAxisString      = "("+unit(varn)+")"
2631                      res@gsnRightString     = " "
2632                      if (xs .EQ. -1) then
2633                        res@trXMinF            = miniwu
2634                      else
2635                        res@trXMinF            = xs
2636                      end if
2637                      if (xe .EQ. -1) then
2638                        res@trXMaxF            = maxiwu
2639                      else
2640                        res@trXMaxF            = xe 
2641                      end if
2642                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2643                  end if
2644                end if
2645
2646                if (vNam(varn) .EQ. "wpvp" .OR. vNam(varn) .EQ. "w"+dq+"v"+dq) then
2647                  miniwpvp=min(data(varn,:,min_z_int:max_z_int))-min_value
2648                  maxiwpvp=max(data(varn,:,min_z_int:max_z_int))+max_value
2649                  if (over .EQ. 1) then
2650                      res@xyDashPattern  = 0
2651                      plot_wpvp = gsn_csm_xy(wks,data(varn,:,:),z,res)
2652                  else
2653                      res@gsnLeftString      = " "
2654                      res@tiXAxisString      = "("+unit(varn)+")"
2655                      res@gsnRightString     = " "
2656                      if (xs .EQ. -1) then
2657                        res@trXMinF            = miniwpvp
2658                      else
2659                        res@trXMinF            = xs
2660                      end if
2661                      if (xe .EQ. -1) then     
2662                        res@trXMaxF            = maxiwpvp
2663                      else
2664                        res@trXMaxF            = xe 
2665                      end if
2666                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2667                  end if
2668                end if
2669                if (vNam(varn) .EQ. "wsvs" .OR. vNam(varn) .EQ. "w*v*") then
2670                  miniwsvs=min(data(varn,:,min_z_int:max_z_int))-min_value
2671                  maxiwsvs=max(data(varn,:,min_z_int:max_z_int))+max_value
2672                  if (over .EQ. 1) then
2673                      res@xyDashPattern  = 1
2674                      plot_wsvs = gsn_csm_xy(wks,data(varn,:,:),z,res)
2675                  else
2676                      res@gsnLeftString      = " "
2677                      res@tiXAxisString      = "("+unit(varn)+")"
2678                      res@gsnRightString     = " "
2679                      if (xs .EQ. -1) then
2680                        res@trXMinF            = miniwsvs
2681                      else
2682                        res@trXMinF            = xs
2683                      end if
2684                      if (xe .EQ. -1) then     
2685                        res@trXMaxF            = maxiwsvs
2686                      else
2687                        res@trXMaxF            = xe 
2688                      end if
2689                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2690                  end if
2691                end if
2692                if (vNam(varn) .EQ. "wv") then
2693                  miniwv=min(data(varn,:,min_z_int:max_z_int))-min_value
2694                  maxiwv=max(data(varn,:,min_z_int:max_z_int))+max_value
2695                  if (over .EQ. 1) then
2696                      res@xyDashPattern  = 2
2697                      plot_wv = gsn_csm_xy(wks,data(varn,:,:),z,res)
2698                  else
2699                      res@gsnLeftString      = " "
2700                      res@tiXAxisString      = "("+unit(varn)+")"
2701                      res@gsnRightString     = " "
2702                      if (xs .EQ. -1) then
2703                        res@trXMinF            = miniwv
2704                      else
2705                        res@trXMinF            = xs
2706                      end if
2707                      if (xe .EQ. -1) then       
2708                        res@trXMaxF            = maxiwv
2709                      else
2710                        res@trXMaxF            = xe 
2711                      end if
2712                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2713                  end if
2714                end if
2715
2716                if (vNam(varn) .EQ. "wpptp" .OR. vNam(varn) \
2717                  .EQ. "w"+dq+"pt"+dq) then
2718                  miniwpptp=min(data(varn,:,min_z_int:max_z_int))-min_value
2719                  maxiwpptp=max(data(varn,:,min_z_int:max_z_int))+max_value
2720                  if (over .EQ. 1) then
2721                      res@xyDashPattern  = 0
2722                      plot_wpptp = gsn_csm_xy(wks,data(varn,:,:),z,res)
2723                  else
2724                      res@gsnLeftString      = " "
2725                      res@tiXAxisString      = "("+unit(varn)+")"
2726                      res@gsnRightString     = " "
2727                      if (xs .EQ. -1) then
2728                        res@trXMinF            = miniwpptp
2729                      else
2730                        res@trXMinF            = xs
2731                      end if
2732                      if (xe .EQ. -1) then       
2733                        res@trXMaxF            = maxiwpptp
2734                      else
2735                        res@trXMaxF            = xe 
2736                      end if
2737                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2738                  end if
2739                end if
2740                if (vNam(varn) .EQ. "wspts" .OR. vNam(varn) .EQ. "w*pt*") then
2741                  miniwspts=min(data(varn,:,min_z_int:max_z_int))-min_value
2742                  maxiwspts=max(data(varn,:,min_z_int:max_z_int))+max_value
2743                  if (over .EQ. 1) then
2744                      res@xyDashPattern  = 1
2745                      plot_wspts = gsn_csm_xy(wks,data(varn,:,:),z,res)
2746                  else
2747                      res@gsnLeftString      = " "
2748                      res@tiXAxisString      = "("+unit(varn)+")"
2749                      res@gsnRightString     = " "
2750                      if (xs .EQ. -1) then
2751                        res@trXMinF            = miniwspts
2752                      else
2753                        res@trXMinF            = xs
2754                      end if
2755                      if (xe .EQ. -1) then       
2756                        res@trXMaxF            = maxiwspts
2757                      else
2758                        res@trXMaxF            = xe 
2759                      end if
2760                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2761                  end if
2762                end if
2763                if (vNam(varn) .EQ. "wpt") then       
2764                  miniwpt=min(data(varn,:,min_z_int:max_z_int))-min_value
2765                  maxiwpt=max(data(varn,:,min_z_int:max_z_int))+max_value
2766                  if (over .EQ. 1) then
2767                      res@xyDashPattern  = 2
2768                      plot_wpt = gsn_csm_xy(wks,data(varn,:,:),z,res)
2769                  else
2770                      res@gsnLeftString      = " "
2771                      res@tiXAxisString      = "("+unit(varn)+")"
2772                      res@gsnRightString     = " "
2773                      if (xs .EQ. -1) then
2774                        res@trXMinF            = miniwpt
2775                      else
2776                        res@trXMinF            = xs
2777                      end if
2778                      if (xe .EQ. -1) then       
2779                        res@trXMaxF            = maxiwpt
2780                      else
2781                        res@trXMaxF            = xe 
2782                      end if
2783                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2784                  end if
2785                end if
2786
2787                if (vNam(varn) .EQ. "wsptsBC".OR. vNam(varn) .EQ. "w*pt*BC" ) then
2788                  miniwsptsBC=min(data(varn,:,min_z_int:max_z_int))-min_value
2789                  maxiwsptsBC=max(data(varn,:,min_z_int:max_z_int))+max_value
2790                  if (over .EQ. 1) then
2791                      res@xyDashPattern  = 0
2792                      plot_wsptsBC = gsn_csm_xy(wks,data(varn,:,:),z,res)
2793                  else
2794                      res@gsnLeftString      = " "
2795                      res@tiXAxisString      = "("+unit(varn)+")"
2796                      res@gsnRightString     = " "
2797                      if (xs .EQ. -1) then
2798                        res@trXMinF            = miniwsptsBC
2799                      else
2800                        res@trXMinF            = xs
2801                      end if
2802                      if (xe .EQ. -1) then       
2803                        res@trXMaxF            = maxiwsptsBC
2804                      else
2805                        res@trXMaxF            = xe 
2806                      end if
2807                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2808                  end if
2809                end if             
2810                if (vNam(varn) .EQ. "wptBC") then
2811                  miniwptBC=min(data(varn,:,min_z_int:max_z_int))-min_value
2812                  maxiwptBC=max(data(varn,:,min_z_int:max_z_int))+max_value
2813                  if (over .EQ. 1) then
2814                      res@xyDashPattern  = 1
2815                      plot_wptBC = gsn_csm_xy(wks,data(varn,:,:),z,res)
2816                  else
2817                      res@gsnLeftString      = " "
2818                      res@tiXAxisString      = "("+unit(varn)+")"
2819                      res@gsnRightString     = " "
2820                      if (xs .EQ. -1) then
2821                        res@trXMinF            = miniwptBC
2822                      else
2823                        res@trXMinF            = xs
2824                      end if
2825                      if (xe .EQ. -1) then       
2826                        res@trXMaxF            = maxiwptBC
2827                      else
2828                        res@trXMaxF            = xe 
2829                      end if
2830                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2831                  end if
2832                end if
2833
2834                if (vNam(varn) .EQ. "wpvptp" .OR. vNam(varn) \
2835                    .EQ. "w"+dq+"vpt"+dq) then
2836                  miniwpvptp=min(data(varn,:,min_z_int:max_z_int))-min_value
2837                  maxiwpvptp=max(data(varn,:,min_z_int:max_z_int))+max_value
2838                  if (over .EQ. 1) then
2839                      res@xyDashPattern  = 0
2840                      plot_wpvptp = gsn_csm_xy(wks,data(varn,:,:),z,res)
2841                  else
2842                      res@gsnLeftString      = " "
2843                      res@tiXAxisString      = "("+unit(varn)+")"
2844                      res@gsnRightString     = " "
2845                      if (xs .EQ. -1) then
2846                        res@trXMinF            = miniwpvptp
2847                      else
2848                        res@trXMinF            = xs
2849                      end if
2850                      if (xe .EQ. -1) then       
2851                        res@trXMaxF            = maxiwpvptp
2852                      else
2853                        res@trXMaxF            = xe 
2854                      end if
2855                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2856                  end if
2857                end if
2858                if (vNam(varn) .EQ. "wsvpts" .OR. vNam(varn) .EQ. "w*vpt*") then
2859                  miniwsvpts=min(data(varn,:,min_z_int:max_z_int))-min_value
2860                  maxiwsvpts=max(data(varn,:,min_z_int:max_z_int))+max_value
2861                  if (over .EQ. 1) then
2862                      res@xyDashPattern  = 1
2863                      plot_wsvpts = gsn_csm_xy(wks,data(varn,:,:),z,res)
2864                  else
2865                      res@gsnLeftString      = " "
2866                      res@tiXAxisString      = "("+unit(varn)+")"
2867                      res@gsnRightString     = " "
2868                      if (xs .EQ. -1) then
2869                        res@trXMinF            = miniwsvpts
2870                      else
2871                        res@trXMinF            = xs
2872                      end if
2873                      if (xe .EQ. -1) then       
2874                        res@trXMaxF            = maxiwsvpts
2875                      else
2876                        res@trXMaxF            = xe 
2877                      end if
2878                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2879                  end if
2880                end if
2881                if (vNam(varn) .EQ. "wvpt") then
2882                  miniwvpt=min(data(varn,:,min_z_int:max_z_int))-min_value
2883                  maxiwvpt=max(data(varn,:,min_z_int:max_z_int))+max_value
2884                  if (over .EQ. 1) then
2885                      res@xyDashPattern  = 2
2886                      plot_wvpt = gsn_csm_xy(wks,data(varn,:,:),z,res)
2887                  else
2888                      res@gsnLeftString      = " "
2889                      res@tiXAxisString      = "("+unit(varn)+")"
2890                      res@gsnRightString     = " "
2891                      if (xs .EQ. -1) then
2892                        res@trXMinF            = miniwvpt
2893                      else
2894                        res@trXMinF            = xs
2895                      end if
2896                      if (xe .EQ. -1) then       
2897                        res@trXMaxF            = maxiwvpt
2898                      else
2899                        res@trXMaxF            = xe 
2900                      end if
2901                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2902                  end if
2903                end if
2904
2905                if (vNam(varn) .EQ. "wpqp" .OR. vNam(varn) .EQ. "w"+dq+"q"+dq) then
2906                  miniwpqp=min(data(varn,:,min_z_int:max_z_int))-min_value
2907                  maxiwpqp=max(data(varn,:,min_z_int:max_z_int))+max_value
2908                  if (over .EQ. 1) then
2909                      res@xyDashPattern  = 0
2910                      plot_wpqp = gsn_csm_xy(wks,data(varn,:,:),z,res)
2911                  else
2912                      res@gsnLeftString      = " "
2913                      res@tiXAxisString      = "("+unit(varn)+")"
2914                      res@gsnRightString     = " "
2915                      if (xs .EQ. -1) then
2916                        res@trXMinF            = miniwpqp
2917                      else
2918                        res@trXMinF            = xs
2919                      end if
2920                      if (xe .EQ. -1) then       
2921                        res@trXMaxF            = maxiwpqp
2922                      else
2923                        res@trXMaxF            = xe 
2924                      end if
2925                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2926                  end if
2927                end if
2928                if (vNam(varn) .EQ. "wsqs".OR. vNam(varn) .EQ. "w*s*" ) then
2929                  miniwsqs=min(data(varn,:,min_z_int:max_z_int))-min_value
2930                  maxiwsqs=max(data(varn,:,min_z_int:max_z_int))+max_value
2931                  if (over .EQ. 1) then
2932                      res@xyDashPattern  = 1
2933                      plot_wsqs = gsn_csm_xy(wks,data(varn,:,:),z,res)
2934                  else
2935                      res@gsnLeftString      = " "
2936                      res@tiXAxisString      = "("+unit(varn)+")"
2937                      res@gsnRightString     = " "
2938                      if (xs .EQ. -1) then
2939                        res@trXMinF            = miniwsqs
2940                      else
2941                        res@trXMinF            = xs
2942                      end if
2943                      if (xe .EQ. -1) then       
2944                        res@trXMaxF            = maxiwsqs
2945                      else
2946                        res@trXMaxF            = xe 
2947                      end if
2948                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2949                  end if
2950                end if
2951                if (vNam(varn) .EQ. "wq") then
2952                  miniwq=min(data(varn,:,min_z_int:max_z_int))-min_value
2953                  maxiwq=max(data(varn,:,min_z_int:max_z_int))+max_value
2954                  if (over .EQ. 1) then
2955                      res@xyDashPattern  = 2
2956                      plot_wq = gsn_csm_xy(wks,data(varn,:,:),z,res)
2957                  else
2958                      res@gsnLeftString      = " "
2959                      res@tiXAxisString      = "("+unit(varn)+")"
2960                      res@gsnRightString     = " "
2961                      if (xs .EQ. -1) then
2962                        res@trXMinF            = miniwq
2963                      else
2964                        res@trXMinF            = xs
2965                      end if
2966                      if (xe .EQ. -1) then       
2967                        res@trXMaxF            = maxiwq
2968                      else
2969                        res@trXMaxF            = xe 
2970                      end if
2971                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2972                  end if
2973                end if
2974
2975                if (vNam(varn) .EQ. "wpqvp" .OR. \
2976                  vNam(varn) .EQ. "w"+dq+"qv"+dq) then
2977                  miniwpqvp=min(data(varn,:,min_z_int:max_z_int))-min_value
2978                  maxiwpqvp=max(data(varn,:,min_z_int:max_z_int))+max_value
2979                  if (over .EQ. 1) then
2980                      res@xyDashPattern  = 0
2981                      plot_wpqvp = gsn_csm_xy(wks,data(varn,:,:),z,res)
2982                  else
2983                      res@gsnLeftString      = " "
2984                      res@tiXAxisString      = "("+unit(varn)+")"
2985                      res@gsnRightString     = " "
2986                      if (xs .EQ. -1) then
2987                        res@trXMinF            = miniwpqvp
2988                      else
2989                        res@trXMinF            = xs
2990                      end if
2991                      if (xe .EQ. -1) then       
2992                        res@trXMaxF            = maxiwpqvp
2993                      else
2994                        res@trXMaxF            = xe 
2995                      end if
2996                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
2997                  end if
2998                end if
2999                if (vNam(varn) .EQ. "wsqvs" .OR. vNam(varn) .EQ. "w*qv*") then
3000                  miniwsqvs=min(data(varn,:,min_z_int:max_z_int))-min_value
3001                  maxiwsqvs=max(data(varn,:,min_z_int:max_z_int))+max_value
3002                  if (over .EQ. 1) then
3003                      res@xyDashPattern  = 1
3004                      plot_wsqvs = gsn_csm_xy(wks,data(varn,:,:),z,res)
3005                  else
3006                      res@gsnLeftString      = " "
3007                      res@tiXAxisString      = "("+unit(varn)+")"
3008                      res@gsnRightString     = " "
3009                      if (xs .EQ. -1) then
3010                        res@trXMinF            = miniwsqvs
3011                      else
3012                        res@trXMinF            = xs
3013                      end if
3014                      if (xe .EQ. -1) then       
3015                        res@trXMaxF            = maxiwsqvs
3016                      else
3017                        res@trXMaxF            = xe 
3018                      end if
3019                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3020                  end if
3021                end if
3022                if (vNam(varn) .EQ. "wqv") then
3023                  miniwqv=min(data(varn,:,min_z_int:max_z_int))-min_value
3024                  maxiwqv=max(data(varn,:,min_z_int:max_z_int))+max_value
3025                  if (over .EQ. 1) then
3026                      res@xyDashPattern  = 2
3027                      plot_wqv = gsn_csm_xy(wks,data(varn,:,:),z,res)
3028                  else
3029                      res@gsnLeftString      = " "
3030                      res@tiXAxisString      = "("+unit(varn)+")"
3031                      res@gsnRightString     = " "
3032                      if (xs .EQ. -1) then
3033                        res@trXMinF            = miniwqv
3034                      else
3035                        res@trXMinF            = xs
3036                      end if
3037                      if (xe .EQ. -1) then       
3038                        res@trXMaxF            = maxiwqv
3039                      else
3040                        res@trXMaxF            = xe 
3041                      end if
3042                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3043                  end if
3044                end if
3045
3046                if (vNam(varn) .EQ. "wpsp" .OR. vNam(varn) .EQ. "w"+dq+"s"+dq) then
3047                  miniwpsp=min(data(varn,:,min_z_int:max_z_int))-min_value
3048                  maxiwpsp=max(data(varn,:,min_z_int:max_z_int))+max_value
3049                  if (over .EQ. 1) then
3050                      res@xyDashPattern  = 0
3051                      plot_wpsp = gsn_csm_xy(wks,data(varn,:,:),z,res)
3052                  else
3053                      res@gsnLeftString      = " "
3054                      res@tiXAxisString      = "("+unit(varn)+")"
3055                      res@gsnRightString     = " "
3056                      if (xs .EQ. -1) then
3057                        res@trXMinF            = miniwpsp
3058                      else
3059                        res@trXMinF            = xs
3060                      end if
3061                      if (xe .EQ. -1) then       
3062                        res@trXMaxF            = maxiwpsp
3063                      else
3064                        res@trXMaxF            = xe 
3065                      end if
3066                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3067                  end if
3068                end if
3069                if (vNam(varn) .EQ. "wsss" .OR. vNam(varn) .EQ. "w*s*" ) then
3070                  miniwsss=min(data(varn,:,min_z_int:max_z_int))-min_value
3071                  maxiwsss=max(data(varn,:,min_z_int:max_z_int))+max_value
3072                  if (over .EQ. 1) then
3073                      res@xyDashPattern  = 1
3074                      plot_wsss = gsn_csm_xy(wks,data(varn,:,:),z,res)
3075                  else
3076                      res@gsnLeftString      = " "
3077                      res@tiXAxisString      = "("+unit(varn)+")"
3078                      res@gsnRightString     = " "
3079                      if (xs .EQ. -1) then
3080                        res@trXMinF            = miniwsss
3081                      else
3082                        res@trXMinF            = xs
3083                      end if
3084                      if (xe .EQ. -1) then       
3085                        res@trXMaxF            = maxiwsss
3086                      else
3087                        res@trXMaxF            = xe 
3088                      end if
3089                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3090                  end if
3091                end if
3092                if (vNam(varn) .EQ. "ws") then
3093                  miniws=min(data(varn,:,min_z_int:max_z_int))-min_value
3094                  maxiws=max(data(varn,:,min_z_int:max_z_int))+max_value
3095                  if (over .EQ. 1) then
3096                      res@xyDashPattern  = 2
3097                      plot_ws = gsn_csm_xy(wks,data(varn,:,:),z,res)
3098                  else
3099                      res@gsnLeftString      = " "
3100                      res@tiXAxisString      = "("+unit(varn)+")"
3101                      res@gsnRightString     = " "
3102                      if (xs .EQ. -1) then
3103                        res@trXMinF            = miniws
3104                      else
3105                        res@trXMinF            = xs
3106                      end if
3107                      if (xe .EQ. -1) then       
3108                        res@trXMaxF            = maxiws
3109                      else
3110                        res@trXMaxF            = xe 
3111                      end if
3112                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3113                  end if
3114                end if
3115
3116                if (vNam(varn) .EQ. "wpsap" .OR. \
3117                    vNam(varn) .EQ. "w"+dq+"sa"+dq) then
3118                  miniwpsap=min(data(varn,:,min_z_int:max_z_int))-min_value
3119                  maxiwpsap=max(data(varn,:,min_z_int:max_z_int))+max_value
3120                  if (over .EQ. 1) then
3121                      res@xyDashPattern  = 0
3122                      plot_wpsap = gsn_csm_xy(wks,data(varn,:,:),z,res)
3123                  else
3124                      res@gsnLeftString      = " "
3125                      res@tiXAxisString      = "("+unit(varn)+")"
3126                      res@gsnRightString     = " "
3127                      if (xs .EQ. -1) then
3128                        res@trXMinF            = miniwpsap
3129                      else
3130                        res@trXMinF            = xs
3131                      end if
3132                      if (xe .EQ. -1) then       
3133                        res@trXMaxF            = maxiwpsap
3134                      else
3135                        res@trXMaxF            = xe 
3136                      end if
3137                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3138                  end if
3139                end if
3140                if (vNam(varn) .EQ. "wssas" .OR. vNam(varn) .EQ. "w*sa*") then
3141                  miniwssas=min(data(varn,:,min_z_int:max_z_int))-min_value
3142                  maxiwssas=max(data(varn,:,min_z_int:max_z_int))+max_value
3143                  if (over .EQ. 1) then
3144                      res@xyDashPattern  = 1
3145                      plot_wssas = gsn_csm_xy(wks,data(varn,:,:),z,res)
3146                  else
3147                      res@gsnLeftString      = " "
3148                      res@tiXAxisString      = "("+unit(varn)+")"
3149                      res@gsnRightString     = " "
3150                      if (xs .EQ. -1) then
3151                        res@trXMinF            = miniwssas
3152                      else
3153                        res@trXMinF            = xs
3154                      end if
3155                      if (xe .EQ. -1) then       
3156                        res@trXMaxF            = maxiwssas
3157                      else
3158                        res@trXMaxF            = xe 
3159                      end if
3160                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3161                  end if
3162                end if
3163                if (vNam(varn) .EQ. "wsa") then
3164                  miniwsa=min(data(varn,:,min_z_int:max_z_int))-min_value
3165                  maxiwsa=max(data(varn,:,min_z_int:max_z_int))+max_value
3166                  if (over .EQ. 1) then
3167                      res@xyDashPattern  = 2
3168                      plot_wsa = gsn_csm_xy(wks,data(varn,:,:),z,res)
3169                  else
3170                      res@gsnLeftString      = " "
3171                      res@tiXAxisString      = "("+unit(varn)+")"
3172                      res@gsnRightString     = " "
3173                      if (xs .EQ. -1) then
3174                        res@trXMinF            = miniwsa
3175                      else
3176                        res@trXMinF            = xs
3177                      end if
3178                      if (xe .EQ. -1) then       
3179                        res@trXMaxF            = maxiwsa
3180                      else
3181                        res@trXMaxF            = xe 
3182                      end if
3183                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3184                  end if
3185                end if
3186
3187                if (vNam(varn) .EQ. "us2" .OR. vNam(varn) .EQ. "u*2") then
3188                  minius2=min(data(varn,:,min_z_int:max_z_int))-min_value
3189                  maxius2=max(data(varn,:,min_z_int:max_z_int))+max_value
3190                  if (over .EQ. 1) then
3191                      res@xyDashPattern  = 0
3192                      plot_us2 = gsn_csm_xy(wks,data(varn,:,:),z,res)
3193                  else
3194                      res@gsnLeftString      = " "
3195                      res@tiXAxisString      = "("+unit(varn)+")"
3196                      res@gsnRightString     = " "
3197                      if (xs .EQ. -1) then
3198                        res@trXMinF            = minius2
3199                      else
3200                        res@trXMinF            = xs
3201                      end if
3202                      if (xe .EQ. -1) then       
3203                        res@trXMaxF            = maxius2
3204                      else
3205                        res@trXMaxF            = xe 
3206                      end if
3207                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3208                  end if
3209                end if
3210                if (vNam(varn) .EQ. "vs2" .OR. vNam(varn) .EQ. "v*2") then
3211                  minivs2=min(data(varn,:,min_z_int:max_z_int))-min_value
3212                  maxivs2=max(data(varn,:,min_z_int:max_z_int))+max_value
3213                  if (over .EQ. 1) then
3214                      res@xyDashPattern  = 1
3215                      plot_vs2 = gsn_csm_xy(wks,data(varn,:,:),z,res)
3216                  else
3217                      res@gsnLeftString      = " "
3218                      res@tiXAxisString      = "("+unit(varn)+")"
3219                      res@gsnRightString     = " "
3220                      if (xs .EQ. -1) then
3221                        res@trXMinF            = minivs2
3222                      else
3223                        res@trXMinF            = xs
3224                      end if
3225                      if (xe .EQ. -1) then       
3226                        res@trXMaxF            = maxivs2
3227                      else
3228                        res@trXMaxF            = xe 
3229                      end if
3230                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3231                  end if
3232                end if
3233                if (vNam(varn) .EQ. "ws2" .OR. vNam(varn) .EQ. "w*2") then
3234                  miniws2=min(data(varn,:,min_z_int:max_z_int))-min_value
3235                  maxiws2=max(data(varn,:,min_z_int:max_z_int))+max_value
3236                  if (over .EQ. 1) then
3237                      res@xyDashPattern  = 2
3238                      plot_ws2 = gsn_csm_xy(wks,data(varn,:,:),z,res)
3239                  else
3240                      res@gsnLeftString      = " "
3241                      res@tiXAxisString      = "("+unit(varn)+")"
3242                      res@gsnRightString     = " "
3243                      if (xs .EQ. -1) then
3244                        res@trXMinF            = miniws2
3245                      else
3246                        res@trXMinF            = xs
3247                      end if
3248                      if (xe .EQ. -1) then       
3249                        res@trXMaxF            = maxiws2
3250                      else
3251                        res@trXMaxF            = xe 
3252                      end if
3253                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3254                  end if
3255                end if
3256
3257                if (vNam(varn) .EQ. "wsususodz" .OR. \
3258                    vNam(varn) .EQ. "w*u*u*:dz") then
3259                  miniwsususodz=min(data(varn,:,min_z_int:max_z_int))-min_value
3260                  maxiwsususodz=max(data(varn,:,min_z_int:max_z_int))+max_value
3261                  if (over .EQ. 1) then
3262                      res@xyDashPattern  = 0
3263                      plot_wsususodz = gsn_csm_xy(wks,data(varn,:,:),z,res)
3264                  else
3265                      res@gsnLeftString      = " "
3266                      res@tiXAxisString      = "("+unit(varn)+")"
3267                      res@gsnRightString     = " "
3268                      if (xs .EQ. -1) then
3269                        res@trXMinF            = miniwsususodz
3270                      else
3271                        res@trXMinF            = xs
3272                      end if
3273                      if (xe .EQ. -1) then       
3274                        res@trXMaxF            = maxiwsususodz
3275                      else
3276                        res@trXMaxF            = xe 
3277                      end if
3278                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3279                  end if 
3280                end if
3281                if (vNam(varn) .EQ. "wspsodz" .OR. vNam(varn) .EQ. "w*p*:dz") then
3282                  miniwspsodz=min(data(varn,:,min_z_int:max_z_int))-min_value
3283                  maxiwspsodz=max(data(varn,:,min_z_int:max_z_int))+max_value
3284                  if (over .EQ. 1) then
3285                      res@xyDashPattern  = 1
3286                      plot_wspsodz = gsn_csm_xy(wks,data(varn,:,:),z,res)
3287                  else
3288                      res@gsnLeftString      = " "
3289                      res@tiXAxisString      = "("+unit(varn)+")"
3290                      res@gsnRightString     = " "
3291                      if (xs .EQ. -1) then
3292                        res@trXMinF            = miniwspsodz
3293                      else
3294                        res@trXMinF            = xs
3295                      end if
3296                      if (xe .EQ. -1) then       
3297                        res@trXMaxF            = maxiwspsodz
3298                      else
3299                        res@trXMaxF            = xe 
3300                      end if
3301                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3302                  end if
3303                end if
3304                if (vNam(varn) .EQ. "wpeodz" .OR. \
3305                    vNam(varn) .EQ. "w"+dq+"p:dz") then
3306                  miniwpeodz=min(data(varn,:,min_z_int:max_z_int))-min_value
3307                  maxiwpeodz=max(data(varn,:,min_z_int:max_z_int))+max_value
3308                  if (over .EQ. 1) then
3309                      res@xyDashPattern  = 2
3310                      plot_wpeodz = gsn_csm_xy(wks,data(varn,:,:),z,res)
3311                  else
3312                      res@gsnLeftString      = " "
3313                      res@tiXAxisString      = "("+unit(varn)+")"
3314                      res@gsnRightString     = " "
3315                      if (xs .EQ. -1) then
3316                        res@trXMinF            = miniwpeodz
3317                      else
3318                        res@trXMinF            = xs
3319                      end if
3320                      if (xe .EQ. -1) then       
3321                        res@trXMaxF            = maxiwpeodz
3322                      else
3323                        res@trXMaxF            = xe 
3324                      end if
3325                      plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res)
3326                  end if
3327                end if
3328
3329                ; ***************************************************
3330                ; legend for each plot
3331                ; ***************************************************
3332               
3333                if (over .EQ. 0) then
3334
3335                   label = vNam(varn)
3336
3337                   lgres1                    = True
3338                   lgMonoDashIndex           = False
3339                   lgres1@lgDashIndexes      = (/0,1,2/)
3340                   lgres1@lgLabelFont        = "helvetica"   
3341                   lgres1@lgLabelFontHeightF = font_size_legend * 10.0 * 0.7
3342                   lgres1@vpWidthF           = 0.12           
3343                   lgres1@vpHeightF          = 0.1 / 3.0 + 0.01         
3344               
3345                   lbid = gsn_create_legend(wks,1,label,lgres1)       
3346
3347                   amres = True
3348                   amres@amParallelPosF   = -0.5                   
3349                   amres@amOrthogonalPosF = -0.5   
3350                   amres@amJust           = "TopLeft"
3351
3352                   annoid1 = gsn_add_annotation(plot(n),lbid,amres)
3353     
3354                end if
3355
3356            if (no_files .GT. 1) then
3357                print("nof="+nof+" und n="+n)
3358                multi_plot(nof,n)=plot(n)
3359                max_nof(nof,n)=max(data(varn,:,min_z_int:max_z_int))
3360                min_nof(nof,n)=min(data(varn,:,min_z_int:max_z_int))
3361                name(nof,n)   =vNam(varn)
3362                unit_(nof,n)  =unit(varn)
3363            end if
3364            if (over .EQ. 0) then
3365                n=n+1 
3366            end if 
3367            if (prof3d .EQ. 0)then   
3368                varn=varn+1
3369            end if
3370            delete(temp)
3371            delete(temp_att)
3372          end if         
3373      end do
3374      if (no_files .GT. 1) then
3375          delete(vNam)
3376          delete(files)
3377      end if
3378     
3379      end do
3380
3381      ;#########ENDE DO LOOP FOR no_files GT 1#############
3382
3383      if (isStrSubset(data@long_name," SR " ) .and. over_remind) then
3384          print(" ")
3385          print("If you have outputs of statistic regions you cannot overlay "+\
3386                "variables;")
3387          print("'over' is set to 0" )
3388          print(" ")
3389          over = 0
3390      end if
3391
3392      if (count_var .EQ. 0) then
3393          print(" ")
3394          print("The variables 'var="+var+"'" )
3395          print("do not exist on your input file;")
3396          print("be sure to have one comma before and after each variable")
3397          print(" ")
3398          exit       
3399      end if
3400     
3401      if (no_files .GT. 1) then
3402          multi_legend=new(6,string)
3403          string_len=new(6,integer)
3404          multi_dash=new(no_files,string)
3405          multi_legend(0)="  "+name_legend_1
3406          string_len(0)=strlen(multi_legend(0))
3407          multi_legend(1)="  "+name_legend_2
3408          string_len(1)=strlen(multi_legend(1))
3409          multi_legend(2)="  "+name_legend_3
3410          string_len(2)=strlen(multi_legend(2))
3411          multi_legend(3)="  "+name_legend_4
3412          string_len(3)=strlen(multi_legend(3))
3413          multi_legend(4)="  "+name_legend_5
3414          string_len(4)=strlen(multi_legend(4))
3415          multi_legend(5)="  "+name_legend_6
3416          string_len(5)=strlen(multi_legend(5))
3417          do ml=1,no_files
3418            multi_dash(ml-1)=ml-1
3419          end do
3420          delete(plot)
3421          plot = new(dim,graphic)
3422          do pl=0,n-1
3423            plot0 = new(1,graphic)
3424            res@trXMinF = min(min_nof(:,pl))
3425            res@trXMaxF = max(max_nof(:,pl))
3426            res@gsnLeftString  = " "; name(0,pl)
3427            res@gsnRightString = " ";unit_(0,pl)
3428            res@tiXAxisString  = "(" + unit_(0,pl) + ")"
3429            data_0(:,:) = min(min_nof(:,pl))
3430
3431            if (mod(no_plots - pl + 1, no_rows * no_columns) .EQ. 2) then
3432                res@pmLegendDisplayMode  = "Always"
3433            else
3434                res@pmLegendDisplayMode  = "Never"
3435            end if
3436
3437            plot0 = gsn_csm_xy(wks,data_0(:,:),z_(pl,:),res)
3438
3439            ; ***************************************************
3440            ; legend for combined plot
3441            ; ***************************************************
3442
3443            if (mod(no_plots - pl + 1, no_rows * no_columns) .EQ. 2)then
3444                lgres                    = True
3445                lgMonoDashIndex          = False
3446                lgres@lgLabelFont        = "helvetica"   
3447                lgres@lgLabelFontHeightF = font_size_legend * 10.0       
3448                lgres@vpWidthF           = max(string_len)*0.02;0.015           
3449                lgres@vpHeightF          = 0.03*no_files   
3450                lgres@lgDashIndexes      = multi_dash(no_files-1:0)
3451                lbid = gsn_create_legend(\
3452                                    wks,no_files,multi_legend(no_files-1:0),lgres)
3453
3454                amres = True
3455                amres@amParallelPosF   = -0.5       
3456                amres@amOrthogonalPosF = -0.55
3457                amres@amJust           = "BottomLeft"
3458                annoid1 = gsn_add_annotation(plot0,lbid,amres)
3459            end if
3460
3461            do plo=0,no_files-1
3462                overlay(plot0,multi_plot(plo,pl))
3463                plot(pl)=plot0
3464            end do
3465            delete(plot0)
3466          end do
3467      end if
3468
3469      if (count_var .EQ. 0) then
3470          print(" ")
3471          print("Select a variable 'var=' or use the default value")
3472          print(" ")
3473          print("Your selection '"+var+"' does not exist on the input file")
3474          print(" ")
3475          exit
3476      end if
3477
3478      if (over .EQ. 1 ) then
3479
3480          overlay(plot_u,plot_v)
3481          overlay(plot_u,plot_w)
3482          u=0
3483          overlay(plot_pt,plot_vpt)
3484          overlay(plot_pt,plot_lpt)
3485          pt=0
3486          overlay(plot_q,plot_qv)
3487          overlay(plot_q,plot_ql)
3488          q=0
3489          overlay(plot_e,plot_es)
3490          e=0
3491          overlay(plot_km,plot_kh)
3492          km=0
3493          overlay(plot_wpup,plot_wsus)
3494          overlay(plot_wpup,plot_wu)
3495          wpup=0
3496          overlay(plot_wpvp,plot_wsvs)
3497          overlay(plot_wpvp,plot_wv)
3498          wpvp=0
3499          overlay(plot_wpptp,plot_wspts)
3500          overlay(plot_wpptp,plot_wpt)
3501          wpptp=0
3502          overlay(plot_wsptsBC,plot_wptBC)
3503          wsptsBC=0
3504          overlay(plot_wpvptp,plot_wsvpts)
3505          overlay(plot_wpvptp,plot_wvpt)
3506          wpvptp=0
3507          overlay(plot_wpqp,plot_wsqs)
3508          overlay(plot_wpqp,plot_wq)
3509          wpqp=0
3510          overlay(plot_wpqvp,plot_wsqvs)
3511          overlay(plot_wpqvp,plot_wqv)
3512          wpqvp=0
3513          overlay(plot_wpsp,plot_wsss)
3514          overlay(plot_wpsp,plot_ws)
3515          wpsp=0
3516          overlay(plot_wpsap,plot_wssas)
3517          overlay(plot_wpsap,plot_wsa)
3518          wpsap=0
3519          overlay(plot_us2,plot_vs2)
3520          overlay(plot_us2,plot_ws2)
3521          us2=0
3522          overlay(plot_wsususodz,plot_wspsodz)
3523          overlay(plot_wsususodz,plot_wpeodz)
3524          wsususodz=0
3525
3526      end if
3527
3528      if (over .EQ. 1) then
3529     
3530          do varn = 0,dim-1   
3531           
3532            check = True
3533         
3534            if (prof3d .EQ. 0) then
3535                if ( isStrSubset( vNam(varn), "time") .OR. \
3536                    isStrSubset( vNam(varn), "NORM")) then
3537                  check = False
3538                end if
3539            else
3540                if ( isStrSubset( vNam(varn), "time") .OR.  \
3541                    isStrSubset( vNam(varn), "zusi") .OR.  \
3542                    isStrSubset( vNam(varn), "zwwi") .OR.  \
3543                    isStrSubset( vNam(varn), "x") .OR.     \
3544                    isStrSubset( vNam(varn), "xu") .OR.    \
3545                    isStrSubset( vNam(varn), "y") .OR.     \
3546                    isStrSubset( vNam(varn), "yv") .OR.    \
3547                    isStrSubset( vNam(varn), "zu_3d") .OR. \
3548                    isStrSubset( vNam(varn), "zw_3d")) then
3549                  check = False
3550                end if
3551            end if
3552
3553            if (var .NE. "all") then     
3554                check = isStrSubset( var,","+vNam(varn)+"," )
3555            end if     
3556
3557            if (check)then
3558
3559                if (prof3d .EQ. 0) then
3560                  if (log_z .EQ. 1) then
3561                      z = f_att->$vNam(varn+1)$(1:dimz-1)
3562                  else
3563                      z = f_att->$vNam(varn+1)$               
3564                  end if
3565                else
3566                  do i=0,b-1           
3567                      if (isStrSubset( a(i),"zu_3d" ))then
3568                        z_v(varn,:) = z_u
3569                        if (log_z .EQ. 1) then
3570                            z = z_v(varn,1:dimz-1)
3571                        else
3572                            z = z_v(varn,:)
3573                        end if
3574                      else
3575                        if (isStrSubset( a(i),"zw_3d" ))then
3576                            z_v(varn,:) = z_w
3577                            if (log_z .EQ. 1) then
3578                              z = z_v(varn,1:dimz-1)
3579                            else
3580                              z = z_v(varn,:)
3581                            end if
3582                        end if                   
3583                      end if
3584                  end do           
3585                end if
3586
3587                z=z/norm_z
3588               
3589                if (abs(min(data(varn,:,min_z_int:max_z_int))) .GT. 10)then
3590                  min_value = abs(0.01*min(data(varn,:,min_z_int:max_z_int)))
3591                  max_value = abs(0.01*max(data(varn,:,min_z_int:max_z_int)))
3592                else
3593                  if (abs(min(data(varn,:,min_z_int:max_z_int))) .LT. 0.01 .AND. \
3594                      abs(max(data(varn,:,min_z_int:max_z_int))) .GT. 0.01)then
3595                      min_value = abs(0.1*max(data(varn,:,min_z_int:max_z_int)))
3596                      max_value = abs(0.1*max(data(varn,:,min_z_int:max_z_int)))
3597                  else
3598                      if (abs(max(data(varn,:,:))) .LT. 0.01 .AND. \
3599                          abs(min(data(varn,:,min_z_int:max_z_int))) .GT. 0.01)then
3600                        min_value = abs(0.1*min(data(varn,:,min_z_int:max_z_int)))
3601                        max_value = abs(0.1*min(data(varn,:,min_z_int:max_z_int)))
3602                      else
3603                        min_value = abs(0.1*min(data(varn,:,min_z_int:max_z_int)))
3604                        max_value = abs(0.1*max(data(varn,:,min_z_int:max_z_int)))
3605                      end if
3606                  end if
3607                end if
3608                if (min(data(varn,:,min_z_int:max_z_int)) .EQ. 0 .AND.\
3609                    max(data(varn,:,min_z_int:max_z_int)) .EQ. 0)then
3610                  min_value = 0.1
3611                  max_value = 0.1
3612                end if
3613
3614                res@gsnLeftString      = vNam(varn)
3615                res@tiXAxisString      = "("+unit(varn)+")"
3616                res@gsnRightString     = " "
3617                res@trYMinF            = min_z
3618                res@trYMaxF            = max_z 
3619                res@xyDashPattern = 0
3620
3621                if (xs .EQ. -1) then
3622                  res@trXMinF = min(data(varn,:,min_z_int:max_z_int))-min_value
3623                else
3624                  res@trXMinF = xs
3625                end if
3626                if (xe .EQ. -1) then
3627                  res@trXMaxF = max(data(varn,:,min_z_int:max_z_int))+max_value
3628                else
3629                  res@trXMaxF = xe 
3630                end if
3631                plot(n) = gsn_csm_xy(wks,data(varn,:,:),z,res) 
3632               
3633                if (vNam(varn) .EQ. "u" .OR. vNam(varn) .EQ. "v" .OR. \
3634                    vNam(varn) .EQ. "w") then
3635                  if (u .EQ. 0) then
3636                      res@gsnLeftString      = "u, v and w"
3637                      res@tiXAxisString      = "("+unit(varn)+")"
3638                      res@gsnRightString     = " "
3639                      if (xs .EQ. -1) then
3640                        res@trXMinF = min((/miniu,miniv,miniw/))
3641                      else
3642                        res@trXMinF = xs
3643                      end if
3644                      if (xe .EQ. -1) then
3645                        res@trXMaxF = max((/maxiu,maxiv,maxiw/))
3646                      else
3647                        res@trXMaxF = xe 
3648                      end if 
3649                      if (vNam(varn) .EQ. "v")then
3650                        res@xyDashPattern = 1
3651                      end if
3652                      if (vNam(varn) .EQ. "w")then
3653                        res@xyDashPattern = 2
3654                      end if
3655                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
3656
3657                      ; ***************************************************
3658                      ; legend for overlaid plot
3659                      ; ***************************************************
3660         
3661                      lgres                    = True
3662                      lgMonoDashIndex          = False
3663                      lgres@lgLabelFont        = "helvetica"   
3664                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
3665                      lgres@vpWidthF           = 0.07           
3666                      lgres@vpHeightF          = 0.12         
3667                      lgres@lgDashIndexes      = (/0,1,2/)
3668                      lbid = gsn_create_legend(wks,3,(/"u","v","w"/),lgres)       
3669
3670                      amres = True
3671                      amres@amParallelPosF   = 0.88                 
3672                      amres@amOrthogonalPosF = 0.33           
3673                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
3674                      overlay(plot(n),plot_u)
3675                      u=1
3676                  else
3677                      if (prof3d .EQ. 0)then
3678                        varn=varn+1
3679                      end if
3680                      continue
3681                  end if       
3682                end if 
3683       
3684                if (vNam(varn) .EQ. "pt" .OR. vNam(varn) .EQ. "vpt" .OR. \
3685                    vNam(varn) .EQ. "lpt") then
3686                  if (pt .EQ. 0) then
3687                      res@gsnLeftString      = "pt, vpt and lpt"
3688                      res@tiXAxisString      = "("+unit(varn)+")"
3689                      res@gsnRightString     = " "
3690                      if (xs .EQ. -1) then
3691                        res@trXMinF = min((/minipt,minivpt,minilpt/))
3692                      else
3693                        res@trXMinF = xs
3694                      end if
3695                      if (xe .EQ. -1) then
3696                        res@trXMaxF = max((/maxipt,maxivpt,maxilpt/))
3697                      else
3698                        res@trXMaxF = xe 
3699                      end if
3700                      if (vNam(varn) .EQ. "vpt")then
3701                        res@xyDashPattern = 1
3702                      end if
3703                      if (vNam(varn) .EQ. "lpt")then
3704                        res@xyDashPattern = 2
3705                      end if
3706                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
3707
3708                      ; ***************************************************
3709                      ; legend for overlaid plot
3710                      ; ***************************************************
3711         
3712                      lgres                    = True
3713                      lgMonoDashIndex          = False
3714                      lgres@lgLabelFont        = "helvetica"   
3715                      lgres@lgLabelFontHeightF = font_size_legend*10.0         
3716                      lgres@vpWidthF           = 0.07           
3717                      lgres@vpHeightF          = 0.12         
3718                      lgres@lgDashIndexes      = (/0,1,2/)
3719                      lbid = gsn_create_legend(wks,3,(/"pt","vpt","lpt"/),lgres)
3720                      amres = True
3721                      amres@amParallelPosF   = 0.88         
3722                      amres@amOrthogonalPosF = 0.33           
3723                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
3724                      overlay(plot(n),plot_pt)
3725                      pt=1
3726                  else
3727                      if (prof3d .EQ. 0)then
3728                        varn=varn+1
3729                      end if
3730                      continue       
3731                  end if
3732                end if           
3733                if (vNam(varn) .EQ. "q" .OR. vNam(varn) .EQ. "qv" .OR. \
3734                    vNam(varn) .EQ. "ql") then
3735                  if (q .EQ. 0) then
3736                      res@gsnLeftString      = "q, qv and ql"
3737                      res@tiXAxisString      = "("+unit(varn)+")"
3738                      res@gsnRightString     = " "
3739                      if (xs .EQ. -1) then
3740                        res@trXMinF = min((/miniq,miniqv,miniql/))
3741                      else
3742                        res@trXMinF = xs
3743                      end if
3744                      if (xe .EQ. -1) then
3745                        res@trXMaxF = max((/maxiq,maxiqv,maxiql/))
3746                      else
3747                        res@trXMaxF = xe 
3748                      end if
3749                      if (vNam(varn) .EQ. "qv")then
3750                        res@xyDashPattern = 1
3751                      end if
3752                      if (vNam(varn) .EQ. "ql")then
3753                        res@xyDashPattern = 2
3754                      end if
3755                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
3756
3757                      ; ***************************************************
3758                      ; legend for overlaid plot
3759                      ; ***************************************************
3760         
3761                      lgres                    = True
3762                      lgMonoDashIndex          = False
3763                      lgres@lgLabelFont        = "helvetica"   
3764                      lgres@lgLabelFontHeightF = font_size_legend*10.0         
3765                      lgres@vpWidthF           = 0.07         
3766                      lgres@vpHeightF          = 0.12         
3767                      lgres@lgDashIndexes      = (/0,1,2/)
3768                      lbid = gsn_create_legend(wks,3,(/"q","qv","ql"/),lgres)
3769
3770                      amres = True
3771                      amres@amParallelPosF   = 0.88                 
3772                      amres@amOrthogonalPosF = 0.33           
3773                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
3774                      overlay(plot(n),plot_q)
3775                      q=1
3776                  else
3777                      if (prof3d .EQ. 0)then
3778                        varn=varn+1
3779                      end if
3780                      continue   
3781                  end if
3782                end if   
3783             
3784                if (vNam(varn) .EQ. "e" .OR. vNam(varn) .EQ. "es" .OR. \
3785                    vNam(varn) .EQ. "e*" ) then
3786                  if (e .EQ. 0) then
3787                      res@gsnLeftString      = "e and e*"
3788                      res@tiXAxisString      = "("+unit(varn)+")"
3789                      res@gsnRightString     = " "
3790                      if (xs .EQ. -1) then
3791                        res@trXMinF = min((/minie,minies/))
3792                      else
3793                        res@trXMinF = xs
3794                      end if
3795                      if (xe .EQ. -1) then
3796                        res@trXMaxF = max((/maxie,maxies/))
3797                      else
3798                        res@trXMaxF = xe 
3799                      end if
3800                      if (vNam(varn) .EQ. "es")then
3801                        res@xyDashPattern = 1
3802                      end if
3803                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
3804
3805                      ; ***************************************************
3806                      ; legend for overlaid plot
3807                      ; ***************************************************
3808         
3809                      lgres                    = True
3810                      lgMonoDashIndex          = False
3811                      lgres@lgLabelFont        = "helvetica"   
3812                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
3813                      lgres@vpWidthF           = 0.07           
3814                      lgres@vpHeightF          = 0.08         
3815                      lgres@lgDashIndexes      = (/0,1,2/)
3816                      lbid = gsn_create_legend(wks,2,(/"e","e*"/),lgres)       
3817
3818                      amres = True
3819                      amres@amParallelPosF   = 0.88                 
3820                      amres@amOrthogonalPosF = 0.365           
3821                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
3822                      overlay(plot(n),plot_e)
3823                      e=1
3824                  else
3825                      if (prof3d .EQ. 0)then
3826                        varn=varn+1
3827                      end if
3828                      continue   
3829                  end if
3830                end if           
3831                if (vNam(varn) .EQ. "km" .OR. vNam(varn) .EQ. "kh") then
3832                  if (km .EQ. 0) then
3833                      res@gsnLeftString      = "km and kh"
3834                      res@tiXAxisString      = "("+unit(varn)+")"
3835                      res@gsnRightString     = " "
3836                      if (xs .EQ. -1) then
3837                        res@trXMinF = min((/minikm,minikh/))
3838                      else
3839                        res@trXMinF = xs
3840                      end if
3841                      if (xe .EQ. -1) then
3842                        res@trXMaxF = max((/maxikm,maxikh/))
3843                      else
3844                        res@trXMaxF = xe 
3845                      end if
3846                      if (vNam(varn) .EQ. "kh")then
3847                        res@xyDashPattern = 1
3848                      end if
3849                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
3850
3851                      ; ***************************************************
3852                      ; legend for overlaid plot
3853                      ; ***************************************************
3854         
3855                      lgres                    = True
3856                      lgMonoDashIndex          = False
3857                      lgres@lgLabelFont        = "helvetica"   
3858                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
3859                      lgres@vpWidthF           = 0.07           
3860                      lgres@vpHeightF          = 0.08         
3861                      lgres@lgDashIndexes      = (/0,1,2/)
3862                      lbid = gsn_create_legend(wks,2,(/"km","kh"/),lgres)       
3863
3864                      amres = True
3865                      amres@amParallelPosF   = 0.88                 
3866                      amres@amOrthogonalPosF = 0.365           
3867                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
3868                      overlay(plot(n),plot_km)
3869                      km=1
3870                  else
3871                      if (prof3d .EQ. 0)then
3872                        varn=varn+1
3873                      end if
3874                      continue   
3875                  end if
3876                end if           
3877             
3878                if (vNam(varn) .EQ. "wpup" .OR. vNam(varn) .EQ. "wsus" .OR.      \
3879                    vNam(varn) .EQ. "wu" .OR. vNam(varn) .EQ. "w"+dq+"u"+dq .OR. \
3880                    vNam(varn) .EQ. "w*u*") then
3881                  if (wpup .EQ. 0) then
3882                      res@gsnLeftString      = "w"+dq+"u"+dq+", w*u* and wu"
3883                      res@tiXAxisString      = "("+unit(varn)+")"
3884                      res@gsnRightString     = " "
3885                      if (xs .EQ. -1) then
3886                        res@trXMinF = min((/miniwpup,miniwsus,miniwu/))
3887                      else
3888                        res@trXMinF = xs
3889                      end if
3890                      if (xe .EQ. -1) then
3891                        res@trXMaxF = max((/maxiwpup,maxiwsus,maxiwu/))
3892                      else
3893                        res@trXMaxF = xe 
3894                      end if 
3895                      if (vNam(varn) .EQ. "wsus")then
3896                        res@xyDashPattern = 1
3897                      end if
3898                      if (vNam(varn) .EQ. "wu")then
3899                        res@xyDashPattern = 2
3900                      end if
3901                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
3902
3903                      ; ***************************************************
3904                      ; legend for overlaid plot
3905                      ; ***************************************************
3906         
3907                      lgres                    = True
3908                      lgMonoDashIndex          = False
3909                      lgres@lgLabelFont        = "helvetica"   
3910                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
3911                      lgres@vpWidthF           = 0.08           
3912                      lgres@vpHeightF          = 0.12         
3913                      lgres@lgDashIndexes      = (/0,1,2/)
3914                      lbid = gsn_create_legend(\
3915                                    wks,3,(/"w"+dq+"u"+dq,"w*u*","wu"/),lgres)
3916
3917                      amres = True
3918                      amres@amParallelPosF   = 0.88                 
3919                      amres@amOrthogonalPosF = 0.33           
3920                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
3921                      overlay(plot(n),plot_wpup)
3922                      wpup=1
3923                  else
3924                      if (prof3d .EQ. 0)then
3925                        varn=varn+1
3926                      end if
3927                      continue   
3928                  end if
3929                end if
3930                if (vNam(varn) .EQ. "wpvp" .OR. vNam(varn) .EQ. "wsvs" .OR.     \
3931                  vNam(varn) .EQ. "wv" .OR. vNam(varn) .EQ. "w"+dq+"v"+dq .OR. \
3932                  vNam(varn) .EQ. "w*v*") then
3933                  if (wpvp .EQ. 0) then
3934                      res@gsnLeftString      = "w"+dq+"v"+dq+", w*v* and wv"
3935                      res@tiXAxisString      = "("+unit(varn)+")"
3936                      res@gsnRightString     = " "
3937                      if (xs .EQ. -1) then
3938                        res@trXMinF = min((/miniwpvp,miniwsvs,miniwv/))
3939                      else
3940                        res@trXMinF = xs
3941                      end if
3942                      if (xe .EQ. -1) then
3943                        res@trXMaxF = max((/maxiwpvp,maxiwsvs,maxiwv/))
3944                      else
3945                        res@trXMaxF = xe 
3946                      end if
3947                      if (vNam(varn) .EQ. "wsvs")then
3948                        res@xyDashPattern = 1
3949                      end if
3950                      if (vNam(varn) .EQ. "wv")then
3951                        res@xyDashPattern = 2
3952                      end if
3953                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
3954
3955                      ; ***************************************************
3956                      ; legend for overlaid plot
3957                      ; ***************************************************
3958         
3959                      lgres                    = True
3960                      lgMonoDashIndex          = False
3961                      lgres@lgLabelFont        = "helvetica"   
3962                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
3963                      lgres@vpWidthF           = 0.08         
3964                      lgres@vpHeightF          = 0.12         
3965                      lgres@lgDashIndexes      = (/0,1,2/)
3966                      lbid = gsn_create_legend(\
3967                                wks,3,(/"w"+dq+"v"+dq,"w*v*","wv"/),lgres)       
3968
3969                      amres = True
3970                      amres@amParallelPosF   = 0.88                 
3971                      amres@amOrthogonalPosF = 0.33           
3972                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
3973                      overlay(plot(n),plot_wpvp)
3974                      wpvp=1
3975                  else
3976                      if (prof3d .EQ. 0)then
3977                        varn=varn+1
3978                      end if
3979                      continue   
3980                  end if
3981                end if
3982                if (vNam(varn) .EQ. "wpptp" .OR. vNam(varn) .EQ. "wspts" .OR.     \
3983                    vNam(varn) .EQ. "wpt" .OR. vNam(varn) .EQ. "w"+dq+"pt"+dq .OR.\
3984                    vNam(varn) .EQ. "w*pt*") then
3985                  if (wpptp .EQ. 0) then                 
3986                      res@gsnLeftString      = "w"+dq+"pt"+dq+", w*pt* and wpt"
3987                      res@tiXAxisString      = "("+unit(varn)+")"
3988                      res@gsnRightString     = " "
3989                      if (xs .EQ. -1) then
3990                        res@trXMinF = min((/miniwpptp,miniwspts,miniwpt/))
3991                      else
3992                        res@trXMinF = xs
3993                      end if
3994                      if (xe .EQ. -1) then
3995                        res@trXMaxF = max((/maxiwpptp,maxiwspts,maxiwpt/))
3996                      else
3997                        res@trXMaxF = xe 
3998                      end if
3999                      if (vNam(varn) .EQ. "wspts" .OR. vNam(varn) .EQ. "w*pt*")then
4000                        res@xyDashPattern = 1
4001                      end if
4002                      if (vNam(varn) .EQ. "wpt")then
4003                        res@xyDashPattern = 2
4004                      end if
4005                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
4006
4007                      ; ***************************************************
4008                      ; legend for overlaid plot
4009                      ; ***************************************************
4010         
4011                      lgres                    = True
4012                      lgMonoDashIndex          = False
4013                      lgres@lgLabelFont        = "helvetica"   
4014                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
4015                      lgres@vpWidthF           = 0.09           
4016                      lgres@vpHeightF          = 0.12         
4017                      lgres@lgDashIndexes      = (/0,1,2/)             
4018                      lbid = gsn_create_legend(\
4019                                  wks,3,(/"w"+dq+"pt"+dq,"w*pt*","wpt"/),lgres)
4020
4021                      amres = True
4022                      amres@amParallelPosF   = 0.88                 
4023                      amres@amOrthogonalPosF = 0.33           
4024                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
4025                      overlay(plot(n),plot_wpptp)
4026                      wpptp=1
4027                  else
4028                      if (prof3d .EQ. 0)then
4029                        varn=varn+1
4030                      end if
4031                      continue   
4032                  end if
4033                end if
4034                if (vNam(varn) .EQ. "wsptsBC" .OR. vNam(varn) .EQ. "wptBC" .OR.\
4035                    vNam(varn) .EQ. "w*pt*BC") then
4036                  if (wsptsBC .EQ. 0) then
4037                      res@gsnLeftString      = "w*pt*BC and wptBC"
4038                      res@tiXAxisString      = "("+unit(varn)+")"
4039                      res@gsnRightString     = " "
4040                      if (xs .EQ. -1) then
4041                        res@trXMinF = min((/miniwsptsBC,miniwptBC/))
4042                      else
4043                        res@trXMinF = xs
4044                      end if
4045                      if (xe .EQ. -1) then
4046                        res@trXMaxF = max((/maxiwsptsBC,maxiwptBC/))
4047                      else
4048                        res@trXMaxF = xe 
4049                      end if
4050                      if (vNam(varn) .EQ. "wptBC")then
4051                        res@xyDashPattern = 1
4052                      end if
4053                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
4054
4055                      ; ***************************************************
4056                      ; legend for overlaid plot
4057                      ; ***************************************************
4058         
4059                      lgres                    = True
4060                      lgMonoDashIndex          = False
4061                      lgres@lgLabelFont        = "helvetica"   
4062                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
4063                      lgres@vpWidthF           = 0.1           
4064                      lgres@vpHeightF          = 0.12         
4065                      lgres@lgDashIndexes      = (/0,1,2/)
4066                      lbid = gsn_create_legend(\
4067                                            wks,3,(/"w*pt*BC","wptBC"/),lgres)
4068
4069                      amres = True
4070                      amres@amParallelPosF   = 0.88                 
4071                      amres@amOrthogonalPosF = 0.33           
4072                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
4073                      overlay(plot(n),plot_wsptsBC)
4074                      wsptsBC=1
4075                  else
4076                      if (prof3d .EQ. 0)then
4077                        varn=varn+1
4078                      end if
4079                      continue   
4080                  end if 
4081                end if             
4082                if (vNam(varn) .EQ. "wpvptp" .OR. vNam(varn) .EQ. "wsvpts" .OR. \
4083                    vNam(varn) .EQ. "wvpt" .OR. vNam(varn) .EQ. \
4084                    "w"+dq+"vpt"+dq .OR. vNam(varn) .EQ. "w*vpt*") then
4085                  if (wpvptp .EQ. 0) then
4086                      res@gsnLeftString      = "w"+dq+"vpt"+dq+", w*vpt* and wvpt"
4087                      res@tiXAxisString      = "("+unit(varn)+")"
4088                      res@gsnRightString     = " "
4089                      if (xs .EQ. -1) then
4090                        res@trXMinF = min((/miniwpvptp,miniwsvpts,miniwvpt/))
4091                      else
4092                        res@trXMinF = xs
4093                      end if
4094                      if (xe .EQ. -1) then
4095                        res@trXMaxF = max((/maxiwpvptp,maxiwsvpts,maxiwvpt/))
4096                      else
4097                        res@trXMaxF = xe 
4098                      end if
4099                      if (vNam(varn) .EQ. "wsvpts")then
4100                        res@xyDashPattern = 1
4101                      end if
4102                      if (vNam(varn) .EQ. "wvpt")then
4103                        res@xyDashPattern = 2
4104                      end if
4105                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
4106
4107                      ; ***************************************************
4108                      ; legend for overlaid plot
4109                      ; ***************************************************
4110         
4111                      lgres                    = True
4112                      lgMonoDashIndex          = False
4113                      lgres@lgLabelFont        = "helvetica"   
4114                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
4115                      lgres@vpWidthF           = 0.1           
4116                      lgres@vpHeightF          = 0.12         
4117                      lgres@lgDashIndexes      = (/0,1,2/)
4118                      lbid = gsn_create_legend(\
4119                                wks,3,(/"w"+dq+"vpt"+dq,"w*vpt*","wvpt"/),lgres)
4120                      amres = True
4121                      amres@amParallelPosF   = 0.88                 
4122                      amres@amOrthogonalPosF = 0.33           
4123                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
4124                      overlay(plot(n),plot_wpvptp)
4125                      wpvptp=1
4126                  else
4127                      if (prof3d .EQ. 0)then
4128                        varn=varn+1
4129                      end if
4130                      continue   
4131                  end if
4132                end if
4133                if (vNam(varn) .EQ. "wpqp" .OR. vNam(varn) .EQ. "wsqs" .OR.      \
4134                    vNam(varn) .EQ. "wq" .OR. vNam(varn) .EQ. "w"+dq+"q"+dq .OR. \
4135                    vNam(varn) .EQ. "w*q*") then
4136                  if (wpqp .EQ. 0) then
4137                      res@gsnLeftString      = "w"+dq+"q"+dq+", w*q* and wq"
4138                      res@tiXAxisString      = "("+unit(varn)+")"
4139                      res@gsnRightString     = " "
4140                      if (xs .EQ. -1) then
4141                        res@trXMinF = min((/miniwpqp,miniwsqs,miniwq/))
4142                      else
4143                        res@trXMinF = xs
4144                      end if
4145                      if (xe .EQ. -1) then
4146                        res@trXMaxF = max((/maxiwpqp,maxiwsqs,maxiwq/))
4147                      else
4148                        res@trXMaxF = xe 
4149                      end if 
4150                      if (vNam(varn) .EQ. "wsqs")then
4151                        res@xyDashPattern = 1
4152                      end if
4153                      if (vNam(varn) .EQ. "wq")then
4154                        res@xyDashPattern = 2
4155                      end if
4156                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
4157
4158                      ; ***************************************************
4159                      ; legend for overlaid plot
4160                      ; ***************************************************
4161         
4162                      lgres                    = True
4163                      lgMonoDashIndex          = False
4164                      lgres@lgLabelFont        = "helvetica"   
4165                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
4166                      lgres@vpWidthF           = 0.08           
4167                      lgres@vpHeightF          = 0.12         
4168                      lgres@lgDashIndexes      = (/0,1,2/)
4169                      lbid = gsn_create_legend(\
4170                                wks,3,(/"w"+dq+"q"+dq,"w*q*","wq"/),lgres)       
4171
4172                      amres = True
4173                      amres@amParallelPosF   = 0.88                 
4174                      amres@amOrthogonalPosF = 0.33           
4175                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
4176                      overlay(plot(n),plot_wpqp)
4177                      wpqp=1
4178                  else
4179                      if (prof3d .EQ. 0)then
4180                        varn=varn+1
4181                      end if
4182                      continue   
4183                  end if
4184                end if
4185                if (vNam(varn) .EQ. "wpqvp" .OR. vNam(varn) .EQ. "wsqvs" .OR.     \
4186                    vNam(varn) .EQ. "wqv" .OR. vNam(varn) .EQ. "w"+dq+"qv"+dq .OR.\
4187                    vNam(varn) .EQ. "w*qv*") then
4188                  if (wpqvp .EQ. 0) then
4189                      res@gsnLeftString      ="w"+dq+"qv"+dq+" , w*qv* and wqv"
4190                      res@tiXAxisString      = "("+unit(varn)+")"
4191                      res@gsnRightString     = " "
4192                      if (xs .EQ. -1) then
4193                        res@trXMinF = min((/miniwpqp,miniwsqvs,miniwqv/))
4194                      else
4195                        res@trXMinF = xs
4196                      end if
4197                      if (xe .EQ. -1) then
4198                        res@trXMaxF = max((/maxiwpqp,maxiwsqvs,maxiwqv/))
4199                      else
4200                        res@trXMaxF = xe 
4201                      end if
4202                      if (vNam(varn) .EQ. "wsqvs")then
4203                        res@xyDashPattern = 1
4204                      end if
4205                      if (vNam(varn) .EQ. "wqv")then
4206                        res@xyDashPattern = 2
4207                      end if
4208                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
4209
4210                      ; ***************************************************
4211                      ; legend for overlaid plot
4212                      ; ***************************************************
4213         
4214                      lgres                    = True
4215                      lgMonoDashIndex          = False
4216                      lgres@lgLabelFont        = "helvetica"   
4217                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
4218                      lgres@vpWidthF           = 0.09           
4219                      lgres@vpHeightF          = 0.12         
4220                      lgres@lgDashIndexes      = (/0,1,2/)
4221                      lbid = gsn_create_legend(\
4222                                    wks,3,(/"w"+dq+"qv"+dq,"w*qv*","wqv"/),lgres)
4223
4224                      amres = True
4225                      amres@amParallelPosF   = 0.88                 
4226                      amres@amOrthogonalPosF = 0.33           
4227                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
4228                      overlay(plot(n),plot_wpqvp)
4229                      wpqvp=1
4230                  else
4231                      if (prof3d .EQ. 0)then
4232                        varn=varn+1
4233                      end if
4234                      continue   
4235                  end if
4236                end if
4237                if (vNam(varn) .EQ. "wpsp" .OR. vNam(varn) .EQ. "wsss" .OR.     \
4238                    vNam(varn) .EQ. "ws" .OR. vNam(varn) .EQ. "w"+dq+"s"+dq .OR.\
4239                    vNam(varn) .EQ. "w*s*") then
4240                  if (wpsp .EQ. 0) then
4241                      res@gsnLeftString      = "w"+dq+"s"+dq+", w*s* and ws"
4242                      res@tiXAxisString      = "("+unit(varn)+")"
4243                      res@gsnRightString     = " "
4244                      if (xs .EQ. -1) then
4245                        res@trXMinF = min((/miniwpsp,miniwsss,miniws/))
4246                      else
4247                        res@trXMinF = xs
4248                      end if
4249                      if (xe .EQ. -1) then
4250                        res@trXMaxF = max((/maxiwpsp,maxiwsss,maxiws/))
4251                      else
4252                        res@trXMaxF = xe 
4253                      end if
4254                      if (vNam(varn) .EQ. "wsss")then
4255                        res@xyDashPattern = 1
4256                      end if
4257                      if (vNam(varn) .EQ. "ws")then
4258                        res@xyDashPattern = 2
4259                      end if
4260                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
4261
4262                      ; ***************************************************
4263                      ; legend for overlaid plot
4264                      ; ***************************************************
4265         
4266                      lgres                    = True
4267                      lgMonoDashIndex          = False
4268                      lgres@lgLabelFont        = "helvetica"   
4269                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
4270                      lgres@vpWidthF           = 0.08           
4271                      lgres@vpHeightF          = 0.12         
4272                      lgres@lgDashIndexes      = (/0,1,2/)
4273                      lbid = gsn_create_legend(\
4274                              wks,3,(/"w"+dq+"s"+dq,"w*s*","ws"/),lgres)       
4275
4276                      amres = True
4277                      amres@amParallelPosF   = 0.88                 
4278                      amres@amOrthogonalPosF = 0.33           
4279                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
4280                      overlay(plot(n),plot_wpsp)
4281                      wpsp=1
4282                  else
4283                      if (prof3d .EQ. 0)then
4284                        varn=varn+1
4285                      end if
4286                      continue   
4287                  end if
4288                end if
4289                if (vNam(varn) .EQ. "wpsap" .OR.vNam(varn) .EQ. "wssas" .OR.      \
4290                    vNam(varn) .EQ. "wsa" .OR. vNam(varn) .EQ. "w"+dq+"sa"+dq .OR.\
4291                    vNam(varn) .EQ. "w*sa*") then
4292                  if (wpsap .EQ. 0) then
4293                      res@gsnLeftString      = "w"+dq+"sa"+dq+", w*sa* and wsa"
4294                      res@tiXAxisString      = "("+unit(varn)+")"
4295                      res@gsnRightString     = " "
4296                      if (xs .EQ. -1) then
4297                        res@trXMinF = min((/miniwpsap,miniwssas,miniwsa/))
4298                      else
4299                        res@trXMinF = xs
4300                      end if
4301                      if (xe .EQ. -1) then
4302                        res@trXMaxF = max((/maxiwpsap,maxiwssas,maxiwsa/))
4303                      else
4304                        res@trXMaxF = xe 
4305                      end if
4306                      if (vNam(varn) .EQ. "wssas")then
4307                        res@xyDashPattern = 1
4308                      end if
4309                      if (vNam(varn) .EQ. "wsa")then
4310                        res@xyDashPattern = 2
4311                      end if
4312                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
4313
4314                      ; ***************************************************
4315                      ; legend for overlaid plot
4316                      ; ***************************************************
4317         
4318                      lgres                    = True
4319                      lgMonoDashIndex          = False
4320                      lgres@lgLabelFont        = "helvetica"   
4321                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
4322                      lgres@vpWidthF           = 0.09           
4323                      lgres@vpHeightF          = 0.12         
4324                      lgres@lgDashIndexes      = (/0,1,2/)
4325                      lbid = gsn_create_legend(\
4326                                wks,3,(/"w"+dq+"sa"+dq,"w*sa*","wsa"/),lgres)
4327                      amres = True
4328                      amres@amParallelPosF   = 0.88                 
4329                      amres@amOrthogonalPosF = 0.33           
4330                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
4331                      overlay(plot(n),plot_wpsap)
4332                      wpsap=1
4333                  else
4334                      if (prof3d .EQ. 0)then
4335                        varn=varn+1
4336                      end if
4337                      continue   
4338                  end if
4339                end if
4340             
4341                if (vNam(varn) .EQ. "us2" .OR. vNam(varn) .EQ. "vs2" .OR. \
4342                    vNam(varn) .EQ. "ws2" .OR. vNam(varn) .EQ. "u*2" .OR. \
4343                    vNam(varn) .EQ. "v*2" .OR. vNam(varn) .EQ. "w*2" ) then
4344                  if (us2 .EQ. 0) then
4345                      res@gsnLeftString      = "u*2, v*2 and w*2"
4346                      res@tiXAxisString      = "("+unit(varn)+")"
4347                      res@gsnRightString     = " "
4348                      if (xs .EQ. -1) then
4349                        res@trXMinF = min((/minius2,minivs2,miniws2/))
4350                      else
4351                        res@trXMinF = xs
4352                      end if
4353                      if (xe .EQ. -1) then
4354                        res@trXMaxF = max((/maxius2,maxivs2,maxiws2/))
4355                      else
4356                        res@trXMaxF = xe 
4357                      end if
4358                      if (vNam(varn) .EQ. "vs2")then
4359                        res@xyDashPattern = 1
4360                      end if
4361                      if (vNam(varn) .EQ. "ws2")then
4362                        res@xyDashPattern = 2
4363                      end if
4364                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
4365
4366                      ; ***************************************************
4367                      ; legend for overlaid plot
4368                      ; ***************************************************
4369         
4370                      lgres                    = True
4371                      lgMonoDashIndex          = False
4372                      lgres@lgLabelFont        = "helvetica"   
4373                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
4374                      lgres@vpWidthF           = 0.07           
4375                      lgres@vpHeightF          = 0.12         
4376                      lgres@lgDashIndexes      = (/0,1,2/)
4377                      lbid = gsn_create_legend(wks,3,(/"u*2","v*2","w*2"/),lgres)
4378                      amres = True
4379                      amres@amParallelPosF   = 0.88                 
4380                      amres@amOrthogonalPosF = 0.33           
4381                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
4382                      overlay(plot(n),plot_us2)
4383                      us2=1
4384                  else
4385                      if (prof3d .EQ. 0)then
4386                        varn=varn+1
4387                      end if
4388                      continue   
4389                  end if
4390                end if
4391             
4392                if (vNam(varn) .EQ. "wsususodz" .OR. \
4393                    vNam(varn) .EQ. "wspsodz" .OR.   \
4394                    vNam(varn) .EQ. "wpeodz" .OR.    \
4395                    vNam(varn) .EQ. "w*u*u*:dz" .OR. \
4396                    vNam(varn) .EQ. "w*p*:dz" .OR.   \
4397                    vNam(varn) .EQ. "w"+dq+"e:dz") then
4398                  if (wsususodz .EQ. 0) then
4399                      res@gsnLeftString      = "w*u*u*:dz, w*p*:dz and w"+dq+"e:dz"
4400                      res@tiXAxisString      = "("+unit(varn)+")"
4401                      res@gsnRightString     = " "
4402                      if (xs .EQ. -1) then
4403                        res@trXMinF = min((/miniwsususodz,\
4404                                                      miniwspsodz,miniwpeodz/))
4405                      else
4406                        res@trXMinF = xs
4407                      end if
4408                      if (xe .EQ. -1) then
4409                        res@trXMaxF = max((/maxiwsususodz,maxiwspsodz,\
4410                                                                  maxiwpeodz/))
4411                      else
4412                        res@trXMaxF = xe 
4413                      end if
4414                      if (vNam(varn) .EQ. "wspsodz")then
4415                        res@xyDashPattern = 1
4416                      end if
4417                      if (vNam(varn) .EQ. "wpeodz")then
4418                        res@xyDashPattern = 2
4419                      end if
4420                      plot(n) =  gsn_csm_xy(wks,data(varn,:,:),z,res)
4421
4422                      ; ***************************************************
4423                      ; legend for overlaid plot
4424                      ; ***************************************************
4425         
4426                      lgres                    = True
4427                      lgMonoDashIndex          = False
4428                      lgres@lgLabelFont        = "helvetica"   
4429                      lgres@lgLabelFontHeightF = font_size_legend*10.0           
4430                      lgres@vpWidthF           = 0.12           
4431                      lgres@vpHeightF          = 0.12         
4432                      lgres@lgDashIndexes      = (/0,1,2/)
4433                      lbid = gsn_create_legend(\
4434                              wks,3,(/"w*u*u*:dz","w*p*:dz","w"+dq+"e:dz"/),lgres)
4435                      amres = True
4436                      amres@amParallelPosF   = 0.88                 
4437                      amres@amOrthogonalPosF = 0.33           
4438                      annoid1 = gsn_add_annotation(plot(n),lbid,amres)
4439                      overlay(plot(n),plot_wsususodz)
4440                      wsususodz=1
4441                  else
4442                      if (prof3d .EQ. 0)then
4443                        varn=varn+1
4444                      end if
4445                      continue   
4446                  end if
4447                end if     
4448                n=n+1
4449                if (prof3d .EQ. 0)then
4450                  varn=varn+1
4451                end if   
4452            end if
4453          end do
4454      end if
4455     
4456      if (com_i .NE. 0) then
4457          delete(com_var_avail)
4458      end if
4459      com_var_avail=new(count_var,string)   
4460
4461      if (combine .EQ. 1) then
4462          co=0     
4463          n_o=0
4464          do varn = 0,dim-1
4465   
4466            check = True
4467         
4468            if (prof3d .EQ. 0) then
4469                if ( isStrSubset( vNam(varn), "time") .OR. \
4470                    isStrSubset( vNam(varn), "NORM")) then
4471                  check = False
4472                end if
4473            else
4474                if ( isStrSubset( vNam(varn), "time") .OR.  \
4475                    isStrSubset( vNam(varn), "zusi") .OR.  \
4476                    isStrSubset( vNam(varn), "zwwi") .OR.  \
4477                    isStrSubset( vNam(varn), "x") .OR.     \
4478                    isStrSubset( vNam(varn), "xu") .OR.    \
4479                    isStrSubset( vNam(varn), "y") .OR.     \
4480                    isStrSubset( vNam(varn), "yv") .OR.    \
4481                    isStrSubset( vNam(varn), "zu_3d") .OR. \
4482                    isStrSubset( vNam(varn), "zw_3d")) then
4483                  check = False
4484                end if
4485            end if
4486
4487            if (var .NE. "all") then         
4488                check = isStrSubset( var,","+vNam(varn)+"," )
4489            end if     
4490
4491            if (check)then
4492               
4493                if (prof3d .EQ. 0) then
4494                  if (log_z .EQ. 1) then
4495                      z = f_att->$vNam(varn+1)$(1:dimz-1)
4496                  else
4497                      z = f_att->$vNam(varn+1)$               
4498                  end if
4499                else
4500                  do i=0,b-1           
4501                      if (isStrSubset( a(i),"zu_3d" ))then
4502                        z_v(varn,:) = z_u
4503                        if (log_z .EQ. 1) then
4504                            z = z_v(varn,1:dimz-1)
4505                        else
4506                            z = z_v(varn,:)
4507                        end if
4508                      else
4509                        if (isStrSubset( a(i),"zw_3d" ))then
4510                            z_v(varn,:) = z_w
4511                            if (log_z .EQ. 1) then
4512                              z = z_v(varn,1:dimz-1)
4513                            else
4514                              z = z_v(varn,:)
4515                            end if
4516                        end if                   
4517                      end if
4518                  end do           
4519                end if
4520
4521                z=z/norm_z
4522               
4523                com_var_avail(n_o)=vNam(varn)
4524               
4525                com=isStrSubset( c_var,","+vNam(varn)+"," )
4526                if (com)then
4527                  co = co+1           
4528                  if (n_o .EQ. 1) then
4529                      res@xyDashPattern  = 1                                   
4530                  else           
4531                      if (n_o .EQ. 2) then
4532                        res@xyDashPattern  = 2
4533                      else
4534                        res@xyDashPattern  = 0
4535                        res@gsnLeftString  = " " ; "Combined Plot of "+c_var_name
4536                        res@tiXAxisString  = "("+unit(varn)+")"
4537                        res@gsnRightString = " "
4538                        if (xs .EQ. -1) then
4539                            res@trXMinF = min(mini)
4540                        else
4541                            res@trXMinF = xs
4542                        end if
4543                        if (xe .EQ. -1) then
4544                            res@trXMaxF = max(maxi)
4545                        else
4546                            res@trXMaxF = xe 
4547                        end if
4548                      end if
4549                  end if
4550                  label(n_o)=vNam(varn)
4551                  color_o(n_o)=237
4552                  plot_o(n_o)=gsn_csm_xy(wks,data(varn,:,:),z,res)
4553                  n_o=n_o+1
4554                end if
4555                if (prof3d .EQ. 0)then
4556                  varn=varn+1
4557                end if           
4558            end if
4559          end do
4560     
4561          if(number_comb .EQ. 2)then
4562            if (co .EQ. 2)then
4563                overlay(plot_o(0),plot_o(1))
4564            else
4565                print(" ")
4566                print("combining is not possible,")
4567                print("'c_var'(= "+c_var+") must include two variables of "+\
4568                      "the general plots = ")
4569                print("- "+com_var_avail)
4570                print("be sure to have one comma before and after the variable")
4571                print(" ")
4572                exit 
4573            end if
4574          end if
4575          if(number_comb .EQ. 3)then
4576            if (co .EQ. 3)then
4577                overlay(plot_o(0),plot_o(1))
4578                overlay(plot_o(0),plot_o(2))
4579            else
4580                print(" ")
4581                print("combining is not possible,")
4582                print("'c_var'(= "+c_var+") must include three variables of "+\
4583                      "the general plots = ")
4584                print("- "+com_var_avail)
4585                print("be sure to have one comma before and after the variable")
4586                print(" ")
4587                exit
4588            end if
4589          end if
4590
4591          ; ***************************************************
4592          ; legend for combined plot
4593          ; ***************************************************
4594         
4595          lgres                    = True
4596          lgMonoDashIndex          = False
4597          lgres@lgDashIndexes      = (/0,1,2/)
4598          lgres@lgLabelFont        = "helvetica"   
4599          lgres@vpWidthF           = 0.12           
4600          lgres@vpHeightF          = 0.1 * number_comb / 3.0 + 0.01       
4601          lgres@lgLabelFontHeightF = 10.0 * font_size_legend * ((number_comb \
4602                                     - 2.0) * 0.35 + 0.65)
4603
4604          lbid = gsn_create_legend(wks,number_comb,label,lgres)       
4605
4606          amres = True
4607          amres@amParallelPosF   = -0.5             
4608          amres@amOrthogonalPosF = -0.5   
4609          amres@amJust           = "TopLeft"
4610
4611          annoid1 = gsn_add_annotation(plot_o(0),lbid,amres)
4612       
4613          plot(0) = plot_o(0)
4614
4615      end if
4616
4617      if (c_var_log .EQ. 1) then
4618          c_var_plot(com_i) = plot(0)
4619      end if
4620
4621   end do
4622
4623   ; ***************************************************
4624   ; merge plots onto one page
4625   ; ***************************************************
4626
4627   if (c_var_log .EQ. 1) then
4628      delete(plot_)
4629      plot_=new(max_com_i+1,graphic)
4630      plot_=c_var_plot
4631      n = max_com_i + 1
4632   else
4633      do m=0,n-1
4634         plot_(m)=plot(n-1-m)
4635      end do
4636   end if
4637
4638   no_frames = 0
4639
4640   if ((format_out .EQ. "eps" .OR. format_out .EQ. "epsi") .AND. \
4641       n .gt. no_rows*no_columns) then
4642      gsn_panel(wks,plot_,(/n,1/),resP)
4643      print(" ")
4644      print("Outputs to .eps or .epsi have only one frame")
4645      print(" ")
4646   else   
4647      do i = 0,n-1, no_rows*no_columns
4648         if( (i+no_rows*no_columns) .gt. (n-1)) then
4649            gsn_panel(wks,plot_(i:n-1),(/no_rows,no_columns/),resP)
4650            no_frames = no_frames + 1 
4651         else
4652            gsn_panel(wks,plot_(i:i+no_rows*no_columns-1),\
4653                                                (/no_rows,no_columns/),resP)
4654            no_frames = no_frames + 1 
4655         end if
4656      end do
4657   end if
4658
4659
4660   if (format_out .EQ. "png" ) then
4661     png_output = new((/no_frames/), string)
4662     j = 0
4663
4664     if (no_frames .eq. 1) then
4665        if (ncl_version .GE. 6 ) then
4666           png_output(0) = file_out+".png"
4667        else
4668           png_output(0) = file_out+".00000"+1+".png"
4669        end if
4670        ;using imagemagick's convert for reducing the white
4671        ;space around the plot
4672        cmd = "convert -geometry 1000x1000 -density 300 -trim " +  \
4673              png_output(0) + " " + png_output(0)
4674       system(cmd)
4675     else
4676
4677       do i=0, no_frames-1
4678         j = i + 1
4679         if (j .LE. 9) then
4680           png_output(i) = file_out+".00000"+j+".png"
4681         end if
4682         if (j .GT. 9 .AND. j .LE. 99) then
4683           png_output(i) = file_out+".0000"+j+".png"
4684         end if
4685         if (j .GT. 99 .AND. j .LE. 999) then
4686           png_output(i) = file_out+".000"+j+".png"
4687         end if
4688         if (j .GT. 999) then
4689           png_output(i) = file_out+".00"+j+".png"
4690         end if
4691
4692         ;using imagemagick's convert for reducing the white
4693         ;space around the plot
4694         cmd = "convert -geometry 1000x1000 -density 300 -trim " +  \
4695                png_output(i) + " " + png_output(i)
4696         system(cmd)
4697       end do
4698     end if
4699
4700     print(" ")
4701     print("Output to: "+ png_output)
4702     print(" ")
4703   else
4704     print(" ")
4705     print("Output to: " + file_out +"."+ format_out)
4706     print(" ")
4707   end if
4708
4709end
Note: See TracBrowser for help on using the repository browser.