# Data Reconciliation

## Overview

Plant measurements are noisy, and taken at face value they contradict the model: nominations do not close, balances do not balance, pressure drops do not match flows. **Data reconciliation** finds the smallest statistically weighted adjustment to the measurements that makes them satisfy the model equations. Because the model carries information the sensors do not, the reconciled estimates are *more* precise than the raw measurements — often by a factor of two or three.

`difflow.reconciliation` provides:

- **Constrained weighted least squares** against any differentiable residual function $F(x, \theta) = 0$ — a flowsheet, a gas network, or an equation set you write yourself.
- **Joint parameter estimation**: mark a variable unmeasured and it is *estimated* rather than reconciled, with a standard error, in the same solve.
- **Covariance of the estimates**, from the inverse KKT matrix, and the sensitivity $\partial \hat x/\partial y$ by automatic differentiation.
- **Gross error detection**: the global $\chi^2$ test, the measurement test, and serial elimination.
- **Observability and redundancy classification**, run *before* the solve so an ill-posed problem raises a named error instead of returning `NaN`.
- **Sensor placement**: what a proposed meter would buy, before anyone buys it.

Everything that makes this work — the constraint Jacobian, the covariance, the sensitivities — is derivative information. A conventional flowsheet package hand-derives it per unit operation; a differentiable one gets it from `jax.jacobian`.

The module is domain-agnostic. `difflow_gas.residuals` shows how a plugin supplies its equation set; see [`examples/28_data_reconciliation.ipynb`](../examples/28_data_reconciliation.ipynb) for the worked gas-network case.

---

## Mathematical Formulation

### The problem

$$\min_x \; (x - y)^T W (x - y) \quad \text{subject to} \quad F(x, \theta) = 0$$

where $y$ are the measurements, $W = \operatorname{diag}(1/\sigma_i^2)$, and $F$ are the model equations. Entries with $\sigma_i = \infty$ get $W_{ii} = 0$: they are **estimated**, not reconciled. A *finite* $\sigma$ on an unmetered variable acts as a Bayesian prior, which is the graceful way to handle a weakly identified parameter.

### The KKT system

The first-order conditions are

$$\begin{bmatrix} W & A^T \\ A & 0 \end{bmatrix} \begin{bmatrix} \Delta x \\ \lambda \end{bmatrix} = \begin{bmatrix} -W(x - y) \\ -F(x) \end{bmatrix}, \qquad A = \frac{\partial F}{\partial x},$$

iterated to convergence. This is Gauss–Newton: it drops the $\sum_k \lambda_k \nabla^2 F_k$ term of the exact Newton Jacobian, which is where a pipe law's non-smooth $q|q|$ second derivative would enter. Exact gradients are recovered afterwards from one implicit-function-theorem correction that *does* use the full Jacobian, so `jax.grad` through `reconcile` is exact with respect to $y$, $\sigma$ and $\theta$.

### Solvability is observability

For $W \succeq 0$, the KKT matrix is nonsingular **iff** $A$ has full row rank and $Z^T W Z \succ 0$ on a basis $Z$ of $\ker A$. With $W$ diagonal and zero exactly on the unmeasured entries, the second condition collapses to

> the unmeasured columns $A_U$ must have full column rank

which is precisely the classical observability condition. `classify` runs this test first and raises `ReconciliationStructureError` naming the culprits.

Two consequences worth knowing:

- **Boundary flows must be state variables, not fixed parameters.** With them fixed, the node-balance block of $A$ is the incidence matrix, whose rank is only $n_{\text{nodes}} - 1$, and the KKT matrix is singular for purely structural reasons.
- **Ranks come from `svd(A)`, never from the eigenvalues of $A^T A$.** Squaring the matrix squares its condition number, and a structurally zero singular value reappears near $\sqrt{\varepsilon}\,\sigma_{\max}$ — an unobservable system reported as full rank.

### Covariance

$$\Sigma_{\hat x} = [K^{-1}]_{11}$$

For a fully measured problem this equals the textbook projection $\Sigma - \Sigma A^T (A\Sigma A^T)^{-1} A \Sigma$, but unlike it, it stays well defined when a variable is unmeasured — so the standard error of an estimated parameter comes from the same expression as the reconciled variance of a metered flow. The adjustment covariance is $\Sigma_{\text{adj}} = \Sigma - \Sigma_{\hat x}$.

