Xarray Interface
A wrapper to manage multiple oopinterface.QGField
objects, with
xarray-based in- and output.
Note
The xarray interface is neither memory- nor dask-optimized. For large datasets it might be necessary to manually split the data and process them chunk-by-chunk to avoid running out of main memory. See Issue #50 on GitHub where progress on memory optimization is tracked.
-
class
xarrayinterface.
QGDataset
(da_u, da_v=None, da_t=None, *, var_names=None, qgfield=falwa.oopinterface.QGFieldNH18, qgfield_args=None, qgfield_kwargs=None)[source] A wrapper for multiple QGField objects with xarray in- and output.
For each combination of timestep, ensemble member, etc. in the input data, a
oopinterface.QGField
object is instanciated. The constructor will automatically flip latitude and pressure dimensions of the input data if necessary to meet the requirements of QGField.This wrapper class imitates the methods of QGField (but not the properties/attributes) and collects and re-organizes output data in xarray Datasets for convenience. All calculations are performed by the QGField routines.
New in version 0.6.1.
- Parameters
da_u (xarray.DataArray | xarray.Dataset) – Input 3D fields of zonal wind. The 3D fields’s dimensions must end with height, latitude and longitude. Other dimensions (e.g. time, ensemble member id) are preserved in the output datasets. Alternatively, a dataset can be given, from which u, v and T fields are then extracted. The da_v and da_t arguments can then be omitted or used as an override.
da_v (xarray.DataArray, optional) – Input 3D fields of meridional wind. The 3D fields’s dimensions must end with height, latitude and longitude. Other dimensions (e.g. time, ensemble member id) are preserved in the output datasets.
da_t (xarray.DataArray, optional) – Input 3D fields of temperature. The 3D fields’s dimensions must end with height, latitude and longitude. Other dimensions (e.g. time, ensemble member id) are preserved in the output datasets.
var_names (dict, optional) – If the auto-detection of variable or coordinate names fails, provide a lookup table that maps plev, ylat, xlon, u, v and/or t to the names used in the dataset.
qgfield (QGField class, optional) – The QGField class to use in the computation. Default:
oopinterface.QGFieldNH18
.qgfield_args (tuple, optional) – Positional arguments given to the QGField constructor.
qgfield_kwargs (dict, optional) – Keyword arguments given to the QGField constructor.
Examples
>>> data = xarray.open_dataset("path/to/some/uvt-data.nc") >>> qgds = QGDataset(data) >>> qgds.interpolate_fields() <xarray.Dataset> ...
Demo script for the analyses done in Nakamura and Huang (2018, Science): using xarray
>>> data_u = xarray.load_dataset("path/to/some/u-data.nc") >>> data_v = xarray.load_dataset("path/to/some/v-data.nc") >>> data_t = xarray.load_dataset("path/to/some/t-data.nc") >>> qgds = QGDataset(data_u, data_v, data_t)
-
property
adv_flux_f1
See
oopinterface.QGFieldBase.adv_flux_f1
.- Returns
- Return type
xarray.DataArray
-
property
adv_flux_f2
See
oopinterface.QGFieldBase.adv_flux_f2
.- Returns
- Return type
xarray.DataArray
-
property
adv_flux_f3
See
oopinterface.QGFieldBase.adv_flux_f3
.- Returns
- Return type
xarray.DataArray
-
property
attrs
Metadata dictionary that is attached to output datasets.
-
compute_lwa_and_barotropic_fluxes
(return_dataset=True)[source] Call compute_lwa_and_barotropic_fluxes on all contained fields.
See
oopinterface.QGFieldBase.compute_lwa_and_barotropic_fluxes()
.- Parameters
return_dataset (bool) – Whether to return the computed fields as a dataset.
- Returns
- Return type
xarray.Dataset or None
-
compute_reference_states
(return_dataset=True)[source] Call compute_reference_states on all contained fields.
See
oopinterface.QGFieldBase.compute_reference_states()
.- Parameters
return_dataset (bool) – Whether to return the computed fields as a dataset.
- Returns
- Return type
xarray.Dataset or None
-
property
convergence_zonal_advective_flux
See
oopinterface.QGFieldBase.convergence_zonal_advective_flux
.- Returns
- Return type
xarray.DataArray
-
property
divergence_eddy_momentum_flux
See
oopinterface.QGFieldBase.divergence_eddy_momentum_flux
.- Returns
- Return type
xarray.DataArray
-
property
fields
Access to the QGField objects created by the QGDataset.
The
oopinterface.QGField
objects are stored in a flattened list.
-
interpolate_fields
(return_dataset=True)[source] Call interpolate_fields on all contained fields.
See
oopinterface.QGFieldBase.interpolate_fields()
.Note
A QGField class may define static stability globally or hemispherically on each level. The output dataset contains a single variable for static stability in case of a global definition and two variables for static stability for a hemispheric definition (suffix
_n
for the northern hemisphere and_s
for the southern hemisphere).- Parameters
return_dataset (bool) – Whether to return the computed fields as a dataset.
- Returns
- Return type
xarray.Dataset or None
-
property
interpolated_theta
See
oopinterface.QGFieldBase.interpolated_theta
.- Returns
- Return type
xarray.DataArray
-
property
interpolated_u
See
oopinterface.QGFieldBase.interpolated_u
.- Returns
- Return type
xarray.DataArray
-
property
interpolated_v
See
oopinterface.QGFieldBase.interpolated_v
.- Returns
- Return type
xarray.DataArray
-
property
lwa
See
oopinterface.QGFieldBase.lwa
.- Returns
- Return type
xarray.DataArray
-
property
lwa_baro
See
oopinterface.QGFieldBase.lwa_baro
.- Returns
- Return type
xarray.DataArray
-
property
meridional_heat_flux
See
oopinterface.QGFieldBase.meridional_heat_flux
.- Returns
- Return type
xarray.DataArray
-
property
ptref
See
oopinterface.QGFieldBase.ptref
.- Returns
- Return type
xarray.DataArray
-
property
qgpv
See
oopinterface.QGFieldBase.qgpv
.- Returns
- Return type
xarray.DataArray
-
property
qref
See
oopinterface.QGFieldBase.qref
.- Returns
- Return type
xarray.DataArray
-
property
static_stability
See
oopinterface.QGFieldBase.static_stability
.- Returns
- Return type
xr.Dataset | Tuple[xr.Dataset, xr.Dataset]
-
property
u_baro
See
oopinterface.QGFieldBase.u_baro
.- Returns
- Return type
xarray.DataArray
-
property
uref
See
oopinterface.QGFieldBase.uref
.- Returns
- Return type
xarray.DataArray
-
xarrayinterface.
hemisphere_to_globe
(ds, var_names=None)[source] Create a global dataset from a hemispheric one.
Takes data from the given hemisphere, mirrors it to the other hemisphere and combines both hemispheres into a global dataset.
If the meridional wind component is found in the dataset, its values will be negated on the created hemisphere. This results in identical fields of local wave activity on both hemispheres (since absolute vorticity is also the same except for the sign), making it possible to use northern_hemisphere_only in the methods of
QGDataset
even if only southern hemisphere data is available. Discontinuities in the meridional wind and derived fields arise due to this at the equator but they generally have only a small effect on the outputs.New in version 0.6.1.
- Parameters
ds (xarray.Dataset) – Input dataset. Must contain the equator (0° latitude).
var_names (dict, optional) – The names of the latitude and meridional wind fields are automatically detected. If the auto-detection of the latitude coordinate and/or the meridional wind component fails, provide a lookup table that maps ylat, and/or v to the names used in the dataset.
- Returns
- Return type
xarray.Dataset
-
xarrayinterface.
integrate_budget
(ds, var_names=None)[source] Compute the integrated LWA budget terms for the given data.
Integrates the LWA tendencies from equation (2) of NH18 in time (over the time interval covered by the input data). The residual (term IV) is determined by subtracting terms (I), (II) and (III) from the LWA difference between the last and first time step in the data. Uses
xarray.DataArray.integrate()
for the time integration of the tendencies.See
QGDataset.compute_lwa_and_barotropic_fluxes()
, which computes all required tendency terms as well as the LWA fields.New in version 0.6.1.
- Parameters
ds (xarray.Dataset) – Input dataset containing the budget tendencies for the time integration interval.
var_names (dict, optional) – The names of LWA and the tendency term variables are automatically detected. If the auto-detection fails, provide a lookup table that maps time, lwa_baro, convergence_zonal_advective_flux, divergence_eddy_momentum_flux, and/or meridional_heat_flux to the names used in the input dataset.
- Returns
- Return type
xarray.Dataset
Examples
>>> qgds = QGDataset(data) >>> ... >>> terms = qgds.compute_lwa_and_barotropic_fluxes() >>> integrate_budget(terms.isel({ "time": slice(5, 10) }))