API reference

Metric sub-components

What are the sub-components?

All statistics and algorithms in this package are built ontop of functions which we refer to as sub-components. These sub-components (should) have one role (e.g. to calculate atmospheric mass at a given kPa level), and can be used to build future metrics implemented into the package upon. They are listed in full here.

Jet statistic sub-components

Algorithms and calculations for the metrics used to identify or classify jet-stream in the literature.

This file is in order of publish year of the metrics (see details_for_all_metrics.py)

Functions and Classes ordered by the year of the paper that uses each metric component first

jsmetrics.metrics.jet_statistics_components.calc_atmospheric_mass_at_kPa(pressure, gravity=9.81, atmospheric_area=510000000.0)[source]

Will calculate the atmospheric mass at a given pressure level.

Radius of earth (R) = 6.372E3 km; Surface area of earth = 4 Pi R^2 = 5.1E8 km^2

Component of method from Archer & Caldiera (2008) https://doi.org/10.1029/2008GL033614 & Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

Parameters:
pressurefloat

in kPa

gravityfloat

m/s^2

jsmetrics.metrics.jet_statistics_components.get_atm_mass_at_one_hPa(hPa)[source]

Get atmospheric mass at one pressure level.

Component of method from Archer & Caldiera (2008) https://doi.org/10.1029/2008GL033614 & Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

Parameters:
hPaint or float

One pressure level in hPa

Returns:
atm_massint or float

Atmospheric mass at given hPa pressure level

jsmetrics.metrics.jet_statistics_components.get_weighted_average_at_one_pressure_level(data, pressure_level, atm_mass, ws_col)[source]

Get weighted average wind speed at one pressure level.

Component of method from Archer & Caldiera (2008) https://doi.org/10.1029/2008GL033614 & Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

Parameters:
dataxarray.Dataset

Data containing windspeed (‘ws’) or u-, v-component wind

pressure_levelint or float

Pressure level to subset

atm_massint or float

Atmospheric mass at given hPa pressure level

ws_colstring

Name of column to calculate weighted average from (e.g. ‘ws’, ‘ua’, ‘va’)

Returns:
outputxarray.Dataset

Data with weighted average at a single pressure level

jsmetrics.metrics.jet_statistics_components.get_mass_weighted_average_wind(data, ws_col, plev_flux=False)[source]

Get mass-weighted average wind-speed from ‘u’, ‘v’ component wind or ‘wind vector’.

Component of method from Archer & Caldiera (2008) https://doi.org/10.1029/2008GL033614 & Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

Parameters:
dataxarray.Dataset

Data containing u- and v-component wind

ws_colstring

Name of column to calculate weighted average from (e.g. ‘ws’, ‘ua’, ‘va’)

Returns:
outputxarray.Dataset

Data containing mass weighted average ws

jsmetrics.metrics.jet_statistics_components.get_sum_atm_mass(data)[source]

Sum atmospheric mass for given pressure levels.

Component of method from Archer & Caldiera (2008) https://doi.org/10.1029/2008GL033614 & Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

Parameters:
dataxarray.Dataset

Data containing plev

Returns:
sum_atm_massint or float

Sum of atmospheric mass

jsmetrics.metrics.jet_statistics_components.get_hPa_or_mbar_and_Pa_or_bar_multiplying_factors_for_plev(data)[source]

Will get the multiplying factor for plev units. E.g. if data contains plev with the unit: Pa the two multiplying factors will be 0.01 and 1 E.g. if data contains plev with the unit: hPa the two multiplying factors will be 1 and 100

Component of method from Archer & Caldiera (2008) https://doi.org/10.1029/2008GL033614 & Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

Parameters:
dataxr.DataArray or xr.Dataset

Data containing plev coord and units for plev

Returns:
hPa_multiplying_factorint

Factor to multiply plev units to convert them to hPa

Pa_multiplying_factorint

Factor to multiply plev units to convert them to Pa

jsmetrics.metrics.jet_statistics_components.calc_mass_weighted_average(data, ws_col)[source]

Calculate mass weighted average wind-speed. Component of method from Archer & Caldiera (2008) https://doi.org/10.1029/2008GL033614 & Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

Parameters:
dataxarray.Dataset

Data containing plev

ws_colstring

Name of column to calculate weighted average from (e.g. ‘ws’, ‘ua’, ‘va’)

Returns:
weighted_averagexr.DataArray

Data with weighted average windspeed based on sum atmospheric mass

jsmetrics.metrics.jet_statistics_components.calc_mass_flux_weighted_pressure(data, ws_col)[source]

Calculate mass-flux weighted pressure.

Component of method from Archer & Caldiera (2008) https://doi.org/10.1029/2008GL033614

Parameters:
dataxarray.Dataset

Data containing windspeed and plev

ws_colstring

Name of column to calculate weighted average from (e.g. ‘ws’, ‘ua’, ‘va’)

Returns:
mass_flux_weighted_pressurexr.DataArray

Data with mass flux weighted pressure

jsmetrics.metrics.jet_statistics_components.calc_mass_flux_weighted_latitude(data, lat_min, lat_max, ws_col)[source]

Calculate mass-flux weighted latitude.

Component of method from Archer & Caldiera (2008) https://doi.org/10.1029/2008GL033614 NOTE: Problem with including 1000 hPa

Parameters:
dataxarray.Dataset

Data containing windspeed and lat

lat_minint or float

Minimum latitude to consider for weighted latitude

lat_maxint or float

Maximum latitude to consider for weighted latitude

ws_colstring

Name of column to calculate weighted average from (e.g. ‘ws’, ‘ua’, ‘va’)

Returns:
mass_flux_weighted_latitudexr.DataArray

Data with mass flux weighted latitude

jsmetrics.metrics.jet_statistics_components.get_climatology(data, freq)[source]

Makes a climatology at given interval (i.e. days, months, season)

Parameters:
dataxarray.Dataset

data with regular time stamp

freqstr

‘day’, ‘month’ or ‘season’

Returns:
climatologyxarray.Dataset

Climatology of a given frequency

jsmetrics.metrics.jet_statistics_components.calc_low_pass_weights(window, cutoff)[source]

Calculate weights for a low pass Lanczos filter. A low-pass filter removes short-term random fluctations in a time series

Component of method from Woollings et al (2010) http://dx.doi.org/10.1002/qj.625 & Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

TAKEN FROM: https://scitools.org.uk/iris/docs/v1.2/examples/graphics/SOI_filtering.html

Parameters:
windowint

The length of the filter window.

cutofffloat

The cutoff frequency in inverse time steps.

Returns:
outputnumpy.array

filtered weights

jsmetrics.metrics.jet_statistics_components.apply_lanczos_filter(dataarray, filter_freq, window_size)[source]

Will carry out Lanczos low-pass filter

Component of method from Woollings et al (2010) http://dx.doi.org/10.1002/qj.625 & Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

Parameters:
datarrayxarray.DataArray

Data to apply filter to containing zonal mean (u-component) wind

filter_freqint

number of days in filter

window_sizeint

number of days in window for Lancoz filter

Returns:
window_consxarray.DataArray

Filtered zonal mean data

jsmetrics.metrics.jet_statistics_components.get_latitude_and_speed_where_max_ws(data_row)[source]

Will return the latitude and windspeed at the index of maximum wind speed from a row of data TODO: expects only one cell to be max

Component of method from Woollings et al (2010) http://dx.doi.org/10.1002/qj.625 & Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1 Barnes & Simpson (2017) https://doi.org/10.1175/JCLI-D-17-0299.1 & Grise & Polvani (2016) https://doi.org/10.1002/2015JD024687

Parameters:
data_rowxarray.DataArray

Data of single time unit containing u-component wind

Returns:
lat_at_maxint or float

Latitude of maximum windspeed

speed_at_maxint or float

Speed at latitude of maximum windspeed

jsmetrics.metrics.jet_statistics_components.calc_latitude_and_speed_where_max_ws(data, var_col='ua')[source]

