import numpy as np
import pandas as pd
from ......utils.datastores.empty import init_empty
[docs]
def make_obs(self, ddi, datastore, runsubdir, mode, tracer, input_type, do_simu=True):
"""Dumps observation locations and time steps to obs.bin to let LMDZ know
where to extract concentrations.
"""
# If empty datastore, do nothing
if (isinstance(datastore, dict) and not datastore) or datastore.size == 0:
return
# Re-initialize auxiliary_indexes
# if first species processed for the current period
if not hasattr(self, "auxiliary_indexes"):
self.auxiliary_indexes = {ddi: {
outcomp: {spec: {}
for spec in self.chemistry.acspecies.attributes}
for outcomp in self.output_components}
}
self.chunk_indexes = {ddi: {
outcomp: {spec: {}
for spec in self.chemistry.acspecies.attributes}
for outcomp in self.output_components}
}
if self.reset_obs[ddi]:
self.auxiliary_indexes[ddi] = {
outcomp: {spec: {} for spec in self.chemistry.acspecies.attributes}
for outcomp in self.output_components
}
self.chunk_indexes[ddi] = {
outcomp: {spec: {} for spec in self.chemistry.acspecies.attributes}
for outcomp in self.output_components
}
self.reset_obs[ddi] = False
# If input type is not "concs", i.e., auxiliary data such as pressure,
# Keep in memory the length of the datastore, but not more
if input_type == "concs":
self.auxiliary_indexes[ddi][input_type][tracer] = [len(datastore)]
# If values to extract come from several transforms, keep in memory
# chunks
if "itransform" in datastore["metadata"]:
self.chunk_indexes[ddi][input_type][tracer] = \
datastore["metadata"]["itransform"].values
else:
return
# If do not need to do LMDz simulation, just update obs datastore
if not do_simu:
if not self.iniobs[ddi]:
self.nbobs_prior[ddi] = 0
self.nbdatatot_prior[ddi] = 0
# Keep in memory that observations were already dumped for that period
self.iniobs[ddi] = True
self.auxiliary_indexes[ddi][input_type][tracer].extend(
[self.nbobs_prior[ddi], self.nbdatatot_prior[ddi]]
)
self.nbobs_prior[ddi] += len(datastore)
self.nbdatatot_prior[ddi] += len(datastore)
return
col2extract = "adj_out" if mode == "adj" else "obs"
# Include only part of the datastore
mask = (
datastore["metadata"]["parameter"]
.str.replace("__sample#", "_").str.upper()
.isin([s.upper() for s in self.chemistry.acspecies.attributes])
)
metadata = datastore["metadata"].loc[mask]
maindata = datastore["maindata"].loc[mask]
data = np.concatenate([
metadata.loc[:, ["tstep", "dtstep", "i", "j"]].values,
maindata.loc[:, [col2extract]].values,
metadata.loc[:, ["level"]].values],
axis=1)
# Converts level, tstep, i, j to fortran levels
data[:, [0, 2, 3, 5]] += 1
# Assumes that stations with no level are in first level
# TODO: make it general
data[np.isnan(data[:, -1]), -1] = 1
data[data[:, -1] <= 0, -1] = 1
data[:, -1] /= 100.0
# Attribute ID to each species
for k, s in enumerate(self.chemistry.acspecies.attributes):
mask = (metadata["parameter"].str.lower().str.replace("__sample#", "_")
== s.lower()).values
data[mask, -1] += k + 1
# Update binary file if necessary
if self.iniobs[ddi]:
nobs_orig = \
int(open(f"{runsubdir}/obs.txt", "r").readlines()[0])
obs_orig = (
np.fromfile(f"{runsubdir}/obs.bin", dtype="float")
.reshape((-1, 6), order="F")
)
data = np.append(data, obs_orig, axis=0)
# Write in a binary
obs_file = f"{runsubdir}/obs.bin"
data.T.tofile(obs_file)
with open(f"{runsubdir}/obs.txt", "w") as f:
f.write(str(len(data)))
# Keep in mind that obs.bin already exists and should be updated if
# another tracer is added
self.iniobs[ddi] = True