|
forrest@0
|
1 |
; ***********************************************
|
|
forrest@0
|
2 |
; combine all scatter
|
|
forrest@0
|
3 |
; ***********************************************
|
|
forrest@0
|
4 |
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
|
|
forrest@0
|
5 |
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
|
|
forrest@0
|
6 |
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
|
|
forrest@0
|
7 |
;************************************************
|
|
forrest@0
|
8 |
begin
|
|
forrest@0
|
9 |
|
|
forrest@0
|
10 |
plot_type = "ps"
|
|
forrest@0
|
11 |
plot_type_new = "png"
|
|
forrest@0
|
12 |
|
|
forrest@0
|
13 |
;************************************************
|
|
forrest@0
|
14 |
; read in ob data
|
|
forrest@0
|
15 |
;************************************************
|
|
forrest@0
|
16 |
|
|
forrest@0
|
17 |
;(A) plot ob scatter_81
|
|
forrest@0
|
18 |
|
|
forrest@0
|
19 |
dir81 = "/fis/cgd/cseg/people/jeff/clamp_data/npp/ob/"
|
|
forrest@0
|
20 |
fil81 = "data.81.nc"
|
|
forrest@0
|
21 |
f81 = addfile (dir81+fil81,"r")
|
|
forrest@0
|
22 |
|
|
forrest@0
|
23 |
x81 = f81->SITE_ID
|
|
forrest@0
|
24 |
y81 = f81->TNPP_C
|
|
forrest@0
|
25 |
xo81 = f81->LONG_DD
|
|
forrest@0
|
26 |
yo81 = f81->LAT_DD
|
|
forrest@0
|
27 |
|
|
forrest@0
|
28 |
x81@long_name = "SITE_ID"
|
|
forrest@0
|
29 |
y81@long_name = "NPP (gC/m2/year)"
|
|
forrest@0
|
30 |
|
|
forrest@0
|
31 |
plot_name = "scatter_ob_81"
|
|
forrest@0
|
32 |
title = "Observed 81 sites"
|
|
forrest@0
|
33 |
|
|
forrest@0
|
34 |
wks = gsn_open_wks (plot_type,plot_name) ; open workstation
|
|
forrest@0
|
35 |
res = True ; plot mods desired
|
|
forrest@0
|
36 |
res@tiMainString = title ; add title
|
|
forrest@0
|
37 |
res@xyMarkLineModes = "Markers" ; choose which have markers
|
|
forrest@0
|
38 |
res@xyMarkers = 16 ; choose type of marker
|
|
forrest@0
|
39 |
res@xyMarkerColor = "red" ; Marker color
|
|
forrest@0
|
40 |
res@xyMarkerSizeF = 0.01 ; Marker size (default 0.01)
|
|
forrest@0
|
41 |
res@tmLabelAutoStride = True ; nice tick mark labels
|
|
forrest@0
|
42 |
|
|
forrest@0
|
43 |
plot = gsn_csm_xy (wks,x81,y81,res) ; create plot
|
|
forrest@0
|
44 |
frame(wks)
|
|
forrest@0
|
45 |
|
|
forrest@0
|
46 |
system("convert "+plot_name+"."+plot_type+" "+plot_name+"."+plot_type_new)
|
|
forrest@0
|
47 |
system("rm "+plot_name+"."+plot_type)
|
|
forrest@0
|
48 |
system("rm "+plot_name+"-1."+plot_type_new)
|
|
forrest@0
|
49 |
system("mv "+plot_name+"-0."+plot_type_new+" "+plot_name+"."+plot_type_new)
|
|
forrest@0
|
50 |
|
|
forrest@0
|
51 |
clear (wks)
|
|
forrest@0
|
52 |
|
|
forrest@0
|
53 |
;(B) plot ob scatter_933
|
|
forrest@0
|
54 |
|
|
forrest@0
|
55 |
dir933 = "/fis/cgd/cseg/people/jeff/clamp_data/npp/ob/"
|
|
forrest@0
|
56 |
fil933 = "data.933.nc"
|
|
forrest@0
|
57 |
f933 = addfile (dir933+fil933,"r")
|
|
forrest@0
|
58 |
|
|
forrest@0
|
59 |
x933 = f933->SITE_ID
|
|
forrest@0
|
60 |
y933 = f933->TNPP_C
|
|
forrest@0
|
61 |
xo933 = f933->LONG_DD
|
|
forrest@0
|
62 |
yo933 = f933->LAT_DD
|
|
forrest@0
|
63 |
|
|
forrest@0
|
64 |
x933@long_name = "SITE_ID"
|
|
forrest@0
|
65 |
y933@long_name = "NPP (gC/m2/year)"
|
|
forrest@0
|
66 |
|
|
forrest@0
|
67 |
plot_name = "scatter_ob_933"
|
|
forrest@0
|
68 |
title = "Observed 933 sites"
|
|
forrest@0
|
69 |
|
|
forrest@0
|
70 |
wks = gsn_open_wks (plot_type,plot_name) ; open workstation
|
|
forrest@0
|
71 |
res = True ; plot mods desired
|
|
forrest@0
|
72 |
res@tiMainString = title ; add title
|
|
forrest@0
|
73 |
res@xyMarkLineModes = "Markers" ; choose which have markers
|
|
forrest@0
|
74 |
res@xyMarkers = 16 ; choose type of marker
|
|
forrest@0
|
75 |
res@xyMarkerColor = "red" ; Marker color
|
|
forrest@0
|
76 |
res@xyMarkerSizeF = 0.01 ; Marker size (default 0.01)
|
|
forrest@0
|
77 |
res@tmLabelAutoStride = True ; nice tick mark labels
|
|
forrest@0
|
78 |
|
|
forrest@0
|
79 |
plot = gsn_csm_xy (wks,x933,y933,res) ; create plot
|
|
forrest@0
|
80 |
frame(wks)
|
|
forrest@0
|
81 |
|
|
forrest@0
|
82 |
system("convert "+plot_name+"."+plot_type+" "+plot_name+"."+plot_type_new)
|
|
forrest@0
|
83 |
system("rm "+plot_name+"."+plot_type)
|
|
forrest@0
|
84 |
system("rm "+plot_name+"-1."+plot_type_new)
|
|
forrest@0
|
85 |
system("mv "+plot_name+"-0."+plot_type_new+" "+plot_name+"."+plot_type_new)
|
|
forrest@0
|
86 |
|
|
forrest@0
|
87 |
clear (wks)
|
|
forrest@0
|
88 |
;************************************************
|
|
forrest@0
|
89 |
; read in model data
|
|
forrest@0
|
90 |
;************************************************
|
|
forrest@0
|
91 |
|
|
forrest@0
|
92 |
;(C) plot model scatter_81
|
|
forrest@0
|
93 |
|
|
forrest@0
|
94 |
model_name = "b30.061n"
|
|
forrest@0
|
95 |
|
|
forrest@0
|
96 |
dirm = "/fis/cgd/cseg/people/jeff/clamp_data/model/"
|
|
forrest@0
|
97 |
film = "b30.061n_1995-2004_ANN_climo_lnd.nc"
|
|
forrest@0
|
98 |
fm = addfile (dirm+film,"r")
|
|
forrest@0
|
99 |
|
|
forrest@0
|
100 |
ymod = fm->NPP
|
|
forrest@0
|
101 |
xm = fm->lon
|
|
forrest@0
|
102 |
ym = fm->lat
|
|
forrest@0
|
103 |
|
|
forrest@0
|
104 |
nx = dimsizes(xo81)
|
|
forrest@0
|
105 |
do i= 0,nx-1
|
|
forrest@0
|
106 |
if (xo81(i) .lt. 0.) then
|
|
forrest@0
|
107 |
xo81(i) = xo81(i)+ 360.
|
|
forrest@0
|
108 |
end if
|
|
forrest@0
|
109 |
end do
|
|
forrest@0
|
110 |
print (xo81)
|
|
forrest@0
|
111 |
|
|
forrest@0
|
112 |
sec_to_year = 86400.*365.
|
|
forrest@0
|
113 |
|
|
forrest@0
|
114 |
ymod81 = linint2_points(xm,ym,ymod(0,:,:),True,xo81,yo81,0) * sec_to_year
|
|
forrest@0
|
115 |
xmod81 = x81
|
|
forrest@0
|
116 |
|
|
forrest@0
|
117 |
ymod81@long_name = "NPP (gC/m2/year)"
|
|
forrest@0
|
118 |
xmod81@long_name = "SITE_ID"
|
|
forrest@0
|
119 |
print (ymod81)
|
|
forrest@0
|
120 |
print (xmod81)
|
|
forrest@0
|
121 |
|
|
forrest@0
|
122 |
plot_name = "scatter_model_81"
|
|
forrest@0
|
123 |
title = "Model "+ model_name +" 81 sites"
|
|
forrest@0
|
124 |
|
|
forrest@0
|
125 |
wks = gsn_open_wks (plot_type,plot_name) ; open workstation
|
|
forrest@0
|
126 |
res = True ; plot mods desired
|
|
forrest@0
|
127 |
res@tiMainString = title ; add title
|
|
forrest@0
|
128 |
res@xyMarkLineModes = "Markers" ; choose which have markers
|
|
forrest@0
|
129 |
res@xyMarkers = 16 ; choose type of marker
|
|
forrest@0
|
130 |
res@xyMarkerColor = "red" ; Marker color
|
|
forrest@0
|
131 |
res@xyMarkerSizeF = 0.01 ; Marker size (default 0.01)
|
|
forrest@0
|
132 |
res@tmLabelAutoStride = True ; nice tick mark labels
|
|
forrest@0
|
133 |
|
|
forrest@0
|
134 |
plot = gsn_csm_xy (wks,xmod81,ymod81,res) ; create plot
|
|
forrest@0
|
135 |
frame(wks)
|
|
forrest@0
|
136 |
|
|
forrest@0
|
137 |
system("convert "+plot_name+"."+plot_type+" "+plot_name+"."+plot_type_new)
|
|
forrest@0
|
138 |
system("rm "+plot_name+"."+plot_type)
|
|
forrest@0
|
139 |
system("rm "+plot_name+"-1."+plot_type_new)
|
|
forrest@0
|
140 |
system("mv "+plot_name+"-0."+plot_type_new+" "+plot_name+"."+plot_type_new)
|
|
forrest@0
|
141 |
|
|
forrest@0
|
142 |
clear (wks)
|
|
forrest@0
|
143 |
|
|
forrest@0
|
144 |
;(D) plot model scatter_933
|
|
forrest@0
|
145 |
|
|
forrest@0
|
146 |
model_name = "b30.061n"
|
|
forrest@0
|
147 |
|
|
forrest@0
|
148 |
dirm = "/fis/cgd/cseg/people/jeff/clamp_data/model/"
|
|
forrest@0
|
149 |
film = "b30.061n_1995-2004_ANN_climo_lnd.nc"
|
|
forrest@0
|
150 |
fm = addfile (dirm+film,"r")
|
|
forrest@0
|
151 |
|
|
forrest@0
|
152 |
ymod = fm->NPP
|
|
forrest@0
|
153 |
xm = fm->lon
|
|
forrest@0
|
154 |
ym = fm->lat
|
|
forrest@0
|
155 |
|
|
forrest@0
|
156 |
nx = dimsizes(xo933)
|
|
forrest@0
|
157 |
do i= 0,nx-1
|
|
forrest@0
|
158 |
if (xo933(i) .lt. 0.) then
|
|
forrest@0
|
159 |
xo933(i) = xo933(i)+ 360.
|
|
forrest@0
|
160 |
end if
|
|
forrest@0
|
161 |
end do
|
|
forrest@0
|
162 |
print (xo933)
|
|
forrest@0
|
163 |
|
|
forrest@0
|
164 |
sec_to_year = 86400.*365.
|
|
forrest@0
|
165 |
|
|
forrest@0
|
166 |
ymod933 = linint2_points(xm,ym,ymod(0,:,:),True,xo933,yo933,0) * sec_to_year
|
|
forrest@0
|
167 |
xmod933 = x933
|
|
forrest@0
|
168 |
|
|
forrest@0
|
169 |
ymod933@long_name = "NPP (gC/m2/year)"
|
|
forrest@0
|
170 |
xmod933@long_name = "SITE_ID"
|
|
forrest@0
|
171 |
;print (ymod933)
|
|
forrest@0
|
172 |
;print (xmod933)
|
|
forrest@0
|
173 |
|
|
forrest@0
|
174 |
plot_name = "scatter_model_933"
|
|
forrest@0
|
175 |
title = "Model "+ model_name +" 933 sites"
|
|
forrest@0
|
176 |
|
|
forrest@0
|
177 |
wks = gsn_open_wks (plot_type,plot_name) ; open workstation
|
|
forrest@0
|
178 |
res = True ; plot mods desired
|
|
forrest@0
|
179 |
res@tiMainString = title ; add title
|
|
forrest@0
|
180 |
res@xyMarkLineModes = "Markers" ; choose which have markers
|
|
forrest@0
|
181 |
res@xyMarkers = 16 ; choose type of marker
|
|
forrest@0
|
182 |
res@xyMarkerColor = "red" ; Marker color
|
|
forrest@0
|
183 |
res@xyMarkerSizeF = 0.01 ; Marker size (default 0.01)
|
|
forrest@0
|
184 |
res@tmLabelAutoStride = True ; nice tick mark labels
|
|
forrest@0
|
185 |
|
|
forrest@0
|
186 |
plot = gsn_csm_xy (wks,xmod933,ymod933,res) ; create plot
|
|
forrest@0
|
187 |
frame(wks)
|
|
forrest@0
|
188 |
|
|
forrest@0
|
189 |
system("convert "+plot_name+"."+plot_type+" "+plot_name+"."+plot_type_new)
|
|
forrest@0
|
190 |
system("rm "+plot_name+"."+plot_type)
|
|
forrest@0
|
191 |
system("rm "+plot_name+"-1."+plot_type_new)
|
|
forrest@0
|
192 |
system("mv "+plot_name+"-0."+plot_type_new+" "+plot_name+"."+plot_type_new)
|
|
forrest@0
|
193 |
|
|
forrest@0
|
194 |
clear (wks)
|
|
forrest@0
|
195 |
|
|
forrest@0
|
196 |
;(E) scatter model vs ob 81
|
|
forrest@0
|
197 |
|
|
forrest@0
|
198 |
ymod81@long_name = "NPP (gC/m2/year)"
|
|
forrest@0
|
199 |
y81@long_name = "NPP (gC/m2/year)"
|
|
forrest@0
|
200 |
;print (ymod81)
|
|
forrest@0
|
201 |
;print (y81)
|
|
forrest@0
|
202 |
|
|
forrest@0
|
203 |
plot_name = "scatter_model_vs_ob_81"
|
|
forrest@0
|
204 |
title = model_name +" vs ob 81 sites"
|
|
forrest@0
|
205 |
|
|
forrest@0
|
206 |
wks = gsn_open_wks (plot_type,plot_name) ; open workstation
|
|
forrest@0
|
207 |
res = True ; plot mods desired
|
|
forrest@0
|
208 |
res@tiMainString = title ; add title
|
|
forrest@0
|
209 |
res@xyMarkLineModes = "Markers" ; choose which have markers
|
|
forrest@0
|
210 |
res@xyMarkers = 16 ; choose type of marker
|
|
forrest@0
|
211 |
res@xyMarkerColor = "red" ; Marker color
|
|
forrest@0
|
212 |
res@xyMarkerSizeF = 0.01 ; Marker size (default 0.01)
|
|
forrest@0
|
213 |
res@tmLabelAutoStride = True ; nice tick mark labels
|
|
forrest@0
|
214 |
;-------------------------------
|
|
forrest@0
|
215 |
;compute correlation coef. and M
|
|
forrest@0
|
216 |
ccr81 = esccr(ymod81,y81,0)
|
|
forrest@0
|
217 |
print (ccr81)
|
|
forrest@0
|
218 |
|
|
forrest@0
|
219 |
;new eq
|
|
forrest@0
|
220 |
bias = sum(abs(ymod81-y81)/(ymod81+y81))
|
|
forrest@0
|
221 |
M81 = (1. - (bias/dimsizes(y81)))*5.
|
|
forrest@0
|
222 |
print (M81)
|
|
forrest@0
|
223 |
|
|
forrest@0
|
224 |
tRes = True
|
|
forrest@0
|
225 |
tRes@txFontHeightF = 0.025
|
|
forrest@0
|
226 |
|
|
forrest@0
|
227 |
correlation_text = "(correlation coef = "+sprintf("%5.2f", ccr81)+")"
|
|
forrest@0
|
228 |
|
|
forrest@0
|
229 |
gsn_text_ndc(wks,correlation_text,0.5,0.95,tRes)
|
|
forrest@0
|
230 |
;--------------------------------
|
|
forrest@0
|
231 |
plot = gsn_csm_xy (wks,y81,ymod81,res) ; create plot
|
|
forrest@0
|
232 |
frame(wks)
|
|
forrest@0
|
233 |
|
|
forrest@0
|
234 |
system("convert "+plot_name+"."+plot_type+" "+plot_name+"."+plot_type_new)
|
|
forrest@0
|
235 |
system("rm "+plot_name+"."+plot_type)
|
|
forrest@0
|
236 |
system("rm "+plot_name+"-1."+plot_type_new)
|
|
forrest@0
|
237 |
system("mv "+plot_name+"-0."+plot_type_new+" "+plot_name+"."+plot_type_new)
|
|
forrest@0
|
238 |
|
|
forrest@0
|
239 |
clear (wks)
|
|
forrest@0
|
240 |
|
|
forrest@0
|
241 |
;(F) scatter model vs ob 933
|
|
forrest@0
|
242 |
|
|
forrest@0
|
243 |
ymod933@long_name = "NPP (gC/m2/year)"
|
|
forrest@0
|
244 |
y933@long_name = "NPP (gC/m2/year)"
|
|
forrest@0
|
245 |
;print (ymod33)
|
|
forrest@0
|
246 |
;print (y933)
|
|
forrest@0
|
247 |
|
|
forrest@0
|
248 |
plot_name = "scatter_model_vs_ob_933"
|
|
forrest@0
|
249 |
title = model_name +" vs ob 933 sites"
|
|
forrest@0
|
250 |
|
|
forrest@0
|
251 |
wks = gsn_open_wks (plot_type,plot_name) ; open workstation
|
|
forrest@0
|
252 |
res = True ; plot mods desired
|
|
forrest@0
|
253 |
res@tiMainString = title ; add title
|
|
forrest@0
|
254 |
res@xyMarkLineModes = "Markers" ; choose which have markers
|
|
forrest@0
|
255 |
res@xyMarkers = 16 ; choose type of marker
|
|
forrest@0
|
256 |
res@xyMarkerColor = "red" ; Marker color
|
|
forrest@0
|
257 |
res@xyMarkerSizeF = 0.01 ; Marker size (default 0.01)
|
|
forrest@0
|
258 |
res@tmLabelAutoStride = True ; nice tick mark labels
|
|
forrest@0
|
259 |
;-------------------------------
|
|
forrest@0
|
260 |
;compute correlation coef. and M
|
|
forrest@0
|
261 |
ccr933 = esccr(ymod933,y933,0)
|
|
forrest@0
|
262 |
print (ccr933)
|
|
forrest@0
|
263 |
|
|
forrest@0
|
264 |
;new eq
|
|
forrest@0
|
265 |
bias = sum(abs(ymod933-y933)/(ymod933+y933))
|
|
forrest@0
|
266 |
M933 = (1. - (bias/dimsizes(y933)))*5.
|
|
forrest@0
|
267 |
print (M933)
|
|
forrest@0
|
268 |
|
|
forrest@0
|
269 |
tRes = True
|
|
forrest@0
|
270 |
tRes@txFontHeightF = 0.025
|
|
forrest@0
|
271 |
|
|
forrest@0
|
272 |
correlation_text = "(correlation coef = "+sprintf("%5.2f", ccr933)+")"
|
|
forrest@0
|
273 |
|
|
forrest@0
|
274 |
gsn_text_ndc(wks,correlation_text,0.5,0.95,tRes)
|
|
forrest@0
|
275 |
;--------------------------------
|
|
forrest@0
|
276 |
plot = gsn_csm_xy (wks,y933,ymod933,res) ; create plot
|
|
forrest@0
|
277 |
frame(wks)
|
|
forrest@0
|
278 |
|
|
forrest@0
|
279 |
system("convert "+plot_name+"."+plot_type+" "+plot_name+"."+plot_type_new)
|
|
forrest@0
|
280 |
system("rm "+plot_name+"."+plot_type)
|
|
forrest@0
|
281 |
system("rm "+plot_name+"-1."+plot_type_new)
|
|
forrest@0
|
282 |
system("mv "+plot_name+"-0."+plot_type_new+" "+plot_name+"."+plot_type_new)
|
|
forrest@0
|
283 |
|
|
forrest@0
|
284 |
clear (wks)
|
|
forrest@0
|
285 |
end
|