"""Calculates responses for each species and scenario."""
import logging
import numpy as np
import xarray as xr
from .calc_ch4 import calc_pmo_rf
from .calc_swv import calc_swv_mass_conc, calc_swv_rf
from .config_model import OUT_TO_INV_REQUIRED
from .interpolate_space import calc_weights, unique_lat_plev_vals
logger = logging.getLogger(__name__)
# CONSTANTS
#
# Correction factor for NO inventory emissions (instead NO2)
# Note: NO2 is the default for emission inventories and response surfaces (O3 and CH4).
CORR_NO = 46.0 / 30.0
[docs]
def calc_resp(spec: str, inv: xr.Dataset, weights: xr.Dataset) -> np.ndarray:
"""Calculate response from response surfaces and emission inventories.
Uses pre-computed weighting parameters.
Args:
spec (str): Name of response species
inv (xarray.Dataset): Emission inventory data
weights (xarray.Dataset): Dataset with weighting parameters
Raises:
KeyError: if species not valid
Returns:
numpy.ndarray: Response array
"""
inv_spec = OUT_TO_INV_REQUIRED[spec]
inv_arr = inv[inv_spec].values
weights_arr = weights["weights"].values
if spec in ["H2O", "O3", "CH4"]:
pass
else:
raise KeyError("calculating response: species not valid")
# Elememt-wise multiplication of inventory emissions and weights
out_arr = (np.multiply(inv_arr.T, weights_arr.T)).T
# Sum over index axis (all steps in emission inventory)
out_arr = np.sum(out_arr, axis=0)
return out_arr
[docs]
def calc_resp_all(config, resp_dict, inv_dict):
"""Loop calc_response function over elements in response dictionary.
Args:
config (dict): Configuration dictionary from config
resp_dict (dict): Dictionary of response xarray Datasets, keys are species
inv_dict (dict): Dictionary of inventory xarray Datasets, keys are years
Returns:
dict: Dictionary of dictionary of numpy arrays of computed responses,
keys are species and inventory years
"""
# "NO" or "NO2" in emission inventory
nox = config["species"]["nox"]
if nox == "NO2":
corr_nox = 1.0
elif nox == "NO":
corr_nox = CORR_NO
else:
raise KeyError("Invalid NOx assumption in config['species']['nox'].")
# unique (lat, plev) locations are shared by all species of an inventory, so these
# values need to only be calculated once. "key" corresponds to the inventory year
unique_lat_plev_dict = {
key: unique_lat_plev_vals(inv) for key, inv in inv_dict.items()
}
out_dict = {}
for spec, resp in resp_dict.items():
# resp_type (str): "conc" or "rf"
resp_type = resp.attrs["resp_type"]
if resp_type not in ["conc", "rf", "tau", "delta"]:
raise ValueError("resp_type not valid")
# Scaling factor of entire response surface
try:
resp_scale = resp.resp_scale.item()
except AttributeError as exc:
msg = "No scaling factor found in " + spec + " response file"
raise AttributeError(msg) from exc
corr = resp_scale
if spec in ["O3", "CH4"]:
corr = corr * corr_nox
out_inv_dict = {}
for key, inv in inv_dict.items():
year = inv.attrs["Inventory_Year"]
weights = calc_weights(spec, resp, inv, unique_lat_plev_dict[key])
# weights = find_weights(spec, resp, inv)
out_arr = corr * calc_resp(spec, inv, weights)
# conc = np.sum(conc_arr)
out_inv_dict[year] = out_arr
out_dict[spec] = out_inv_dict
return out_dict
[docs]
def calc_resp_sub(
species_sub: list[str], config: dict, output_dict: dict, ac: str
) -> tuple[dict, dict]:
"""Calculates responses for specified sub-species.
The calculation of sub-species responses depends on the results
of main species which must be calculated and written to output beforehand.
Args:
species_sub (list[str]): List of sub-species names, such as 'PMO'
config (dict): Configuration dictionary from config
output_dict (dict): Dictionary with computed responses for main
species, keyed by aircraft identifier
ac (str): Aircraft identifier
Returns:
tuple[dict, dict]: ``rf_sub_dict``, dictionary with computed RF
responses, and ``conc_sub_dict``, dictionary with computed
concentration responses; keys are sub-species
Raises:
KeyError: If no method defined for the sub-species
"""
# Get results computed for other species
rf_sub_dict: dict = {}
conc_sub_dict: dict = {}
for spec in species_sub:
if spec == "PMO":
rf_pmo_dict = calc_pmo_rf(output_dict[ac])
rf_sub_dict = rf_sub_dict | rf_pmo_dict
elif spec == "SWV":
if "conc_CH4" in output_dict[ac]:
mass_swv_dict = {}
conc_swv_dict = {}
mass_swv_dict["SWV"], conc_swv_dict["SWV"], _ = calc_swv_mass_conc(
output_dict[ac]["conc_CH4"],
config,
)
rf_swv_dict = calc_swv_rf(mass_swv_dict)
rf_sub_dict = rf_sub_dict | rf_swv_dict
conc_sub_dict = conc_sub_dict | conc_swv_dict
else:
raise KeyError("SWV RF response requires a CH4 concentration")
else:
msg = "No method defined for sub species " + spec
raise KeyError(msg)
return rf_sub_dict, conc_sub_dict