source: palm/tags/release-3.5/SCRIPTS/NCL/timeseries.ncl @ 4109

Last change on this file since 4109 was 194, checked in by letzel, 14 years ago
  • NCL scripts in trunk/SCRIPTS/NCL updated
File size: 47.0 KB
Line 
1load "$NCARG_ROOT/lib/ncarg/nclex/gsun/gsn_code.ncl"
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 
6;***************************************************
7; load ncl_preferences.ncl
8;***************************************************
9
10function which_script()
11local script
12begin
13   script="timeseries"
14   return(script)
15end
16   
17if (isfilepresent("~/ncl_preferences.ncl")) then
18   loadscript("~/ncl_preferences.ncl")
19else
20  if (isfilepresent("~/palm/current_version/trunk/SCRIPTS/NCL/ncl_preferences.ncl")) then
21     loadscript( "~/palm/current_version/trunk/SCRIPTS/NCL/ncl_preferences.ncl")
22  else
23      print(" ")
24      print("'ncl_preferences.ncl' does not exist in $home or $home/palm/current_version/trunk/SCRIPTS/NCL/")
25      print(" ")
26      exit
27   end if
28end if
29   
30begin
31
32   ;***************************************************
33   ; set up default parameter values and strings
34   ;***************************************************
35 
36   if (file_1 .EQ. "File in") then
37      print(" ")
38      print("Please provide 1st input file 'file_1=' in 'ncl_preferences.ncl'")
39      print(" ")
40      exit
41   else
42      file_in = file_1   
43   end if
44   if (.not. isfilepresent(file_in)) then
45      print(" ")
46      print("1st input file: '"+file_in+"' does not exist")
47      print(" ")
48      exit
49   end if
50
51   if (format_out .NE. "x11" .AND. format_out .NE. "pdf" .AND. format_out .NE. "eps" .AND. format_out .NE. "ps" .AND. format_out .NE. "epsi" .AND. format_out .NE. "ncgm")then
52      print(" ")
53      print("'format_out = "+format_out+"' is invalid and set to'x11'")
54      print(" ")
55      format_out="x11"
56   end if 
57
58   if (over .NE. 0 .AND. over .NE. 1) then
59      print(" ")
60      print("'over'= "+over+" is invalid and set to 0")
61      print(" ")
62      over = 0
63   end if   
64 
65
66   ;***************************************************
67   ; open input file
68   ;***************************************************
69
70   f  = addfile(file_in , "r" )
71
72   vNam  = getfilevarnames(f)
73   print(" ")
74   print("Variables in input file:")
75   print("- "+ vNam)
76   print(" ")
77   dim = dimsizes(vNam)
78   if (dim .EQ. 0) then
79      print(" ")
80      print("There are no data on file")
81      print(" ")
82   end if
83
84   t_all = f->time
85   nt  = dimsizes(t_all)
86   delta_t=t_all(nt-1)/nt
87   
88   ;****************************************************       
89   ; start of time step and different types of mistakes that could be done
90   ;****************************************************
91
92   if (start_time_step .EQ. -1.) then           
93      start_time_step=t_all(0)/3600     
94   else
95      if (start_time_step .GE. t_all(nt-1)/3600)
96         print(" ")
97         print("'start_time_step' = "+ start_time_step +"h is equal or greater than last time step = " + t_all(nt-1)+"s = "+t_all(nt-1)/3600+"h")
98         print(" ")
99         print("Please select another 'start_time_step'")
100         print(" ")
101         exit
102      end if
103      if (start_time_step .LT. t_all(0)/3600)
104         print(" ")
105         print("'start_time_step' = "+ start_time_step +"h is lower than first time step = " + t_all(0)+"s = "+t_all(0)/3600+"h")
106         print(" ")
107         print("Please select another 'start_time_step'")
108         print(" ")
109         exit
110      end if
111   end if
112   do i=0,nt-2     
113      if (start_time_step .GE. (t_all(i)-delta_t/2)/3600 .AND. start_time_step .LT. (t_all(i)+delta_t/2)/3600)then
114         st=i
115         break
116      end if
117   end do
118   if (start_time_step .GE. t_all(nt-1)-delta_t/2 .AND. start_time_step .LT. t_all(nt-1)) then
119      st=nt-2   
120   end if
121   if (.not. isvar("st"))then
122      print(" ")
123      print("'start_time_step' = "+ start_time_step +"h is invalid")
124      print(" ")
125      print("Please select another 'start_time_step'")
126      print(" ")
127      exit
128   end if
129     
130   ;****************************************************
131   ; end of time step and different types of mistakes that could be done
132   ;****************************************************
133
134   if (end_time_step .EQ. -1.) then             
135      end_time_step = t_all(nt-1)/3600 
136   else
137      if (end_time_step .LE. t_all(0)/3600)
138         print(" ")
139         print("'end_time_step' = "+end_time_step+ "h is lower or equal than first time step = " + t_all(0)+"s = "+t_all(0)/3600+"h")
140         print(" ")
141         print("Please select another 'end_time_step'")
142         print(" ")
143         exit
144      end if
145      if (end_time_step .GT. t_all(nt-1)/3600)
146         print(" ")
147         print("'end_time_step' = "+ end_time_step +"h is greater than last time step = " + t_all(nt-1)+"s = "+t_all(nt-1)/3600+"h")
148         print(" ")
149         print("Please select another 'end_time_step'") 
150         print(" ")
151         exit
152      end if
153      if (end_time_step .LE. start_time_step/3600)
154         print(" ")
155         print("'end_time_step' = "+ end_time_step +"h is equal or lower than 'start_time_step' = "+start_time_step+"h")
156         print(" ")
157         print("Please select another 'start_time_step' or 'end_time_step'")
158         print(" ")
159         exit
160      end if
161   end if
162   do i=0,nt-1     
163      if (end_time_step .GE. (t_all(i)-delta_t/2)/3600 .AND. end_time_step .LT. (t_all(i)+delta_t/2)/3600)then
164         et=i
165         break
166      end if
167   end do
168   
169   if (.not. isvar("et"))then
170      print(" ")
171      print("'end_time_step' = "+ end_time_step +"h is invalid")
172      print(" ")
173      print("Please select another 'end_time_step'")
174      print(" ")
175      exit
176   end if
177 
178   delete(start_time_step)
179   start_time_step=round(st,3)
180   delete(end_time_step)
181   end_time_step=round(et,3)
182
183   print(" ")
184   print("Output of time steps from "+t_all(start_time_step)/3600+" h = "+t_all(start_time_step)+" s => index = "+start_time_step)
185   print("                     till "+t_all(end_time_step)/3600+" h = "+t_all(end_time_step)+" s => index = "+end_time_step)
186   print(" ")
187
188   t = f->time(start_time_step:end_time_step)
189
190 
191   ; ***************************************************
192   ; set up recourses
193   ; ***************************************************
194
195   res                         = True
196   res@gsnDraw                 = False
197   res@gsnFrame                = False
198   res@gsnPaperOrientation     = "portrait"
199   res@gsnPaperWidth           = 8.27
200   res@gsnPaperHeight          = 11.69
201   res@gsnPaperMargin          = 0.79
202   res@tmXBMode                = True
203   res@tmYLMode                = True
204   res@txFont                  = "helvetica"
205   res@tiMainFont              = "helvetica"
206   res@tiXAxisFont             = "helvetica"
207   res@tiYAxisFont             = "helvetica"
208   res@tmXBLabelFont           = "helvetica"
209   res@tmYLLabelFont           = "helvetica"
210   res@xyLineColors            = (/237/)
211   
212   res@lgLabelFontHeightF     = .02
213
214   resP                        = True
215   resP@txFont                 = "helvetica"
216   resP@txString               = f@title+" time series "
217   resP@txFuncCode             = "~"
218   resP@txFontHeightF          = 0.015
219
220   res@vpWidthF=4
221
222   txres = True
223
224   ; ***************************************************
225   ; read data and create plots
226   ; ***************************************************
227   
228   wks_ps = gsn_open_wks(format_out,file_out)           
229   gsn_define_colormap(wks_ps,"rainbow+white")
230   plot_ps=new(dim,graphic)
231                                 
232   n=0
233   minE=1.E27
234   maxE=-1.E27
235   minus=1.E27
236   maxus=-1.E27
237   minu=1.E27
238   maxu=-1.E27
239   minz=1.E27
240   maxz=-1.E27
241   minw=1.E27
242   maxw=-1.E27
243   minp=1.E27
244   maxp=-1.E27
245   mins=1.E27
246   maxs=-1.E27
247
248   data   = new((/dim,(end_time_step-start_time_step)+1/),float)
249   unit   = new(dim,string)
250   data_0 = new((end_time_step-start_time_step)+1,float)
251   data_0 = 0.0
252   mini   = new(dim,float)
253   maxi   = new(dim,float)
254   
255   if (over .EQ. 1) then
256      plot_E       = gsn_csm_xy(wks_ps,t,data_0(:),res)
257      plot_Es      = gsn_csm_xy(wks_ps,t,data_0(:),res)
258      plot_us      = gsn_csm_xy(wks_ps,t,data_0(:),res)
259      plot_ws      = gsn_csm_xy(wks_ps,t,data_0(:),res)
260      plot_umax    = gsn_csm_xy(wks_ps,t,data_0(:),res)
261      plot_vmax    = gsn_csm_xy(wks_ps,t,data_0(:),res)
262      plot_wmax    = gsn_csm_xy(wks_ps,t,data_0(:),res)
263      plot_z_i_wpt = gsn_csm_xy(wks_ps,t,data_0(:),res)
264      plot_z_i_pt  = gsn_csm_xy(wks_ps,t,data_0(:),res)
265      plot_wpptp0  = gsn_csm_xy(wks_ps,t,data_0(:),res)
266      plot_wpptp   = gsn_csm_xy(wks_ps,t,data_0(:),res)
267      plot_wpt     = gsn_csm_xy(wks_ps,t,data_0(:),res)
268      plot_pt_0_   = gsn_csm_xy(wks_ps,t,data_0(:),res)
269      plot_pt_zp_  = gsn_csm_xy(wks_ps,t,data_0(:),res)
270      plot_splptx  = gsn_csm_xy(wks_ps,t,data_0(:),res)
271      plot_splpty  = gsn_csm_xy(wks_ps,t,data_0(:),res)
272      plot_splptz  = gsn_csm_xy(wks_ps,t,data_0(:),res)
273   end if
274
275   count_var=0
276   do varn = dim-1,0,1
277
278      if( isStrSubset (vNam(varn), "time") )
279         check = False
280      else
281         check = True
282      end if
283 
284      if (var .NE. "all") then
285         check = isStrSubset( var,","+vNam(varn)+"," )
286      end if
287   
288     
289      if(check) then
290         count_var=count_var+1
291
292         data_all = f ->$vNam(varn)$
293         unit(varn) = data_all@units
294               
295         data(varn,:)=data_all(start_time_step:end_time_step)
296         
297         if (over .EQ. 1) then
298
299            mini(varn) = min(data(varn,:))
300            maxi(varn) = max(data(varn,:))
301           
302            if (vNam(varn) .EQ. "E" .OR. vNam(varn) .EQ. "Es") then
303               if (mini(varn) .EQ. maxi(varn)) then
304                  if (min(data(varn,:)) .EQ. 0)then
305                     mini(varn)= mini(varn)-1.
306                     maxi(varn)= maxi(varn)+1.
307                  end if
308                  if (min(data(varn,:)) .LT. 0)then
309                     mini(varn)= mini(varn)-1.+(mini(varn))/2
310                     maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
311                  end if
312                  if (min(data(varn,:)) .GT. 0)then
313                     mini(varn)= mini(varn)-1.-(mini(varn))/2
314                     maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
315                  end if
316               end if
317               minE=min((/minE,mini(varn)/))
318               maxE=max((/maxE,maxi(varn)/))
319            end if
320
321            if (vNam(varn) .EQ. "us" .OR. vNam(varn) .EQ. "ws") then
322               if (mini(varn) .EQ. maxi(varn)) then
323                  if (min(data(varn,:)) .EQ. 0)then
324                     mini(varn)= mini(varn)-1.
325                     maxi(varn)= maxi(varn)+1.
326                  end if
327                  if (min(data(varn,:)) .LT. 0)then
328                     mini(varn)= mini(varn)-1.+(mini(varn))/2
329                     maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
330                  end if
331                  if (min(data(varn,:)) .GT. 0)then
332                     mini(varn)= mini(varn)-1.-(mini(varn))/2
333                     maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
334                  end if
335               end if
336               minus=min((/minus,mini(varn)/))
337               maxus=max((/maxus,maxi(varn)/))
338            end if
339
340            if (vNam(varn) .EQ. "umax" .OR. vNam(varn) .EQ. "vmax" .OR. vNam(varn) .EQ. "wmax") then
341               if (mini(varn) .EQ. maxi(varn)) then
342                  if (mini(varn) .EQ. 0)then
343                     mini(varn)= mini(varn)-1.
344                     maxi(varn)= maxi(varn)+1.
345                  end if
346                  if (mini(varn) .LT. 0)then
347                     mini(varn)= mini(varn)-1.+(mini(varn))/2
348                     maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
349                  end if
350                  if (mini(varn) .GT. 0)then
351                     mini(varn)= mini(varn)-1.-(mini(varn))/2
352                     maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
353                  end if
354               end if
355               minu=min((/minu,mini(varn)/))
356               maxu=max((/maxu,maxi(varn)/))
357            end if
358
359            if (vNam(varn) .EQ. "z_i_wpt" .OR. vNam(varn) .EQ. "z_i_pt") then
360               if (mini(varn) .EQ. maxi(varn)) then
361                  if (min(data(varn,:)) .EQ. 0)then
362                     mini(varn)= mini(varn)-1.
363                     maxi(varn)= maxi(varn)+1.
364                  end if
365                  if (min(data(varn,:)) .LT. 0)then
366                     mini(varn)= mini(varn)-1.+(mini(varn))/2
367                     maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
368                  end if
369                  if (min(data(varn,:)) .GT. 0)then
370                     mini(varn)= mini(varn)-1.-(mini(varn))/2
371                     maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
372                  end if
373               end if
374               minz=min((/minz,mini(varn)/))
375               maxz=max((/maxz,maxi(varn)/))
376            end if
377
378            if (vNam(varn) .EQ. "wpptp0" .OR. vNam(varn) .EQ. "wpptp" .OR. vNam(varn) .EQ. "wpt") then
379               if (mini(varn) .EQ. maxi(varn)) then
380                  if (min(data(varn,:)) .EQ. 0)then
381                     mini(varn)= mini(varn)-1.
382                     maxi(varn)= maxi(varn)+1.
383                  end if
384                  if (min(data(varn,:)) .LT. 0)then
385                     mini(varn)= mini(varn)-1.+(mini(varn))/2
386                     maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
387                  end if
388                  if (min(data(varn,:)) .GT. 0)then
389                     mini(varn)= mini(varn)-1.-(mini(varn))/2
390                     maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
391                  end if
392               end if
393               minw=min((/minw,mini(varn)/))
394               maxw=max((/maxw,maxi(varn)/))
395            end if
396
397            if (vNam(varn) .EQ. "pt_0_" .OR. vNam(varn) .EQ. "pt_zp_") then
398               if (mini(varn) .EQ. maxi(varn)) then
399                  if (min(data(varn,:)) .EQ. 0)then
400                     mini(varn)= mini(varn)-1.
401                     maxi(varn)= maxi(varn)+1.
402                  end if
403                  if (min(data(varn,:)) .LT. 0)then
404                     mini(varn)= mini(varn)-1.+(mini(varn))/2
405                     maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
406                  end if
407                  if (min(data(varn,:)) .GT. 0)then
408                     mini(varn)= mini(varn)-1.-(mini(varn))/2
409                     maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
410                  end if
411               end if
412               minp=min((/minp,mini(varn)/))
413               maxp=max((/maxp,maxi(varn)/))
414            end if
415
416            if (vNam(varn) .EQ. "splptx" .OR. vNam(varn) .EQ. "splpty" .OR. vNam(varn) .EQ. "splptz") then
417               if (mini(varn) .EQ. maxi(varn)) then
418                  if (min(data(varn,:)) .EQ. 0)then
419                     mini(varn)= mini(varn)-1.
420                     maxi(varn)= maxi(varn)+1.
421                  end if
422                  if (min(data(varn,:)) .LT. 0)then
423                     mini(varn)= mini(varn)-1.+(mini(varn))/2
424                     maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
425                  end if
426                  if (min(data(varn,:)) .GT. 0)then
427                     mini(varn)= mini(varn)-1.-(mini(varn))/2
428                     maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
429                  end if
430               end if
431               mins=min((/mins,mini(varn)/))
432               maxs=max((/maxs,maxi(varn)/))
433            end if
434           
435         end if
436      end if
437   end do
438
439   if (count_var .EQ. 0) then
440      print(" ")
441      print("The variables 'var="+var+"' do not exist on your input file;")
442      print("be sure to have one comma berfore and after each variable")
443      print(" ")
444      exit
445   end if
446
447   do varn = dim-1,0,1
448     
449      if( isStrSubset (vNam(varn), "time") )
450         check = False
451      else
452         check = True
453      end if     
454   
455      if (var .NE. "all") then
456         check = isStrSubset( var,","+vNam(varn)+"," )
457      end if
458 
459      if(check) then
460
461        if (isStrSubset(vNam(varn),"_0" ))then
462            print(" ")
463            print("If you have Outputs of statistic regions you cannot overlay variables; 'over' is set to 0")
464            print(" ")
465            over = 0
466         end if
467
468        if (over .EQ. 1) then
469 
470            res@gsnLeftString       = "overlayed plot"
471            res@gsnRightString      = unit(varn)
472            res@tiXAxisString        = " time [s] "
473            res@tiYAxisString        = " "
474            res@tiXAxisFontHeightF   = 0.07
475            res@txFontHeightF        = 0.07
476            res@tiYAxisFontHeightF   = 0.07           
477
478            if (vNam(varn) .EQ. "E")
479               E=0
480               res@xyLineColors     = (/237/)
481               res@xyLineLabelFontHeightF = 0.05
482               res@xyLineLabelFontColor   = 237
483               res@trYMaxF           = minE
484               res@trYMinF           = maxE
485               plot_E = gsn_csm_xy(wks_ps,t,data(varn,:),res)
486            end if
487            if (vNam(varn) .EQ. "Es")
488               Es=0
489               res@xyLineColors            = (/144/)
490               res@xyLineLabelFontHeightF = 0.05
491               res@xyLineLabelFontColor   = 144
492               plot_Es = gsn_csm_xy(wks_ps,t,data(varn,:),res)
493            end if
494             
495            if (vNam(varn) .EQ. "us")
496               us=0
497               res@xyLineColors            = (/237/)
498               res@xyLineLabelFontHeightF = 0.05
499               res@xyLineLabelFontColor   = 237
500               res@trYMaxF           = minus
501               res@trYMinF           = maxus
502               plot_us = gsn_csm_xy(wks_ps,t,data(varn,:),res)
503            end if             
504            if (vNam(varn) .EQ. "ws")
505               ws=0
506               res@xyLineColors            = (/144/)
507               res@xyLineLabelFontHeightF = 0.05
508               res@xyLineLabelFontColor   = 144
509               plot_ws = gsn_csm_xy(wks_ps,t,data(varn,:),res)
510            end if
511             
512            if (vNam(varn) .EQ. "umax")
513               u=0
514               res@xyLineColors            = (/237/)
515               res@xyLineLabelFontHeightF = 0.05
516               res@xyLineLabelFontColor   = 237
517               res@trYMaxF           = minu
518               res@trYMinF           = maxu
519               plot_umax = gsn_csm_xy(wks_ps,t,data(varn,:),res)
520            end if
521            if (vNam(varn) .EQ. "vmax")
522               v=0
523               res@xyLineColors            = (/144/)
524               res@xyLineLabelFontHeightF = 0.05
525               res@xyLineLabelFontColor   = 144
526               plot_vmax = gsn_csm_xy(wks_ps,t,data(varn,:),res)
527            end if
528            if (vNam(varn) .EQ. "wmax")
529               w=0
530               res@xyLineColors            = (/80/)
531               res@xyLineLabelFontHeightF = 0.05
532               res@xyLineLabelFontColor   = 80
533               plot_wmax = gsn_csm_xy(wks_ps,t,data(varn,:),res)
534            end if
535   
536            if (vNam(varn) .EQ. "z_i_wpt")
537               zw=0
538               res@xyLineColors            = (/237/)
539               res@xyLineLabelFontHeightF = 0.05
540               res@xyLineLabelFontColor   = 237
541               res@trYMaxF           = minz
542               res@trYMinF           = maxz
543               plot_z_i_wpt = gsn_csm_xy(wks_ps,t,data(varn,:),res)
544            end if
545            if (vNam(varn) .EQ. "z_i_pt")
546               z=0
547               res@xyLineColors            = (/144/) 
548               res@xyLineLabelFontHeightF = 0.05
549               res@xyLineLabelFontColor   = 144
550               plot_z_i_pt = gsn_csm_xy(wks_ps,t,data(varn,:),res)
551            end if
552
553            if (vNam(varn) .EQ. "wpptp0")
554               w0=0
555               res@xyLineColors            = (/237/)
556               res@xyLineLabelFontHeightF = 0.05
557               res@xyLineLabelFontColor   = 237
558               res@trYMaxF           = minw
559               res@trYMinF           = maxw
560               plot_wpptp0 = gsn_csm_xy(wks_ps,t,data(varn,:),res)
561            end if
562            if (vNam(varn) .EQ. "wpptp")
563               wp=0
564               res@xyLineColors            = (/144/)
565               res@xyLineLabelFontHeightF = 0.05
566               res@xyLineLabelFontColor   = 144
567               plot_wpptp = gsn_csm_xy(wks_ps,t,data(varn,:),res) 
568            end if
569            if (vNam(varn) .EQ. "wpt")
570               wt=0
571               res@xyLineColors            = (/80/)
572               res@xyLineLabelFontHeightF = 0.05
573               res@xyLineLabelFontColor   = 80
574               plot_wpt = gsn_csm_xy(wks_ps,t,data(varn,:),res)
575            end if
576
577            if (vNam(varn) .EQ. "pt_0_")
578               p=0
579               res@xyLineColors            = (/237/)
580               res@xyLineLabelFontHeightF = 0.05
581               res@xyLineLabelFontColor   = 237
582               res@trYMaxF           = minp
583               res@trYMinF           = maxp
584               plot_pt_0_ = gsn_csm_xy(wks_ps,t,data(varn,:),res) 
585            end if
586            if (vNam(varn) .EQ. "pt_zp_")
587               pz=0
588               res@xyLineColors            = (/144/)
589               res@xyLineLabelFontHeightF = 0.05
590               res@xyLineLabelFontColor   = 144
591               plot_pt_zp_ = gsn_csm_xy(wks_ps,t,data(varn,:),res)
592            end if
593
594            if (vNam(varn) .EQ. "splptx")
595               x=0
596               res@xyLineColors            = (/237/)
597               res@xyLineLabelFontHeightF = 0.05
598               res@xyLineLabelFontColor   = 237
599               res@trYMaxF           = mins
600               res@trYMinF           = maxs
601               plot_splptx = gsn_csm_xy(wks_ps,t,data(varn,:),res)   
602            end if
603            if (vNam(varn) .EQ. "splpty")
604               y=0
605               res@xyLineColors            = (/144/)
606               res@xyLineLabelFontHeightF = 0.05
607               res@xyLineLabelFontColor   = 144
608               plot_splpty = gsn_csm_xy(wks_ps,t,data(varn,:),res)
609            end if
610            if (vNam(varn) .EQ. "splptz")
611               z=0
612               res@xyLineColors            = (/80/)
613               res@xyLineLabelFontHeightF = 0.05
614               res@xyLineLabelFontColor   = 80
615               plot_splptz = gsn_csm_xy(wks_ps,t,data(varn,:),res)         
616            end if
617
618         end if
619      end if
620   end do
621   
622   do varn = dim-1,0,1
623 
624      if( isStrSubset (vNam(varn), "time") )
625         check = False
626      else
627         check = True
628      end if   
629 
630      if (var.NE. "all") then
631         check = isStrSubset( var,","+vNam(varn)+"," )
632      end if
633   
634      if(check) then
635       
636         if (over .EQ. 1) then       
637       
638            if (vNam(varn) .EQ. "E" .AND. Es .NE. 1) then
639               E=1   
640               overlay(plot_E,plot_Es)
641               n=n+1
642               plot_ps(n) = plot_E   
643
644               ; ***************************************************
645               ; legend for combined plot
646               ; ***************************************************
647     
648               lgres                    = True
649               lgMonoDashIndex          = False
650               lgres@lgLabelFont        = "helvetica"   
651               lgres@lgLabelFontHeightF = .1           
652               lgres@vpWidthF           = 0.4           
653               lgres@vpHeightF          = 0.4         
654               lgres@lgDashIndexes      = (/0,0,0/)
655               lgres@lgLineColors       = (/237,144,80/)
656               lbid = gsn_create_legend(wks_ps,2,(/"E","Es"/),lgres)       
657
658               amres = True
659               amres@amParallelPosF   = 0.6                 
660               amres@amOrthogonalPosF = -0.2           
661               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)           
662            end if
663            if (vNam(varn) .EQ. "Es" .AND. E .NE. 1) then
664               Es=1
665               overlay(plot_E,plot_Es)
666               n=n+1
667               plot_ps(n) = plot_E 
668
669               ; ***************************************************
670               ; legend for combined plot
671               ; ***************************************************
672     
673               lgres                    = True
674               lgMonoDashIndex          = False
675               lgres@lgLabelFont        = "helvetica"   
676               lgres@lgLabelFontHeightF = .1           
677               lgres@vpWidthF           = 0.4           
678               lgres@vpHeightF          = 0.4         
679               lgres@lgDashIndexes      = (/0,0,0/)
680               lgres@lgLineColors       = (/237,144,80/)
681               lbid = gsn_create_legend(wks_ps,2,(/"E","Es"/),lgres)       
682
683               amres = True
684               amres@amParallelPosF   = 0.6                 
685               amres@amOrthogonalPosF = -0.2           
686               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)         
687            end if
688
689            if (vNam(varn) .EQ. "us" .AND. ws .NE. 1) then
690               us=1
691               overlay(plot_us,plot_ws)
692               n=n+1
693               plot_ps(n) = plot_us
694
695               ; ***************************************************
696               ; legend for combined plot
697               ; ***************************************************
698     
699               lgres                    = True
700               lgMonoDashIndex          = False
701               lgres@lgLabelFont        = "helvetica"   
702               lgres@lgLabelFontHeightF = .1           
703               lgres@vpWidthF           = 0.4           
704               lgres@vpHeightF          = 0.4         
705               lgres@lgDashIndexes      = (/0,0,0/)
706               lgres@lgLineColors       = (/237,144,80/)
707               lbid = gsn_create_legend(wks_ps,2,(/"us","ws"/),lgres)       
708
709               amres = True
710               amres@amParallelPosF   = 0.6                 
711               amres@amOrthogonalPosF = -0.2           
712               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
713            end if
714            if (vNam(varn) .EQ. "ws" .AND. us .NE. 1) then
715               ws=1
716               overlay(plot_us,plot_ws)
717               n=n+1
718               plot_ps(n) = plot_us
719
720               ; ***************************************************
721               ; legend for combined plot
722               ; ***************************************************
723     
724               lgres                    = True
725               lgMonoDashIndex          = False
726               lgres@lgLabelFont        = "helvetica"   
727               lgres@lgLabelFontHeightF = .1           
728               lgres@vpWidthF           = 0.4           
729               lgres@vpHeightF          = 0.4         
730               lgres@lgDashIndexes      = (/0,0,0/)
731               lgres@lgLineColors       = (/237,144,80/)
732               lbid = gsn_create_legend(wks_ps,2,(/"us","ws"/),lgres)       
733
734               amres = True
735               amres@amParallelPosF   = 0.6                 
736               amres@amOrthogonalPosF = -0.2           
737               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
738            end if
739         
740            if (vNam(varn) .EQ. "umax" .AND. v .NE. 1)
741               if (w .NE. 1) then
742                  u=1         
743                  overlay(plot_umax,plot_vmax)
744                  overlay(plot_umax,plot_wmax)
745                  n=n+1
746                  plot_ps(n) = plot_umax
747
748                  ; ***************************************************
749                  ; legend for combined plot
750                  ; ***************************************************
751     
752                  lgres                    = True
753                  lgMonoDashIndex          = False
754                  lgres@lgLabelFont        = "helvetica"   
755                  lgres@lgLabelFontHeightF = .1           
756                  lgres@vpWidthF           = 0.4           
757                  lgres@vpHeightF          = 0.4         
758                  lgres@lgDashIndexes      = (/0,0,0/)
759                  lgres@lgLineColors       = (/237,144,80/)
760                  lbid = gsn_create_legend(wks_ps,3,(/"umax","vmax","wmax"/),lgres)       
761
762                  amres = True
763                  amres@amParallelPosF   = 0.6             
764                  amres@amOrthogonalPosF = -0.2           
765                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
766               end if
767            end if
768            if (vNam(varn) .EQ. "vmax" .AND. u .NE. 1)
769               if (w .NE. 1) then
770                  v=1 
771                  overlay(plot_umax,plot_vmax)
772                  overlay(plot_umax,plot_wmax)
773                  n=n+1
774                  plot_ps(n) = plot_umax
775
776                  ; ***************************************************
777                  ; legend for combined plot
778                  ; ***************************************************
779     
780                  lgres                    = True
781                  lgMonoDashIndex          = False
782                  lgres@lgLabelFont        = "helvetica"   
783                  lgres@lgLabelFontHeightF = .1         
784                  lgres@vpWidthF           = 0.4           
785                  lgres@vpHeightF          = 0.4         
786                  lgres@lgDashIndexes      = (/0,0,0/)
787                  lgres@lgLineColors       = (/237,144,80/)
788                  lbid = gsn_create_legend(wks_ps,3,(/"umax","vmax","wmax"/),lgres)       
789
790                  amres = True
791                  amres@amParallelPosF   = 0.6             
792                  amres@amOrthogonalPosF = -0.2           
793                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
794               end if
795            end if
796            if (vNam(varn) .EQ. "wmax" .AND. v .NE. 1)
797               if(u .NE. 1) then 
798                  w=1       
799                  overlay(plot_umax,plot_vmax)
800                  overlay(plot_umax,plot_wmax)
801                  n=n+1
802                  plot_ps(n) = plot_umax
803
804                  ; ***************************************************
805                  ; legend for combined plot
806                  ; ***************************************************
807     
808                  lgres                    = True
809                  lgMonoDashIndex          = False
810                  lgres@lgLabelFont        = "helvetica"   
811                  lgres@lgLabelFontHeightF = .1           
812                  lgres@vpWidthF           = 0.4           
813                  lgres@vpHeightF          = 0.4         
814                  lgres@lgDashIndexes      = (/0,0,0/)
815                  lgres@lgLineColors       = (/237,144,80/)
816                  lbid = gsn_create_legend(wks_ps,3,(/"umax","vmax","wmax"/),lgres)       
817
818                  amres = True
819                  amres@amParallelPosF   = 0.6             
820                  amres@amOrthogonalPosF = -0.2           
821                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
822               end if
823            end if
824           
825            if (vNam(varn) .EQ. "z_i_wpt" .AND. z .NE. 1) then
826               zw=1       
827               overlay(plot_z_i_wpt,plot_z_i_pt)
828               n=n+1
829               plot_ps(n) = plot_z_i_wpt
830       
831               ; ***************************************************
832               ; legend for combined plot
833               ; ***************************************************
834     
835               lgres                    = True
836               lgMonoDashIndex          = False
837               lgres@lgLabelFont        = "helvetica"   
838               lgres@lgLabelFontHeightF = .1           
839               lgres@vpWidthF           = 0.4           
840               lgres@vpHeightF          = 0.4         
841               lgres@lgDashIndexes      = (/0,0,0/)
842               lgres@lgLineColors       = (/237,144,80/)
843               lbid = gsn_create_legend(wks_ps,2,(/"z_i_wpt","z_i_pt"/),lgres)       
844
845               amres = True
846               amres@amParallelPosF   = 0.6                 
847               amres@amOrthogonalPosF = -0.2           
848               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)         
849            end if
850            if (vNam(varn) .EQ. "z_i_pt" .AND. zw .NE. 1) then
851               z=1     
852               overlay(plot_z_i_wpt,plot_z_i_pt)
853               n=n+1
854               plot_ps(n) = plot_z_i_wpt
855
856               ; ***************************************************
857               ; legend for combined plot
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/)
868               lbid = gsn_create_legend(wks_ps,2,(/"z_i_wpt","z_i_pt"/),lgres)       
869
870               amres = True
871               amres@amParallelPosF   = 0.6                 
872               amres@amOrthogonalPosF = -0.2           
873               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)           
874            end if   
875         
876            if (vNam(varn) .EQ. "wpptp0" .AND. wp .NE. 1)
877               if (wt .NE. 1) then
878                  w0=1
879                  overlay(plot_wpptp0,plot_wpptp)
880                  overlay(plot_wpptp0,plot_wpt)
881                  n=n+1
882                  plot_ps(n) = plot_wpptp0
883
884                  ; ***************************************************
885                  ; legend for combined plot
886                  ; ***************************************************
887     
888                  lgres                    = True
889                  lgMonoDashIndex          = False
890                  lgres@lgLabelFont        = "helvetica"   
891                  lgres@lgLabelFontHeightF = .1           
892                  lgres@vpWidthF           = 0.4           
893                  lgres@vpHeightF          = 0.4         
894                  lgres@lgDashIndexes      = (/0,0,0/)
895                  lgres@lgLineColors       = (/237,144,80/)
896                  lbid = gsn_create_legend(wks_ps,3,(/"wpptp0","wpptp","wpt"/),lgres)       
897
898                  amres = True
899                  amres@amParallelPosF   = 0.6             
900                  amres@amOrthogonalPosF = -0.2           
901                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
902               end if           
903            end if
904            if (vNam(varn) .EQ. "wpptp" .AND. w0 .NE. 1)
905               if (wt .NE. 1) then
906                 wp=1
907                 overlay(plot_wpptp0,plot_wpptp)
908                 overlay(plot_wpptp0,plot_wpt)
909                 n=n+1
910                 plot_ps(n) = plot_wpptp0
911                 
912                 ; ***************************************************
913                  ; legend for combined plot
914                  ; ***************************************************
915     
916                  lgres                    = True
917                  lgMonoDashIndex          = False
918                  lgres@lgLabelFont        = "helvetica"   
919                  lgres@lgLabelFontHeightF = .1           
920                  lgres@vpWidthF           = 0.4           
921                  lgres@vpHeightF          = 0.4         
922                  lgres@lgDashIndexes      = (/0,0,0/)
923                  lgres@lgLineColors       = (/237,144,80/)
924                  lbid = gsn_create_legend(wks_ps,3,(/"wpptp0","wpptp","wpt"/),lgres)       
925
926                  amres = True
927                  amres@amParallelPosF   = 0.6             
928                  amres@amOrthogonalPosF = -0.2           
929                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres) 
930               end if           
931            end if
932            if (vNam(varn) .EQ. "wpt" .AND. wp .NE. 1)
933               if (w0 .NE. 1) then
934                  wt=1
935                  overlay(plot_wpptp0,plot_wpptp)
936                  overlay(plot_wpptp0,plot_wpt)
937                  n=n+1
938                  plot_ps(n) = plot_wpptp0
939
940                  ; ***************************************************
941                  ; legend for combined plot
942                  ; ***************************************************
943     
944                  lgres                    = True
945                  lgMonoDashIndex          = False
946                  lgres@lgLabelFont        = "helvetica"   
947                  lgres@lgLabelFontHeightF = .1           
948                  lgres@vpWidthF           = 0.4           
949                  lgres@vpHeightF          = 0.4         
950                  lgres@lgDashIndexes      = (/0,0,0/)
951                  lgres@lgLineColors       = (/237,144,80/)
952                  lbid = gsn_create_legend(wks_ps,3,(/"wpptp0","wpptp","wpt"/),lgres)       
953
954                  amres = True
955                  amres@amParallelPosF   = 0.6             
956                  amres@amOrthogonalPosF = -0.2           
957                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
958               end if
959            end if
960
961            if (vNam(varn) .EQ. "pt_0_" .AND. pz .NE. 1) then
962               p=1     
963               overlay(plot_pt_0_,plot_pt_zp_)
964               n=n+1
965               plot_ps(n) = plot_pt_0_ 
966     
967               ; ***************************************************
968               ; legend for combined plot
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/)
979               lbid = gsn_create_legend(wks_ps,2,(/"pt_0_","pt_zp_"/),lgres)       
980
981               amres = True
982               amres@amParallelPosF   = 0.6                 
983               amres@amOrthogonalPosF = -0.2           
984               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
985            end if
986            if (vNam(varn) .EQ. "pt_zp_" .AND. p .NE. 1) then
987               pz=1       
988               overlay(plot_pt_0_,plot_pt_zp_)
989               n=n+1
990               plot_ps(n) = plot_pt_0_   
991
992               ; ***************************************************
993               ; legend for combined plot
994               ; ***************************************************
995     
996               lgres                    = True
997               lgMonoDashIndex          = False
998               lgres@lgLabelFont        = "helvetica"   
999               lgres@lgLabelFontHeightF = .1           
1000               lgres@vpWidthF           = 0.4           
1001               lgres@vpHeightF          = 0.4         
1002               lgres@lgDashIndexes      = (/0,0,0/)
1003               lgres@lgLineColors       = (/237,144,80/)
1004               lbid = gsn_create_legend(wks_ps,2,(/"pt_0_","pt_zp_"/),lgres)       
1005
1006               amres = True
1007               amres@amParallelPosF   = 0.6                 
1008               amres@amOrthogonalPosF = -0.2           
1009               annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
1010            end if
1011           
1012            if (vNam(varn) .EQ. "splptx" .AND. y .NE. 1)
1013               if (z .NE.1 ) then
1014                  x=1     
1015                  overlay(plot_splptx,plot_splpty)
1016                  overlay(plot_splptx,plot_splptz)
1017                  n=n+1
1018                  plot_ps(n) = plot_splptx   
1019
1020                  ; ***************************************************
1021                  ; legend for combined plot
1022                  ; ***************************************************
1023     
1024                  lgres                    = True
1025                  lgMonoDashIndex          = False
1026                  lgres@lgLabelFont        = "helvetica"   
1027                  lgres@lgLabelFontHeightF = .1           
1028                  lgres@vpWidthF           = 0.4           
1029                  lgres@vpHeightF          = 0.4         
1030                  lgres@lgDashIndexes      = (/0,0,0/)
1031                  lgres@lgLineColors       = (/237,144,80/)
1032                  lbid = gsn_create_legend(wks_ps,3,(/"splptx","splpty","splptz"/),lgres)       
1033
1034                  amres = True
1035                  amres@amParallelPosF   = 0.6             
1036                  amres@amOrthogonalPosF = -0.2           
1037                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
1038               end if
1039            end if
1040            if (vNam(varn) .EQ. "splpty" .AND. x .NE. 1)
1041               if(z .NE.1 ) then
1042                  y=1         
1043                  overlay(plot_splptx,plot_splpty)
1044                  overlay(plot_splptx,plot_splptz)
1045                  n=n+1
1046                  plot_ps(n) = plot_splptx
1047
1048                  ; ***************************************************
1049                  ; legend for combined plot
1050                  ; ***************************************************
1051     
1052                  lgres                    = True
1053                  lgMonoDashIndex          = False
1054                  lgres@lgLabelFont        = "helvetica"   
1055                  lgres@lgLabelFontHeightF = .1           
1056                  lgres@vpWidthF           = 0.4           
1057                  lgres@vpHeightF          = 0.4         
1058                  lgres@lgDashIndexes      = (/0,0,0/)
1059                  lgres@lgLineColors       = (/237,144,80/)
1060                  lbid = gsn_create_legend(wks_ps,3,(/"splptx","splpty","splptz"/),lgres)       
1061
1062                  amres = True
1063                  amres@amParallelPosF   = 0.6             
1064                  amres@amOrthogonalPosF = -0.2           
1065                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
1066               end if           
1067            end if
1068            if (vNam(varn) .EQ. "splptz" .AND. y .NE. 1)
1069               if(x .NE.1 ) then
1070                  z=1         
1071                  overlay(plot_splptx,plot_splpty)
1072                  overlay(plot_splptx,plot_splptz)
1073                  n=n+1
1074                  plot_ps(n) = plot_splptx 
1075
1076                  ; ***************************************************
1077                  ; legend for combined plot
1078                  ; ***************************************************
1079     
1080                  lgres                    = True
1081                  lgMonoDashIndex          = False
1082                  lgres@lgLabelFont        = "helvetica"   
1083                  lgres@lgLabelFontHeightF = .1           
1084                  lgres@vpWidthF           = 0.4           
1085                  lgres@vpHeightF          = 0.4         
1086                  lgres@lgDashIndexes      = (/0,0,0/)
1087                  lgres@lgLineColors       = (/237,144,80/)
1088                  lbid = gsn_create_legend(wks_ps,3,(/"splptx","splpty","splptz"/),lgres)       
1089
1090                  amres = True
1091                  amres@amParallelPosF   = 0.6             
1092                  amres@amOrthogonalPosF = -0.2           
1093                  annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
1094               end if       
1095            end if
1096
1097            if(vNam(varn) .NE. "splptz" .AND. vNam(varn) .NE. "splpty" .AND. vNam(varn) .NE. "splptx" .AND. vNam(varn) .NE. "pt_zp_" .AND. vNam(varn) .NE. "pt_0_" .AND. vNam(varn) .NE. "wpt" .AND. vNam(varn) .NE. "wpptp" .AND. vNam(varn) .NE. "wpptp0" .AND. vNam(varn) .NE. "z_i_pt" .AND. vNam(varn) .NE. "z_i_wpt" .AND. vNam(varn) .NE. "wmax" .AND. vNam(varn) .NE. "vmax" .AND. vNam(varn) .NE. "umax" .AND. vNam(varn) .NE. "ws" .AND.  vNam(varn) .NE. "us" .AND. vNam(varn) .NE. "Es" .AND. vNam(varn) .NE. "E") then
1098               n=n+1
1099               res@xyLineColors   = (/237/)
1100               res@xyLabelMode    = False
1101               res@gsnLeftString  = vNam(varn)
1102               res@gsnRightString = unit(varn)
1103               res@trYMaxF        = min(data(varn,:))
1104               res@trYMinF        = max(data(varn,:))
1105               if (min(data(varn,:)) .EQ. max(data(varn,:))) then
1106                  if (min(data(varn,:)) .EQ. 0)then
1107                     res@trYMaxF = min(data(varn,:))-1.
1108                     res@trYMinF = max(data(varn,:))+1.
1109                  end if
1110                  if (min(data(varn,:)) .LT. 0)then
1111                     res@trYMaxF = min(data(varn,:))+(min(data(varn,:)))/2
1112                     res@trYMinF = max(data(varn,:))-(max(data(varn,:)))/2
1113                  end if
1114                  if (min(data(varn,:)) .GT. 0)then
1115                     res@trYMaxF = min(data(varn,:))-(min(data(varn,:)))/2
1116                     res@trYMinF = max(data(varn,:))+(max(data(varn,:)))/2
1117                  end if
1118               end if
1119               plot_ps(n) = gsn_csm_xy(wks_ps,t,data(varn,:),res) 
1120            end if
1121       
1122         else
1123           
1124            print("plot of " + vNam(varn))
1125           
1126            n=n+1
1127            res@xyLineColors            = (/237/)
1128            res@gsnLeftString       = vNam(varn)
1129            res@gsnRightString      = unit(varn)
1130            res@tiXAxisString        = " time [s] "
1131            res@tiYAxisString        = " "
1132            res@tiXAxisFontHeightF   = 0.07
1133            res@txFontHeightF        = 0.07
1134            res@tiYAxisFontHeightF   = 0.07
1135            res@trYMaxF        = min(data(varn,:))
1136            res@trYMinF        = max(data(varn,:))
1137            if (min(data(varn,:)) .EQ. max(data(varn,:))) then
1138               if (min(data(varn,:)) .EQ. 0)then
1139                     res@trYMaxF = min(data(varn,:))-1.
1140                     res@trYMinF = max(data(varn,:))+1.
1141                  end if
1142                  if (min(data(varn,:)) .LT. 0)then
1143                     res@trYMaxF = min(data(varn,:))+(min(data(varn,:)))/2
1144                     res@trYMinF = max(data(varn,:))-(max(data(varn,:)))/2
1145                  end if
1146                  if (min(data(varn,:)) .GT. 0)then
1147                     res@trYMaxF = min(data(varn,:))-(min(data(varn,:)))/2
1148                     res@trYMinF = max(data(varn,:))+(max(data(varn,:)))/2
1149                  end if
1150            end if
1151            plot_ps(n) = gsn_csm_xy(wks_ps,t,data(varn,:),res)
1152           
1153         end if     
1154      end if
1155   end do
1156 
1157   ; ***************************************************
1158   ; merge plots onto one page
1159   ; ***************************************************
1160 
1161   if (format_out .EQ. "eps" .OR. format_out .EQ. "epsi") then
1162      gsn_panel(wks_ps,plot_ps(1:n),(/n,1/),resP)
1163   else
1164      do np = 1,n,no_lines*no_columns   
1165         if ( np + no_lines*no_columns .gt. n) then   
1166            gsn_panel(wks_ps, plot_ps(np:n),(/no_lines,no_columns/),resP)
1167         else
1168            gsn_panel(wks_ps, plot_ps(np:np+no_lines*no_columns-1),(/no_lines,no_columns/),resP)
1169         end if
1170      end do
1171   end if
1172
1173   print(" ")
1174   print("Output to: " + file_out +"."+ format_out)
1175   print(" ")
1176 
1177end
Note: See TracBrowser for help on using the repository browser.