Will return the latitude and windspeed at the index of maximum wind speed from a row of data TODO: expects only one cell to be max

Component of method from Woollings et al (2010) http://dx.doi.org/10.1002/qj.625 & Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1 Barnes & Simpson (2017) https://doi.org/10.1175/JCLI-D-17-0299.1 & Grise & Polvani (2016) https://doi.org/10.1002/2015JD024687

Parameters:
dataarrayxarray.dataset

Dataset containing u-component wind

Returns:
jet_lats

Latitude of maximum windspeed

jet_speedsint or float

Speed at latitude of maximum windspeed

jsmetrics.metrics.jet_statistics_components.assign_jet_lat_and_speed_to_data(data, max_lat_ws, max_lats_col='jet_lat', max_ws_col='jet_speed')[source]

Will return a data array with the maximum windspeed and latitude of that maximum wind speed

Component of method from Woollings et al (2010) http://dx.doi.org/10.1002/qj.625 & Barnes & Simpson 2017 https://doi.org/10.1175/JCLI-D-17-0299.1

Parameters:
dataxarray.DataArray

Data of single time unit containing u-component wind

max_lat_wsnp.array

Array of maximum latitude and windspeed at those maximums for each timeunit

Returns:
data_with_max_lats_wsxarray.Dataset

Data with maximum latitude and windspeed attached

jsmetrics.metrics.jet_statistics_components.apply_low_freq_fourier_filter(data, highest_freq_to_keep)[source]

Carries out a Fourier transform for filtering keeping only low frequencies.

Component of method from Woollings et al (2010) http://dx.doi.org/10.1002/qj.625

ADAPTED FROM: https://scipy-lectures.org/intro/scipy/auto_examples/plot_fftpack.html

Parameters:
datanp.array - 1-d

time series data at regular intervals

highest_freq_to_keepint

highest frequency to keep in the fourier transform expression NOTE: starts at 0, so highest_freq_to_keep=1 will only keep the constant and first expresion

Returns:
filtered_signp.array

Fourier-filtered signal

jsmetrics.metrics.jet_statistics_components.assign_filtered_lats_and_ws_to_data(data, filtered_max_lats, filtered_max_ws, dim)[source]

Assigns the filtered data back to the returned dataset

Component of method from Woollings et al (2010) http://dx.doi.org/10.1002/qj.625

Parameters:
dataxarray.Dataset

Data to have fourier filtered data assigned to

filtered_max_latsnumpy.array

Fourier-filtered maximum latitude by given timeunit

filtered_max_latsnumpy.array

Fourier-filtered maximum speed at maximum latitude by given timeunit

Returns:
filtered_dataxarray.Dataset

Data with fourier-filtered maximum lat and windspeed attached to it

jsmetrics.metrics.jet_statistics_components.calc_jet_width_for_one_day(data_row, jet_lat, jet_speed)[source]

Calculates jet width on a given data array using the calculated jet lat and speed at that lat

Component of method from Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

Parameters:
data_rowxarray.DataArray

Data of single time unit containing u-component wind

jet_latint or float

Latitude of jet

jet_speedint or float

Speed of jet at jet_lat

Returns:
jet_widthsarray-like

Width of jet-stream in latitude. Units for latitude are the same as the input units i.e. likely to be ‘degrees_north’

jsmetrics.metrics.jet_statistics_components.get_possible_surrounding_jet_lat_vals_for_one_day(one_day_wind_data, jet_speed)[source]

Returns array of windspeed values that are within half speed of the jet_speed given.

Component of method from Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

Parameters:
one_day_wind_dataxarray.DataArray

Data of single time unit containing u-component wind

jet_speedint or float

Speed of jet at jet_lat

Returns:
possible_surrounding_jet_valsarray-like

DataArray of True or False as to whether the windspeed value is within half speed of the jet_speed

jsmetrics.metrics.jet_statistics_components.get_width_of_jet_from_surrounding_lat_vals(jet_lat, possible_surrounding_jet_lat_vals, lat_resolution)[source]

Get latitude width of jet around the given jet latitude value. Width given in units of latitude

Component of method from Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

Parameters:
jet_latint or float

Latitude of jet

possible_surrounding_jet_lat_valsarray-like

DataArray of lats where the windspeed value is within half speed of the jet_speed

lat_resolutionint or float

Latitude resolution in degrees

Returns:
jet_widtharray-like

Width of jet-stream in latitude. Unit for latitude are the same as the input units i.e. likely to be ‘degrees_north’

jsmetrics.metrics.jet_statistics_components.get_all_surrounding_jet_lats_around_max_lat(jet_lat, possible_surrounding_jet_lat_vals, lat_resolution)[source]

Returns all the latitudes in a continuous group (i.e. within one unit of lat resolution) around selected jet (jet_lat).

Component of method from Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

Parameters:
jet_latint or float

Latitude of jet

possible_surrounding_jet_lat_valsarray-like

DataArray of lats where the windspeed value is within half speed of the jet_speed

lat_resolutionint or float

Latitude resolution in degrees

Returns:
lats_within_jetlist

All latitudes in continous group (i.e. within lat_resolution) around jet_lat

jsmetrics.metrics.jet_statistics_components.check_lats_within_jet_in_one_direction(jet_lat, possible_surrounding_jet_lat_vals, amount_to_add_to_lat_val)[source]

Creates a list of latitudes that are within a given distance from the jet_lat and only in one direction i.e. increasing or decreasing

Component of method from Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1

Parameters:
jet_latint or float

Latitude of jet

possible_surrounding_jet_lat_valsarray-like

DataArray of lats where the windspeed value is within half speed of the jet_speed

amount_to_add_to_lat_valint or float

Check if lat within this amount

Returns:
lats_within_jetlist

All latitudes in continous group (i.e. within lat_resolution) around jet_lat

jsmetrics.metrics.jet_statistics_components.scale_lat_vals_with_quadratic_func(lats, speeds, scaled_lats)[source]

Will downscale or upscale the resolution of latitude using a quadratic func.

Component of method from Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1 & Grise & Polvani (2016) https://doi.org/10.1002/2015JD024687

Parameters:
latsxr.DataArray or array-like

Array of latitude values

speedsxr.DataArray or array-like

Array of wind-speeds

scaled_latsarray-like

Array of scaled latitude values of a given resolution (see rescale_lat_resolution function)

Returns:
scaled_lat_valsxr.DataArray or array-like

Array of rescaled latitude values based scaled_lats

jsmetrics.metrics.jet_statistics_components.get_latitude_and_speed_where_max_ws_at_reduced_resolution(lats_and_ws, lat_resolution)[source]

Get latitude and speed where max windspeed at a reduced resolution. Makes use of the quadratic func to scale latitude values.

Component of method from Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1 & Grise & Polvani (2016) https://doi.org/10.1002/2015JD024687

Parameters:
lats_and_wsxr.DataArray or array-like

Array of latitudes and windspeeds

lat_resolutionint or float

Latitude resolution in degrees

Returns:
outputnumpy.array

latitude and wind-speed value scaled by quadratic func

jsmetrics.metrics.jet_statistics_components.get_3_latitudes_and_speed_around_max_ws(row)[source]

Will get the latitudes neighbouring to east, at and to west where the max windspeed is found

Component of method from Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1 & Grise & Polvani (2016) https://doi.org/10.1002/2015JD024687

Parameters:
rowxr.DataArray or array-like

Array containing latitude and wind-speed values

Returns:
neighbouring_latsarray-like

array of 3 neighbouring latitude coordinates around where maximum windspeed is found in input (row)

neighbouring_speedsarray-like

array of 3 neighbouring winspeeds values around where maximum windspeed is found in input (row)

jsmetrics.metrics.jet_statistics_components.add_nan_value_to_arr_if_not_len_enough(arr, len_required)[source]

Add NaN value to array if any values in the array so it is a certain length

Component of method from Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1 & Grise & Polvani (2016) https://doi.org/10.1002/2015JD024687

