openghg_inversions.basis#
Functions for creating basis functions and applying them to sensitivity matrices.
- openghg_inversions.basis.basis_functions_wrapper(fp_all: dict, species: str, domain: str, start_date: str, emissions_name: list[str] | None, nbasis: int, use_bc: bool, basis_algorithm: str | None = None, fix_outer_regions: bool = False, fp_basis_case: str | None = None, bc_basis_case: str | None = None, basis_directory: str | None = None, bc_basis_directory: str | None = None, outputname: str | None = None, output_path: str | None = None)#
Wrapper function for selecting basis function algorithm.
- Parameters:
fp_all (dict) – Dictionary object produced from get_data functions
species (str) – Atmospheric trace gas species of interest
domain (str) – Model domain
start_date (str) – Start date of period of inference
emissions_name (str/list) – Emissions dataset key words for retrieving from object store
nbasis (int) – Number of basis function regions to calculated in domain
use_bc (bool) – Option to include/exclude boundary conditions in inversion
basis_algorithm (str, optional) – One of “quadtree” (for using Quadtree algorithm) or “weighted” (for using an algorihtm that splits region by input data). Land-sea separation is not imposed in the quadtree basis functions, but is imposed by default in “weighted” Default None
fixed_outer_region (bool) – When set to True uses InTEM regions to derive basis functions for inner region Default False
fp_basis_case (str) – Name of basis function to use for emissions. Default None
bc_basis_case (str, optional) – Name of basis case type for boundary conditions (NOTE, I don’t think that currently you can do anything apart from scaling NSEW boundary conditions if you want to scale these monthly.) Default None
basis_directory (str, optional) – Directory containing the basis function if not default. Default None
bc_basis_directory (str, optional) – Directory containing the boundary condition basis functions (e.g. files starting with “NESW”) Default None
outputname (str, optional) – File output name Default None
output_path (str, optional) – Passed to outputdir argument of quadtreebasisfunction. Used for testing. Default None
- Returns:
Dictionary object similar to fp_all but with information on basis functions and sensitivities
- Return type:
fp_data (dict)
- openghg_inversions.basis.bucketbasisfunction(fp_all: dict, start_date: str, domain: str, emissions_name: list[str] | None = None, nbasis: int = 100, abs_flux: bool = False, mask: DataArray | None = None) DataArray#
Basis functions calculated using a weighted region approach where each basis function / scaling region contains approximately the same value.
- Parameters:
fp_all (dict) – fp_all dictionary of datasets as produced from get_data functions
start_date (str) – Start date of period of inversion
domain (str) – domain for the basis functions to be calculated over
emissions_name (list) – List of keyword “source” args used for retrieving emissions files from the Object store Default None
nbasis (int) – Desired number of basis function regions Default 100
abs_flux (bool) – When set to True uses absolute values of a flux array Default False
mask (xarray.DataArray) – Boolean mask on lat/lon coordinates. Used to find basis on sub-region Default None
- Returns:
Array with lat/lon dimensions and basis regions encoded by integers.
- Return type:
bucket_basis (xarray.DataArray)
- openghg_inversions.basis.fixed_outer_regions_basis(fp_all: dict, start_date: str, basis_algorithm: str, domain: str, emissions_name: list[str] | None = None, nbasis: int = 100, abs_flux: bool = False) DataArray#
Fix outer region of basis functions to InTEM regions, and fit the inner regions using basis_algorithm.
- Parameters:
fp_all (dict) – fp_all dictionary object as produced from get_data functions
start_date (str) – Start date of period of inference
basis_algorithm (str) – Name of the basis algorithm used. Options are “quadtree”, “weighted”
domain (str) – domain for the basis functions to be calculated over
emissions_name (list) – List of keyword “source” args used for retrieving emissions files from the Object store.
nbasis (int) – Desired number of basis function regions Default 100
abs_flux – When set to True uses absolute values of a flux array Default False
- Returns:
Array with lat/lon dimensions and basis regions encoded by integers.
- Return type:
- openghg_inversions.basis.quadtreebasisfunction(fp_all: dict, start_date: str, domain: str, emissions_name: list[str] | None = None, nbasis: int = 100, abs_flux: bool = False, seed: int | None = None, mask: DataArray | None = None) DataArray#
Creates a basis function with nbasis grid cells using a quadtree algorithm.
The domain is split with smaller grid cells for regions which contribute more to the a priori (above basline) mole fraction. This is based on the average footprint over the inversion period and the a priori emissions field.
The number of basis functions is optimised using dual annealing. Probably not the best or fastest method as there should only be one minima, but doesn’t require the Jacobian or Hessian for optimisation.
- Parameters:
fp_all (dict) – fp_all dictionary of datasets as produced from get_data functions
start_date (str) – Start date of period of inversion
domain (str) – Domain across which to calculate basis functions.
emissions_name (list) – List of keyword “source” args used for retrieving emissions files from the Object store Default None
nbasis (int) – Desired number of basis function regions Default 100
abs_flux (bool) – When set to True uses absolute values of a flux array Default False
seed (int) – Optional seed to pass to scipy.optimize.dual_annealing. Used for testing. Default None
mask (xarray.DataArray) – Boolean mask on lat/lon coordinates. Used to find basis on sub-region Default None
- Returns:
Array with lat/lon dimensions and basis regions encoded by integers.
- Return type:
quad_basis (xarray.DataArray)