`measurement_sensitivity` computes $S = \partial \hat x/\partial y$ by differentiating the solver. For **linear** constraints $S \Sigma S^T = \Sigma_{\hat x}$ exactly. For nonlinear ones the two differ by the curvature term the covariance formula drops: they agree to machine precision when the data are consistent (the multipliers vanish) and diverge in proportion to how inconsistent the data are. The classical formula is itself a linearization; the discrepancy is a useful diagnostic of how nonlinear the model is over the range the adjustments span.

### Scaling

The KKT matrix mixes $W$ (units of 1/variable²) with $A$ (units of residual/variable), so its conditioning depends on the unit system — the same gas network posed in Pa and Pa² rather than bar and bar² is many orders worse conditioned, past what float64 resolves. Variables are scaled by $d_i = \sigma_i$, which makes the scaled weight matrix a 0/1 mask (so $1/\sigma^2$ is never evaluated and an infinite $\sigma$ cannot produce a `NaN`), and residual rows are equilibrated to unit 2-norm. This is on by default; pass `scaling=False` only to observe what it prevents.

### Gross error detection

The **global test** statistic is the optimal objective itself, distributed $\chi^2$ on the degrees of redundancy $m - \operatorname{rank}(A_U)$. The often-quoted form $F(y)^T (A\Sigma A^T)^{-1} F(y)$ is a special case that cannot be evaluated at all when a variable is unmeasured, and carries the wrong degrees of freedom.

The **measurement test** standardizes each adjustment, $z_i = (\hat x_i - y_i)/\sqrt{\Sigma_{\text{adj},ii}}$, which is standard normal under the null hypothesis. A measurement nothing checks has $\Sigma_{\text{adj},ii} = 0$ and is reported as untestable rather than given a spurious $z$. Because that and the redundancy classification are read from the same $\Sigma_{\hat x}$, they cannot disagree.

Least squares smears a gross error across neighbouring measurements, so identification is reliable only where redundancy is high; `serial_elimination` discards the prime suspect and re-tests until the data are clean.

---

## API Reference

### `reconcile`

<!-- doc-test: skip: signature listing, not runnable code -->
```python
def reconcile(
    residual_fn,                    # F(x, params) -> (m,) Array
    y: Array,
    sigma: Array,                   # inf entry => unmeasured
    *,
    params=None, names=None, x0=None,
    unmeasured_init=None, unmeasured_scale=None,
    scaling: Scaling | bool = True,
    max_steps: int = 20, tol: float = 1e-9,
    method: str = "gauss_newton",
    check_structure: bool = True, rank_tol=None,
) -> ReconcileResult
```

**Parameters:**
- `residual_fn` — the model equations, JAX-traceable.
- `y` — measurements. Entries with infinite `sigma` are ignored and may be `nan`.
- `sigma` — standard deviations; `inf` marks a variable to estimate.
- `params` — extra argument threaded to `residual_fn`; `jax.grad` with respect to it answers *how would the reconciled state move if this fixed parameter changed* — a different question from estimating it.
- `check_structure` — run the observability test first. Disable only inside a `jit`/`vmap` sweep whose structure you have already validated.

**Returns** a `ReconcileResult` with `x`, `x_named`, `adjustment`, `objective`, `covariance`, `std`, `converged`, `structure` and `scaling`, plus a `summary()` table.

### Other entry points

<!-- doc-test: skip: signature listing, not runnable code -->
```python
classify(residual_fn, x, sigma, *, scaling, names=None) -> StructureReport
reconciled_covariance(residual_fn, x, sigma, *, scaling) -> Array
measurement_sensitivity(residual_fn, y, sigma, *, x0, scaling) -> Array

global_test(result, alpha=0.05) -> GlobalTestResult
measurement_test(result, alpha=0.05, bonferroni=True) -> MeasurementTestResult
serial_elimination(residual_fn, y, sigma, *, alpha=0.05, max_removed=3) -> list

sensor_value(residual_fn, x, sigma, *, target, candidate, candidate_sigma) -> dict
sensor_ranking(residual_fn, x, sigma, *, target, candidates, candidate_sigma) -> list

monitor(residual_fn, measurements, sigma, *, names=None, alpha=0.05) -> MonitorResult
blame_concentration(suspects, window=15) -> (fraction, culprit)
reconcile_multi(residual_fn, measurements, sigma, *, shared, names=None) -> MultiReconcileResult
```

`StructureReport.classes` assigns every variable one of `measured-redundant`, `measured-just-determined`, `unmeasured-observable` or `unmeasured-unobservable`, and `summary()` prints the table.

---

## Worked Example

A three-stream splitter whose meters do not close:

