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

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

Bugfix: enable plot of data if it is of kind double instead of kind float

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