source: palm/trunk/SCRIPTS/NCL/timeseries.ncl @ 771

Last change on this file since 771 was 769, checked in by heinze, 13 years ago

Bugfixes in case of plot of t=0h and plot of topography zusi/zwwi possible

File size: 50.4 KB
RevLine 
[194]1load "$NCARG_ROOT/lib/ncarg/nclex/gsun/gsn_code.ncl"
[154]2load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
3load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
4load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
5 
[190]6;***************************************************
[418]7; load .ncl.config or .ncl.config.default
[190]8;***************************************************
[194]9
10function which_script()
11local script
12begin
13   script="timeseries"
14   return(script)
15end
[418]16
17if (isfilepresent("$PALM_BIN/../../.ncl.config")) then
18   loadscript("$PALM_BIN/../../.ncl.config")
[190]19else
[418]20  if (isfilepresent("$PALM_BIN/NCL/.ncl.config.default")) then
21     loadscript( "$PALM_BIN/NCL/.ncl.config.default")
[190]22  else
[418]23      palm_bin_path = getenv("PALM_BIN")
[190]24      print(" ")
[418]25      print("Neither the personal configuration file '.ncl.config' exists in")
26      print("~/palm/current_version")
[566]27      print("nor the default configuration file '.ncl.config.default' "+\
28            "exists in")
[418]29      print(palm_bin_path + "/NCL")
[190]30      print(" ")
31      exit
[154]32   end if
[418]33end if   
[190]34   
35begin
[154]36
[190]37   ;***************************************************
[532]38   ; Retrieving the NCL version used
39   ;***************************************************
40   
41   ncl_version_ch = systemfunc("ncl -V")
42   ncl_version    = stringtofloat(ncl_version_ch)
43
44   ;***************************************************
[526]45   ; Retrieving the double quote character
46   ;***************************************************
47   
48   dq=str_get_dq()
49
50   ;***************************************************
[190]51   ; set up default parameter values and strings
52   ;***************************************************
[218]53   
[190]54   if (file_1 .EQ. "File in") then
55      print(" ")
[418]56      print("Declare input file 'file_1=' in '.ncl.config' or prompt")
[190]57      print(" ")
58      exit
[175]59   else
60      file_in = file_1   
[154]61   end if
[175]62
[566]63   if (format_out .NE. "x11" .AND. format_out .NE. "pdf" .AND.   \
64       format_out .NE. "eps" .AND. format_out .NE. "ps" .AND.    \
65       format_out .NE. "epsi" .AND. format_out .NE. "ncgm" .AND. \
66       format_out .NE. "png")then
[190]67      print(" ")
68      print("'format_out = "+format_out+"' is invalid and set to'x11'")
69      print(" ")
70      format_out="x11"
[154]71   end if 
[175]72
[532]73   if (ncl_version .LE. 5.1 .AND. format_out .EQ. "png") then
74      print(" ")
75      print("Output of png files not available")
76      print("png output is avaiable with NCL version 5.2.0 and higher ")
77      print("NCL version used: " + ncl_version_ch)
78      print(" ")
79      exit
80   end if
81
[190]82   if (over .NE. 0 .AND. over .NE. 1) then
83      print(" ")
84      print("'over'= "+over+" is invalid and set to 0")
85      print(" ")
[161]86      over = 0
[218]87   end if 
88   
89   if (norm_t .EQ. 0) then
90      print(" ")
91      print("Normalising with 0 is not allowed, 'norm_t' is set to 1.0")
92      print(" ")
93      norm_t = 1.0
94   end if
[190]95 
[175]96
[190]97   ;***************************************************
[154]98   ; open input file
[190]99   ;***************************************************
[218]100   
101   file_in_1 = False
102   if (isStrSubset(file_in, ".nc"))then
103      start_f = -2
104      end_f = -2
105      file_in_1 = True     
106   end if
[154]107
[218]108   if (start_f .EQ. -1)then
109      print(" ")
[566]110      print("'start_f' must be one of the cyclic numbers (at least 0) of "+\
111            "your input file(s)")
[218]112      print(" ") 
113      exit
114   end if
115   if (end_f .EQ. -1)then
116      print(" ")
[566]117      print("'end_f' must be one of the cyclic numbers (at least 0) of "+\
118            "your input file(s)")
[218]119      print(" ") 
120      exit
121   end if
[154]122
[218]123   files=new(end_f-start_f+1,string)
124
125   if (file_in_1)then
126      if (isfilepresent(file_in))then
127         files(0)=file_in
128      else
129         print(" ")
130         print("1st input file: '"+file_in+"' does not exist")
131         print(" ")
132         exit
133      end if
134   else   
135      if (start_f .EQ. 0)then
136         if (isfilepresent(file_in+".nc"))then
137            files(0)=file_in+".nc"
138            do i=1,end_f
139               if (isfilepresent(file_in+"."+i+".nc"))then   
140                  files(i)=file_in+"."+i+".nc"
141               else
142                  print(" ")
143                  print("Input file: '"+file_in+"."+i+".nc' does not exist")
144                  print(" ")
145                  exit 
146               end if       
147            end do         
148         else
149            print(" ")
150            print("Input file: '"+file_in+".nc' does not exist")
151            print(" ")
152            exit
153         end if
154      else
155         do i=start_f,end_f
156            if (isfilepresent(file_in+"."+i+".nc"))then   
157               files(i-start_f)=file_in+"."+i+".nc"
158            else
159               print(" ")
160               print("Input file: '"+file_in+"."+i+".nc' does not exist")
161               print(" ")
162               exit 
163            end if
164         end do
165      end if
166   end if
167
168   f=addfiles(files,"r")
169   f_att=addfile(files(0),"r")
170   ListSetType(f,"cat")
171
[585]172   vNam  =getfilevarnames(f_att)
173   vType = getfilevartypes(f_att,vNam)
174   vNam_static = vNam
[218]175
[585]176   if ((all(vType .eq. "double"))) then ;distinction if data is double or float
177      check_vType = True
178   else
179      check_vType = False
180   end if
181
[154]182   print(" ")
[190]183   print("Variables in input file:")
184   print("- "+ vNam)
[154]185   print(" ")
[161]186   dim = dimsizes(vNam)
187   if (dim .EQ. 0) then
188      print(" ")
[310]189      print("There is no data on file")
[161]190      print(" ")
191   end if
192
[310]193   t_all   = f[:]->time
194   nt      = dimsizes(t_all)
195   delta_t = t_all(nt-1)/nt
[218]196
[769]197   if (nt .LE. 1) then
198      print(" ")
199      print("Input file contains only one time step -> " +\
200            "plot of timeseries is not possible ")
201      print(" ")
202      exit
203   end if
204
[190]205   ;****************************************************       
[154]206   ; start of time step and different types of mistakes that could be done
[190]207   ;****************************************************
[154]208
[190]209   if (start_time_step .EQ. -1.) then           
210      start_time_step=t_all(0)/3600     
[154]211   else
[162]212      if (start_time_step .GE. t_all(nt-1)/3600)
[154]213         print(" ")
[566]214         print("'start_time_step' = "+ start_time_step +\
215               "h is equal or greater than last time step = " + \
216               t_all(nt-1)+"s = "+t_all(nt-1)/3600+"h")
[154]217         print(" ")
[218]218         print("Select another 'start_time_step'")
[162]219         print(" ")
[154]220         exit
221      end if
[162]222      if (start_time_step .LT. t_all(0)/3600)
[154]223         print(" ")
[566]224         print("'start_time_step' = "+ start_time_step +\
225               "h is lower than first time step = " + t_all(0)+"s = "+\
226               t_all(0)/3600+"h")
[154]227         print(" ")
[218]228         print("Select another 'start_time_step'")
[162]229         print(" ")
[154]230         exit
231      end if
232   end if
[162]233   do i=0,nt-2     
[566]234      if (start_time_step .GE. (t_all(i)-delta_t/2)/3600 .AND. \
235          start_time_step .LT. (t_all(i)+delta_t/2)/3600)then
[162]236         st=i
237         break
[769]238       else
239         st=0
[162]240      end if
241   end do
[566]242   if (start_time_step .GE. t_all(nt-1)-delta_t/2 .AND. \
243       start_time_step .LT. t_all(nt-1)) then
[162]244      st=nt-2   
[769]245    else
246      st=0
[162]247   end if
248     
[190]249   ;****************************************************
[154]250   ; end of time step and different types of mistakes that could be done
[190]251   ;****************************************************
[154]252
[190]253   if (end_time_step .EQ. -1.) then             
254      end_time_step = t_all(nt-1)/3600 
[154]255   else
[162]256      if (end_time_step .LE. t_all(0)/3600)
[154]257         print(" ")
[566]258         print("'end_time_step' = "+end_time_step+ \
259              "h is lower or equal than first time step = " + \
260               t_all(0)+"s = "+t_all(0)/3600+"h")
[154]261         print(" ")
[218]262         print("Select another 'end_time_step'")
[162]263         print(" ")
[154]264         exit
265      end if
[162]266      if (end_time_step .GT. t_all(nt-1)/3600)
[154]267         print(" ")
[566]268         print("'end_time_step' = "+ end_time_step +\
269               "h is greater than last time step = " + t_all(nt-1)+"s = "+\
270               t_all(nt-1)/3600+"h")
[154]271         print(" ")
[218]272         print("Select another 'end_time_step'") 
[162]273         print(" ")
[154]274         exit
275      end if
[162]276      if (end_time_step .LE. start_time_step/3600)
[161]277         print(" ")
[566]278         print("'end_time_step' = "+ end_time_step +\
279               "h is equal or lower than 'start_time_step' = "+\
280               start_time_step+"h")
[161]281         print(" ")
[218]282         print("Select another 'start_time_step' or 'end_time_step'")
[162]283         print(" ")
[161]284         exit
[154]285      end if
[162]286   end if
287   do i=0,nt-1     
[566]288      if (end_time_step .GE. (t_all(i)-delta_t/2)/3600 .AND. \
289          end_time_step .LT. (t_all(i)+delta_t/2)/3600)then
[162]290         et=i
291         break
[769]292       else
293         et=0
[162]294      end if
295   end do
296 
297   delete(start_time_step)
298   start_time_step=round(st,3)
299   delete(end_time_step)
300   end_time_step=round(et,3)
[154]301
[175]302   print(" ")
[566]303   print("Output of time steps from "+t_all(start_time_step)/3600+" h = "+\
304         t_all(start_time_step)+" s => index = "+start_time_step)
305   print("                     till "+t_all(end_time_step)/3600+" h = "+\
306         t_all(end_time_step)+" s => index = "+end_time_step)
[175]307   print(" ")
308
[218]309   t = t_all(start_time_step:end_time_step)/norm_t
[154]310 
311   ; ***************************************************
312   ; set up recourses
313   ; ***************************************************
314
315   res                         = True
316   res@gsnDraw                 = False
317   res@gsnFrame                = False
318   res@tmXBMode                = True
319   res@tmYLMode                = True
320   res@txFont                  = "helvetica"
321   res@tiMainFont              = "helvetica"
322   res@tiXAxisFont             = "helvetica"
323   res@tiYAxisFont             = "helvetica"
324   res@tmXBLabelFont           = "helvetica"
[161]325   res@tmYLLabelFont           = "helvetica"
326   res@xyLineColors            = (/237/)
[758]327   res@trXMaxF                 = t(start_time_step)
328   res@trXMinF                 = t(end_time_step)
[161]329   
[218]330   res@lgLabelFontHeightF     = 0.02
[161]331
[154]332   resP                        = True
333   resP@txFont                 = "helvetica"
[218]334   resP@txString               = f_att@title+" time series "
[154]335   resP@txFuncCode             = "~"
336   resP@txFontHeightF          = 0.015
337
[218]338   res@tmXBMinorPerMajor = 4
339   res@tmYLMinorPerMajor = 4
340
[154]341   res@vpWidthF=4
342
[161]343   txres = True
[154]344
345   ; ***************************************************
346   ; read data and create plots
347   ; ***************************************************
348   
[532]349   if ( format_out .EQ. "pdf" .OR. format_out .EQ. "ps" ) then
350      format_out@wkPaperSize = "A4"
351   end if
352   if ( format_out .EQ. "png" ) then
353      format_out@wkWidth  = 1000
354      format_out@wkHeight = 1000
355   end if
356
[154]357   wks_ps = gsn_open_wks(format_out,file_out)           
358   gsn_define_colormap(wks_ps,"rainbow+white")
[161]359   plot_ps=new(dim,graphic)
360                                 
[154]361   n=0
[585]362   if (check_vType) then
363      minE  =  100000.d
364      maxE  = -100000.d
365      minus =  100000.d
366      maxus = -100000.d
367      minu  =  100000.d
368      maxu  = -100000.d
369      minz  =  100000.d
370      maxz  = -100000.d
371      minw  =  100000.d
372      maxw  = -100000.d
373      minp  =  100000.d
374      maxp  = -100000.d
375      mins  =  100000.d
376      maxs  = -100000.d
377   else
378      minE  =  100000.
379      maxE  = -100000.
380      minus =  100000.
381      maxus = -100000.
382      minu  =  100000.
383      maxu  = -100000.
384      minz  =  100000.
385      maxz  = -100000.
386      minw  =  100000.
387      maxw  = -100000.
388      minp  =  100000.
389      maxp  = -100000.
390      mins  =  100000.
391      maxs  = -100000.
392   end if
[161]393
[585]394   if (check_vType) then
395      data   = new((/dim,(end_time_step-start_time_step)+1/),double)
396      data_0 = new((end_time_step-start_time_step)+1,double)
397      mini   = new(dim,double)
398      maxi   = new(dim,double)
399   else
400      data   = new((/dim,(end_time_step-start_time_step)+1/),float)
401      data_0 = new((end_time_step-start_time_step)+1,float)
402      mini   = new(dim,float)
403      maxi   = new(dim,float)
404   end if
405   unit   = new(dim,string)   
[161]406   data_0 = 0.0
[585]407 
[161]408   
409   if (over .EQ. 1) then
410      plot_E       = gsn_csm_xy(wks_ps,t,data_0(:),res)
411      plot_Es      = gsn_csm_xy(wks_ps,t,data_0(:),res)
412      plot_us      = gsn_csm_xy(wks_ps,t,data_0(:),res)
413      plot_ws      = gsn_csm_xy(wks_ps,t,data_0(:),res)
414      plot_umax    = gsn_csm_xy(wks_ps,t,data_0(:),res)
415      plot_vmax    = gsn_csm_xy(wks_ps,t,data_0(:),res)
416      plot_wmax    = gsn_csm_xy(wks_ps,t,data_0(:),res)
417      plot_z_i_wpt = gsn_csm_xy(wks_ps,t,data_0(:),res)
418      plot_z_i_pt  = gsn_csm_xy(wks_ps,t,data_0(:),res)
419      plot_wpptp0  = gsn_csm_xy(wks_ps,t,data_0(:),res)
420      plot_wpptp   = gsn_csm_xy(wks_ps,t,data_0(:),res)
421      plot_wpt     = gsn_csm_xy(wks_ps,t,data_0(:),res)
422      plot_pt_0_   = gsn_csm_xy(wks_ps,t,data_0(:),res)
423      plot_pt_zp_  = gsn_csm_xy(wks_ps,t,data_0(:),res)
424   end if
[162]425
[161]426
[310]427   if (var .NE. "all") then
428      vNam_temp = new((/dim/),string)
429      comma     = 0
430      var_char  = stringtocharacter(var)
431      no_char   = dimsizes(var_char)
432       
433      do j = 0, no_char-1
434         if(var_char(j) .eq. ",")
435            comma = comma + 1
436         end if
437      end do
438
[317]439      if(comma .eq. 0 .or. comma .eq. 1)
440          print(" ")
441          print("The variables 'var="+var+"'" )
442          print("do not exist on your input file;")
443          print("be sure to have one comma before and after each variable")
444          print(" ")
445          exit
446      end if
447
[310]448      indices = new((/comma/), integer)
449      comma   = 0
450       
451      do j = 0, no_char-1
452         if(var_char(j) .eq. ",")
453           indices(comma) = j
454           comma = comma + 1
455         end if
456      end do
457
458      do j = 0, comma -2
[566]459         vNam_temp(j) = charactertostring(\
460                                   var_char(indices(j)+1:indices(j+1)-1))
[310]461      end do
462     
463      delete(vNam)
464      vNam = new((/comma-1/),string)
465      do varn = 0, comma - 2
466         vNam(varn) = vNam_temp(comma - 2 - varn)
467      end do
468     
469      do j = 0, comma-2
470        count_check = 0
471        do varn = 0, dim-1
472           if(vNam(j) .ne. vNam_static(varn))
473             count_check=count_check+1
474           end if
[317]475        end do   
476
[310]477        if(count_check .eq. dim)
478            print(" ")
479            print("The variables 'var="+var+"'" )
480            print("do not exist on your input file;")
481            print("be sure to have one comma before and after each variable")
482            print(" ")
483            exit
484        end if
485      end do
486
487      dim = comma-1
488    end if
489
490
491   do varn = dim-1, 0, 1
492
[161]493      if( isStrSubset (vNam(varn), "time") )
494         check = False
495      else
496         check = True
[310]497      end if
[161]498     
499      if(check) then
[162]500
[218]501         data_all = f[:]->$vNam(varn)$
502         data_att = f_att->$vNam(varn)$
503         unit(varn) = data_att@units
[310]504 
505     
[161]506         data(varn,:)=data_all(start_time_step:end_time_step)
507         
508         if (over .EQ. 1) then
509
510            mini(varn) = min(data(varn,:))
511            maxi(varn) = max(data(varn,:))
512           
[566]513            if (vNam(varn) .EQ. "E" .OR. vNam(varn) .EQ. "Es" .OR. \
514                vNam(varn) .EQ. "E*") then
[161]515               if (mini(varn) .EQ. maxi(varn)) then
516                  if (min(data(varn,:)) .EQ. 0)then
[218]517                     mini(varn)= mini(varn)-0.1
518                     maxi(varn)= maxi(varn)+0.1
[161]519                  end if
520                  if (min(data(varn,:)) .LT. 0)then
521                     mini(varn)= mini(varn)-1.+(mini(varn))/2
522                     maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
523                  end if
524                  if (min(data(varn,:)) .GT. 0)then
525                     mini(varn)= mini(varn)-1.-(mini(varn))/2
526                     maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
527                  end if
528               end if
529               minE=min((/minE,mini(varn)/))
530               maxE=max((/maxE,maxi(varn)/))
531            end if
532
[526]533            if (vNam(varn) .EQ. "us" .OR. vNam(varn) .EQ. "ws"\
534                .OR. vNam(varn) .EQ. "u*" .OR. vNam(varn) .EQ. "w*") then
[161]535               if (mini(varn) .EQ. maxi(varn)) then
536                  if (min(data(varn,:)) .EQ. 0)then
[218]537                     mini(varn)= mini(varn)-0.1
538                     maxi(varn)= maxi(varn)+0.1
[161]539                  end if
540                  if (min(data(varn,:)) .LT. 0)then
541                     mini(varn)= mini(varn)-1.+(mini(varn))/2
542                     maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
543                  end if
544                  if (min(data(varn,:)) .GT. 0)then
545                     mini(varn)= mini(varn)-1.-(mini(varn))/2
546                     maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
547                  end if
548               end if
549               minus=min((/minus,mini(varn)/))
550               maxus=max((/maxus,maxi(varn)/))
551            end if
552
[566]553            if (vNam(varn) .EQ. "umax" .OR. vNam(varn) .EQ. "vmax" .OR. \
554                vNam(varn) .EQ. "wmax") then
[161]555               if (mini(varn) .EQ. maxi(varn)) then
556                  if (mini(varn) .EQ. 0)then
[218]557                     mini(varn)= mini(varn)-0.1
558                     maxi(varn)= maxi(varn)+0.1
[161]559                  end if
560                  if (mini(varn) .LT. 0)then
561                     mini(varn)= mini(varn)-1.+(mini(varn))/2
562                     maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
563                  end if
564                  if (mini(varn) .GT. 0)then
565                     mini(varn)= mini(varn)-1.-(mini(varn))/2
566                     maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
567                  end if
568               end if
569               minu=min((/minu,mini(varn)/))
570               maxu=max((/maxu,maxi(varn)/))
571            end if
572
573            if (vNam(varn) .EQ. "z_i_wpt" .OR. vNam(varn) .EQ. "z_i_pt") then
574               if (mini(varn) .EQ. maxi(varn)) then
575                  if (min(data(varn,:)) .EQ. 0)then
[218]576                     mini(varn)= mini(varn)-0.1
577                     maxi(varn)= maxi(varn)+0.1
[161]578                  end if
579                  if (min(data(varn,:)) .LT. 0)then
580                     mini(varn)= mini(varn)-1.+(mini(varn))/2
581                     maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
582                  end if
583                  if (min(data(varn,:)) .GT. 0)then
584                     mini(varn)= mini(varn)-1.-(mini(varn))/2
585                     maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
586                  end if
587               end if
588               minz=min((/minz,mini(varn)/))
589               maxz=max((/maxz,maxi(varn)/))
590            end if
591
[566]592            if (vNam(varn) .EQ. "wpptp0" .OR. vNam(varn) .EQ. "wpptp" .OR. \
593                vNam(varn) .EQ. "wpt" .OR. vNam(varn) .EQ. \
594                "w"+dq+"pt"+dq+"0" .OR. vNam(varn) .EQ. "w"+dq+"pt"+dq ) then
[161]595               if (mini(varn) .EQ. maxi(varn)) then
596                  if (min(data(varn,:)) .EQ. 0)then
[218]597                     mini(varn)= mini(varn)-0.1
598                     maxi(varn)= maxi(varn)+0.1
[161]599                  end if
600                  if (min(data(varn,:)) .LT. 0)then
601                     mini(varn)= mini(varn)-1.+(mini(varn))/2
602                     maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
603                  end if
604                  if (min(data(varn,:)) .GT. 0)then
605                     mini(varn)= mini(varn)-1.-(mini(varn))/2
606                     maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
607                  end if
608               end if
609               minw=min((/minw,mini(varn)/))
610               maxw=max((/maxw,maxi(varn)/))
611            end if
612
[526]613            if (vNam(varn) .EQ. "pt_0_" .OR. vNam(varn) .EQ. "pt_zp_"\
614               .OR. vNam(varn) .EQ. "pt(0)" .OR. vNam(varn) .EQ. "pt(zp)") then
[161]615               if (mini(varn) .EQ. maxi(varn)) then
616                  if (min(data(varn,:)) .EQ. 0)then
[218]617                     mini(varn)= mini(varn)-0.1
618                     maxi(varn)= maxi(varn)+0.1
[161]619                  end if
620                  if (min(data(varn,:)) .LT. 0)then
621                     mini(varn)= mini(varn)-1.+(mini(varn))/2
622                     maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
623                  end if
624                  if (min(data(varn,:)) .GT. 0)then
625                     mini(varn)= mini(varn)-1.-(mini(varn))/2
626                     maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
627                  end if
628               end if
629               minp=min((/minp,mini(varn)/))
630               maxp=max((/maxp,maxi(varn)/))
631            end if
632           
633         end if
634      end if
635   end do
[310]636 
637  if (isStrSubset(data@long_name," SR " ))then
[162]638      print(" ")
[566]639      print("If you have outputs of statistic regions you cannot overlay "+\
640            "variables;")
[310]641      print("'over' is set to 0" )
[162]642      print(" ")
[310]643      over = 0
[175]644   end if
[162]645
[310]646   do varn = dim-1, 0, 1
[161]647     
[154]648      if( isStrSubset (vNam(varn), "time") )
649         check = False
[161]650      else
651         check = True
652      end if     
[190]653   
654      if (var .NE. "all") then
[161]655         check = isStrSubset( var,","+vNam(varn)+"," )
656      end if
[154]657 
658      if(check) then
[310]659       
[161]660        if (over .EQ. 1) then
[154]661 
[161]662            res@gsnLeftString       = "overlayed plot"
663            res@gsnRightString      = unit(varn)
[218]664            if (norm_t .NE. 1.)then
[534]665               res@tiXAxisString        = "t ("+unit_t+")"
[218]666            else
[534]667               res@tiXAxisString        = "t (s)"
[218]668            end if
669            res@tiYAxisString = " "
670            res@tiXAxisFontHeightF   = font_size
671            res@txFontHeightF        = font_size
672            res@tiYAxisFontHeightF   = font_size           
[154]673
[161]674            if (vNam(varn) .EQ. "E")
675               E=0
676               res@xyLineColors     = (/237/)
677               res@xyLineLabelFontHeightF = 0.05
678               res@xyLineLabelFontColor   = 237
[310]679               res@trYMaxF           = maxE
680               res@trYMinF           = minE
[161]681               plot_E = gsn_csm_xy(wks_ps,t,data(varn,:),res)
682            end if
[526]683            if (vNam(varn) .EQ. "Es" .OR. vNam(varn) .EQ. "E*")
[161]684               Es=0
685               res@xyLineColors            = (/144/)
686               res@xyLineLabelFontHeightF = 0.05
687               res@xyLineLabelFontColor   = 144
688               plot_Es = gsn_csm_xy(wks_ps,t,data(varn,:),res)
689            end if
690             
[526]691            if (vNam(varn) .EQ. "us" .OR. vNam(varn) .EQ. "u*")
[161]692               us=0
693               res@xyLineColors            = (/237/)
694               res@xyLineLabelFontHeightF = 0.05
695               res@xyLineLabelFontColor   = 237
[310]696               res@trYMaxF           = maxus
697               res@trYMinF           = minus
[161]698               plot_us = gsn_csm_xy(wks_ps,t,data(varn,:),res)
699            end if             
[526]700            if (vNam(varn) .EQ. "ws" .OR. vNam(varn) .EQ. "w*")
[161]701               ws=0
702               res@xyLineColors            = (/144/)
703               res@xyLineLabelFontHeightF = 0.05
704               res@xyLineLabelFontColor   = 144
705               plot_ws = gsn_csm_xy(wks_ps,t,data(varn,:),res)
706            end if
707             
708            if (vNam(varn) .EQ. "umax")
709               u=0
710               res@xyLineColors            = (/237/)
711               res@xyLineLabelFontHeightF = 0.05
712               res@xyLineLabelFontColor   = 237
[310]713               res@trYMaxF           = maxu
714               res@trYMinF           = minu
[161]715               plot_umax = gsn_csm_xy(wks_ps,t,data(varn,:),res)
716            end if
717            if (vNam(varn) .EQ. "vmax")
718               v=0
719               res@xyLineColors            = (/144/)
720               res@xyLineLabelFontHeightF = 0.05
721               res@xyLineLabelFontColor   = 144
722               plot_vmax = gsn_csm_xy(wks_ps,t,data(varn,:),res)
723            end if
724            if (vNam(varn) .EQ. "wmax")
725               w=0
726               res@xyLineColors            = (/80/)
727               res@xyLineLabelFontHeightF = 0.05
728               res@xyLineLabelFontColor   = 80
729               plot_wmax = gsn_csm_xy(wks_ps,t,data(varn,:),res)
730            end if
[154]731   
[161]732            if (vNam(varn) .EQ. "z_i_wpt")
733               zw=0
734               res@xyLineColors            = (/237/)
735               res@xyLineLabelFontHeightF = 0.05
736               res@xyLineLabelFontColor   = 237
[310]737               res@trYMaxF           = maxz
738               res@trYMinF           = minz
[161]739               plot_z_i_wpt = gsn_csm_xy(wks_ps,t,data(varn,:),res)
740            end if
741            if (vNam(varn) .EQ. "z_i_pt")
742               z=0
743               res@xyLineColors            = (/144/) 
744               res@xyLineLabelFontHeightF = 0.05
745               res@xyLineLabelFontColor   = 144
746               plot_z_i_pt = gsn_csm_xy(wks_ps,t,data(varn,:),res)
747            end if
748
[566]749            if (vNam(varn) .EQ. "wpptp0" .OR. \
750                vNam(varn) .EQ. "w"+dq+"pt"+dq+"0" )
[161]751               w0=0
752               res@xyLineColors            = (/237/)
753               res@xyLineLabelFontHeightF = 0.05
754               res@xyLineLabelFontColor   = 237
[310]755               res@trYMaxF           = maxw
756               res@trYMinF           = minw
[161]757               plot_wpptp0 = gsn_csm_xy(wks_ps,t,data(varn,:),res)
758            end if
[526]759            if (vNam(varn) .EQ. "wpptp" .OR. vNam(varn) .EQ. "w"+dq+"pt"+dq)
[161]760               wp=0
761               res@xyLineColors            = (/144/)
762               res@xyLineLabelFontHeightF = 0.05
763               res@xyLineLabelFontColor   = 144
764               plot_wpptp = gsn_csm_xy(wks_ps,t,data(varn,:),res) 
765            end if
766            if (vNam(varn) .EQ. "wpt")
767               wt=0
768               res@xyLineColors            = (/80/)
769               res@xyLineLabelFontHeightF = 0.05
770               res@xyLineLabelFontColor   = 80
771               plot_wpt = gsn_csm_xy(wks_ps,t,data(varn,:),res)
772            end if
773
[526]774            if (vNam(varn) .EQ. "pt_0_" .OR. vNam(varn) .EQ. "pt(0)")
[161]775               p=0
776               res@xyLineColors            = (/237/)
777               res@xyLineLabelFontHeightF = 0.05
778               res@xyLineLabelFontColor   = 237
[310]779               res@trYMaxF           = maxp
780               res@trYMinF           = minp
[161]781               plot_pt_0_ = gsn_csm_xy(wks_ps,t,data(varn,:),res) 
782            end if
[526]783            if (vNam(varn) .EQ. "pt_zp_" .OR. vNam(varn) .EQ. "pt(zp)")
[161]784               pz=0
785               res@xyLineColors            = (/144/)
786               res@xyLineLabelFontHeightF = 0.05
787               res@xyLineLabelFontColor   = 144
788               plot_pt_zp_ = gsn_csm_xy(wks_ps,t,data(varn,:),res)
789            end if
790
791         end if
[154]792      end if
[161]793   end do
[310]794
795
796   do varn = dim-1, 0, 1
[161]797   
798      if( isStrSubset (vNam(varn), "time") )
799         check = False
800      else
801         check = True
802      end if   
[190]803 
804      if (var.NE. "all") then
[161]805         check = isStrSubset( var,","+vNam(varn)+"," )
806      end if
807   
808      if(check) then
809       
810         if (over .EQ. 1) then       
811       
812            if (vNam(varn) .EQ. "E" .AND. Es .NE. 1) then
813               E=1   
814               overlay(plot_E,plot_Es)
815               n=n+1
[162]816               plot_ps(n) = plot_E   
817
818               ; ***************************************************
[566]819               ; legend for overlaid plot
[162]820               ; ***************************************************
821     
822               lgres                    = True
823               lgMonoDashIndex          = False
824               lgres@lgLabelFont        = "helvetica"   
825               lgres@lgLabelFontHeightF = .1           
826               lgres@vpWidthF           = 0.4           
827               lgres@vpHeightF          = 0.4         
828               lgres@lgDashIndexes      = (/0,0,0/)
829               lgres@lgLineColors       = (/237,144,80/)
[526]830               lbid = gsn_create_legend(wks_ps,2,(/"E","E*"/),lgres)       
[162]831
832               amres = True
833               amres@amParallelPosF   = 0.6                 
834               amres@amOrthogonalPosF = -0.2           
835               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)           
[161]836            end if
[566]837            if ((vNam(varn) .EQ. "Es" .OR. vNam(varn) .EQ. "E*") .AND. \
838                E .NE. 1) then
[161]839               Es=1
840               overlay(plot_E,plot_Es)
841               n=n+1
[162]842               plot_ps(n) = plot_E 
843
844               ; ***************************************************
[566]845               ; legend for overlaid plot
[162]846               ; ***************************************************
847     
848               lgres                    = True
849               lgMonoDashIndex          = False
850               lgres@lgLabelFont        = "helvetica"   
851               lgres@lgLabelFontHeightF = .1           
852               lgres@vpWidthF           = 0.4           
853               lgres@vpHeightF          = 0.4         
854               lgres@lgDashIndexes      = (/0,0,0/)
855               lgres@lgLineColors       = (/237,144,80/)
[526]856               lbid = gsn_create_legend(wks_ps,2,(/"E","E*"/),lgres)       
[162]857
858               amres = True
859               amres@amParallelPosF   = 0.6                 
860               amres@amOrthogonalPosF = -0.2           
861               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)         
[161]862            end if
863
[566]864            if ((vNam(varn) .EQ. "us" .OR. vNam(varn) .EQ. "u*") .AND. \
865                 ws .NE. 1) then
[161]866               us=1
867               overlay(plot_us,plot_ws)
868               n=n+1
869               plot_ps(n) = plot_us
[162]870
871               ; ***************************************************
[566]872               ; legend for overlaid plot
[162]873               ; ***************************************************
874     
875               lgres                    = True
876               lgMonoDashIndex          = False
877               lgres@lgLabelFont        = "helvetica"   
878               lgres@lgLabelFontHeightF = .1           
879               lgres@vpWidthF           = 0.4           
880               lgres@vpHeightF          = 0.4         
881               lgres@lgDashIndexes      = (/0,0,0/)
882               lgres@lgLineColors       = (/237,144,80/)
[526]883               lbid = gsn_create_legend(wks_ps,2,(/"u*","w*"/),lgres)       
[162]884
885               amres = True
886               amres@amParallelPosF   = 0.6                 
887               amres@amOrthogonalPosF = -0.2           
888               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]889            end if
[566]890            if ((vNam(varn) .EQ. "ws" .OR. vNam(varn) .EQ. "w*") .AND. \
891                 us .NE. 1) then
[161]892               ws=1
893               overlay(plot_us,plot_ws)
894               n=n+1
895               plot_ps(n) = plot_us
[162]896
897               ; ***************************************************
[566]898               ; legend for overlaid  plot
[162]899               ; ***************************************************
900     
901               lgres                    = True
902               lgMonoDashIndex          = False
903               lgres@lgLabelFont        = "helvetica"   
904               lgres@lgLabelFontHeightF = .1           
905               lgres@vpWidthF           = 0.4           
906               lgres@vpHeightF          = 0.4         
907               lgres@lgDashIndexes      = (/0,0,0/)
908               lgres@lgLineColors       = (/237,144,80/)
[526]909               lbid = gsn_create_legend(wks_ps,2,(/"u*","w*"/),lgres)       
[162]910
911               amres = True
912               amres@amParallelPosF   = 0.6                 
913               amres@amOrthogonalPosF = -0.2           
914               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]915            end if
916         
917            if (vNam(varn) .EQ. "umax" .AND. v .NE. 1)
918               if (w .NE. 1) then
919                  u=1         
920                  overlay(plot_umax,plot_vmax)
921                  overlay(plot_umax,plot_wmax)
[310]922               
[161]923                  n=n+1
924                  plot_ps(n) = plot_umax
[162]925
926                  ; ***************************************************
[566]927                  ; legend for overlaid plot
[162]928                  ; ***************************************************
929     
930                  lgres                    = True
931                  lgMonoDashIndex          = False
932                  lgres@lgLabelFont        = "helvetica"   
933                  lgres@lgLabelFontHeightF = .1           
934                  lgres@vpWidthF           = 0.4           
935                  lgres@vpHeightF          = 0.4         
936                  lgres@lgDashIndexes      = (/0,0,0/)
937                  lgres@lgLineColors       = (/237,144,80/)
[566]938                  lbid = gsn_create_legend(\
939                                  wks_ps,3,(/"umax","vmax","wmax"/),lgres)
[162]940
941                  amres = True
942                  amres@amParallelPosF   = 0.6             
943                  amres@amOrthogonalPosF = -0.2           
944                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]945               end if
946            end if
947            if (vNam(varn) .EQ. "vmax" .AND. u .NE. 1)
948               if (w .NE. 1) then
949                  v=1 
950                  overlay(plot_umax,plot_vmax)
951                  overlay(plot_umax,plot_wmax)
952                  n=n+1
953                  plot_ps(n) = plot_umax
[162]954
955                  ; ***************************************************
[566]956                  ; legend for overlaid plot
[162]957                  ; ***************************************************
958     
959                  lgres                    = True
960                  lgMonoDashIndex          = False
961                  lgres@lgLabelFont        = "helvetica"   
962                  lgres@lgLabelFontHeightF = .1         
963                  lgres@vpWidthF           = 0.4           
964                  lgres@vpHeightF          = 0.4         
965                  lgres@lgDashIndexes      = (/0,0,0/)
966                  lgres@lgLineColors       = (/237,144,80/)
[566]967                  lbid = gsn_create_legend(\
968                                      wks_ps,3,(/"umax","vmax","wmax"/),lgres)
[162]969                  amres = True
970                  amres@amParallelPosF   = 0.6             
971                  amres@amOrthogonalPosF = -0.2           
972                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]973               end if
974            end if
975            if (vNam(varn) .EQ. "wmax" .AND. v .NE. 1)
976               if(u .NE. 1) then 
977                  w=1       
978                  overlay(plot_umax,plot_vmax)
979                  overlay(plot_umax,plot_wmax)
980                  n=n+1
981                  plot_ps(n) = plot_umax
[162]982
983                  ; ***************************************************
[566]984                  ; legend for overlaid plot
[162]985                  ; ***************************************************
986     
987                  lgres                    = True
988                  lgMonoDashIndex          = False
989                  lgres@lgLabelFont        = "helvetica"   
990                  lgres@lgLabelFontHeightF = .1           
991                  lgres@vpWidthF           = 0.4           
992                  lgres@vpHeightF          = 0.4         
993                  lgres@lgDashIndexes      = (/0,0,0/)
994                  lgres@lgLineColors       = (/237,144,80/)
[566]995                  lbid = gsn_create_legend(\
996                                       wks_ps,3,(/"umax","vmax","wmax"/),lgres)
[162]997                  amres = True
998                  amres@amParallelPosF   = 0.6             
999                  amres@amOrthogonalPosF = -0.2           
1000                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]1001               end if
1002            end if
1003           
1004            if (vNam(varn) .EQ. "z_i_wpt" .AND. z .NE. 1) then
1005               zw=1       
1006               overlay(plot_z_i_wpt,plot_z_i_pt)
1007               n=n+1
[162]1008               plot_ps(n) = plot_z_i_wpt
1009       
1010               ; ***************************************************
[566]1011               ; legend for overlaid plot
[162]1012               ; ***************************************************
1013     
1014               lgres                    = True
1015               lgMonoDashIndex          = False
1016               lgres@lgLabelFont        = "helvetica"   
1017               lgres@lgLabelFontHeightF = .1           
1018               lgres@vpWidthF           = 0.4           
1019               lgres@vpHeightF          = 0.4         
1020               lgres@lgDashIndexes      = (/0,0,0/)
1021               lgres@lgLineColors       = (/237,144,80/)
[566]1022               lbid = gsn_create_legend(wks_ps,2,(/"z_i_wpt","z_i_pt"/),lgres)
[162]1023
1024               amres = True
1025               amres@amParallelPosF   = 0.6                 
1026               amres@amOrthogonalPosF = -0.2           
1027               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)         
[161]1028            end if
1029            if (vNam(varn) .EQ. "z_i_pt" .AND. zw .NE. 1) then
1030               z=1     
1031               overlay(plot_z_i_wpt,plot_z_i_pt)
1032               n=n+1
[162]1033               plot_ps(n) = plot_z_i_wpt
1034
1035               ; ***************************************************
[566]1036               ; legend for overlaid plot
[162]1037               ; ***************************************************
1038     
1039               lgres                    = True
1040               lgMonoDashIndex          = False
1041               lgres@lgLabelFont        = "helvetica"   
1042               lgres@lgLabelFontHeightF = .1           
1043               lgres@vpWidthF           = 0.4           
1044               lgres@vpHeightF          = 0.4         
1045               lgres@lgDashIndexes      = (/0,0,0/)
1046               lgres@lgLineColors       = (/237,144,80/)
[566]1047               lbid = gsn_create_legend(wks_ps,2,(/"z_i_wpt","z_i_pt"/),lgres)
[162]1048
1049               amres = True
1050               amres@amParallelPosF   = 0.6                 
1051               amres@amOrthogonalPosF = -0.2           
1052               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)           
[161]1053            end if   
1054         
[526]1055            if ((vNam(varn) .EQ. "wpptp0" .OR. vNam(varn) .EQ. "w"+dq+"pt"+dq+"0") .AND. wp .NE. 1)
[161]1056               if (wt .NE. 1) then
1057                  w0=1
1058                  overlay(plot_wpptp0,plot_wpptp)
1059                  overlay(plot_wpptp0,plot_wpt)
1060                  n=n+1
1061                  plot_ps(n) = plot_wpptp0
[162]1062
1063                  ; ***************************************************
[566]1064                  ; legend for overlaid plot
[162]1065                  ; ***************************************************
1066     
1067                  lgres                    = True
1068                  lgMonoDashIndex          = False
1069                  lgres@lgLabelFont        = "helvetica"   
1070                  lgres@lgLabelFontHeightF = .1           
1071                  lgres@vpWidthF           = 0.4           
1072                  lgres@vpHeightF          = 0.4         
1073                  lgres@lgDashIndexes      = (/0,0,0/)
1074                  lgres@lgLineColors       = (/237,144,80/)
[566]1075                  lbid = gsn_create_legend(\
1076                                    wks_ps,3,(/"wpptp0","wpptp","wpt"/),lgres)
[162]1077                  amres = True
1078                  amres@amParallelPosF   = 0.6             
1079                  amres@amOrthogonalPosF = -0.2           
1080                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]1081               end if           
1082            end if
[566]1083            if ((vNam(varn) .EQ. "wpptp" .OR. \
1084                 vNam(varn) .EQ. "w"+dq+"pt"+dq) .AND. w0 .NE. 1)
[161]1085               if (wt .NE. 1) then
1086                 wp=1
1087                 overlay(plot_wpptp0,plot_wpptp)
1088                 overlay(plot_wpptp0,plot_wpt)
1089                 n=n+1
[162]1090                 plot_ps(n) = plot_wpptp0
1091                 
1092                 ; ***************************************************
[566]1093                  ; legend for overlaid plot
[162]1094                  ; ***************************************************
1095     
1096                  lgres                    = True
1097                  lgMonoDashIndex          = False
1098                  lgres@lgLabelFont        = "helvetica"   
1099                  lgres@lgLabelFontHeightF = .1           
1100                  lgres@vpWidthF           = 0.4           
1101                  lgres@vpHeightF          = 0.4         
1102                  lgres@lgDashIndexes      = (/0,0,0/)
1103                  lgres@lgLineColors       = (/237,144,80/)
[566]1104                  lbid = gsn_create_legend(\
1105                                     wks_ps,3,(/"wpptp0","wpptp","wpt"/),lgres)
[162]1106
1107                  amres = True
1108                  amres@amParallelPosF   = 0.6             
1109                  amres@amOrthogonalPosF = -0.2           
1110                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres) 
[161]1111               end if           
1112            end if
1113            if (vNam(varn) .EQ. "wpt" .AND. wp .NE. 1)
1114               if (w0 .NE. 1) then
1115                  wt=1
1116                  overlay(plot_wpptp0,plot_wpptp)
1117                  overlay(plot_wpptp0,plot_wpt)
1118                  n=n+1
[162]1119                  plot_ps(n) = plot_wpptp0
1120
1121                  ; ***************************************************
[566]1122                  ; legend for overlaid plot
[162]1123                  ; ***************************************************
1124     
1125                  lgres                    = True
1126                  lgMonoDashIndex          = False
1127                  lgres@lgLabelFont        = "helvetica"   
1128                  lgres@lgLabelFontHeightF = .1           
1129                  lgres@vpWidthF           = 0.4           
1130                  lgres@vpHeightF          = 0.4         
1131                  lgres@lgDashIndexes      = (/0,0,0/)
1132                  lgres@lgLineColors       = (/237,144,80/)
[566]1133                  lbid = gsn_create_legend(\
1134                                  wks_ps,3,(/"wpptp0","wpptp","wpt"/),lgres)
[162]1135                  amres = True
1136                  amres@amParallelPosF   = 0.6             
1137                  amres@amOrthogonalPosF = -0.2           
1138                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]1139               end if
1140            end if
1141
[566]1142            if ((vNam(varn) .EQ. "pt_0_" .OR. vNam(varn) .EQ. "pt(0)") .AND. \
1143                pz .NE. 1) then
[161]1144               p=1     
1145               overlay(plot_pt_0_,plot_pt_zp_)
1146               n=n+1
[162]1147               plot_ps(n) = plot_pt_0_ 
1148     
1149               ; ***************************************************
[566]1150               ; legend for overlaid plot
[162]1151               ; ***************************************************
1152     
1153               lgres                    = True
1154               lgMonoDashIndex          = False
1155               lgres@lgLabelFont        = "helvetica"   
1156               lgres@lgLabelFontHeightF = .1           
1157               lgres@vpWidthF           = 0.4           
1158               lgres@vpHeightF          = 0.4         
1159               lgres@lgDashIndexes      = (/0,0,0/)
1160               lgres@lgLineColors       = (/237,144,80/)
[566]1161               lbid = gsn_create_legend(wks_ps,2,(/"pt_0_","pt_zp_"/),lgres)
[162]1162
1163               amres = True
1164               amres@amParallelPosF   = 0.6                 
1165               amres@amOrthogonalPosF = -0.2           
1166               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]1167            end if
[566]1168            if ((vNam(varn) .EQ. "pt_zp_" .OR. vNam(varn) .EQ. "pt(zp)") .AND.\
1169                 p .NE. 1) then
[161]1170               pz=1       
1171               overlay(plot_pt_0_,plot_pt_zp_)
1172               n=n+1
[162]1173               plot_ps(n) = plot_pt_0_   
1174
1175               ; ***************************************************
[566]1176               ; legend for overlaid plot
[162]1177               ; ***************************************************
1178     
1179               lgres                    = True
1180               lgMonoDashIndex          = False
1181               lgres@lgLabelFont        = "helvetica"   
1182               lgres@lgLabelFontHeightF = .1           
1183               lgres@vpWidthF           = 0.4           
1184               lgres@vpHeightF          = 0.4         
1185               lgres@lgDashIndexes      = (/0,0,0/)
1186               lgres@lgLineColors       = (/237,144,80/)
[566]1187               lbid = gsn_create_legend(wks_ps,2,(/"pt_0_","pt_zp_"/),lgres)
[162]1188
1189               amres = True
1190               amres@amParallelPosF   = 0.6                 
1191               amres@amOrthogonalPosF = -0.2           
1192               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]1193            end if
[162]1194
[526]1195            if(vNam(varn) .NE. "pt_zp_" .AND. vNam(varn) .NE. "pt(zp)" .AND. \
[566]1196               vNam(varn) .NE. "pt_0_" .AND. vNam(varn) .NE. "pt(0)" .AND.   \
1197               vNam(varn) .NE. "wpt" .AND. vNam(varn) .NE. "wpptp" .AND.     \
1198               vNam(varn) .NE. "w"+dq+"pt"+dq .AND.                          \
1199               vNam(varn) .NE. "wpptp0" .AND.                                \
1200               vNam(varn) .NE. "w"+dq+"pt"+dq+"0" .AND.                      \
1201               vNam(varn) .NE. "z_i_pt" .AND.                                \
1202               vNam(varn) .NE. "z_i_wpt" .AND. vNam(varn) .NE. "wmax" .AND.  \
1203               vNam(varn) .NE. "vmax" .AND. vNam(varn) .NE. "umax" .AND.     \
1204               vNam(varn) .NE. "ws" .AND. vNam(varn) .NE. "w*" .AND.         \
1205               vNam(varn) .NE. "us" .AND. vNam(varn) .NE. "u*" .AND.         \
1206               vNam(varn) .NE. "Es" .AND. vNam(varn) .NE. "E*" .AND.         \
[526]1207               vNam(varn) .NE. "E") then
[162]1208
[161]1209               n=n+1
1210               res@xyLineColors   = (/237/)
1211               res@xyLabelMode    = False
1212               res@gsnLeftString  = vNam(varn)
1213               res@gsnRightString = unit(varn)
[310]1214               res@trYMaxF        = max(data(varn,:))
1215               res@trYMinF        = min(data(varn,:))
[161]1216               if (min(data(varn,:)) .EQ. max(data(varn,:))) then
1217                  if (min(data(varn,:)) .EQ. 0)then
[310]1218                     res@trYMaxF = max(data(varn,:))-0.1
1219                     res@trYMinF = min(data(varn,:))+0.1
[161]1220                  end if
1221                  if (min(data(varn,:)) .LT. 0)then
[310]1222                     res@trYMaxF = max(data(varn,:))+(min(data(varn,:)))/2
1223                     res@trYMinF = min(data(varn,:))-(max(data(varn,:)))/2
[161]1224                  end if
1225                  if (min(data(varn,:)) .GT. 0)then
[310]1226                     res@trYMaxF = max(data(varn,:))-(min(data(varn,:)))/2
1227                     res@trYMinF = min(data(varn,:))+(max(data(varn,:)))/2
[161]1228                  end if
1229               end if
1230               plot_ps(n) = gsn_csm_xy(wks_ps,t,data(varn,:),res) 
1231            end if
[310]1232
1233            ; ***************************************************
1234            ; transparent legend as placeholder for single plots
1235            ; ***************************************************
1236            lgresmono                    = True
1237            lgMonoDashIndex              = False
1238            lgresmono@lgMonoLineColor    = True
1239            lgresmono@lgLabelFont        = "helvetica"   
1240            lgresmono@lgLabelFontHeightF = .1           
1241            lgresmono@vpWidthF           = 0.4           
1242            lgresmono@vpHeightF          = 0.4         
1243            lgresmono@lgLineColor        = 0
1244            lgresmono@lgPerimOn          = False
1245            lbid = gsn_create_legend(wks_ps,1,(/""/),lgresmono)       
1246
1247            amres = True
1248            amres@amParallelPosF   = 0.6                   
1249            amres@amOrthogonalPosF = -0.2           
1250            annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]1251       
1252         else
1253           
1254            print("plot of " + vNam(varn))
1255           
1256            n=n+1
1257            res@xyLineColors            = (/237/)
1258            res@gsnLeftString       = vNam(varn)
1259            res@gsnRightString      = unit(varn)
[310]1260       
[218]1261            if (norm_t .NE. 1.)then
[534]1262               res@tiXAxisString        = "t ("+unit_t+")"
[218]1263            else
[534]1264               res@tiXAxisString        = "t (s)"
[218]1265            end if
1266            res@tiYAxisString = " "
1267            res@tiXAxisFontHeightF   = font_size
1268            res@txFontHeightF        = font_size
1269            res@tiYAxisFontHeightF   = font_size
1270            res@tmXBLabelFontHeightF = font_size
1271            res@tmYLLabelFontHeightF = font_size
[310]1272            res@trYMaxF        = max(data(varn,:))
1273            res@trYMinF        = min(data(varn,:))
[161]1274            if (min(data(varn,:)) .EQ. max(data(varn,:))) then
1275               if (min(data(varn,:)) .EQ. 0)then
[310]1276                     res@trYMaxF = max(data(varn,:))-0.1
1277                     res@trYMinF = min(data(varn,:))+0.1
[161]1278                  end if
1279                  if (min(data(varn,:)) .LT. 0)then
[310]1280                     res@trYMaxF = max(data(varn,:))+(min(data(varn,:)))/2
1281                     res@trYMinF = min(data(varn,:))-(max(data(varn,:)))/2
[161]1282                  end if
1283                  if (min(data(varn,:)) .GT. 0)then
[310]1284                     res@trYMaxF = max(data(varn,:))-(min(data(varn,:)))/2
1285                     res@trYMinF = min(data(varn,:))+(max(data(varn,:)))/2
[161]1286                  end if
1287            end if
1288            plot_ps(n) = gsn_csm_xy(wks_ps,t,data(varn,:),res)
1289           
1290         end if     
1291      end if
[154]1292   end do
[161]1293 
[154]1294   ; ***************************************************
1295   ; merge plots onto one page
1296   ; ***************************************************
[532]1297
1298   no_frames = 0
[162]1299 
[566]1300   if ((format_out .EQ. "eps" .OR. format_out .EQ. "epsi") .AND. \
1301        n .gt. no_rows*no_columns) then
[161]1302      gsn_panel(wks_ps,plot_ps(1:n),(/n,1/),resP)
[250]1303      print(" ")
1304      print("Outputs to .eps or .epsi have only one frame")
1305      print(" ")
[161]1306   else
[218]1307      do np = 1,n,no_rows*no_columns   
1308         if ( np + no_rows*no_columns .gt. n) then   
1309            gsn_panel(wks_ps, plot_ps(np:n),(/no_rows,no_columns/),resP)
[532]1310            no_frames = no_frames + 1           
[161]1311         else
[566]1312            gsn_panel(wks_ps, plot_ps(np:np+no_rows*no_columns-1),\
1313                                                   (/no_rows,no_columns/),resP)
[532]1314            no_frames = no_frames + 1
[161]1315         end if
1316      end do
1317   end if
[154]1318
[532]1319   if (format_out .EQ. "png" ) then
1320     png_output = new((/no_frames/), string)
1321     j = 0
1322     do i=0, no_frames-1
1323       j = i + 1
1324       if (j .LE. 9) then
1325         png_output(i) = file_out+".00000"+j+".png"
1326       end if
1327       if (j .GT. 9 .AND. j .LE. 99) then
1328         png_output(i) = file_out+".0000"+j+".png"
1329       end if
1330       if (j .GT. 99 .AND. j .LE. 999) then
1331         png_output(i) = file_out+".000"+j+".png"
1332       end if
1333       if (j .GT. 999) then
1334         png_output(i) = file_out+".00"+j+".png"
1335       end if
1336
1337       ;using imagemagick's convert for reducing the white
1338       ;space around the plot
1339       cmd = "convert -geometry 1000x1000 -density 300 -trim " +  \
1340              png_output(i) + " " + png_output(i)
1341       system(cmd)
1342     end do
1343
1344     print(" ")
1345     print("Output to: "+ png_output)
1346     print(" ")
1347   else
1348     print(" ")
1349     print("Output to: " + file_out +"."+ format_out)
1350     print(" ")
1351   end if
[154]1352 
[532]1353 
[154]1354end
Note: See TracBrowser for help on using the repository browser.