```python
import jax.numpy as jnp
from difflow.reconciliation import reconcile, global_test, measurement_test

def balance(x, params=None):
    return jnp.array([x[0] - x[1] - x[2]])       # feed = top + bottom

y     = jnp.array([100.0, 62.0, 40.0])           # out by 2 units
sigma = jnp.array([2.0, 1.0, 1.0])               # the feed meter is worst

res = reconcile(balance, y, sigma, names=["feed", "top", "bottom"])
print(res.summary())
print(global_test(res))
```

The feed meter, being the least trusted, absorbs most of the adjustment; the reconciled standard deviations are all below the raw `sigma`. With one equation and three measurements the degree of redundancy is 1.

To *estimate* a quantity instead of reconciling it, give it an infinite sigma:

```python
sigma = jnp.array([2.0, 1.0, jnp.inf])           # the bottoms meter failed
res = reconcile(balance, y, sigma, names=["feed", "top", "bottom"])
res.x_named["bottom"]                            # inferred from the balance
res.std["bottom"]                                # and its standard error
```

Ask for one unknown too many and the problem is diagnosed rather than silently failing:

```python
from difflow.reconciliation import ReconciliationStructureError

sigma = jnp.array([2.0, jnp.inf, jnp.inf])
try:
    reconcile(balance, y, sigma, names=["feed", "top", "bottom"])
except ReconciliationStructureError as err:
    print(err)
# ReconciliationStructureError: reconciliation problem is not solvable:
# 1 unmeasured variable(s) cannot be determined from the constraints.
# Unobservable: top, bottom. ...
```

---

## Gas Networks

`difflow_gas` supplies the equation set for transmission networks. `difflow_gas.residuals.network_residuals` is the single definition of that set — nodal balances, resistance laws, valve relations, **and** the compressor relation $p_{\text{to}} = r\, p_{\text{from}}$. [`difflow_gas.verify`](unit-operations-gas.md) is the reporting layer over it, unflattening the same residual vector into labelled dicts of floats.

`verify` reports every block except the compressor relation, which a sequential solve satisfies by construction — there is nothing to check. A reconciliation cannot drop it, and on the example network it turns out to be exactly what makes the loop observable.

```python
import jax
import difflow_gas as dg

# A small looped network: src -> a -> (compressor) -> b -> c/d, with a cycle
RATIOS = {"cs1": 1.2}
net = dg.GasNetwork(
    arcs={"p1": ("src", "a", "pipe"), "cs1": ("a", "b", "compressor"),
          "p2": ("b", "c", "pipe"), "p3": ("b", "d", "pipe"),
          "p4": ("c", "d", "pipe")},
    beta={aid: dg.weymouth_beta(length_m=L, diameter_m=0.6, roughness_m=1e-4)
          for aid, L in [("p1", 20e3), ("p2", 40e3), ("p3", 60e3), ("p4", 80e3)]},
    supply_kg_s={"src": 120.0, "c": -50.0, "d": -70.0},
    pressure_bounds_bar={n: (30.0, 80.0) for n in ["src", "a", "b", "c", "d"]},
)
fs, dec = dg.build_network_flowsheet(net, root="src", p_slack_pa=60.0e5, ratios=RATIOS)
streams = fs.solve(tol=1e-12, max_iter=500)
p_true = dg.verify.node_pressures_bar(streams, dec)       # the simulated truth
q_true = dg.verify.arc_flows_kg_s(streams, dec)

layout = dg.gas_state_layout(net, efficiency_arcs=["p3"])   # estimate fouling
x_true = layout.pack(p_true, q_true, net.supply_kg_s, {"eta_p3": 1.0})
sigma  = dg.measurement_sigma(layout)                       # meter accuracies
y      = dg.perturb(x_true, sigma, jax.random.PRNGKey(0))   # simulated data

res = dg.reconcile_network(net, y, sigma, layout, ratios={"cs1": 1.2})
p_bar, q_kg_s, supply = dg.reconciled_values(res, layout)
# the reconciliation was free to move the nominations, so check against those
net_rec = dg.GasNetwork(arcs=net.arcs, beta=net.beta, supply_kg_s=supply)
dg.verify.residuals_from_values(p_bar, q_kg_s, net_rec).ok  # True
```

`dg.monitor_network` and `dg.reconcile_network_multi` are the campaign-scale versions of the same call, and take the same `(network, ..., layout, ratios=...)` arguments. All three fill in `names` and `unmeasured_scale` from the layout, so a gas reconciliation never has to restate them.

