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

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