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

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

Calculation of q*2 added

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