Measured nominations are deliberately kept out of `GasNetwork`, which rejects supplies that do not sum to zero — real nominations do not, and making them close is what the reconciliation is for.

See [`examples/28_data_reconciliation.ipynb`](../examples/28_data_reconciliation.ipynb) for the full treatment: variance reduction, a biased flow meter found and eliminated, an unmeasured pipe-fouling factor estimated with its standard error, the observability boundary, and sensor placement.

## Reconciling data vs. updating the model

Both are the same optimisation — a variable with a finite `sigma` is a measurement you may move at a cost, one with `sigma = inf` is a parameter you may move for free — so the interesting question is not *how* to update a model but **when you are entitled to**. Letting a parameter float absorbs whatever is wrong, including a broken sensor, and hands back a confident wrong number.

The practical discipline is two clocks: reconcile routinely with parameters **fixed**, so the global test stays a genuine instrument-health monitor; re-estimate parameters only as a deliberate campaign. The trigger is the pair of tests read together — persistent rejection with the *same* sensor blamed every day is an instrument fault, while persistent rejection with a *wandering* suspect is model drift.

### The routine clock: `monitor`

`monitor` runs the first clock. It reconciles a sequence of data sets against one fixed model and keeps the diagnostics, which is what makes the resulting statistic series readable: the model and the sigmas never move, so every wobble in it comes from the data.

```python
import difflow_gas as dg

layout0 = dg.gas_state_layout(net)               # the plain layout, no eta
sigma0 = dg.measurement_sigma(layout0)
x0_true = layout0.pack(p_true, q_true, net.supply_kg_s)
daily_measurements = [dg.perturb(x0_true, sigma0, jax.random.PRNGKey(k))
                      for k in range(6)]         # six days of simulated data

mon = dg.monitor_network(net, daily_measurements, sigma0, layout0, ratios=RATIOS)
mon.statistic          # the chi-squared series, shape (n_days,)
mon.suspects           # who the measurement test blamed each day
mon.diagnose(window=15)
# model drift: 93% of the last 15 steps reject, blame concentration 40%
```

`monitor(residual_fn, measurements, sigma, names=...)` is the domain-agnostic form; `dg.monitor_network` is the same call with the network's residual closure and the layout's names and scales filled in.

`diagnose` applies the rule above. `blame_concentration` is the measurement behind it — the fraction of the window blaming the single most-blamed sensor, counted over the *whole* window including the quiet days, so a campaign that rejects rarely does not read as a concentrated fault just because its few rejections agreed. A verdict of `instrument fault` names the culprit and means go and calibrate it; `model drift` is the one verdict that entitles you to the second clock. Both thresholds are arguments, so the rule can be tuned to a plant's own noise.

A data set whose problem cannot be posed at all records `failed` on its step rather than aborting the campaign.

### The campaign clock: `reconcile_multi`

One day's data gives a noisy parameter estimate, so pool a window. `reconcile_multi` gives each data set its own copy of the plant state while the variables in `shared` appear **once**, estimated from all of them in a single solve:

```python
from difflow.reconciliation import global_test

layout = dg.gas_state_layout(net, efficiency_arcs=["p3"])   # carries eta
sigma = dg.measurement_sigma(layout)
plain_layout, week_of_data = layout0, daily_measurements
window = [layout.embed(y, plain_layout) for y in week_of_data]

res = dg.reconcile_network_multi(
    net, window, sigma, layout, shared=["eta_p3"], ratios=RATIOS,
)
res.shared["eta_p3"], res.shared_std["eta_p3"]
res.states[0]                       # day 0's reconciled state
global_test(res)                    # on the pooled problem
```

Again `reconcile_multi(residual_fn, measurements, sigma, shared=...)` is the general form. `GasStateLayout.embed` re-packs measurements taken against a plainer layout into the one carrying the parameter, by name — the added entry is unmeasured, so it is filled with `nan`.

This is not the same as reconciling each day separately and averaging the estimates, in two ways that matter. The standard error it reports is that *of the pooled estimate*, roughly $\sqrt{K}$ tighter, whereas averaging point estimates leaves you holding one day's error bar for a quantity $K$ of them informed. And because the parameter is one unknown rather than $K$ private copies, the degrees of redundancy count it once — pooling recovers the $K-1$ that separate estimations throw away, and a parameter too weakly identified to be recovered from one data set can still be observable from several. The structure check runs on the stacked problem and reports this.

