In this notebook I'll show you how to use windng to perform a shear analysis on wind data and extrapolate wind speed beyond measurement height. First we need to import some stuff.
%pylab inline
import windng.io as wio
import windng.plot.windplot as wplot
import windng.weibull as weibull
import windng.wind as wind
Welcome to pylab, a matplotlib-based Python environment [backend: module://IPython.zmq.pylab.backend_inline]. For more information, type 'help(pylab)'.
I'll use a sample data file with one year of hourly data at 80 m and 50 m a.g.l.
with open('sample01.txt') as f:
for i in range(10):
print f.readline(),
date time s1a s1s s2a s2s d1a d2a 1998-10-01 0:00:00 10.3 2 9.5 1.8 290 294 1998-10-01 1:00:00 12.3 1.6 11.6 1.4 295 297 1998-10-01 2:00:00 10.1 1.5 9.3 1.4 292 294 1998-10-01 3:00:00 8.5 1.5 7.9 1.3 291 295 1998-10-01 4:00:00 9.3 1.3 8.7 1.2 291 292 1998-10-01 5:00:00 8.3 1.4 8.2 1.2 287 287 1998-10-01 6:00:00 9.9 1.7 9.2 1.4 292 294 1998-10-01 7:00:00 8.7 1.4 8.9 1.4 290 290 1998-10-01 8:00:00 10.2 1.7 9.9 1.6 294 296
The meaning of the column headers is:
s1a
/ s2a
-> Average speed at 80 m / 50 m
s1s
/ s2s
-> Standard deviation of speed at 80 m / 50 m
d1a
/ d2a
-> Average direction at 80 m / 50 m
Let's load the data into a pandas DataFrame and change the column names for clarity.
df = wio.read_wind('sample01.txt', freq='H', index_name=None)
df.columns = ['s80a', 's80s', 's50a', 's50s', 'd80a', 'd50a']
Now let's print some tables and graphics. This is not necessary for shear analysis but it's always good to know what's in your data.
df.head()
s80a | s80s | s50a | s50s | d80a | d50a | |
---|---|---|---|---|---|---|
1998-10-01 00:00:00 | 10.3 | 2.0 | 9.5 | 1.8 | 290 | 294 |
1998-10-01 01:00:00 | 12.3 | 1.6 | 11.6 | 1.4 | 295 | 297 |
1998-10-01 02:00:00 | 10.1 | 1.5 | 9.3 | 1.4 | 292 | 294 |
1998-10-01 03:00:00 | 8.5 | 1.5 | 7.9 | 1.3 | 291 | 295 |
1998-10-01 04:00:00 | 9.3 | 1.3 | 8.7 | 1.2 | 291 | 292 |
df.tail()
s80a | s80s | s50a | s50s | d80a | d50a | |
---|---|---|---|---|---|---|
1999-09-30 19:00:00 | 10.4 | 1.2 | 9.6 | 1.2 | 323 | 323 |
1999-09-30 20:00:00 | 10.3 | 1.3 | 9.6 | 1.3 | 341 | 342 |
1999-09-30 21:00:00 | 9.8 | 1.4 | 8.7 | 1.4 | 335 | 335 |
1999-09-30 22:00:00 | 3.9 | 0.7 | 3.6 | 0.7 | 342 | 343 |
1999-09-30 23:00:00 | 4.8 | 0.4 | 4.6 | 0.5 | 341 | 339 |
from IPython.core.display import HTML
display(HTML(df.describe().to_html()))
s80a | s80s | s50a | s50s | d80a | d50a | |
---|---|---|---|---|---|---|
count | 8760.000000 | 8760.000000 | 8760.000000 | 8760.000000 | 8760.000000 | 8760.000000 |
mean | 7.374726 | 0.879361 | 7.018299 | 0.886256 | 239.064612 | 238.924087 |
std | 3.545199 | 0.578978 | 3.375945 | 0.533688 | 89.723788 | 90.247493 |
min | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 |
25% | 4.700000 | 0.400000 | 4.500000 | 0.500000 | 164.000000 | 163.000000 |
50% | 6.900000 | 0.800000 | 6.600000 | 0.800000 | 275.000000 | 275.000000 |
75% | 9.600000 | 1.300000 | 9.100000 | 1.200000 | 316.000000 | 317.000000 |
max | 23.600000 | 3.800000 | 22.200000 | 3.800000 | 360.000000 | 360.000000 |
df['1-may-1999':'31-may-1999'][['s80a', 's50a']].plot(figsize=(14, 4))
plt.grid(True, which='both')
monthly80 = wind.downsample(df['s80a'], freq='M')
monthly50 = wind.downsample(df['s50a'], freq='M')
monthly80
mean | count | possible | rec_rate | |
---|---|---|---|---|
1998-10 | 7.471237 | 744 | 744 | 1 |
1998-11 | 7.956250 | 720 | 720 | 1 |
1998-12 | 7.252554 | 744 | 744 | 1 |
1999-01 | 8.721774 | 744 | 744 | 1 |
1999-02 | 8.612054 | 672 | 672 | 1 |
1999-03 | 8.145968 | 744 | 744 | 1 |
1999-04 | 7.568333 | 720 | 720 | 1 |
1999-05 | 7.238575 | 744 | 744 | 1 |
1999-06 | 6.074167 | 720 | 720 | 1 |
1999-07 | 6.368952 | 744 | 744 | 1 |
1999-08 | 6.544355 | 744 | 744 | 1 |
1999-09 | 6.620972 | 720 | 720 | 1 |
monthly50
mean | count | possible | rec_rate | |
---|---|---|---|---|
1998-10 | 7.066532 | 744 | 744 | 1 |
1998-11 | 7.485139 | 720 | 720 | 1 |
1998-12 | 6.832392 | 744 | 744 | 1 |
1999-01 | 8.177688 | 744 | 744 | 1 |
1999-02 | 8.128274 | 672 | 672 | 1 |
1999-03 | 7.854167 | 744 | 744 | 1 |
1999-04 | 7.230417 | 720 | 720 | 1 |
1999-05 | 6.934946 | 744 | 744 | 1 |
1999-06 | 5.812639 | 720 | 720 | 1 |
1999-07 | 6.151210 | 744 | 744 | 1 |
1999-08 | 6.279570 | 744 | 744 | 1 |
1999-09 | 6.335000 | 720 | 720 | 1 |
monthly80.describe()
mean | count | possible | rec_rate | |
---|---|---|---|---|
count | 12.000000 | 12.000000 | 12.000000 | 12 |
mean | 7.381266 | 730.000000 | 730.000000 | 1 |
std | 0.869637 | 21.608079 | 21.608079 | 0 |
min | 6.074167 | 672.000000 | 672.000000 | 1 |
25% | 6.601818 | 720.000000 | 720.000000 | 1 |
50% | 7.361895 | 744.000000 | 744.000000 | 1 |
75% | 8.003679 | 744.000000 | 744.000000 | 1 |
max | 8.721774 | 744.000000 | 744.000000 | 1 |
monthly50.describe()
mean | count | possible | rec_rate | |
---|---|---|---|---|
count | 12.000000 | 12.000000 | 12.000000 | 12 |
mean | 7.023998 | 730.000000 | 730.000000 | 1 |
std | 0.786208 | 21.608079 | 21.608079 | 0 |
min | 5.812639 | 672.000000 | 672.000000 | 1 |
25% | 6.321142 | 720.000000 | 720.000000 | 1 |
50% | 7.000739 | 744.000000 | 744.000000 | 1 |
75% | 7.577396 | 744.000000 | 744.000000 | 1 |
max | 8.177688 | 744.000000 | 744.000000 | 1 |
wplot.plot_means(monthly80['mean'], style='-o', label='80 m')
wplot.plot_means(monthly50['mean'], style='-o', label='50 m')
plt.legend(loc='best')
<matplotlib.legend.Legend at 0x3b1e150>
Now in order to compute a global shear profile the function shear_analysis()
can be used. It takes a sequence of speed series and a sequence of measurement heights. Then it computes the mean speed of each series over the common time steps and fits a power function of the form $s=\beta h^\alpha \space$, where $s$ is the speed and $h$ is the height. It returns a Series containing the exponent alpha
, the coefficient beta
, the number of common time steps
and the computed means at each height.
shear = wind.shear_analysis([df['s80a'], df['s50a']], [80, 50])
shear
alpha 0.105399 beta 4.646889 steps 8760.000000 80 7.374726 50 7.018299 dtype: float64
Use the function plot_shear()
to compare the fitted power function with the observed means.
wplot.plot_shear(shear)
([<matplotlib.lines.Line2D at 0x3f4f410>], [<matplotlib.lines.Line2D at 0x3f4f310>])
But a global shear is not enough. Wind shear and mean speed vary across directions so you usually want to estimate the wind shear for each direction sector. The function shear_analysis()
can do this too if you provide an additional dir
argument with directions.
dshear = wind.shear_analysis(series=(df['s80a'], df['s50a']), heights=(80, 50), dir=df['d80a'])
dshear
alpha | beta | steps | 80 | 50 | |
---|---|---|---|---|---|
0 | 0.161174 | 2.583762 | 355 | 5.235775 | 4.853803 |
1 | 0.078343 | 3.153772 | 112 | 4.445536 | 4.284821 |
2 | -0.137823 | 7.863582 | 73 | 4.298630 | 4.586301 |
3 | 0.003041 | 7.484595 | 120 | 7.585000 | 7.574167 |
4 | 0.198817 | 2.946425 | 181 | 7.041436 | 6.413260 |
5 | 0.160196 | 3.533805 | 379 | 7.130343 | 6.613193 |
6 | 0.089558 | 4.547162 | 538 | 6.732528 | 6.455019 |
7 | 0.066892 | 4.677869 | 847 | 6.271192 | 6.077096 |
8 | 0.081763 | 3.938632 | 706 | 5.635694 | 5.423229 |
9 | 0.143505 | 2.567121 | 428 | 4.814486 | 4.500467 |
10 | 0.173924 | 1.971267 | 377 | 4.224138 | 3.892573 |
11 | 0.159336 | 2.123808 | 224 | 4.269196 | 3.961161 |
12 | 0.106811 | 5.109761 | 414 | 8.159662 | 7.760145 |
13 | 0.085239 | 6.667894 | 981 | 9.687360 | 9.306932 |
14 | 0.096104 | 6.268679 | 1983 | 9.551488 | 9.129652 |
15 | 0.146228 | 3.894038 | 1042 | 7.390691 | 6.899808 |
For a polar plot of the power law exponents use plot_rose()
or bar_rose()
.
wplot.bar_rose(dshear['alpha'], ec='w')
wplot.bar_rose(-dshear['alpha'], fc='r', ec='w')
plt.title('Shear rose')
<matplotlib.text.Text at 0x4044e50>
The function plot_shear_by_sector()
creates a subplot for each direction sector.
wplot.plot_shear_by_sector(dshear, figsize=(12, 8))
In order to extrapolate your data from measurement height to another height you can apply the corresponding power law exponent to each speed value. This is what shear_extrapolation()
does. It takes speed, direction, measurement height, destination height and the exponents for each direction sector as arguments and returns the extrapolated speed.
df['s120a'] = wind.shear_extrapolation(df['s80a'], df['d80a'], 80, 120, dshear['alpha'])
df['1-may-1999':'31-may-1999'][['s120a', 's80a', 's50a']].plot(figsize=(14, 4))
plt.grid(True, which='both')
monthly120 = wind.downsample(df['s120a'], freq='M')
monthly120
mean | count | possible | rec_rate | |
---|---|---|---|---|
1998-10 | 7.807622 | 744 | 744 | 1 |
1998-11 | 8.323802 | 720 | 720 | 1 |
1998-12 | 7.572116 | 744 | 744 | 1 |
1999-01 | 9.125926 | 744 | 744 | 1 |
1999-02 | 9.028491 | 672 | 672 | 1 |
1999-03 | 8.466868 | 744 | 744 | 1 |
1999-04 | 7.910519 | 720 | 720 | 1 |
1999-05 | 7.540650 | 744 | 744 | 1 |
1999-06 | 6.337107 | 720 | 720 | 1 |
1999-07 | 6.633770 | 744 | 744 | 1 |
1999-08 | 6.822522 | 744 | 744 | 1 |
1999-09 | 6.901698 | 720 | 720 | 1 |
monthly120.describe()
mean | count | possible | rec_rate | |
---|---|---|---|---|
count | 12.000000 | 12.000000 | 12.000000 | 12 |
mean | 7.705924 | 730.000000 | 730.000000 | 1 |
std | 0.918208 | 21.608079 | 21.608079 | 0 |
min | 6.337107 | 672.000000 | 672.000000 | 1 |
25% | 6.881904 | 720.000000 | 720.000000 | 1 |
50% | 7.689869 | 744.000000 | 744.000000 | 1 |
75% | 8.359568 | 744.000000 | 744.000000 | 1 |
max | 9.125926 | 744.000000 | 744.000000 | 1 |
wplot.plot_means(monthly120['mean'], style='-o', label='120 m')
wplot.plot_means(monthly80['mean'], style='-o', label='80 m')
wplot.plot_means(monthly50['mean'], style='-o', label='50 m')
plt.legend(loc='best')
<matplotlib.legend.Legend at 0x5489b90>