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

Last change on this file since 515 was 513, checked in by heinze, 15 years ago

Using the NCL scripts by means of the shell script palmplot

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