The catch is lag: pooling estimates the *average* truth over its window, so choose a window short enough that the parameter is genuinely constant across it. A finite `sigma` on a shared variable is a prior, and it is applied once — $K$ copies of one prior would count it $K$ times.

[`examples/29_model_updating.ipynb`](../examples/29_model_updating.ipynb) works this through on a pipe that fouls over a 45-day campaign, including the case where a free parameter manufactures a fouling estimate out of a biased flow meter.

## Tracking a drifting parameter: `track_parameters`

The two clocks are a discipline a person applies. A digital twin — a model kept current against the plant so that what you optimise, price or plan against is the plant as it *is* — has nobody to apply it, and the discipline does not survive being automated naively. Two things have to change.

**The verdict has to gate the update, not advise it.** `diagnose()` already decides whether a sensor or the model is at fault. Left as a report, nothing stops a scheduled re-estimation from running on a day when the honest answer was "go calibrate the dp meter", and section 6 of [`examples/29_model_updating.ipynb`](../examples/29_model_updating.ipynb) shows what that costs: the free parameter absorbs part of the bias, the twin reports a confident fouling estimate for a pipe nobody inspected, **and the $\chi^2$ statistic falls while it happens**. The model looks healthier as it gets wronger. No filter gain prevents this — a slow filter reaches the wrong answer gracefully. Only a gate prevents it.

**The update has to have a memory.** Pooling a window is a rectangular filter: hard edges, a ten-day-old period weighted like this morning's, recomputed from scratch each time. Run it on a schedule and the estimate lurches as periods fall off the back of the window.

`track_parameters` is the two clocks written as an update law:

```python
from difflow.reconciliation import (
    TrackerState, drift_std_from_time_constant, track_parameters,
)

def F(x, params):
    """The network residuals, with the tracked parameter passed through `params`."""
    return dg.residuals.network_residuals(
        x, net, layout0, ratios=RATIOS, efficiencies={"p3": params["eta"]},
    )

run = track_parameters(
    F, daily_measurements, sigma0,
    state=TrackerState.initial(["eta"], [1.0], std=[0.02]),
    drift_std=drift_std_from_time_constant(0.05, 30.0),   # 5% over a month
    names=layout0.names, unmeasured_scale=layout0.default_scale,
)

run.final.as_params()      # {'eta': ...} -- feeds a model, a Block, an LP
run.final.std              # {'eta': ...}
print(run.summary())
```

Each period does four things:

1. reconcile with the parameters **frozen at the current estimate** and record both gross-error tests — the routine clock, and the only reason the tests mean anything;
2. draw a verdict from the campaign so far and put it to `update_gate`;
3. `time_update` the tracker, **whatever the verdict said**;
4. only if the gate opened, estimate the parameters from this period and fold the result in with `measurement_update`.

### The filter

The parameter is given the random walk

$$\theta_{k+1} = \theta_k + w_k, \qquad \operatorname{cov}(w_k) = Q\,\Delta t,$$

which is the same model `difflow.mhe.augment_parameters` puts on a drifting parameter, and `drift_std` = $\sqrt{\operatorname{diag} Q}$ is the same knob as `process_std` there, carrying the same warning: **too large and the parameter absorbs sensor noise, too small and a genuine drift is rejected**. It is a required argument, never a default, because it is the bandwidth of the twin rather than a nuisance.

A rate is hard to have an opinion about; the question an engineer can answer is *how far does this move, and over how long?* `drift_std_from_time_constant(spread, tau)` converts one to the other — a random walk accumulates $\sqrt{Q}\sqrt{t}$, so fouling that costs five percent of duty in a month is `drift_std_from_time_constant(0.05, 30.0)` on a daily clock.

No new estimator is involved. A single period's reconciliation already returns both halves of a Kalman measurement update — the estimate, and the block of `reconciled_covariance` belonging to it — so the update is the combination of two Gaussians and all the physics stays inside `reconcile`:

$$S = P^- + R, \quad K = P^- S^{-1}, \quad \theta^+ = \theta^- + K(\hat\theta - \theta^-).$$

Four properties of that combination are worth stating, because they are what make it better than a rolling refit:

