source: trunk/UTIL/NCL/visu-utils.ncl @ 823

Last change on this file since 823 was 808, checked in by slebonnois, 12 years ago

SL: Many changes for VENUS (related to newstart) and TITAN (related to clouds). Please read DOC/chantiers/commit_importants.log (cf v808).

File size: 16.0 KB
Line 
1;;;;;;;;;;;; positionDims ;;;;;;;;;;;;;;;;
2
3undef("positionDims")
4function positionDims(var,nbdim,names)
5;============================
6
7; var   : variable
8; nbdim : number of dimensions for var
9; names : names of the 4 dimensions: time,levs,lat,lon
10
11begin
12
13  flag=(/10,10,10,10/)
14  do i=0,nbdim-1
15     if (var!i.eq.names(0)) then
16        flag(0)=i
17     end if
18     if (var!i.eq.names(1)) then
19        flag(1)=i
20     end if
21     if (var!i.eq.names(2)) then
22        flag(2)=i
23     end if
24     if (var!i.eq.names(3)) then
25        flag(3)=i
26     end if
27  end do
28
29  return(flag)
30
31end
32
33;;;;;;;;;;;; trimVar ;;;;;;;;;;;;;;;;
34
35undef("trimVar")
36function trimVar(var,nbdim,dimname,mindimval,maxdimval)
37;============================
38
39; var      : variable
40; nbdim    : number of dimensions for var
41; dimname  : names of the dimensions for this var
42; mindimval: min values for the dimensions
43; maxdimval: max values for the dimensions
44
45begin
46
47  if (nbdim.eq.1) then
48     vartrim=var({$dimname(0)$|mindimval(0):maxdimval(0)})
49  end if
50  if (nbdim.eq.2) then
51       vartrim=var({$dimname(0)$|mindimval(0):maxdimval(0)},{$dimname(1)$|mindimval(1):maxdimval(1)})
52  end if
53  if (nbdim.eq.3) then
54       vartrim=var({$dimname(0)$|mindimval(0):maxdimval(0)},{$dimname(1)$|mindimval(1):maxdimval(1)},{$dimname(2)$|mindimval(2):maxdimval(2)})
55  end if
56  if (nbdim.eq.4) then
57       vartrim=var({$dimname(0)$|mindimval(0):maxdimval(0)},{$dimname(1)$|mindimval(1):maxdimval(1)},{$dimname(2)$|mindimval(2):maxdimval(2)},{$dimname(3)$|mindimval(3):maxdimval(3)})
58;       vartrim=var($dimname(0)$|0:0,{$dimname(1)$|mindimval(1):maxdimval(1)},{$dimname(2)$|mindimval(2):maxdimval(2)},{$dimname(3)$|mindimval(3):maxdimval(3)})
59  end if
60
61  return(vartrim)
62
63end
64
65;;;;;;;;;;;; avgVar ;;;;;;;;;;;;;;;;
66
67undef("avgVar")
68function avgVar(vartrim,nbdim,dimavg)
69;============================
70
71; vartrim : variable
72; nbdim   : number of dimensions for var
73; dimavg  : logicals for averaging the dimensions
74
75begin
76
77  if ((nbdim.eq.4).and.(dimavg(3))) then
78       vartrim_3=dim_avg_n_Wrap(vartrim,3)
79  else
80       vartrim_3=vartrim
81  end if
82  if ((nbdim.ge.3).and.(dimavg(2))) then
83       vartrim_2=dim_avg_n_Wrap(vartrim_3,2)
84  else
85       vartrim_2=vartrim_3
86  end if
87  if ((nbdim.ge.2).and.(dimavg(1))) then
88       vartrim_1=dim_avg_n_Wrap(vartrim_2,1)
89  else
90       vartrim_1=vartrim_2
91  end if
92  if ((nbdim.ge.1).and.(dimavg(0))) then
93       varavg=dim_avg_n_Wrap(vartrim_1,0)
94  else
95       varavg=vartrim_1
96  end if
97
98  return(varavg)
99
100end
101
102;;;;;;;;;;;; deflimits_2d ;;;;;;;;;;;;;;;;
103
104undef("deflimits_2d")
105function deflimits_2d(var,res)
106;============================
107
108; var : 2D variable
109; res : ressource
110
111begin
112
113  dims=dimsizes(var)
114  dimname=getvardims(var)
115  res@trXMinF = tofloat(var&$dimname(1)$(0))
116  res@trXMaxF = tofloat(var&$dimname(1)$(dims(1)-1))
117  res@trYMinF = tofloat(var&$dimname(0)$(0))
118  res@trYMaxF = tofloat(var&$dimname(0)$(dims(0)-1))
119
120  return(res)
121
122end
123
124;;;;;;;;;;;; contours_2d ;;;;;;;;;;;;;;;;
125
126undef("contours_2d")
127function contours_2d(res,limits)
128;============================
129
130; res    : modified ressource for the plot
131; limits : min, max and values for contours (-888 => automatic)
132
133begin
134
135  if (limits(0).ne.-888) then
136     res@cnLevelSelectionMode = "ManualLevels"
137     res@cnMinLevelValF       = limits(0)
138     res@cnMaxLevelValF       = limits(1)
139     if (limits(2).ne.-888) then
140        res@cnLevelSpacingF   = limits(2)
141     end if
142  end if
143
144  return(res)
145
146end
147
148;;;;;;;;;;;; optionsOrtho ;;;;;;;;;;;;;;;;
149
150undef("optionsOrtho")
151function optionsOrtho(res,center,grid,spacing,gridcol)
152;============================
153
154; res    : modified ressource for the plot
155; center : lon,lat coordinates for center of projection
156; grid   : logical for printing grid
157; spacing: lon,lat spacing for grid
158; gridcol: color for grid
159
160begin
161
162  res@mpCenterLonF      = center(0)  ; choose center lon
163  res@mpCenterLatF      = center(1)  ; choose center lat
164  res@mpGridAndLimbOn   = grid
165  res@mpGridLonSpacingF = spacing(0)
166  res@mpGridLatSpacingF = spacing(1)
167  res@mpGridLineColor   = gridcol
168
169  return(res)
170
171end
172
173;;;;;;;;;;;; optionsStereo ;;;;;;;;;;;;;;;;
174
175undef("optionsStereo")
176function optionsStereo(res,hm,grid,spacing,gridcol)
177;============================
178
179; res    : modified ressource for the plot
180; hm     : hemisphere to be plotted
181; grid   : logical for printing grid
182; spacing: lon,lat spacing for grid
183; gridcol: color for grid
184
185begin
186
187  res@gsnPolar          = hm
188
189  res@mpGridAndLimbOn   = grid
190  res@mpGridLonSpacingF = spacing(0)
191  res@mpGridLatSpacingF = spacing(1)
192  res@mpGridLineColor   = gridcol
193
194  return(res)
195
196end
197
198;;;;;;;;;;;; optionsOverPlot_2d ;;;;;;;;;;;;;;;;
199
200undef("optionsOverPlot_2d")
201function optionsOverPlot_2d(res)
202;============================
203
204; res    : modified ressource for the plot
205
206begin
207
208  res@cnLineLabelFormat            = "@5.1f"
209  res@cnLineLabelsOn               = True
210  res@cnLabelMasking               = True               
211  res@cnLineLabelBackgroundColor   = -1
212  res@cnInfoLabelOn                = False
213  res@gsnContourNegLineDashPattern = 1 
214
215  return(res)
216
217end
218
219;;;;;;;;;;;; axePress_2d ;;;;;;;;;;;;;;;;
220
221undef("axePress_2d")
222function axePress_2d(res)
223;============================
224
225; res : modified ressource for the plot
226
227begin
228
229  res@tiYAxisString  = "Pressure"
230  res@trYReverse     = True
231  res@gsnYAxisIrregular2Log = True
232
233  return(res)
234
235end
236
237;;;;;;;;;;;; axePress_1d ;;;;;;;;;;;;;;;;
238
239undef("axePress_1d")
240function axePress_1d(res,axe,limits)
241;============================
242
243; res    : modified ressource for the plot
244; axe    : Y axis for plot (pressure)
245; limits : min and max values for variable axis
246
247begin
248
249  res@tiYAxisString  = "Pressure"
250  res@trYReverse     = True
251  res@xyYStyle       = "Log"
252  res@trYMinF         = axe(0)
253  res@trYMaxF         = axe(dimsizes(axe)-1)
254  if (limits(0).ne.-888) then
255    res@trXMinF       = limits(0)
256    res@trXMaxF       = limits(1)
257  end if
258
259  return(res)
260
261end
262
263
264undef("axeNoPress_1d")
265function axeNoPress_1d(res,axe,limits)
266;============================
267
268; res    : modified ressource for the plot
269; axe    : X axis for plot (not pressure)
270; limits : min and max values for variable axis
271
272begin
273
274  res@trXMinF         = axe(0)
275  res@trXMaxF         = axe(dimsizes(axe)-1)
276  if (limits(0).ne.-888) then
277    res@trYMinF       = limits(0)
278    res@trYMaxF       = limits(1)
279  end if
280
281  return(res)
282
283end
284
285;;;;;;;;;;;; axeTimeX_2d ;;;;;;;;;;;;;;;;
286
287undef("axeTimeX_2d")
288function axeTimeX_2d(res,axe)
289;============================
290
291; res : modified ressource for the plot
292; axe : X axis is time
293
294begin
295
296  res@tiXAxisString        = "Time (local days)"
297  res@sfXArray             = axe
298  res@tmXBLabelFontHeightF = 0.015
299  res@trXMinF              = axe(0)
300  res@trXMaxF              = axe(dimsizes(axe)-1)
301
302  return(res)
303
304end
305
306;;;;;;;;;;;; axeTimeXls_2d ;;;;;;;;;;;;;;;;
307
308undef("axeTimeXls_2d")
309function axeTimeXls_2d(res,axe)
310;============================
311
312; res : modified ressource for the plot
313; axe : X axis is Ls
314
315begin
316
317  res@tiXAxisString            = "Solar longitude"
318  res@sfXArray                 = axe
319  res@gsnXAxisIrregular2Linear = True
320  res@tmXBMode                 = "Manual"       
321  res@tmXBTickStartF           = 0
322  res@tmXBTickEndF             = 360
323  res@tmXBTickSpacingF         = 30
324  res@tmXBMinorPerMajor        = 2
325  res@tmXBLabelFontHeightF     = 0.015
326  res@trXMinF                  = axe(0)
327  res@trXMaxF                  = axe(dimsizes(axe)-1)
328
329  return(res)
330
331end
332
333;;;;;;;;;;;; axeTimeY_2d ;;;;;;;;;;;;;;;;
334
335undef("axeTimeY_2d")
336function axeTimeY_2d(res,axe)
337;============================
338
339; res : modified ressource for the plot
340; axe : Y axis is time
341
342begin
343
344  res@tiYAxisString        = "Time (local days)"
345  res@sfYArray             = axe
346  res@tmYLLabelFontHeightF = 0.015
347  res@trYMinF              = axe(0)
348  res@trYMaxF              = axe(dimsizes(axe)-1)
349
350  return(res)
351
352end
353
354;;;;;;;;;;;; axeTimeYls_2d ;;;;;;;;;;;;;;;;
355
356undef("axeTimeYls_2d")
357function axeTimeYls_2d(res,axe)
358;============================
359
360; res : modified ressource for the plot
361; axe : Y axis is Ls
362
363begin
364
365  res@tiYAxisString            = "Solar longitude"
366  res@sfYArray                 = axe
367  res@gsnYAxisIrregular2Linear = True
368  res@tmYLMode                 = "Manual"       
369  res@tmYLTickStartF           = 0
370  res@tmYLTickEndF             = 360
371  res@tmYLTickSpacingF         = 30
372  res@tmYLMinorPerMajor        = 2
373  res@tmYLLabelFontHeightF     = 0.015
374  res@trYMinF                  = axe(0)
375  res@trYMaxF                  = axe(dimsizes(axe)-1)
376
377  return(res)
378
379end
380
381;;;;;;;;;;;; axeTime_1d ;;;;;;;;;;;;;;;;
382
383undef("axeTime_1d")
384function axeTime_1d(res,axe)
385;============================
386
387; res : modified ressource for the plot
388; axe : X axis is time
389
390begin
391
392  res@tiXAxisString        = "Time (local days)"
393  res@tmXBLabelFontHeightF = 0.015
394
395  return(res)
396
397end
398
399;;;;;;;;;;;; axeTimels_1d ;;;;;;;;;;;;;;;;
400
401undef("axeTimels_1d")
402function axeTimels_1d(res,axe)
403;============================
404
405; res : modified ressource for the plot
406; axe : X axis is Ls
407
408begin
409
410  res@tiXAxisString            = "Solar longitude"
411  res@tmXBMode                 = "Manual"       
412  res@tmXBTickStartF           = 0
413  res@tmXBTickEndF             = 360
414  res@tmXBTickSpacingF         = 30
415  res@tmXBMinorPerMajor        = 2
416  res@tmXBLabelFontHeightF     = 0.015
417
418  return(res)
419
420end
421
422;;;;;;;;;;;; inversDim_2d ;;;;;;;;;;;;;;;;
423
424undef("inversDim_2d")
425function inversDim_2d(var)
426;============================
427
428; var is a 2d variable
429
430begin
431
432  dimname=getvardims(var)
433  tmpvar = var($dimname(1)$|:,$dimname(0)$|:)
434
435  return(tmpvar)
436
437end
438
439;;;;;;;;;;;; revertLon_2dmap ;;;;;;;;;;;;;;;;
440
441undef("revertLon_2dmap")
442function revertLon_2dmap(var)
443;============================
444
445; var is a 2d lat-lon variable
446
447begin
448
449  tmpvar  = var
450  dimname = getvardims(var)
451  Lname   = dimname(1)
452  nbL     = dimsizes(var&$Lname$)
453  tmpaxL  = new(nbL,float)
454
455;revert lon axis:
456  do i=0,nbL-1
457    tmpaxL(i)=var&$Lname$(nbL-1-i)
458  end do
459  tmpvar&$Lname$ = tmpaxL
460
461;revert values along lon axis:
462  do i=0,nbL-1
463    tmpvar(:,i) = var(:,nbL-1-i)
464  end do
465
466  return(tmpvar)
467
468end
469
470;;;;;;;;;;;; Map_2d ;;;;;;;;;;;;;;;;
471
472undef("Map_2d")
473function Map_2d(var,wks,res)
474;============================
475
476; var : 2D variable
477; wks : workstation
478; res : ressource
479
480begin
481
482; main plot characteristics
483
484  res@gsnSpreadColors = True
485  res@tmXTOn          = False
486  res@cnFillOn        = True
487  res@cnLineLabelsOn  = True
488  res@cnLabelMasking  = True
489  res@cnLineLabelBackgroundColor = "white"
490  res@lbOrientation   = "vertical"
491
492; Create plot
493
494  plotvar  = gsn_csm_contour(wks, var, res )
495
496  return(plotvar)
497
498end
499
500;;;;;;;;;;;; OPlot_2d ;;;;;;;;;;;;;;;;
501
502undef("OPlot_2d")
503function OPlot_2d(var,wks,res)
504;============================
505
506; var : 2D variable
507; wks : workstation
508; res : ressource
509
510begin
511
512; main plot characteristics
513
514  res@tmXTOn          = False
515  res@tmYROn          = False
516  res@cnFillOn        = False
517
518  res@gsnLeftString   = ""
519  res@gsnRightString  = ""
520
521; Create plot
522
523  plotvar  = gsn_csm_contour(wks, var, res )
524
525  return(plotvar)
526
527end
528
529;;;;;;;;;;;; Plot_1d ;;;;;;;;;;;;;;;;
530
531undef("Plot_1d")
532function Plot_1d(var,wks,res,axe,axep,limits)
533;============================
534
535; var            : 1D variable
536; wks            : workstation
537; res            : ressource
538; axe            : axis for plot
539; axep           : logical indication for pressure axis
540; limits         : min and max values for variable axis
541
542begin
543
544; main plot characteristics
545
546; Create plot
547; pressure axis ?
548  if (axep) then
549    res=axePress_1d(res,axe,limits)
550    plotvar  = gsn_csm_xy(wks, var, axe, res )
551  else
552    res1d=axeNoPress_1d(res,axe,limits)
553    plotvar  = gsn_csm_xy(wks, axe, var, res )
554  end if
555
556  return(plotvar)
557
558end
559
560;;;;;;;;;;;; Ortho_2d ;;;;;;;;;;;;;;;;
561
562undef("Ortho_2d")
563function Ortho_2d(var,wks,res)
564;============================
565
566; var : 2D variable
567; wks : workstation
568; res : ressource
569
570begin
571
572; main plot characteristics
573
574  res@gsnSpreadColors = True
575  res@tmXTOn          = False
576  res@cnFillOn        = True
577  res@cnLineLabelsOn  = True
578  res@cnLabelMasking  = True
579  res@cnLineLabelBackgroundColor = "white"
580  res@lbOrientation   = "vertical"
581
582  res@mpProjection      = "Orthographic"       ; choose projection
583  res@mpOutlineOn       = False
584  res@mpPerimOn         = False             ; turn off box around plot
585  res@mpFillOn          = False
586  res@mpGridLineDashPattern = 2
587 
588; Create plot
589
590  plotvar  = gsn_csm_contour_map(wks, var, res )
591
592  return(plotvar)
593
594end
595
596;;;;;;;;;;;; Stereo_2d ;;;;;;;;;;;;;;;;
597
598undef("Stereo_2d")
599function Stereo_2d(var,wks,res)
600;============================
601
602; var : 2D variable
603; wks : workstation
604; res : ressource
605
606begin
607
608; main plot characteristics
609
610  res@gsnSpreadColors = True
611  res@tmXTOn          = False
612  res@cnFillOn        = True
613  res@cnLineLabelsOn  = True
614  res@cnLabelMasking  = True
615  res@cnLineLabelBackgroundColor = "white"
616  res@lbOrientation   = "vertical"
617
618  res@mpProjection      = "Stereographic"       ; choose projection
619  res@mpOutlineOn       = False
620  res@mpPerimOn         = False             ; turn off box around plot
621  res@mpFillOn          = False
622 
623; Create plot
624
625  plotvar  = gsn_csm_contour_map_polar(wks, var, res )
626
627  return(plotvar)
628
629end
630
631;;;;;;;;;;;; customVar ;;;;;;;;;;;;;;;;
632
633undef("customVar")
634function customVar(infile,labelvar)
635;============================
636
637; infile   : input file
638; labelvar : input label for variable
639
640; THIS IS TO BE CUSTOMIZED EACH TIME ACCORDING TO YOUR NEEDS
641; because it is impossible to make it automatic...
642
643begin
644
645prepared = "dummy"    ; DEFAULT
646
647if (labelvar.eq."dTdyn") then
648  if (isfilevar(infile,"vitu")) then
649    dtdyn = infile->dtdyn
650  else
651    dtdyn = infile->DTCORE
652  end if
653  myvar = dtdyn*1e7
654  copy_VarMeta(dtdyn,myvar)
655  myvar@units    ="K/Vd"
656  myvar@long_name="Dynamics heating rate"
657  prepared = labelvar
658end if
659
660if (labelvar.eq."dTadj") then
661  if (isfilevar(infile,"vitu")) then
662    dtadj = infile->dtajs
663  else
664    dtadj = infile->DADJ
665  end if
666  myvar = dtadj*1e7
667  copy_VarMeta(dtadj,myvar)
668  myvar@units    ="K/Vd"
669  myvar@long_name="Dry convection heating rate"
670  prepared = labelvar
671end if
672
673if (labelvar.eq."dTvdf") then
674  if (isfilevar(infile,"vitu")) then
675    dtvdf = infile->dtvdf
676  else
677    dtvdf = infile->DTV
678  end if
679  myvar = dtvdf*1e7
680  copy_VarMeta(dtvdf,myvar)
681  myvar@units    ="K/Vd"
682  myvar@long_name="PBL heating rate"
683  prepared = labelvar
684end if
685
686if (labelvar.eq."dTrad") then
687  if (isfilevar(infile,"vitu")) then
688    dtlwr = infile->dtlwr
689    dtswr = infile->dtswr
690  else
691    dtlwr = infile->QRL
692    dtswr = infile->QRS
693  end if
694  myvar = (dtlwr+dtswr)*1e7
695  copy_VarMeta(dtlwr,myvar)
696  myvar@units    ="K/Vd"
697  myvar@long_name="Radiative balance"
698  prepared = labelvar
699end if
700
701if (labelvar.eq."Tbal") then
702  if (isfilevar(infile,"vitu")) then
703    dtlwr = infile->dtlwr
704    dtswr = infile->dtswr
705    dtajs = infile->dtajs
706    dtdyn = infile->dtdyn
707    dtvdf = infile->dtvdf
708  else
709    dtlwr = infile->QRL
710    dtswr = infile->QRS
711    dtajs = infile->DADJ
712    dtdyn = infile->DTCORE
713    dtvdf = infile->DTV
714  end if
715  myvar = (dtlwr+dtswr+dtajs+dtdyn+dtvdf)*1e7
716  copy_VarMeta(dtlwr,myvar)
717  myvar@units    ="K/Vd"
718  myvar@long_name="Thermal balance"
719  prepared = labelvar
720end if
721
722if (labelvar.eq."Tforc") then
723  if (isfilevar(infile,"vitu")) then
724    dtlwr = infile->dtlwr
725    dtswr = infile->dtswr
726    dtajs = infile->dtajs
727    if (isfilevar(infile,"dtvdf")) then
728          dtvdf = infile->dtvdf
729    else
730          dtvdf = dtajs*0.
731    end if
732  else
733    dtlwr = infile->QRL
734    dtswr = infile->QRS
735    dtajs = infile->DADJ
736    if (isfilevar(infile,"DTV")) then
737          dtvdf = infile->DTV
738    else
739          dtvdf = dtajs*0.
740    end if
741  end if
742  myvar = (dtlwr+dtswr+dtajs+dtvdf)*1e7
743  copy_VarMeta(dtlwr,myvar)
744  myvar@units    ="K/Vd"
745  myvar@long_name="Thermal physical forcing"
746  prepared = labelvar
747end if
748
749if (labelvar.eq."ske") then
750  if (isfilevar(infile,"vitu")) then
751    u = infile->vitu
752    v = infile->vitv
753  else
754    u = infile->U
755    v = infile->V
756  end if
757  myvar = u*u+v*v
758  copy_VarMeta(u,myvar)
759  myvar@units    ="m2/s2"
760  myvar@long_name="Specific kinetic energy"
761  prepared = labelvar
762end if
763
764if (prepared.eq."dummy") then   ; DEFAULT: LEAVE AS IS
765  print("You chose a customized variable !")
766  print("Modify the function customVar in visu-utils.ncl to prepare this variable")
767  exit
768end if
769
770  return(myvar)
771
772end
773
Note: See TracBrowser for help on using the repository browser.