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

Last change on this file since 771 was 769, checked in by heinze, 13 years ago

Bugfixes in case of plot of t=0h and plot of topography zusi/zwwi possible

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