Parameters:
arrarray-like

Array to check.

len_requiredint

length of array required

Returns:
arrarray-like

new array

jsmetrics.metrics.jet_statistics_components.get_3_neighbouring_coord_values(coord_val, coord_resolution)[source]

Get the three neighbouring values to a given input coord (coord_val)

Component of method from Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1 & Grise & Polvani (2016) https://doi.org/10.1002/2015JD024687

Parameters:
coord_valint or float

Central coord value to get neighbours from

coord_resolutionint or float

in degrees

Returns:
outputarray-like

array of 3 neighbouring latitude coordinates

jsmetrics.metrics.jet_statistics_components.quadratic_func(x, y)[source]

Quadratic function.

Component of method from Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1 & Grise & Polvani (2016) https://doi.org/10.1002/2015JD024687

Parameters:
xxr.DataArray or array-like

Array 1

yxr.DataArray or array-like

Array 2

jsmetrics.metrics.jet_statistics_components.apply_quadratic_func(x, y, vals)[source]

Apply quadratic function to an array of values (vals).

Component of method from Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1 & Grise & Polvani (2016) https://doi.org/10.1002/2015JD024687

Parameters:
xxr.DataArray or array-like

Array 1

yxr.DataArray or array-like

Array 2

valsarray-like

Values to apply function to

Returns:
outputarray-like

quadratic function output

jsmetrics.metrics.jet_statistics_components.get_jet_lat_and_speed_using_parabola_by_day(data_row)[source]

Will get jet latitude and speed by fitting a parabola and taking maximum

Component of method from Barnes & Polvani (2015) http://journals.ametsoc.org/doi/10.1175/JCLI-D-14-00589.1

Parameters:
data_rowxarray.DataArray

Data of single time unit containing u-component wind

Returns:
data_rowxarray.DataArray

Data of single time unit containing jet_lat and jet_speed variables

jsmetrics.metrics.jet_statistics_components.fit_parabola(x, y)[source]

Fits a parabola.

Component of method from Barnes & Polvani (2015) http://journals.ametsoc.org/doi/10.1175/JCLI-D-14-00589.1

jsmetrics.metrics.jet_statistics_components.parabola(x, lat, speed, a)[source]

Parabola.

Component of method from Barnes & Polvani (2015) http://journals.ametsoc.org/doi/10.1175/JCLI-D-14-00589.1

jsmetrics.metrics.jet_statistics_components.assign_jet_lat_speed_to_ten_day_mean_data(ten_day_mean)[source]

Joins the values for latitude and windspeed of the point of maximum windspeed for a given sector to the ten_day_mean data.

Component of method from Barnes & Simpson 2017 https://doi.org/10.1175/JCLI-D-17-0299.1

Parameters:
ten_day_meanxarray.Dataset

Data containing u-component wind and resampled by 10 days

Returns:
ten_day_meanxarray.Dataset

Data with jet latitude and windspeed

jsmetrics.metrics.jet_statistics_components.assign_ten_day_average_jet_lat_speed_to_data(data, ten_day_mean)[source]

Joins the values for latitude and windspeed of the point of maximum windspeed to the main data

Component of method from Barnes & Simpson 2017 https://doi.org/10.1175/JCLI-D-17-0299.1

Parameters:
dataxarray.Dataset

Data containing u-component wind

Returns:
outputxarray.Dataset

Data with jet latitude and windspeed

jsmetrics.metrics.jet_statistics_components.cubic_spline_interpolation(x, y)[source]

Runs a cubic spline interpolation using 2 equal-length arrays.

Component of method from Bracegirdle et al (2018) https://doi.org/10.1175/JCLI-D-17-0320.1

Parameters:
xxarray.DataArray or array-like

array 1 to run cubic spline interpolation

yxarray.DataArray or array-like

array 2 to run cubic spline interpolation

Returns:
cubic_interpolationnp.array

Cubic spline interpolation

jsmetrics.metrics.jet_statistics_components.run_cubic_spline_interpolation_to_get_max_lat_and_ws(data, lat_resolution, ua_col='ua')[source]

Runs a cubic spline interpolation to find maximum latitude and maximum windspeed at a given resolution of latitude.

Component of method from Bracegirdle et al (2018) https://doi.org/10.1175/JCLI-D-17-0320.1

Parameters:
dataxr.Dataset

must contain coords lat

lat_resolutionint or float

Latitude resolution in degrees for cubic spline interpolation

ua_colstr

u-component windspeed column (default=’ua’)

Returns:
max_latfloat

latitude with the maximum wind-speed for time period as deterimined by cubic spline interpolation to a given lat resolution

max_wsfloat

maximum wind-speed at the given latitude for time period

jsmetrics.metrics.jet_statistics_components.run_cubic_spline_interpolation_for_each_unit_of_climatology_to_get_max_lat_and_ws(data, lat_resolution, time_col)[source]

Runs a cubic spline interpolation to find maximum latitude and maximum windspeed at a given resolution of latitude for each.

Component of method from Bracegirdle et al (2018) https://doi.org/10.1175/JCLI-D-17-0320.1

Parameters:
dataxarray.Dataset

Data containing u-component wind-speed. In the case of Bracegirdle et al. 2018, uses a seasonal and annual mean/climatology

resolutionint or float

Latitude resolution in degrees for cubic spline interpolation

time_colstr

Column name containing period (i.e. season, year, month, etc.)

Returns:
max_latslist

list of latitudes with the maximum wind-speed

max_wslist

list of maximum wind-speed at the given latitudes

jsmetrics.metrics.jet_statistics_components.calc_centroid_jet_lat_from_zonal_mean(zonal_mean, area_by_lat)[source]

Will get the centroid latitude of the u-component wind using: jet_lat = integral(60deg, 30deg)(zonal_u**2*lat) dlat / integral(60deg, 30deg)(zonal_u**2) dlat

Component of method from Ceppi et al (2018) https://doi.org/10.1175/JCLI-D-17-0323.1 & Zappa et al. 2018 https://doi.org/10.1029/2019GL083653

Parameters:
zonal_meanxarray.Dataset

zonally-average data containing u-component wind data

area_by_latxarray.DataArray

Information on area in each latitude (i.e. m2 per latitude)

jsmetrics.metrics.jet_statistics_components.get_latitude_value_in_data_row(data_row, lat_val)[source]

Will return the latitude value of a given data row. Needed to loop through data and external array

Component of method from Ceppi et al (2018) https://doi.org/10.1175/JCLI-D-17-0323.1 & Zappa et al. 2018 https://doi.org/10.1029/2019GL083653

Parameters:
data_rowxarray.DataArray

Data row to extract variable from using given latitude

lat_valfloat or int

latitude value to extract from data row

jsmetrics.metrics.jet_statistics_components.get_moving_averaged_smoothed_jet_lats_for_one_day(data_row, width_of_pulse=10)[source]

Get moving average of smoothed_jet_lats for one day. Used in combination with a groupby mapping.

Component of method from Kerr et al. (2020) https://onlinelibrary.wiley.com/doi/10.1029/2020JD032735

Parameters:
data_rowxarray.DataArray

Data containing u-component windspeed of one time unit

Returns:
smoothed_jet_latxarray.Dataset

Data containing jet-stream position

jsmetrics.metrics.jet_statistics_components.smooth_jet_lat_across_lon_with_rectangular_pulse(jet_lat_data, width_of_pulse)[source]

Smooth jet position (jet latitude) by carrying out a convolution with a rectangular pulse.

Component of method from Kerr et al. (2020) https://onlinelibrary.wiley.com/doi/10.1029/2020JD032735

Parameters:
jet_lat_dataxarray.DataArray

Data detailing the latitude of jet-stream at each longitude of one time unit

width_of_pulsefloat or int

Width to make rectangular pulse (likely to be in ‘degrees_east’)

Returns:
filtered_dataxarray.DataArray

Data detailing the latitude of jet-stream that has been smoothed using rectangular pulse

