3
3
Calculating Seasonal Averages from Timeseries of Monthly Means
4
4
==============================================================
5
5
6
- Author: `Joe Hamman <http ://uw-hydro. github.io/current_member/joe_hamman/ >`_
6
+ Author: `Joe Hamman <https ://github.com/jhamman/ >`__
7
7
8
- The data for this example can be found in the `xray-data <https://github.com/xray/xray-data >`_ repository. This example is also available in an IPython Notebook that is available `here <https://github.com/xray/xray/tree/master/examples/xray_seasonal_means.ipynb >`_.
8
+ The data used for this example can be found in the
9
+ `xarray-data <https://github.com/pydata/xarray-data >`__ repository.
9
10
10
11
Suppose we have a netCDF or xray Dataset of monthly mean data and we
11
12
want to calculate the seasonal average. To do this properly, we need to
12
13
calculate the weighted average considering that each month has a
13
14
different number of days.
14
15
16
+ Suppose we have a netCDF or ``xarray.Dataset `` of monthly mean data and
17
+ we want to calculate the seasonal average. To do this properly, we need
18
+ to calculate the weighted average considering that each month has a
19
+ different number of days.
20
+
15
21
.. code :: python
16
22
17
23
% matplotlib inline
18
24
import numpy as np
19
25
import pandas as pd
20
- import xray
26
+ import xarray as xr
21
27
from netCDF4 import num2date
22
28
import matplotlib.pyplot as plt
23
29
24
30
print (" numpy version : " , np.__version__ )
25
- print (" pandas version : " , pd.version.version)
26
- print (" xray version : " , xray.version.version)
31
+ print (" pandas version : " , pd.__version__ )
32
+ print (" xarray version : " , xr.__version__ )
33
+
27
34
28
35
.. parsed-literal ::
29
36
30
- numpy version : 1.9.2
31
- pandas version : 0.16.2
32
- xray version : 0.5.1
37
+ numpy version : 1.11.1
38
+ pandas version : 0.18.1
39
+ xarray version : 0.8.2
33
40
34
41
35
42
Some calendar information so we can support any netCDF calendar.
@@ -90,21 +97,22 @@ Open the ``Dataset``
90
97
91
98
.. code :: python
92
99
93
- monthly_mean_file = ' RASM_example_data.nc'
94
- ds = xray.open_dataset(monthly_mean_file, decode_coords = False )
100
+ ds = xr.tutorial.load_dataset(' rasm' )
95
101
print (ds)
96
102
97
103
98
104
.. parsed-literal ::
99
105
100
- <xray .Dataset>
106
+ <xarray .Dataset>
101
107
Dimensions: (time: 36, x: 275, y: 205)
102
108
Coordinates:
103
109
* time (time) datetime64[ns] 1980-09-16T12:00:00 1980-10-17 ...
104
- * x (x) int64 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 ...
105
110
* y (y) int64 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 ...
111
+ * x (x) int64 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 ...
106
112
Data variables:
107
113
Tair (time, y, x) float64 nan nan nan nan nan nan nan nan nan nan ...
114
+ yc (y, x) float64 16.53 16.78 17.02 17.27 17.51 17.76 18.0 18.25 ...
115
+ xc (y, x) float64 189.2 189.4 189.6 189.7 189.9 190.1 190.2 190.4 ...
108
116
Attributes:
109
117
title: /workspace/jhamman/processed/R1002RBRxaaa01a/lnd/temp/R1002RBRxaaa01a.vic.ha.1979-09-01.nc
110
118
institution: U.W.
@@ -132,9 +140,8 @@ allong the time dimension.
132
140
.. code :: python
133
141
134
142
# Make a DataArray with the number of days in each month, size = len(time)
135
- month_length = xray.DataArray(get_dpm(ds.time.to_index(),
136
- calendar = ' noleap' ),
137
- coords = [ds.time], name = ' month_length' )
143
+ month_length = xr.DataArray(get_dpm(ds.time.to_index(), calendar = ' noleap' ),
144
+ coords = [ds.time], name = ' month_length' )
138
145
139
146
# Calculate the weights by grouping by 'time.season'.
140
147
# Conversion to float type ('astype(float)') only necessary for Python 2.x
@@ -153,14 +160,16 @@ allong the time dimension.
153
160
154
161
.. parsed-literal ::
155
162
156
- <xray .Dataset>
163
+ <xarray .Dataset>
157
164
Dimensions: (season: 4, x: 275, y: 205)
158
165
Coordinates:
159
166
* y (y) int64 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 ...
160
167
* x (x) int64 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 ...
161
168
* season (season) object 'DJF' 'JJA' 'MAM' 'SON'
162
169
Data variables:
163
170
Tair (season, y, x) float64 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ...
171
+ xc (season, y, x) float64 189.2 189.4 189.6 189.7 189.9 190.1 ...
172
+ yc (season, y, x) float64 16.53 16.78 17.02 17.27 17.51 17.76 18.0 ...
164
173
165
174
166
175
.. code :: python
@@ -172,30 +181,31 @@ allong the time dimension.
172
181
.. code :: python
173
182
174
183
# Quick plot to show the results
175
- is_null = np.isnan (ds_unweighted[' Tair' ][0 ].values )
184
+ notnull = pd.notnull (ds_unweighted[' Tair' ][0 ])
176
185
177
186
fig, axes = plt.subplots(nrows = 4 , ncols = 3 , figsize = (14 ,12 ))
178
187
for i, season in enumerate ((' DJF' , ' MAM' , ' JJA' , ' SON' )):
179
- plt.sca(axes[i, 0 ])
180
- plt.pcolormesh(np.ma.masked_where(is_null, ds_weighted[' Tair' ].sel(season = season).values),
181
- vmin = - 30 , vmax = 30 , cmap = ' Spectral_r' )
182
- plt.colorbar(extend = ' both' )
183
-
184
- plt.sca(axes[i, 1 ])
185
- plt.pcolormesh(np.ma.masked_where(is_null, ds_unweighted[' Tair' ].sel(season = season).values),
186
- vmin = - 30 , vmax = 30 , cmap = ' Spectral_r' )
187
- plt.colorbar(extend = ' both' )
188
-
189
- plt.sca(axes[i, 2 ])
190
- plt.pcolormesh(np.ma.masked_where(is_null, ds_diff[' Tair' ].sel(season = season).values),
191
- vmin = - 0.1 , vmax = .1 , cmap = ' RdBu_r' )
192
- plt.colorbar(extend = ' both' )
193
- for j in range (3 ):
194
- axes[i, j].axes.get_xaxis().set_ticklabels([])
195
- axes[i, j].axes.get_yaxis().set_ticklabels([])
196
- axes[i, j].axes.axis(' tight' )
188
+ ds_weighted[' Tair' ].sel(season = season).where(notnull).plot.pcolormesh(
189
+ ax = axes[i, 0 ], vmin = - 30 , vmax = 30 , cmap = ' Spectral_r' ,
190
+ add_colorbar = True , extend = ' both' )
191
+
192
+ ds_unweighted[' Tair' ].sel(season = season).where(notnull).plot.pcolormesh(
193
+ ax = axes[i, 1 ], vmin = - 30 , vmax = 30 , cmap = ' Spectral_r' ,
194
+ add_colorbar = True , extend = ' both' )
195
+
196
+ ds_diff[' Tair' ].sel(season = season).where(notnull).plot.pcolormesh(
197
+ ax = axes[i, 2 ], vmin = - 0.1 , vmax = .1 , cmap = ' RdBu_r' ,
198
+ add_colorbar = True , extend = ' both' )
197
199
198
200
axes[i, 0 ].set_ylabel(season)
201
+ axes[i, 1 ].set_ylabel(' ' )
202
+ axes[i, 2 ].set_ylabel(' ' )
203
+
204
+ for ax in axes.flat:
205
+ ax.axes.get_xaxis().set_ticklabels([])
206
+ ax.axes.get_yaxis().set_ticklabels([])
207
+ ax.axes.axis(' tight' )
208
+ ax.set_xlabel(' ' )
199
209
200
210
axes[0 , 0 ].set_title(' Weighted by DPM' )
201
211
axes[0 , 1 ].set_title(' Equal Weighting' )
@@ -206,6 +216,15 @@ allong the time dimension.
206
216
fig.suptitle(' Seasonal Surface Air Temperature' , fontsize = 16 , y = 1.02 )
207
217
208
218
219
+
220
+
221
+ .. parsed-literal ::
222
+
223
+ <matplotlib.text.Text at 0x117c18048>
224
+
225
+
226
+
227
+
209
228
.. image :: monthly_means_output.png
210
229
211
230
@@ -214,12 +233,12 @@ allong the time dimension.
214
233
# Wrap it into a simple function
215
234
def season_mean (ds , calendar = ' standard' ):
216
235
# Make a DataArray of season/year groups
217
- year_season = xray .DataArray(ds.time.to_index().to_period(freq = ' Q-NOV' ).to_timestamp(how = ' E' ),
218
- coords = [ds.time], name = ' year_season' )
236
+ year_season = xr .DataArray(ds.time.to_index().to_period(freq = ' Q-NOV' ).to_timestamp(how = ' E' ),
237
+ coords = [ds.time], name = ' year_season' )
219
238
220
239
# Make a DataArray with the number of days in each month, size = len(time)
221
- month_length = xray .DataArray(get_dpm(ds.time.to_index(), calendar = calendar),
222
- coords = [ds.time], name = ' month_length' )
240
+ month_length = xr .DataArray(get_dpm(ds.time.to_index(), calendar = calendar),
241
+ coords = [ds.time], name = ' month_length' )
223
242
# Calculate the weights by grouping by 'time.season'
224
243
weights = month_length.groupby(' time.season' ) / month_length.groupby(' time.season' ).sum()
225
244
0 commit comments