improver.wind_calculations.wind_downscaling module

Module containing wind downscaling plugins.

class improver.wind_calculations.wind_downscaling.FrictionVelocity(u_href, h_ref, z_0, mask)[source]

Bases: improver.BasePlugin

“Class to calculate the friction velocity.

This holds the function to calculate the friction velocity u_star, given a reference height h_ref, the velocity at the reference height u_href and the surface roughness z_0.

__init__(u_href, h_ref, z_0, mask)[source]

Initialise the class.

Parameters
  • u_href (numpy.ndarray) – A 2D array of float32 for the wind speed at h_ref

  • h_ref (numpy.ndarray) – A 2D array of float32 for the reference heights

  • z_0 (numpy.ndarray) – A 2D array of float32 for the vegetative roughness lengths

  • mask (numpy.ndarray) – A 2D array of booleans where True indicates calculate u*

Notes

  • z_0 and h_ref need to have identical units.

  • the calculated friction velocity will have the units of the

    supplied velocity u_href.

_abc_cache = <_weakrefset.WeakSet object>
_abc_negative_cache = <_weakrefset.WeakSet object>
_abc_negative_cache_version = 213
_abc_registry = <_weakrefset.WeakSet object>
process()[source]

Function to calculate the friction velocity.

ustar = K * (u_href / ln(h_ref / z_0))

where ustar is the friction velocity, K is Von Karman’s constant, u_ref is the wind speed at the reference height, h_ref is the reference height and z_0 is the vegetative roughness length.

Returns

A 2D array of float32 friction velocities

Return type

numpy.ndarray

class improver.wind_calculations.wind_downscaling.RoughnessCorrection(a_over_s_cube, sigma_cube, pporo_cube, modoro_cube, modres, z0_cube=None, height_levels_cube=None)[source]

Bases: improver.BasePlugin

Plugin to orographically-correct 3d wind speeds.

__init__(a_over_s_cube, sigma_cube, pporo_cube, modoro_cube, modres, z0_cube=None, height_levels_cube=None)[source]

Initialise the RoughnessCorrection instance.

Parameters
  • a_over_s_cube (iris.cube.Cube) – 2D - model silhouette roughness on pp grid. dimensionless

  • sigma_cube (iris.cube.Cube) – 2D - standard deviation of model orography height on pp grid. In m.

  • pporo_cube (iris.cube.Cube) – 2D - pp orography. In m

  • modoro_cube (iris.cube.Cube) – 2D - model orography interpolated on pp grid. In m

  • modres (float) – original average model resolution in m

  • height_levels_cube (iris.cube.Cube) – 3D or 1D - height of input velocity field. Can be position dependent

  • z0_cube (iris.cube.Cube) – 2D - vegetative roughness length in m. If not given, do not do any RC

_abc_cache = <_weakrefset.WeakSet object>
_abc_negative_cache = <_weakrefset.WeakSet object>
_abc_negative_cache_version = 213
_abc_registry = <_weakrefset.WeakSet object>
calc_av_ppgrid_res(a_cube)[source]

Calculate average grid resolution from a cube.

Parameters

a_cube (iris.cube.Cube) – Cube to calculate average resolution of

Returns

Average grid resolution.

Return type

np.float32

static check_ancils(a_over_s_cube, sigma_cube, z0_cube, pp_oro_cube, model_oro_cube)[source]

Check ancils grid and units.

Check if ancil cubes are on the same grid and if they have the expected units. The testing for “same grid” might be replaced if there is a general utils function made for it or so.

Parameters
  • a_over_s_cube (iris.cube.Cube) – holding the silhouette roughness field

  • sigma_cube (iris.cube.Cube) – holding the standard deviation of height in a grid cell

  • z0_cube (iris.cube.Cube or None) – holding the vegetative roughness field

  • pp_oro_cube (iris.cube.Cube) – holding the post processing grid orography

  • model_oro_cube (iris.cube.Cube) – holding the model orography on post processing grid

Returns

Containing bools describing whether or not the tests passed

Return type

numpy.ndarray

check_wind_ancil(xwp, ywp)[source]

Check wind vs ancillary file grids.

Check if wind and ancillary files are on the same grid and if they have the same ordering.

Parameters
  • xwp (int) – representing the position of the x-axis in the wind cube

  • ywp (int) – representing the position of the y-axis of the wind cube

find_coord_names(cube)[source]

Extract x, y, z, and time coordinate names.

Parameters

cube (iris.cube.Cube) – some iris cube to find coordinate names from

Returns

tuple containing:
xname (str):

name of the axis name in x-direction

yname (str):

name of the axis name in y-direction

zname (str):

name of the axis name in z-direction

tname (str):

name of the axis name in t-direction

Return type

(tuple)

find_coord_order(mcube)[source]

Extract coordinate ordering within a cube.

Use coord_dims to assess the dimension associated with a particular dimension coordinate. If a coordinate is not a dimension coordinate, then a NaN value will be returned for that coordinate.