jsmetrics.metrics.jet_statistics_components.get_jet_lat_by_lon(data_row)[source]

Gets all the latitudes of maximum wind-speeds by longitude.

Component of method from Kerr et al. (2020) https://onlinelibrary.wiley.com/doi/10.1029/2020JD032735

Parameters:
data_rowxarray.DataArray

Data containing u-component windspeed of one time unit

Returns:
outputxarray.DataArray

Data detailing the latitude of jet-stream at each longitude of one time unit

Jet core algorithm sub-components

Algorithms and calculations for the metrics used to identify or classify jet-stream in the literature.

This file is in order of publish year of the metrics (see details_for_all_metrics.py)

Functions and Classes ordered by the year of the paper that uses each algorithm first.

jsmetrics.metrics.jet_core_algorithms_components.get_sum_weighted_ws(data, all_plevs_hPa)[source]

Get sum of weighted windspeed. sum weighted windspeed is calculated as follows:

\[\int_{p1}^{p2} (u^2+v^2)^{1/2} \,dp\]

where p1, p2 is min, max pressure level

Component of method from Koch et al (2006) https://doi.org/10.1002/joc.1255

Parameters:
dataxarray.Dataset

Data which should containing the variables: ‘ua’ and ‘va’, and the coordinates: ‘lon’, ‘lat’, ‘plev’ and ‘time’.

all_plevs_hPaarray-like

list of hPa unit pressure levels

Returns:
sum_weighted_wsxarray.Dataset

Data containing sum weighted windspeed values

jsmetrics.metrics.jet_core_algorithms_components.get_weighted_average_ws(sum_weighted_ws, all_plevs_hPa)[source]

Calculates weighted average windspeed as follows:

\[\alpha vel = \frac{1}{(p2-p1)} \int_{p1}^{p2} (u^2+v^2)^{1/2} \,dp\]

where p1, p2 is min, max pressure level.

Component of method from Koch et al (2006) https://doi.org/10.1002/joc.1255

Parameters:
sum_weighted_wsxarray.Dataset

Data containing sum weighted windspeed values

all_plevs_hPaarray-like

list of hPa unit pressure levels

Returns:
weighted_average_wsxarray.Dataset

Data containing weighted average windspeed values

jsmetrics.metrics.jet_core_algorithms_components.get_all_hPa_list(data)[source]

Will get a list of all the pressure levels in the data in hPa/mbar

Parameters:
dataxarray.Dataset

data with plev coord

Returns:
plevnp.array

arr containing pressure level list in data in hPa/mbar

jsmetrics.metrics.jet_core_algorithms_components.get_local_jet_occurence(row, ws_threshold, u_threshold)[source]

Each jet occurence is detected based on three rules applied to inputted wind speed (V = [u, v]):

  1. |V| is a local maxima in latitude and altitude plane

  2. |V| ≥ 30 m/s

  3. |u| ≥ 0 m/s.

Component of method from Schiemann et al 2009 https://doi.org/10.1175/2008JCLI2625.1

NOTE: will only work if 1 day is the resolution

Parameters:
rowxarray.Dataset

Data of a single time unit containing windspeed (ws), plev, lat, lon

ws_thresholdint or float

Windspeed threshold used to extract jet events (default: 30 ms-1)

u_thresholdint or float

Windspeed threshold used to extract u-component wind speed (default: 0 ms^{-1})

Returns:
rowxarray.Dataset

Data of a single time unit with value for jet-maxima (1 == maxima, 0 == none)

jsmetrics.metrics.jet_core_algorithms_components.run_jet_core_and_region_algorithm_on_one_day(row, jet_core_plev_limit, jet_core_ws_threshold, jet_boundary_ws_threshold, ws_drop_threshold, jet_core_lat_distance, check_diagonals)[source]

This method detects jet cores and defines a boundary region beside those cores based on two windspeed thresholds. Two additional checks are applied after initial detection of cores to check whether cores within the same windspeed region are part of the same feature (see ‘jet_boundary_ws_threshold’). This is achieved by checking whether regions with multiple jet cores are more than a certain distance apart (see ‘jet_core_lat_distance’) and the windspeed between two cores does not drop below a threshold (see ‘ws_drop_threshold’). This function runs this method on a single time unit.

This method returns four outputs
  1. ‘jet_core_mask’ – Regions within each latitude/altitude that are local maxima have windspeeds above the ‘jet_core_ws_threshold’

  2. ‘jet_region_mask’ – Regions above, below, left and right of the jet core with windspeed above the ‘jet_boundary_ws_threshold’

  3. ‘jet_region_above_ws_threshold_mask’ – All contigious regions of windspeeds emcompassing a jet core above the ‘jet_boundary_ws_threshold’ (i.e. not just above, below, left and right)

  4. ‘ws’ – Wind speed calculated from ‘ua’, ‘va’ inputs.

Component of method of algorithm originally introduced in Manney et al. (2011) https://doi.org/10.5194/acp-11-6115-2011

Parameters:
rowxarray.Dataset

Data of single time unit containing the variables: ‘ua’ and ‘va’, and the coordinates: ‘lon’, ‘lat’, ‘plev’.

jet_core_plev_limit: tuple or array

Sequence of two values relating to the pressure level limit of the jet cores (original paper uses (100, 400) hPa).

jet_core_ws_thresholdint or float

Threshold used for jet-stream core point.

jet_boundary_ws_thresholdint or float

Threshold for jet-stream boundary point.

ws_drop_thresholdint or float

Threshold for drop in windspeed along the line between cores.

jet_core_lat_distanceint or float

Threshold for maximum distance between cores to be counted the same.

check_diagonalsbool

Whether to check the diagonal edges of each maxima in the latitude-altitude plane.

Returns:
rowxarray.Dataset

Data for one time unit containing four new variables (ws, jet_core_mask, jet_region_mask, jet_region_above_ws_threshold_mask)

jsmetrics.metrics.jet_core_algorithms_components.get_local_maximas_at_each_longitude(row, check_diagonals, var_name)[source]

Runs ‘find_local_maxima_in_2d_dataarray’ on each longitude (see docs for find_local_maxima_in_2d_dataarray for more information)

Parameters:
rowxr.DataArray

Data of single time unit containing the variables: var_name, and the coordinates: ‘lon’, ‘lat’, ‘plev’

check_diagonalsbool

Whether to check the diagonal edges of each maxima in the latitude-altitude plane.

var_namestr

Name of data array variable to find local maxima from e.g. ‘potential_jet_cores’

Returns:
local_maximas_dictdict

Keys are longitude coordinates, values are local maxima locations

jsmetrics.metrics.jet_core_algorithms_components.find_local_maxima_in_2d_dataarray(arr)[source]

Find indices of local maximas within a 2-D array. Should return two values which relate to the position of local maximas (if any).

Component of method of algorithm originally introduced in Manney et al. (2011) https://doi.org/10.5194/acp-11-6115-2011

Parameters:
arrxr.DataArray

A 2-D array with numeric values (i.e. dtypes float or int)

Returns:
local_maximaxr.DataArray

A 2-D array with values relating to the index of the local maximas

Examples

import numpy as np
import xarray as xr

# Example array
data = np.array([
    [1, 2, 3, 4, 5],
    [2, 3, 4, 5, 4],
    [3, 4, 5, 6, 7],
    [4, 5, 6, 7, 6],
    [5, 6, 7, 8, 9]
])

# Convert NumPy array to xarray DataArray
data_array = xr.DataArray(data)

local_maxima_indices = find_local_maxima_in_2d_dataarray(data_array)

for i, j in local_maxima_indices:
    print(f"Local maximum at ({i}, {j}): {data_array[i, j]}")
jsmetrics.metrics.jet_core_algorithms_components.find_local_maxima_in_2d_dataarray_with_diagonals(arr, threshold=10)[source]

Find indices of local maximas within a 2-D array with check for the diagonals. Should return two values which relate to the position of local maximas (if any).

