Source code for pycif.plugins.transforms.complex.isotopes.native2obsvect

def native2obsvect(
    transf,
    y0,
    mapper,
    mod_input,
    ddi,
    ddf,
    mode,
    runsubdir,
    workdir,
    **kwargs
):
    """Convert model isotopologue outputs to observation-vector δ-values and totals.

    Only acts when ``mod_input == 'obs'``.  From the two simulated
    isotopologue columns, computes the δ-value and the total (or
    mass-corrected) concentration:

    * **δ-value** = :math:`([iso_2] / [iso_1] / R_{std} - 1) \\times 1000`
    * **total** = :math:`[iso_1] + [iso_2]` (volume) or
      mass-corrected equivalent (``unit = 'mass'``)

    In TL mode the tangent-linear of both quantities is computed.

    Forward data (isotopologue and reference concentrations) are saved on
    ``transf.fwd_data[ddi]`` for use by the adjoint in ``obsvect2native``.

    Args:
        transf (Plugin): isotopes plugin instance (carries ``parameters_in``,
            ``parameters_out``, ``unit``, ``fwd_data``).
        y0 (pd.DataFrame): obs-vector DataFrame with ``'sim'`` and ``'sim_tl'``
            columns; modified in-place.
        mapper (dict): transform mapper.
        mod_input (str): component type; the function returns ``y0`` unchanged
            unless ``mod_input == 'obs'``.
        ddi (datetime): sub-simulation start date.
        ddf (datetime): sub-simulation end date.
        mode (str): ``'fwd'`` or ``'tl'``.
        runsubdir (str): unused.
        workdir (str): unused.
        **kwargs: unused.

    Returns:
        pd.DataFrame: updated obs-vector DataFrame.
    """

    inputs = transf.parameters_in.names
    outputs = transf.parameters_out.names
    unit = transf.unit
    r_std = transf.parameters_in.standard
    iso_mass = transf.parameters_out.iso_mass
    spec_mass = transf.parameters_out.spec_mass

    if not mod_input == "obs":
        return y0

    dfs = [y0.loc[y0["parameter"] == s.lower()] for s in outputs]

    mask_iso = dfs[0]["parameter_ref"] == inputs[0].lower()
    mask_tot = dfs[0]["parameter_ref"] == inputs[1].lower()

    isotopologue_iso = dfs[1].loc[mask_iso, "sim"]
    spec_ref_iso = dfs[0].loc[mask_iso, "sim"]
    isotopologue_tot = dfs[1].loc[mask_tot, "sim"]
    spec_ref_tot = dfs[0].loc[mask_tot, "sim"]

    # Saving forward simulations for later use by adjoint
    transf.fwd_data[ddi] = {
        "isotopologue_iso": isotopologue_iso,
        "spec_ref_iso": spec_ref_iso,
        "isotopologue_tot": isotopologue_tot,
        "spec_ref_tot": spec_ref_tot,
    }

    iso = isotopologue_iso / spec_ref_iso
    if unit == "volume":
        tot = isotopologue_tot + spec_ref_tot

    else:
        tot = (
            spec_ref_tot / iso_mass[0] * spec_mass
            + isotopologue_tot / iso_mass[1] * spec_mass
        )
        iso *= iso_mass[0] / iso_mass[1]

    iso /= r_std
    iso = (iso - 1) * 1000

    # Filling the original datastore with correct values
    # Changing names as well

    ind_iso = dfs[0].loc[mask_iso, "indorig"]
    ind_tot = dfs[0].loc[mask_tot, "indorig"]

    dfout = y0.drop_duplicates(subset='indorig', keep='first')
    dfout.iloc[ind_iso, dfout.columns.get_loc("sim")] = iso
    dfout.iloc[ind_tot, dfout.columns.get_loc("sim")] = tot

    dfout.iloc[ind_iso, dfout.columns.get_loc("parameter")] = inputs[0].lower()
    dfout.iloc[ind_tot, dfout.columns.get_loc("parameter")] = inputs[1].lower()

    # Applying tangent-linear
    if mode == "tl":
        isotopologue_iso_tl = dfs[1].loc[mask_iso, "sim_tl"]
        spec_ref_iso_tl = dfs[0].loc[mask_iso, "sim_tl"]
        isotopologue_tot_tl = dfs[1].loc[mask_tot, "sim_tl"]
        spec_ref_tot_tl = dfs[0].loc[mask_tot, "sim_tl"]

        iso_tl = (
            isotopologue_iso_tl / spec_ref_iso
            - isotopologue_iso / spec_ref_iso ** 2 * spec_ref_iso_tl
        )
        if unit == "volume":
            tot_tl = isotopologue_tot_tl + spec_ref_tot_tl

        else:
            tot_tl = (
                isotopologue_iso_tl * iso_mass[0] / spec_mass
                + spec_ref_tot_tl * iso_mass[1] / spec_mass
            )
            iso_tl *= iso_mass[0] / iso_mass[1]

        iso_tl /= r_std / 1000

        dfout.iloc[ind_iso, dfout.columns.get_loc("sim_tl")] = iso_tl
        dfout.iloc[ind_tot, dfout.columns.get_loc("sim_tl")] = tot_tl

    return dfout