Parameters

mcube (iris.cube.Cube) – cube to check the order of coordinate axis

Returns

tuple containing:
xpos (int):

position of x axis.

ypos (int):

position of y axis.

zpos (int):

position of z axis.

tpos (int):

position of t axis.

Return type

(tuple)

find_heightgrid(wind)[source]

Setup the height grid.

Setup the height grid either from the 1D or 3D height grid that was supplied to the plugin or from the z-axis information from the wind grid.

Parameters

wind (iris.cube.Cube) – 3D or 4D - representing the wind data.

Returns

1D or 3D array - representing the height grid.

Return type

numpy.ndarray

process(input_cube)[source]

Adjust the 4d wind field - cube - (x, y, z including times).

Parameters

input_cube (iris.cube.Cube) – The wind cube to be operated upon. Should be wind speed on height_levels for all desired forecast times.

Returns

The 4d wind field with roughness and height correction applied in the same order as the input cube.

Return type

iris.cube.Cube

Raises

TypeError – If input_cube is not a cube.:

tcoordnames = ['time', 'forecast_time']
zcoordnames = ['height', 'model_level_number']
class improver.wind_calculations.wind_downscaling.RoughnessCorrectionUtilities(a_over_s, sigma, z_0, pporo, modoro, ppres, modres)[source]

Bases: object

Class to calculate the height and roughness wind corrections.

This holds functions to calculate the roughness and height corrections given the ancil files:

  • standard deviation of height in grid cell as sigma (model grid on pp grid)

  • Silhouette roughness as a_over_s (model grid on pp grid)

  • vegetative roughness length z_0 (model grid on pp grid)

  • post-processing grid orography pporo

  • model grid orography interpolated on post-processing grid modoro

  • height level 3D/ 1D grid

  • windspeed 3D field on height level 3D grid (from above).

__init__(a_over_s, sigma, z_0, pporo, modoro, ppres, modres)[source]

Set up roughness and height correction.

This sets up the parameters used for roughness and height correction given the ancillary file inputs:

Parameters
  • a_over_s (numpy.ndarray) – 2D array float32 - Silhouette roughness field, dimensionless ancillary data, calculated according to Robinson, D. (2008) - Ancillary file creation for the UM, Unified Model Documentation Paper 73.

  • sigma (numpy.ndarray) – 2D array float32 - Standard deviation field of height in the grid cell, units of length

  • z_0 (numpy.ndarray) – 2D array float32 - Vegetative roughness height field, units of length

  • pporo (numpy.ndarray) – 2D array float32 - Post processing grid orography field

  • modoro (numpy.ndarray) – 2D array float32 - Model orography field interpolated to post processing grid

  • ppres (float) – Float - Grid cell length of post processing grid

  • modres (float) – Float - Grid cell length of model grid

_calc_h_ref()[source]

Calculate the reference height for roughness correction.

The reference height below which the flow is in equilibrium with the vegetative roughness is proportional to 1/wavenum (Howard & Clark, 2007).

Vosper (2009) and Clark (2009) argue that at the reference height, the perturbation should have decayed to a fraction epsilon (ABSOLUTE_CORRECTION_TOL). The factor alpha implements eq. 1.3 in Clark (2009): UK Climatology - Wind Screening Tool. See also Vosper (2009) for a motivation. For a freely available external reference, see the Virtual Met Mast Version 1 Methodology and Verification paper under www.thecrownestate.co.uk.

alpha is the log of scale parameter to determine reference height which is currently set to 0.04 (this corresponds to epsilon in both Vosper and Clark)

Returns

2D array float32 - reference height for roughness correction

Return type

numpy.ndarray

_calc_height_corr(u_a, heightg, mask, onemfrac)[source]

Function to calculate the additive height correction.

Parameters
  • u_a (numpy.ndarray) – 2D array float32 - outer velocity, e.g. velocity at h_ref_orig

  • heightg (numpy.ndarray) – 1D or 3D array float32 - heights above orography

  • mask (numpy.ndarray) – 3D array of bools - Masks the hc_add result

  • onemfrac (float or numpy.ndarray) – Currently, scalar = 1. But can be a function of position and height, e.g. a 3D array (float32)

Returns

3D array float32 - additive height correction to wind speed

Return type

numpy.ndarray

Comments:

The height correction is a disturbance of the flow that decays exponentially with height. The larger the vertical offset h_at0 (the higher the unresolved hill), the larger is the disturbance.

The more smooth the disturbance (the larger the horizontal scale of the disturbance), the smaller the height correction (hence, a larger wavenumber results in a larger disturbance). hc_add = exp(-height*wavenumber)*u(href)*h_at_0*wavenumber

A final factor of 1 is assumed and omitted for the Bessel function term.

_calc_u_at_h(u_in, h_in, hhere, mask, dolog=False)[source]

Function to interpolate u_in on h_in at hhere.

