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

Last change on this file since 827 was 827, checked in by heinze, 12 years ago

allow plotting of data with very small time increments

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