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

Last change on this file since 957 was 827, checked in by heinze, 13 years ago

allow plotting of data with very small time increments

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