Component of method from Manney et al. (2011) https://doi.org/10.5194/acp-11-6115-2011 & Kuang et al (2014) (https://doi.org/10.1007/s00704-013-0994-x)

Parameters:
arrxr.DataArray

A 2-D array with numeric values (i.e. dtypes float or int)

Returns:
local_maximaxr.DataArray

A 1-D array with values relating to the index of the local maximas

Examples

import numpy as np
import xarray as xr

# Example array
data = np.array([
    [1, 2, 3, 4, 5],
    [2, 3, 4, 5, 4],
    [3, 4, 5, 6, 7],
    [4, 5, 6, 7, 6],
    [5, 6, 7, 8, 9]
])

# Convert NumPy array to xarray DataArray
data_array = xr.DataArray(data)

local_maxima_indices = find_local_maxima_in_2d_dataarray_with_diagonals(data_array)

for i, j in local_maxima_indices:
    print(f"Local maximum at ({i}, {j}): {data_array[i, j]}")
jsmetrics.metrics.jet_core_algorithms_components.get_all_jet_core_mask(local_maximas_dict, mask_shape)[source]

Runs ‘create_mask_using_local_maximas’ using each

Component of method of algorithm originally introduced in Manney et al. (2011) https://doi.org/10.5194/acp-11-6115-2011

Parameters:
local_maximas_dictdict

Keys are longitude coordinates, values are local maxima locations

mask_shapetuple

Values to create a mask of given shape from .

Returns
———-
all_jet_core_masknumpy.array

Mask of the jet cores

jsmetrics.metrics.jet_core_algorithms_components.create_mask_using_local_maximas(local_maximas, mask_shape)[source]

Will create a mask with the same dimensions as the inputted mask_shape and only the local maxima as value 1. All other values will be 0.

Component of method of algorithm originally introduced in Manney et al. (2011) https://doi.org/10.5194/acp-11-6115-2011

Parameters:
local_maximasnp.array

An array containing values relating to the index of the local maximas in a mask with shape: mask_shape

mask_shapetuple

Values to create a mask of given shape from.

Returns:
masknp.array

An array with shape: mask_shape with only the indexes of the local maximas as 1 All other values will 0.

jsmetrics.metrics.jet_core_algorithms_components.get_all_jet_region_contour_mask(row, local_maximas_dict)[source]

Runs ‘get_jet_region_contour_mask’ on each longitude in data with one time step

Component of method of algorithm originally introduced in Manney et al. (2011) https://doi.org/10.5194/acp-11-6115-2011

Parameters:
rowxr.DataArray

Data of single time unit containing the variables: ‘potential_jet_cores’, and the coordinates: ‘lon’, ‘lat’, ‘plev’

local_maximasnp.array

An array containing values relating to the index of the local maximas in a mask with shape: mask_shape

Returns:
all_jet_region_contour_masknp.array

A 3-D array the indexes (lat, plev) of the local maximas for each longitude as 1 All other values will 0.

jsmetrics.metrics.jet_core_algorithms_components.get_jet_region_contour_mask(potential_jet_regions, local_maximas)[source]

Will create a mask based on regions (above ws_threshold) around the jet cores that are contiguous. The mask will contain categorical values (e.g. 1, 2) for each cluster of jet regions contain a jet core. All other values will be 0 (i.e not a jet core or jet region)

Component of method of algorithm originally introduced in Manney et al. (2011) https://doi.org/10.5194/acp-11-6115-2011

Parameters:
potential_jet_regionsxr.DataArray

Windspeeds of regions above the jet region threshold

local_maximasnp.array

An array containing values relating to the index of the local maximas (should be the same shape as input: potential_jet_regions)

Returns:
jet_region_contour_masknp.array

A 2-D array with shape: mask_shape with only the indexes of the local maximas as 1 All other values will 0.

jsmetrics.metrics.jet_core_algorithms_components.get_jet_region_numbers(local_maximas, potential_jet_regions_mask)[source]

Will only return clusters ID numbers of jet regions with an actual jet core in them.

Component of method of algorithm originally introduced in Manney et al. (2011) https://doi.org/10.5194/acp-11-6115-2011

Parameters:
local_maximasnp.array

An array containing values relating to the index of the local maximas (should be the same shape as input: potential_jet_regions)

potential_jet_regions_masknumpy.array

Array of jet regions clusters as returned by scipy.ndimage.label

Returns:
actual_jet_region_numslist

A list of cluster ID numbers with local maximas in them

jsmetrics.metrics.jet_core_algorithms_components.subset_jet_region_mask_to_regions_with_cores(potential_jet_regions_mask, actual_jet_region_nums)[source]

Will subset jet region mask to only those with an actual jet core in them.

Component of method of algorithm originally introduced in Manney et al. (2011) https://doi.org/10.5194/acp-11-6115-2011

Parameters:
potential_jet_regions_masknumpy.array

Array of jet regions clusters as returned by scipy.ndimage.label

actual_jet_region_numslist

A list of cluster ID numbers with local maximas in them

Returns:
potential_jet_regions_masknumpy.array

Subset version of jet regions mask containing only jet region clusters with local maxima in them

jsmetrics.metrics.jet_core_algorithms_components.get_all_jet_regions_mask(row, all_jet_region_contour_mask, local_maximas_dict)[source]

Get all jet regions mask by looping over each longitude in row. See also docs of ‘refine_jet_region_to_leftright_and_abovebelow’

Component of method of algorithm originally introduced in Manney et al. (2011) https://doi.org/10.5194/acp-11-6115-2011

Parameters:
rowxr.DataArray

Data of single time unit containing the coordinates: ‘lon’

all_jet_region_contour_masknp.array

A 3-D array the indexes (lat, plev) of the local maximas for each longitude as 1 All other values will 0.

local_maximas_dictdict

Keys are longitude coordinates, values are local maxima locations

Returns:
all_jet_regions_masknp.array

A refined version of the inputted ‘all_jet_region_contour_mask’ of the local maximas for each longitude as 1 All other values will 0.

jsmetrics.metrics.jet_core_algorithms_components.refine_jet_region_to_leftright_and_abovebelow(array, x, y)[source]

This method will remove all values not left, right, above or below the input x, y coordinates in a 2-D array.

Component of method of algorithm originally introduced in Manney et al. (2011) https://doi.org/10.5194/acp-11-6115-2011

Parameters:
arraynp.array

A 2-D array with numeric values relating to the jet region (i.e. dtypes float or int)

xint

x-coordinate of a local maxima to subset jet region by

y: int

y-coordinate of a local maxima to subset input array by

Returns:
refined_arraynp.array

New 2-D array with values only around x and y coordinates

jsmetrics.metrics.jet_core_algorithms_components.get_values_along_a_line_between_two_coordinates(data, start_point, end_point)[source]

Get all values along a shortest path between two coordinates in a 2-D numpy array.

Component of method of algorithm originally introduced in Manney et al. (2011) https://doi.org/10.5194/acp-11-6115-2011

Parameters:
arrnp.array

A 2-D array with numeric values (i.e. dtypes float or int)

Returns:
local_maximanp.array

A 2-D array with values relating to the index of the local maximas

Examples

import numpy as np

# Example array
data = np.array([
    [1, 2, 3, 4, 5],
    [2, 3, 4, 5, 4],
    [3, 4, 5, 6, 7],
    [4, 5, 6, 7, 6],
    [5, 6, 7, 8, 9]
])

local_maxima_indices = find_local_maxima(data)

for i, j in local_maxima_indices:
    print(f"Local maximum at ({i}, {j}): {data[i, j]}")
jsmetrics.metrics.jet_core_algorithms_components.has_ws_drop_between_cores(ws_between_cores, ws_drop_threshold)[source]

Will check for a windspeed drop of a given threshold between cores

Parameters:
ws_between_coresnp.array

A 1-D array with windspeed between two cores as determined by ‘get_values_along_a_line_between_two_coordinates’

