source: palm/tags/release-3.6/SCRIPTS/NCL/profiles.ncl @ 2172

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