improver.wind_calculations.wind_downscaling module

Module containing wind downscaling plugins.

class FrictionVelocity(u_href, h_ref, z_0, mask)[source]

Bases: 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 (ndarray) – A 2D array of float32 for the wind speed at h_ref

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

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

  • mask (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_impl = <_abc_data 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.

Return type:

ndarray

Returns:

A 2D array of float32 friction velocities

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

Bases: PostProcessingPlugin

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 (Cube) – 2D - model silhouette roughness on pp grid. dimensionless

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

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

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

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

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

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

_abc_impl = <_abc_data object>
calc_av_ppgrid_res(a_cube)[source]

Calculate average grid resolution from a cube.

Parameters:

a_cube (Cube) – Cube to calculate average resolution of

Return type:

float

Returns:

Average grid resolution.

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 (Cube) – holding the silhouette roughness field

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

  • z0_cube (Optional[Cube]) – holding the vegetative roughness field

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

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

Return type:

ndarray

Returns:

Containing bools describing whether or not the tests passed

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

Return type:

None

find_coord_names(cube)[source]

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

Parameters:

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

Return type:

Tuple[str, str, str, str]

Returns:

  • name of the axis name in x-direction

  • name of the axis name in y-direction

  • name of the axis name in z-direction

  • name of the axis name in t-direction

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 (Cube) – cube to check the order of coordinate axis

Return type:

Tuple[int, int, int, int]

Returns:

  • position of x axis.

  • position of y axis.

  • position of z axis.

  • position of t axis.

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 (Cube) – 3D or 4D - representing the wind data.

Return type:

ndarray

Returns:

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

process(input_cube)[source]

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

Parameters:

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

Return type:

Cube

Returns:

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

Raises:

TypeError – If input_cube is not a cube.

tcoordnames = ['time', 'forecast_time']
zcoordnames = ['height', 'model_level_number']
class 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 (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 (ndarray) – 2D array float32 - Standard deviation field of height in the grid cell, units of length

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

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

  • modoro (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)

Return type:

ndarray

Returns:

2D array float32 - reference height for roughness correction

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

Function to calculate the additive height correction.

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

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

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

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

Return type:

ndarray

Returns:

3D array float32 - additive height correction to wind speed

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 (ndarray) – 3D array float32 - velocity on h_in layer, last dim is height

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

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

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

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

Return type:

ndarray

Returns:

2D array float32 - velocity interpolated at h

_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}\]
Return type:

ndarray

Returns:

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

_delta_height()[source]

Function to calculate pp-grid diff from model grid.

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

Return type:

ndarray

Returns:

2D array float32 - height difference, ppgrid-model

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

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

Parameters:
  • xup (ndarray) – 2D array float32 - upper x-bins

  • xlow (ndarray) – 2D array float32 - lower x-bins

  • at_x (ndarray) – 2D array float32 - x values to interpolate y at

  • yup (ndarray) – 2D array float32 - y(xup)

  • ylow (ndarray) – 2D array float32 - y(xlow)

Return type:

ndarray

Returns:

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

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

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

Parameters:
  • xup (ndarray) – 2D array float32 - upper x-bins

  • xlow (ndarray) – 2D array float32 - lower x-bins

  • at_x (ndarray) – 2D array float32 - x values to interpolate y at

  • yup (ndarray) – 2D array float32 - y(xup)

  • ylow (ndarray) – 2D array float32 -y(xlow)

Return type:

ndarray

Returns:

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

_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.

Return type:

None

_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.

Return type:

Tuple[ndarray, ndarray]

Returns:

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

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

calc_roughness_correction(hgrid, uold, mask)[source]

Function to perform the roughness correction.

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

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

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

Return type:

ndarray

Returns:

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

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 (ndarray) – 1D or 3D array float32 - height grid of wind input

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

Return type:

ndarray

Returns:

  • RC corrected windspeed

  • HC additional part

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 (ndarray) – 2D array of float32 - standard deviation of height in grid cell.

Return type:

ndarray

Returns:

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

Comments:

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