ws_drop_thresholdint or float

Windspeed threshold to check

Returns:
outboolean

True if windspeeds input have a drop greater than ‘ws_drop_threshold’

jsmetrics.metrics.jet_core_algorithms_components.get_index_of_cores_to_drop_based_on_multicore_regions(local_maxima_ind_slices, current_local_maximas, multi_core_region_ws)[source]

Get index of cores to remove in same region based on whether they contain the first occurence of the maximum windspeed or not.

Parameters:
local_maxima_ind_slicesnp.array of np.arrays

Collection of index slices for current local maximas

current_local_maximasnp.array

Collection of indexes of current local maxima

multi_core_region_wsxr.DataArray

Data containing region windspeeds

Returns:
index_of_cores_to_dropnp.array

Mask of jet cores in current longitude plev/lat slice with multi core regions formatted

jsmetrics.metrics.jet_core_algorithms_components.get_current_region_ws_maxima_lat_and_plev_ind(current_local_maximas, local_maxima_ind_slice, multi_core_region_ws)[source]

Get the plev and lat index of current region windspeed maxima. Will select the first occurence of maxima, if there is more than one (rare).

Parameters:
current_local_maximasnp.array

Collection of indexes of current local maxima

local_maxima_ind_slicenp.array

Index slice for current local maxima

multi_core_region_wsxr.DataArray

Data containing region windspeeds

Returns:
outputnp.array

Array of the plev and latitude

jsmetrics.metrics.jet_core_algorithms_components.run_checks_on_jet_cores_and_return_jet_cores(row, initial_jet_core_masks, local_maximas_dict, jet_core_lat_distance, ws_drop_threshold)[source]

This method runs two checks on the jet cores to check whether there are regions with multiple jet cores. Firstly, it checks whether regions with multiple jet cores are more than a certain distance apart (default is 15 degrees, see ‘jet_core_lat_distance’), and hence seperate cores. Secondly, it will check whether the windspeed between two cores drops below a threshold (default is 25 m/s, see ‘ws_drop_threshold’), if so it will remove the latter core.

Component of method of algorithm originally introduced in Manney et al. (2011) https://doi.org/10.5194/acp-11-6115-2011

Parameters:
rowxarray.Dataset

Data of single time unit containing the variables: ‘jet_region_contour_mask’, and the coordinates: ‘lon’, ‘lat’, ‘plev’

initial_jet_core_masks

Initial mask of jet cores to check.

ws_drop_thresholdint or float

Threshold for drop in windspeed along the line between cores (default: 25 m/s)

jet_core_lat_distanceint or float

Threshold for maximum distance between cores to be counted the same (default: 15 degrees)

Returns:
jet_core_masksnumpy.array

Final jet cores mask of shape of ‘initial_jet_core_masks’.

jsmetrics.metrics.jet_core_algorithms_components.get_empty_local_wind_maxima_data(data, expected_dims=('time', 'plev', 'lat', 'lon'))[source]

Will add a new data var of zeros for local wind maxima

Component of method from Pena-Ortiz (2013) https://doi.org/10.1002/jgrd.50305

Parameters:
dataxarray.Dataset

Data which should containing the variables: ‘ua’ and ‘va’, and the expected dims e.g.: ‘lon’, ‘lat’, ‘plev’ and ‘time’.

expected_dimstuple

Expected dimensions of input data (default: “time”, “plev”, “lat”, “lon”)

Returns:
dataxarray.Dataset

Data containing zeros array of (time, plev, lat, lon) dimensions

jsmetrics.metrics.jet_core_algorithms_components.get_potential_local_wind_maximas_by_ws_threshold(ws_slice, ws_threshold=30)[source]

Will return a 2-d array of potential local windspeed maximas

Component of method from Pena-Ortiz (2013) https://doi.org/10.1002/jgrd.50305

Parameters:
ws_slicexarray.Dataset

Data slice of windspeed that has only ‘lat’ and ‘lon’ dims

ws_thresholdint or float

windspeed threshold to apply (default=30 ms-1)

Returns
———-
ws_slicexarray.Dataset

Data slice of windspeed (lat, lon only) with ws_threshold applied

jsmetrics.metrics.jet_core_algorithms_components.get_local_wind_maxima_by_timeunit(row, ws_threshold)[source]

Get local wind maxima by timeunit (i.e. day)

Component of method from Pena-Ortiz (2013) https://doi.org/10.1002/jgrd.50305 who originally use 30 m/s as their ws threshold

Parameters:
rowxarray.Dataset

Data of single time unit containing the variables: ‘ua’ and ‘va’, and the coordinates: ‘lon’, ‘lat’, ‘plev’

ws_thresholdint or float

windspeed threshold to apply

Returns:
rowxarray.Dataset

Data of a single time unit containing 0 or 1 value (local wind maxima) for that time unit

jsmetrics.metrics.jet_core_algorithms_components.get_number_of_timeunits_per_monthyear_with_local_wind_maxima(data)[source]

Will resample by each month and return number of timeunits (i.e. day) with local wind maxima.

Component of method from Pena-Ortiz (2013) https://doi.org/10.1002/jgrd.50305

Parameters:
dataxarray.Dataset

Input data with local wind maxima values by time unit (i.e. day)

Returns:
dataxarray.Dataset

Data containing zeros array of (time, plev, lat, lon) dimensions

jsmetrics.metrics.jet_core_algorithms_components.subdivide_local_wind_maxima_into_stj_pfj(data, local_wind_column_name='local_wind_maxima_by_monthyear')[source]

Subdivide the local_wind_maxima values into the Subtropical Jet (STJ) and Polar Front Jet (PFJ) based on Table 1 pg. 2709 from Pena-Ortiz (2013) https://doi.org/10.1002/jgrd.50305. After the method in that paper, categorisation for the Northern Hemisphere STJ is only possible in DJF, and the latitude bands used are not based on December, March, June or September in Table 1 of that study.

Component of method from Pena-Ortiz (2013) https://doi.org/10.1002/jgrd.50305

Parameters:
dataxarray.Dataset

Input data with local wind maxima values

Returns:
dataxarray.Dataset

data with polar_front_jet and subtropical_jet subdivisions

jsmetrics.metrics.jet_core_algorithms_components.run_jet_occurence_and_centre_alg_on_one_day(row, occurence_ws_threshold)[source]

Runs JetStreamCoreIdentificationAlgorithm method on a single day.

Component of method from Kuang et al (2014) https://doi.org/10.1007/s00704-013-0994-x

Parameters:
rowxarray.Dataset

Data of single time unit containing the variables: ‘ua’ and ‘va’, and the coordinates: ‘lon’, ‘lat’, ‘plev’

occurence_ws_thresholdint or float

Threshold used to identify a jet-stream occurence point

Returns:
occ_alg.output_dataxarray.Dataset

Data with jet occurence and centre points (1 for occurence, 2 for centre)

class jsmetrics.metrics.jet_core_algorithms_components.JetStreamOccurenceAndCentreAlgorithm(data, occurence_ws_threshold=30)[source]

Bases: object

Jet-stream occurence and centre algorithm.

Component of method from Kuang et al (2014) https://doi.org/10.1007/s00704-013-0994-x

Methods

run()

Run the algorithm.

run_algorithm(data)

Class method for running the algorithm.

classmethod run_algorithm(data)[source]

Class method for running the algorithm.

run()[source]

Run the algorithm.

jsmetrics.metrics.jet_core_algorithms_components.run_jet_core_algorithm_on_one_day(row, ws_core_threshold, ws_boundary_threshold)[source]

Runs JetStreamCoreIdentificationAlgorithm method on a single time unit.