- **The covariance is kept full.** Correlated parameters have a *difference* variance a diagonal covariance gets wrong by a factor of a few, so the prior is a matrix. Pass one to `TrackerState.initial(..., covariance=...)` — `difflow.estimation.predicted_covariance` and `reconciled_covariance` both return the right shape.
- **The covariance is propagated in Joseph form**, $P^+ = (I-K)P^-(I-K)^T + KRK^T$. The short form $(I-K)P^-$ is algebraically equal and numerically worse: it loses symmetry over a long run and can go indefinite. A twin is a long run.
- **A weakly informative period needs no special case.** Its $R$ is large, $K$ goes to zero, the estimate does not move.
- **The time update runs on held periods too.** Holding is not knowing: while the twin refuses to move a parameter, its error bar widens at exactly the rate the drift model claims, and the next permitted update takes a correspondingly larger step. `max_std` bounds that growth if a quiet year would otherwise make the next step a jump.

`Innovation.nis` is the filter's own global test — $v^T S^{-1} v$, which is $\chi^2$ on $p$ degrees of freedom when the drift model and the sigmas are right. A series running well above $p$ says the parameter is moving faster than `drift_std` admits.

### The gate

`update_gate` reads a `MonitorDiagnosis` and permits an update only on `model drift`:

| verdict | gate | why |
|---|---|---|
| `model drift` | **update** | persistent rejection, blame wanders — the signature of a model fault |
| `instrument fault` | hold | blame is concentrated; go calibrate the named sensor |
| `consistent` | hold | the data gave the model nothing to correct |
| `undiagnosed` | hold | persistent rejection with nothing testable to blame |

Holding on `consistent` makes the loop **event-triggered**: the parameter sits still until the evidence is strong enough to reject, then moves. That deadband is deliberate. A parameter re-estimated every period tracks whatever that period's noise favoured, and the twin stops being a model.

The policy is an argument (`allow=`), so it can be widened — deliberately, and at a known cost. `tests/test_tracking.py::TestGate::test_an_ungated_loop_would_have_invented_a_fouling_factor` is the regression: the same filter on the same biased-meter data, with the gate opened, reports a pipe several percent off clean.

The gate is a **statistical rule, not a guarantee**, and the rule has two knobs. `track_parameters` forwards `window`, `rejection_threshold` and `concentration_threshold` to `diagnose`, because the defaults are not right for every plant: on a small network a few days early in a sensor bias can read as diffuse, slip through as `model drift`, and move the parameter before the rule settles. When that happens the fix is to tune the rule to the plant's own noise — lengthen the window, lower the concentration threshold — never to widen `allow`, which removes the rule instead of sharpening it. Section 7 of [`examples/29_model_updating.ipynb`](../examples/29_model_updating.ipynb) shows both the leak and the tuning on a real gas network.

### Where the parameter lives

The tracked parameters are threaded through `params`, so the frozen and free problems are the **same `residual_fn`** — a twin whose two clocks run different code drifts apart in a second, less interesting way. `parameter_measurement` appends them to the state vector with `sigma = inf` and hands the augmented problem to `reconcile` unchanged.

They are left *free* there rather than having the prior passed in as a finite `sigma`, on two counts. `sigma` is a vector, so a prior smuggled through it would be diagonal and would discard exactly the correlations the filter keeps. And with the prior outside, the reconciliation's objective stays a test of data against model, uninflated by how confident the twin already was.

The cost is that each period must identify the parameters on its own. If it cannot, the structure check raises `ReconciliationStructureError` naming them — loudly, rather than returning a NaN. Pool several periods with `reconcile_multi` and hand its `shared` estimate and covariance straight to `measurement_update` instead.

### What it does not fix

A filter on parameters assumes the model *form* is right. Under structural mismatch the parameter converges to something that is not the physical quantity and that shifts with operating point, so the twin's **gradients** are wrong even where its values match — which is fatal, since everything downstream of a differentiable flowsheet is a derivative. The symptom is visible in `run.monitor.statistic`: the gate opens, the parameter moves, and the statistic does not come back down, because no value of the parameter fits. The answer then is `difflow.planning.modifiers.update_modifiers`, which corrects values *and* gradients against the plant, not a faster filter.

`track_parameters` is the offline driver, and it is also a backtest — run it over a recorded campaign to choose `drift_std` before trusting the loop live. The four steps are public and stateless (`update_gate`, `time_update`, `parameter_measurement`, `measurement_update`), so an online loop is the same four calls with a `TrackerState` carried between them.

Section 7 of [`examples/29_model_updating.ipynb`](../examples/29_model_updating.ipynb) runs the whole loop on the gas network the notebook builds, over both campaigns. The gated loop tracks the fouling pipe to 1.257 against a truth of 1.300 and never moves at all on the biased-meter campaign; wiring the gate open buys 0.03 of that lag back and reports 1.165 — a confident 16% fouling claim — on the pipe that is clean.