from logging import info
import xarray as xr
import numpy as np
import os
import glob
import re
from .process_output import process_output
from logging import info
from ......utils import path
from ......utils.hdf5 import _hdf5_lock
from ..utils import OUTPUT_DIRNAME, OUTPUT_PATTERN
from ......utils.check.errclass import CifFileNotFoundError, CifValueError
[docs]
def fetch_sim(self, runsubdir, mode, ddi):
"""Read ICON-ART output files and interpolate to observation locations.
Scans the ``OUTPUT/`` sub-directory for NetCDF output files, calls
:func:`apply_interpolation` to interpolate each field to the observation
metadata (station coordinates and levels), and stores results in
``self.sim_data``.
Args:
self: ICON-ART model plugin instance.
runsubdir (str): path to the period run directory.
mode (str): ``'fwd'``, ``'tl'``, or ``'adj'``.
ddi (datetime): sub-simulation period start.
Raises:
FileNotFoundError: if no ICON-ART output files are found.
"""
# File pattern as specified in the icon nml for the output data
output_path = os.path.join(runsubdir, OUTPUT_DIRNAME)
info(f'Processing output files and performing interpolation...')
if len([f for f in glob.glob(os.path.join(output_path, OUTPUT_PATTERN))]) == 0:
raise CifFileNotFoundError(
f"No output files of format {OUTPUT_PATTERN} were produced by ICON-ART."
)
# Create the directory to process the outputs
path.init_dir(os.path.join(output_path, 'reduced'))
path.init_dir(os.path.join(output_path, 'concatenated_byday'))
path.init_dir(os.path.join(output_path, 'interpolation'))
path.init_dir(os.path.join(output_path, 'dataout'))
# Reduce files and then apply interpolation
process_output(self, runsubdir, ddi)
return
[docs]
def read_sim(self, data2load, runsubdir, mode, ddi, ddf):
"""Extract simulated concentrations from pre-fetched ICON-ART output data.
For each tracer in *data2load*, reads the corresponding interpolated
concentration values from ``self.sim_data`` and writes them into the
CIF data-store (``'spec'`` column for forward, ``'incr'`` for TL).
Args:
self: ICON-ART model plugin instance (carries ``sim_data``).
data2load (dict): tracer-ID-keyed CIF data-store entries to fill.
runsubdir (str): path to the period run directory (unused).
mode (str): ``'fwd'``, ``'tl'``, or ``'adj'``.
ddi (datetime): period start date.
ddf (datetime): period end date.
Returns:
dict: updated data-store with simulated concentrations.
"""
dataout = {}
# Loop over species in data2load
for trid in data2load:
ref_spec = spec = spec_icon = trid[1]
info(f'Reading concentrations of {spec}')
if not len(data2load[trid][ddi]) :
info(f'Empty datastore for {spec}')
continue
dataout[trid] = data2load[trid][ddi].copy()
# Fetch pre-computed dataout
if not hasattr(self, "icon_dataout"):
filepath_dataout = os.path.join(runsubdir, 'outputs/dataout/dataout.nc')
with _hdf5_lock:
self.icon_dataout = xr.open_dataset(filepath_dataout).to_dataframe()
# Change the name if ensemble
if "__sample#" in spec:
ref_spec = spec.split("__sample#")[0]
sample_id = spec.split("__sample#")[1]
new_sample_id = f"-{int(sample_id) + 1:03d}"
spec_icon = ref_spec + new_sample_id
# Filter values for the current tracer only (if multiple species)
mask_trcr = self.icon_dataout['parameter'] == ref_spec
# Put simulated value into correct column
# Different case if concs, or other parameters such as pressure
# Put pressure and other auxiliary data into spec column for later
# interpolation
if trid[0] == 'concs':
M_DRYAIR = 28.97 # mol. weight dry air [kg/mol]
# Get molecular mass
if hasattr(self.chemistry.acspecies, spec):
m_spec = getattr(self.chemistry.acspecies, spec).mass
else:
matches = re.findall(r"([a-zA-Z0-9]+)(?=[-_])", spec)
if hasattr(self.chemistry.acspecies, matches[0]):
m_spec = getattr(self.chemistry.acspecies, matches[0]).mass
else:
try:
m_spec = getattr(self.chemistry.acspecies, matches[0] + "__sample#000").mass
except AttributeError:
raise CifValueError(f"Molecular mass for {spec} not found.")
# Load specific humidity
dataout_qv = self.icon_dataout["qv"]
# Convert concentrations from (moist air mmr) to (dry air vmr in ppbv)
self.icon_dataout[spec_icon] *= (M_DRYAIR / m_spec) * 1e9 / (1 - dataout_qv)
# Fill dataout with new data
dataout[trid][('maindata', 'spec')] = \
self.icon_dataout.loc[mask_trcr, spec_icon].values
dataout[trid][('maindata', 'incr')] = np.full_like(
dataout[trid][('maindata', 'spec')], np.nan)
else:
col = trid[0]
dataout[trid][('maindata', 'spec')] = \
self.icon_dataout.loc[mask_trcr, col].values
dataout[trid][('maindata', 'incr')] = np.full_like(
dataout[trid][('maindata', 'spec')], np.nan)
delattr(self, "icon_dataout")
return dataout