Component of method of jet_core_identification_algorithm in jsmetrics and is inspired by the method from Manney et al. (2011) (https://doi.org/10.5194/acp-11-6115-2011), which is described in Section 3.1 of that study.

Parameters:
rowxarray.Dataset

Data of single time unit containing the variables: ‘ua’ and ‘va’, and the coordinates: ‘lon’, ‘lat’, ‘plev’

ws_core_thresholdint or float

Threshold used for jet-stream core point

ws_boundary_thresholdint or float

Threshold for jet-stream boundary point

Returns:
rowxarray.Dataset

Data for one time unit containing jet-cores (ID number relates to each unique core)

class jsmetrics.metrics.jet_core_algorithms_components.JetStreamCoreIdentificationAlgorithm(data, ws_core_threshold=40, ws_boundary_threshold=30)[source]

Bases: object

Jet-stream core identification algorithm.

Component of method of jet_core_identification_algorithm in jsmetrics and is inspired by the method from Manney et al. (2011) (https://doi.org/10.5194/acp-11-6115-2011), which is described in Section 3.1 of that study.

Methods

run:

run algorithm

run_algorithm:

class method for running algorithm

classmethod run_algorithm(data, ws_core_threshold=40, ws_boundary_threshold=30)[source]

Class method for running algorithm

run()[source]

Runs algorithm.

Utility functions

All functions used in package that provide generic utility (such as for xarray compatibility) are listed here, and are seperate from metric sub-components.

Data Utils

Various utility functions needed for the jet-stream metrics and algorithms not belonging to windspeed or spatial utils. The module is built from xarray data-structures, so this file contains all the stuff that helps the library handle xarray dataset/dataarrays

Classes and Functions ordered alphabetically

jsmetrics.utils.data_utils.add_num_of_days_to_360Datetime(datetime_360day, num_of_days_to_add)[source]

Adds a number of days to cftime.Datetime360Day format

Parameters:
datetime_360daycftime.Datetime360Day

Date to add days to

num_of_days_to_addint or float

Number of days to add to 360 day format

Returns:
new_360day_datecftime.Datetime360Day

Date with added number of days

Raises:
AssertionError

If not ‘360_day’ format or not a cftime.Datetime object

ValueError

If number of days not more than 0

jsmetrics.utils.data_utils.add_num_of_days_to_NoLeapDatetime(datetime_noleap, num_of_days_to_add)[source]

Adds a number of days to cftime.DatetimeNoLeap format

Parameters:
datetime_noleapcftime.DatetimeNoLeap

Date to add days to

num_of_days_to_addint or float

Number of days to add to noleap day format

Returns:
new_noleap_datecftime.DatetimeNoLeap

Date with added number of days

Raises:
AssertionError

If not ‘noleap’ format or not a cftime.Datetime object

ValueError

If number of days not more than 0

jsmetrics.utils.data_utils.add_pad_to_array(arr, pad_width=1, constant_values=0)[source]

Add a edge of constant values around a numpy array

Parameters:
arrarray-like

A 2-D array with numeric values (i.e. dtypes float or int)

pad_widthint

Width of values to add to edges (default: 1)

constant_valuesint

Value to add to edges (default: 0)

Returns:
padded_arrarray-like

A 2-D array with new dimensions as a pad has been added at the edges

jsmetrics.utils.data_utils.check_at_least_n_plevs_in_data(data, n_plevs)[source]

Checks there are at least two pressure-levels (plevs) in xarray dataset

Parameters:
dataxarray.Dataset

Data to check

Raises:
ValueError

If not two pressure levels

jsmetrics.utils.data_utils.check_coords_in_data(data, req_coords)[source]

Will check if a given coordinates are in the data and raises a key error if any of them are not. This function is needed to check before algorithm continues and uses too much memory. Built from xarray.

Parameters:
dataxarray.Dataset

Data to check

req_coordsarray-like or tuple

Coordinates to check if in data

Raises:
KeyError

If any given coord is not in data

jsmetrics.utils.data_utils.check_if_data_is_xarray_datatype(data)[source]

What it says on the tin.

Parameters:
dataxarray.Dataset

Data to check

Raises:
TypeError

If input is not an xarray data type

jsmetrics.utils.data_utils.check_plev_units(data, expected_plev_units)[source]

Will check that pressure level units are in data plev coord (i.e. Pa, hPa, mbar, millibars).

Parameters:
dataxr.DataArray or xr.Dataset

Data containing plev coord and units for plev

expected_plev_unitslist

List of names of plev units (i.e. mbar, Pa, hPa, etc.)

Returns:
unitsstr

Units of plev coord

Raises:
KeyError

If plev is not a coord in input ‘data’. And if the units of plev are not in expected plev units

jsmetrics.utils.data_utils.check_variables_in_data(data, req_variables)[source]

What it says on the tin. Built from xarray

Parameters:
dataxarray.Dataset

Data to check

req_variablesarray-like

Variables needed in data

Raises:
KeyError

If variables not in data

jsmetrics.utils.data_utils.find_nearest_value_to_array(array, value)[source]

Will find the nearest value to a given array Built for use with xarray adapted from: https://stackoverflow.com/questions/2566412/find-nearest-value-in-numpy-array

Parameters:
arrayarray-like

Array to reference

valuenumeric-like

Number to check nearest value in input array

Returns:
outputnumeric

Closest value to input value in input array

jsmetrics.utils.data_utils.get_local_maxima(arr, axis=0)[source]

Uses scipy.signal.argrelextrema to get index location of minimum value in array

Taken from https://stackoverflow.com/questions/4624970/finding-local-maxima-minima-with-numpy-in-a-1d-numpy-array

Parameters:
arrarray-like

array to find index location of minima value from

axisint

axis for scipy.signal.argrelextrema

jsmetrics.utils.data_utils.get_local_minima(arr, axis=0)[source]

Uses scipy.signal.argrelextrema to get index location of maximum value in array

Taken from https://stackoverflow.com/questions/4624970/finding-local-maxima-minima-with-numpy-in-a-1d-numpy-array

Parameters:
arrarray-like

array to find index location of maxima value from

axisint

axis for scipy.signal.argrelextrema

jsmetrics.utils.data_utils.filter_local_extremes_to_min_distance(local_extrema, min_distance_threshold=2)[source]

Filter local extremes (i.e. outputs of minima or maxima from scipy.signal.argrelextrema) so that no neighbours remain in array.

Parameters:
local_extremaarray-like

likely (stacked) outputs of scipy.signal.argrelextrema

min_distance_thresholdint

Minimum distances between indexes in array (default: 2 indexes)

Returns:
filtered_extremaarray-like

(stacked) outputs of scipy.signal.argrelextrema with values less than a min_distance_threshold away removed

Examples

maxima_indices = np.column_stack(jsmetrics.utils.data_utils.get_local_maxima(current[‘ws’].data)) filtered_indices = filter_local_extremes_so_no_neighbours(maxima_indices, min_distance_threshold=2)

jsmetrics.utils.data_utils.find_intersection_between_two_array_of_arrays(array1, array2)[source]

Find the intersection between two arrays of arrays. See examples for example.

Parameters:
array1np.array

First array of arrays to compare with array2

array2np.array

Second array of arrays to compare with array1

Returns:
intersectionnp.array

Intesection of the two arrays returning only the arrays within the original two that are in both.

Examples

array1 = [[1, 2], [3, 4]] array2 = [[1, 2], [4, 3]] find_intersection_between_two_array_of_arrays(array1, array2) # returns [[1, 2]]

jsmetrics.utils.data_utils.get_num_of_decimal_places(num)[source]

Gets number of decimal places in a float

Parameters:
numfloat or int

input number to get decimal places from

Returns:
decimal_placesint

number of decimal places

jsmetrics.utils.data_utils.remove_duplicates(arr)[source]

Removes duplicates from array. From: https://stackoverflow.com/questions/2213923/removing-duplicates-from-a-list-of-lists

Parameters:
arrarray-like

arr to remove duplicates from

Returns:
arrarray-like

arr with no duplicates

jsmetrics.utils.data_utils.remove_unwanted_coords_from_data(data, wanted_coords, unwanted_coords=(), show_error=False)[source]

What it says on the tin. Built from xarray

Parameters:
dataxarray.Dataset

Data to check

wanted_coordsarray-like or tuple

Coords to retain in data

unwanted_coordsarray-like or tuple

Coords to remove from data

Raises:
ValueError

If coord cannot be removed from data

jsmetrics.utils.data_utils.rescale_lat_resolution(lats, lat_resolution)[source]

Rescale latitude resolution to input lat resolution in degrees.

Component of method from Barnes & Polvani (2013) https://doi.org/10.1175/JCLI-D-12-00536.1 & Grise & Polvani 2014 https://doi.org/10.1002/2013GL058466 & Bracegirdle et al (2018) https://doi.org/10.1175/JCLI-D-17-0320.1

Parameters:
latsxr.DataArray or array-like

Array of latitude values

lat_resolutionint or float

Latitude resolution in degrees

Returns:
outputnumpy.array

Rescaled array of latitude values

jsmetrics.utils.data_utils.slice_array_by_index_breaks(array_to_slice, index_breaks)[source]

Will break an array down into segments based on a list of indexes containing information about where to create slices

Parameters:
array_to_slicearray-like

Array to break down based on slice breaks

index_breaksarray-like

Indexes from which to slice array

Returns:
output: list

Broken down slices of original array_to_slice

Spatial Utils

All spatial operations needed for the jet-stream metrics and algorithms

Classes and Functions ordered alphabetically.

jsmetrics.utils.spatial_utils.calc_spatial_integral(xr_da, lon_name='longitude', lat_name='latitude', radius=6371000.0)[source]

Calculate spatial integral of xarray.DataArray with grid cell weighting.

Author: Denis Sergev (https://github.com/dennissergeev) https://gist.github.com/dennissergeev/60bf7b03443f1b2c8eb96ce0b1880150

Parameters:
xr_da: xarray.DataArray

Data to average

lon_name: str, optional

Name of x-coordinate

lat_name: str, optional

Name of y-coordinate

radius: float

Radius of the planet [metres], currently assumed spherical (not important anyway)

Returns:
Spatially averaged xarray.DataArray.
jsmetrics.utils.spatial_utils.calc_spatial_mean(xr_da, lon_name='longitude', lat_name='latitude', radius=6371000.0)[source]

Calculate spatial mean of xarray.DataArray with grid cell weighting.

Author: Denis Sergev (https://github.com/dennissergeev) https://gist.github.com/dennissergeev/60bf7b03443f1b2c8eb96ce0b1880150

Parameters:
xr_da: xarray.DataArray

Data to average

lon_name: str, optional

Name of x-coordinate

lat_name: str, optional

Name of y-coordinate

radius: float

Radius of the planet [metres], currently assumed spherical (not important anyway)

Returns:
Spatially averaged xarray.DataArray.
jsmetrics.utils.spatial_utils.calc_total_great_circle_distance_along_line(line)[source]

Returns the total great circle (haversine) distance along a linestring or multilinestring

Component of method from Cattiaux et al (2016) https://doi.org/10.1002/2016GL070309

Parameters:
lineshapely.geometry.LineString or shapely.geometry.MultiLineString

Line to calculate great circle distance from

Returns:
total_distanceint or float

Total distance in degrees or m

jsmetrics.utils.spatial_utils.get_great_circle_distance_along_linestring(line)[source]

Calculate great circle distance along the length of linestring

Parameters:
lineshapely.geometry.LineString

Line to calculate great circle (haversine) distance along

Returns:
distancefloat

Great circle (haversine) distance along input line

jsmetrics.utils.spatial_utils.get_latitude_circle_linestring(latitude, lon_min, lon_max)[source]

Will return a linestring of a latitude circle.

Component of method from Cattiaux et al (2016) https://doi.org/10.1002/2016GL070309

Parameters:
latitudeint or float

given latitude to calculate circle from

lon_minint or float

Minimum longitude for circle to extend to

lon_maxint or float

Maximum longitude for circle to extend to

Returns:
circleshapely.geometry.LineString

Linestring of latitude circle around a hemisphere

jsmetrics.utils.spatial_utils.get_one_contour_linestring(dataarray, contour_level)[source]

Returns a linestring or multi-linestring of a given contour.

Component of method from Cattiaux et al (2016) https://doi.org/10.1002/2016GL070309

Parameters:
dataarrayxarray.DataArray

Array of Geopotential height (zg) values to calculate contour from

contour_level

Value with which to calculate a contour from geopotential height

Returns:
contour_lineshapely.geometry.LineString or shapely.geometry.MultiLineString

Contour line of geopotential height (zg) a given contour

jsmetrics.utils.spatial_utils.grid_cell_areas(lon1d, lat1d, radius=6371000.0)[source]

Calculate grid cell areas given 1D arrays of longitudes and latitudes for a planet with the given radius.

Author: Denis Sergev (https://github.com/dennissergeev) https://gist.github.com/dennissergeev/60bf7b03443f1b2c8eb96ce0b1880150

Parameters:
lon1d: numpy.array

Array of longitude points [degrees] of shape (M,)

lat1d: numpy.array

Array of latitude points [degrees] of shape (M,)

radius: float, optional

Radius of the planet [metres] (currently assumed spherical)

Returns:
Array of grid cell areas [metres**2] of shape (M, N).
jsmetrics.utils.spatial_utils.haversine(lon1, lat1, lon2, lat2)[source]

Calculate the great circle distance between two points on the earth (specified in decimal degrees)

jsmetrics.utils.spatial_utils.seperate_one_contour_into_line_segments(one_contour)[source]

Seperates a list of vertices of a given contour into multiple segments for eventual conversion to multi-linestring.

Needed with new Matplotlib depreciation to collections and allsegs of the contour plots.

Component of method from Cattiaux et al (2016) https://doi.org/10.1002/2016GL070309

Parameters:
one_contourmatplotlib.path.Path

Path object of one contour from a contour plot

Returns:
all_segmentslist of lists

List of all segments that can be converted to a multilinestring

Windspeed Utils

Operations needed for the jet-stream metrics and jet-stream algorithms that specifically operate on windspeed data. Includes the base classes and function for dealing with windspeed vectors and lat/lon or lat/plev slices of windspeed (so called: WindSpeedSlice class)

Classes and Functions ordered alphabetically.

jsmetrics.utils.windspeed_utils.get_resultant_wind(u, v)[source]

Gets wind vector from u-wind and v-wind

jsmetrics.utils.windspeed_utils.get_wind_direction_in_degrees(u, v)[source]

Gets wind direction from u-wind and v-wind In degrees (0-360)

jsmetrics.utils.windspeed_utils.get_zonal_mean(data)[source]

Will get the zonal mean either by pressure level (plev) or for one layer

Parameters:
dataxarray.Dataset

Data containing lon and plev coords

Returns:
zonal_meanxarray.DataSet

zonal mean data

Raises:
KeyError

when ‘lon’ not discovered as coord

class jsmetrics.utils.windspeed_utils.WindSpeedSlice(data, req_variables=('ua', 'va'))[source]

Bases: object

Base class for windspeed slice.

Methods

label_slice:

Label the windspeed slice using where condtion

get_values:

Get values using xarray

label_slice(condition, label)[source]

Returns labels where condition is met.

get_values()[source]

Get values.

class jsmetrics.utils.windspeed_utils.PressureLevelWindSpeedSlice(data, req_variables=('ua', 'va'))[source]

Bases: WindSpeedSlice

Data will be lon*lat

Methods

label_slice:

Label the windspeed slice using where condtion

get_values:

Get values using xarray

req_coords = ('lat', 'lon')
class jsmetrics.utils.windspeed_utils.LatitudeWindSpeedSlice(data, req_variables=('ua', 'va'))[source]

Bases: WindSpeedSlice

Data will be lon*plev

Methods

label_slice:

Label the windspeed slice using where condtion

get_values:

Get values using xarray

req_coords = ('lat', 'plev')