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

Last change on this file since 560 was 534, checked in by heinze, 15 years ago

Bugfix in spectra.ncl concerning output in png and small changes in layout

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