Supplementary material 1/6¶
This notebook serves as supplementary material to the manuscript on low-level jets submitted to Wind Energy Science. The supplementary material consists of 6 notebooks:
- A very quick overview of the measurement data from the LiDARs.
- Overview of the ERA5 data and especially the procedure of aligning it with the observations.
- Some error statistics on wind speed from ERA5 and observations
- Main results on low-level jets
- Reconstruction of seasonal cycle using a machine learning package
- Spatial figures from ERA5 (data-intensive).
In this notebook (1/6), we provide a brief introduction into the observation data.
Getting up to speed with the observation data¶
The observation data come from 7 sites. LiDARs were installed at each site, and three sites were equipped with two LiDARs. These are the abbreviations.
- BWF - Borssele Wind Farm (Zuid)
- ELP - Europlatform
- HKN - Hollandse Kust Noord
- HKZ - Hollandse Kust Zuid
- K13 - K13a oil platform
- LEG - Lichteiland Goeree
- MMIJ - Met Mast IJmuiden
A website dedicated to these measurement is available at https://windopzee.net/en/home/. The data are stored as NetCDF format, 1 file per platform:
ls ../../data/external/ECN/*.nc
Inspecting one of the files¶
The postprocessing of the platform data is described in the appendix of the paper. The resulting output files all have the following format.
import xarray as xr
ds = xr.open_dataset('../../data/external/ECN/BWFZwinds_lot2.nc')
ds
Some notes/log¶
The output files were produced using matlab, while the further analysis in done in python. To achieve good compatibility, I wrote a preprocessing function, which addresses some remaining inconsistencies.
def clean(ds):
""" This function modifies the platform datasets on load to addresses the following:
- Convert the units of the time arrays in order to read them properly
- Change variable names (remove spaces etc.)
- Discard UTC time and model reference time (not necessary).
- Set time and height variables as coordinates rather than variables
- Re-index time array to have uniform spacing (convert gaps to NaN values)
"""
import xarray as xr
from pandas import date_range
# Convert units of time arrays (convert capitals to lower case)
ds['Model Reference Time'].attrs['units'] = ds['Model Reference Time'].attrs['units'].lower()
ds['File Reference Time'].attrs['units'] = ds['File Reference Time'].attrs['units'].lower()
# Rename variables
ds.rename({
"File Reference Time": "time",
"Measurement Levels": "height",
"Site Latitude": "latitude",
"Site Longitude": "longitude",
"WS": "wspd",
"WD": "wdir"}, inplace=True)
# Discard redundant time arrays
ds = ds.drop(['Model Reference Time', 'UTC Time'])
# This last function correctly parses the time arrays and makes the time equidistant (filling up with nan)
ds = xr.decode_cf(ds) #.resample(time='10min').asfreq() # Equivalent function, but slower (I think)
return ds.reindex(time=date_range(ds.time.values[0], ds.time.values[-1], freq='10min'))
ds = xr.open_dataset('../../data/external/ECN/BWFZwinds_lot1.nc')
clean(ds)
Load all stations at once¶
We now load all the datasets into one dictionary. This code is copied to the code base for later re-use (code/datasets/ecn.py)
data_path = '../../data/external/ECN/'
files = {
'BWF1': data_path+'BWFZwinds_lot1.nc',
'BWF2': data_path+'BWFZwinds_lot2.nc',
'EPL': data_path+'EPLwinds.nc',
'HKNA': data_path+'HKNAwinds.nc',
'HKNB': data_path+'HKNBwinds.nc',
'HKZA': data_path+'HKZAwinds.nc',
'HKZB': data_path+'HKZBwinds.nc',
'K13': data_path+'K13winds.nc',
'LEG': data_path+'LEGwinds.nc',
'MMIJ': data_path+'MMIJwinds.nc'}
# Save the platform data in a dictionary
obs = {}
for name, file in files.items():
print(name)
obs[name] = clean(xr.open_dataset(file))
Quickly inspecting the data¶
import matplotlib.pyplot as plt
plt.xkcd()
fig, ax = plt.subplots()
for name, ds in obs.items():
ax.plot(ds.wspd.mean(dim='time').values, ds.height.values, 'o-', label=name)
ax.legend(bbox_to_anchor=(1.05, 1))
ax.set_xlabel('mean wind speed (m/s)')
ax.set_ylabel('altitude (m)')
ax.set_title('Vertical wind profiles for each of the LiDARs')
plt.show()
Notes¶
- LEG shows a lot of shear in the upper layer
- Large variability between different platforms
- Some platforms observe up to higher altitudes than others.
fig, ax = plt.subplots()
for name, ds in obs.items():
ax.plot(range(24), ds.wspd.mean(dim='height').groupby('time.hour').mean(dim='time'), label=name)
ax.legend(bbox_to_anchor=(1.05, 1))
ax.set_ylabel('height-/time-averaged wind speed (m/s)')
ax.set_xlabel('hour of the day')
ax.set_title('diurnal cycle of wind speed')
plt.show()
Notes¶
- K13, LEG and MMIJ are further offshore.
- Consequently, these have larger wind speeds
- Borssele have a stronger diurnal cycle --> more impact of the surrounding land areas.
Barcode plots¶
The figure below shows the time/height evolution of wind speed for each dataset. It illustrates the temporal extent and overlap of the datasets, as well as the missing data.
fig, axs = plt.subplots(11, 1, figsize=(12, 16), sharex=True, sharey=True)
axs = axs.flatten()
for i, (name, ds) in enumerate(obs.items()):
mesh = axs[i].pcolormesh(ds.time, ds.height, ds.wspd, vmin=0, vmax=30)
axs[i].set_title(name)
plt.tight_layout()
axs[-2].xaxis.set_tick_params(labelbottom=True)
cax = fig.add_axes(axs[-1].get_position())
axs[-1].set_axis_off()
plt.colorbar(mesh, cax=cax, orientation='horizontal', label='wind speed (m/s)')
plt.show()
Notes¶
- MMIJ is by far the longest measurement record.
- HKN does not cover a complete season
- BWF2 very poor time extent
Wind directions¶
For wind direction, I want a cyclic colormap The code below is not relevant for understanding the data, but is necessary for the subsequent visualization.
# For wind direction, I want a cyclic colormap
# Use function from great answer: https://stackoverflow.com/a/34557535/6012085
import numpy as np
import matplotlib.colors as col
!pip install hsluv
import hsluv
##### generate custom colormaps
def make_segmented_cmap():
white = '#ffffff'
black = '#000000'
red = '#ff0000'
blue = '#0000ff'
anglemap = col.LinearSegmentedColormap.from_list(
'anglemap', [black, red, white, blue, black], N=256, gamma=1)
return anglemap
def make_anglemap( N = 256, use_hpl = True ):
h = np.ones(N) # hue
h[:N//2] = 11.6 # red
h[N//2:] = 258.6 # blue
s = 100 # saturation
l = np.linspace(0, 100, N//2) # luminosity
l = np.hstack( (l,l[::-1] ) )
colorlist = np.zeros((N,3))
for ii in range(N):
if use_hpl:
colorlist[ii,:] = hsluv.hpluv_to_rgb( (h[ii], s, l[ii]) )
else:
colorlist[ii,:] = hsluv.hsluv_to_rgb( (h[ii], s, l[ii]) )
colorlist[colorlist > 1] = 1 # correct numeric errors
colorlist[colorlist < 0] = 0
return col.ListedColormap( colorlist )
N = 256
hpluv_anglemap = make_anglemap( use_hpl = True )
Wind direction bar codes¶
The plot below is similar to the previous one for wind speed, but now for wind direction. This is not used in the paper so it is less relevant. In this case, the time axes are not shared between the axis, so it is in fact a 'zoomed in' version of the previous plot.
fig, axs = plt.subplots(11, 1, figsize=(12, 16), sharex=False, sharey=False)
axs = axs.flatten()
for i, (name, ds) in enumerate(obs.items()):
mesh = axs[i].pcolormesh(ds.time, ds.height, ds.wdir, vmin=0, vmax=360, cmap=hpluv_anglemap)
axs[i].set_title(name)
plt.tight_layout()
axs[-2].xaxis.set_tick_params(labelbottom=True)
cax = fig.add_axes(axs[-1].get_position())
axs[-1].set_axis_off()
plt.colorbar(mesh, cax=cax, orientation='horizontal', label='wind direction')
plt.show()
Notes¶
- Collocated datasets show similar time evoluation as expected.
- Borssele has more North-Easterly winds, because other directions have been masked due to wind farm wake effects.
Save cleaned data in intermediate files¶
for name, ds in obs.items():
ds.to_netcdf('../../data/interim/ECNplatforms/'+name+'.nc')
Conclusion¶
So far so good! Some of the scripts used in this notebook have been refactored to a shared code base in order to load the data more quickly in other notebooks.