Parameters
  • u_in (numpy.ndarray) – 3D array float32 - velocity on h_in layer, last dim is height

  • h_in (numpy.ndarray) – 3D or 1D array float32 - height layer array

  • hhere (numpy.ndarray) – 2D array float32 - height grid to interpolate at

  • mask (numpy.ndarray) – 2D array of bools - mask the final result for uath

  • dolog (bool) – if True, log interpolation, default False

Returns

2D array float32 - velocity interpolated at h

Return type

numpy.ndarray

_calc_wav()[source]

Calculate wavenumber k of typical orographic lengthscale.

Function to calculate wavenumber k of typical orographic lengthscale L:

(1)\[ k = 2 * \pi / L\]

L is approximated from half the peak-to-trough height h_over_2 and the silhoutte roughness a_over_s (average of up-slopes per unit length over several cross-sections through a grid cell) as:

(2)\[ L = 2 * \rm{h\_over\_2} / \rm{a\_over\_s}\]

a_over_s is dimensionless since it is the sum of up-slopes measured in the same unit lengths as it is calculated over.

h_over_2 is calculated from the standard deviation of height in a grid cell, sigma, as:

(3)\[ \rm{h\_over\_2} = \sqrt{2} * \rm{sigma}\]

which is based on the assumptions of sine waves, see sigma2hover2.

From eq. (1) and (2) it follows that:

(4)\[ k = 2*\pi / (2 * \rm{h\_over\_2} / \rm{a\_over\_s)} = \rm{a\_over\_s} * \pi / \rm{h\_over\_2}\]
Returns

2D array float32 - wavenumber in units of inverse units of supplied h_over_2.

Return type

numpy.ndarray

_delta_height()[source]

Function to calculate pp-grid diff from model grid.

Calculate the difference between pp-grid height and model grid height.

Returns

2D array float32 - height difference, ppgrid-model

Return type

numpy.ndarray

static _interpolate_1d(xup, xlow, at_x, yup, ylow)[source]

Simple 1D linear interpolation for 2D grid inputs level.

Parameters
Returns

2D array float32 - y(at_x) assuming a lin function between xlow and xup

Return type

numpy.ndarray

static _interpolate_log(xup, xlow, at_x, yup, ylow)[source]

Simple 1D log interpolation y(x), except if lowest layer is ground level.

Parameters
Returns

2D array float32 - y(at_x) assuming a log function between xlow and xup

Return type

numpy.ndarray

_refinemask()[source]

Remask over RMDI and NaN orography.

The mask for HC needs to be False where either of the orographies (model or pp) has an invalid number. This cannot be done before because the mask is used to calculate the wavenumber which can and should be calculated for all points where h_over_2 and a_over_s is a valid number.

_setmask()[source]

Create a ~land-sea mask.

Create a mask that is basically a land-sea mask: Both, the standard deviation and the silouette roughness, are 0 over the sea. A standard deviation of 0 results in a RMDI for h_over_2.

Returns

tuple containing:
hcmask (numpy.ndarray):

2D array of booleans- True for land-points, false for Sea (HC)

rcmask (numpy.ndarray):

2D array of booleans- additionally False for invalid z_0 (RC)

Return type

(tuple)

calc_roughness_correction(hgrid, uold, mask)[source]

Function to perform the roughness correction.

Parameters
  • hgrid (numpy.ndarray) – 3D or 1D array float32 - height above orography

  • uold (numpy.ndarray) – 3D array float32 - original velocities at hgrid.

  • mask (numpy.ndarray) – 2D array of bools that is True for land-points, False for Sea and False for invalid z_0.

Returns

3D np.array float32 - Corrected wind speed on hgrid. Above href, this is equal to uold.

Return type

numpy.ndarray

Comments:

Replace the windspeed profile below the reference height with one that increases logarithmically with height, bound by the original velocity uhref at the reference height h_ref and by a 0 velocity at the vegetative roughness height z_0

do_rc_hc_all(hgrid, uorig)[source]

Function to call HC and RC (height and roughness corrections).

Parameters
  • hgrid (numpy.ndarray) – 1D or 3D array float32 - height grid of wind input

  • uorig (numpy.ndarray) – 3D array float32 - wind speed on these levels

Returns

sum of unew: 3D array float32 - RC corrected windspeed on levels HC: 3D array float32 - HC additional part

Return type

numpy.ndarray

Friedrich, M. M., 2016 Wind Downscaling Program (Internal Met Office Report)

static sigma2hover2(sigma)[source]

Calculate the half peak-to-trough height.

The ancillary data used to estimate the peak to trough height contains the standard deviation of height in a cell. For sine-waves, this relates to the amplitude of the wave as:

Amplitude = sigma * sqrt(2)

The amplitude would correspond to half the peak-to-trough height (h_o_2).

Parameters

sigma (numpy.ndarray) – 2D array of float32 - standard deviation of height in grid cell.

Returns

2D array of float32 - half peak-to-trough height.

Return type

numpy.ndarray

Comments:

Points that had sigma = 0 (i.e. sea points) are set to RMDI.