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

Last change on this file since 566 was 566, checked in by heinze, 14 years ago

Formatting of all NCL scripts.
Items of .ncl.config are re-sorted and xyc, xzc and yzc are no parameters any more.
Parameters start_f_1/end_f_1 are renamed to start_f/end_f in profiles.ncl.
Bugfix in cross_sections: enable vector plot if var is set explicitly.
Deletion of NCL user guide because guide is no available online.

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