-
Notifications
You must be signed in to change notification settings - Fork 2
/
plot-sequence.ncl
142 lines (93 loc) · 3.83 KB
/
plot-sequence.ncl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
;;; NCL script for creating contour plots of 2D time-varying NARCCAP RCM data
;;; Creates a sequence of plots for multiple timesteps. Set variables
;;; 'start' and 'finish' on command line
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
;; Uncomment the following variables to hard-code them here, or pass values in via command-line. Example:
;; ncl -x varname=\"prtot\" infile=\"table6/prtot_CRCM_climatology.nc\" outfile=\"plots/prtot_CRCM_s0\" title=\"CRCM\" timestep=\"0\" narccap-generic-plot.ncl
; varname = "prtot"
; infile = "table6/prtot_CRCM_climatology.nc"
; outfile = "plots/prtot_CRCM_s0"
; title = "CRCM+NCEP, s0 seasonal ave, 1979-2004"
; timestep = 0 ;; defaults to last step in file
;; open file
fin = addfile(infile, "r")
;; read in coordinate variables first
lat = fin->lat
lon = fin->lon
time = fin->time
nx = dimsizes(fin->xc)
ny = dimsizes(fin->yc)
;;;;; default to last timestep in file for QC convenience
;;;if (.not. isvar("timestep")) then
;;; timestep = dimsizes(time) - 1
;;;end if
;;;
;;; ;; default title for convenience
;;; if (.not. isvar("title")) then
;;; date = ut_calendar(time(timestep),0)
;;; timestamp = ""+date(0,0)+"/"+date(0,1)+"/"+date(0,2)+" "+sprintf("%02.0f", date(0,3))+":"+sprintf("%02.0f", date(0,4))
;;; title = systemfunc("basename "+infile)+", "+timestamp
;;; end if
;; read in data.
data = fin->$varname$(start:finish,:,:)
;; hook up coordinate variables
data@lat2d = lat
data@lon2d = lon
;; read in map projection info
pname = data@grid_mapping
proj = fin->$pname$
;; set plotting parameters
res = True
;; map projection stuff.
if (lower_case(proj@grid_mapping_name) .eq. "lambert_conformal_conic") then
res@mpProjection = "LambertConformal"
res@mpLambertMeridianF = proj@longitude_of_central_meridian
res@mpLambertParallel1F = proj@standard_parallel(0)
res@mpLambertParallel2F = proj@standard_parallel(1)
end if
if (lower_case(proj@grid_mapping_name) .eq. "transverse_mercator") then
res@mpProjection = "Mercator"
res@mpCenterLatF = proj@latitude_of_projection_origin
res@mpCenterLonF = proj@longitude_of_central_meridian
end if
if (lower_case(pname) .eq. "polar_stereographic") then
res@mpProjection= "Stereographic"
res@mpCenterLonF= proj@straight_vertical_longitude_from_pole - 360
res@mpCenterLatF = proj@latitude_of_projection_origin
end if
if (lower_case(proj@grid_mapping_name) .eq. "rotated_latitude_longitude") then
res@mpProjection = "Mercator"
res@mpCenterLonF = proj@grid_north_pole_longitude - 180
res@mpCenterLatF = 90 - proj@grid_north_pole_latitude
end if
;; /map projection stuff
;; map boundaries
;; expand slighltly to be sure we plot things on the borders
res@mpLimitMode = "Corners"
res@mpLeftCornerLatF = lat(0,0) -1 ;; SW corner
res@mpLeftCornerLonF = lon(0,0) -1 ;; SW corner
res@mpRightCornerLatF = lat(yc|ny-1,xc|nx-1) + 1 ;; NE corner
res@mpRightCornerLonF = lon(yc|ny-1,xc|nx-1) + 1 ;; NE corner
;; general plotting stuff
res@mpFillOn = False
res@gsnMaximize = True
res@cnFillOn = True
res@cnLinesOn = False
res@cnFillMode = "RasterFill"
res@lbLabelAutoStride = True
res@lbBoxLinesOn = False
res@gsnSpreadColors = True
;; res@tiMainString = title
do timestep = start,finish
res@tiMainString = timestep+" - "+time(timestep)
;; create the plot
print(""+timestep)
wks = gsn_open_wks("ps", outfile+".t"+timestep)
gsn_define_colormap(wks,"nrl_sirkes") ;; 21 colors + fg/bg
plot = gsn_csm_contour_map(wks, data(timestep-start,:,:), res)
end do
exit
;; Copyright 2009-2012 Univ. Corp. for Atmos. Research
;; Author: Seth McGinnis, mcginnis@ucar.edu