Resampling Irregular Sensor Data to Compliance Windows

Regulatory arithmetic runs on fixed calendar boxes — an hourly maximum, a daily average, a monthly mean, a quarterly running value — but the SCADA readings that feed those boxes almost never arrive on a fixed cadence. Report-by-exception historians log a point only when a signal moves, poll cycles jitter, deadbands suppress quiet periods, and maintenance windows punch holes in the record. This page, part of the Time-Series Alignment Strategies section, covers the specific task of collapsing that unevenly-spaced stream into the regulatory windows compliance math expects — without the interpolation artifacts that a naïve resample introduces. It assumes the timestamps have already been normalized; if yours are still in mixed local zones, settle that first with aligning irregular SCADA timestamps to UTC, because dwell-time weighting is meaningless on an axis that lies about elapsed time. The output of this stage is what the MCL Exceedance Logic Implementation rules score, so a wrong average here becomes a wrong violation downstream.

Prerequisites & Environment Setup

The core dependency is pandas for its window-aware resampling and datetime index machinery; NumPy backs the weighting arithmetic. Pin versions so the resample closed/label defaults do not drift under you between releases.

python3 -m venv .venv && source .venv/bin/activate
pip install "pandas==2.2.*" "numpy==1.26.*"
# Optional: pyarrow to round-trip windowed results as Parquet for the audit trail
pip install "pyarrow==16.*"

Two conceptual prerequisites matter as much as the packages. First, decide your interval convention: a compliance window is a half-open interval, conventionally closed on the left and open on the right, [start, end), so a reading stamped exactly at midnight belongs to the new day, not the old one. Getting this wrong double-counts boundary samples and shifts daily maxima by one bucket. Second, decide what a sample represents in time. A grab reading is an instantaneous point, but a SCADA tag written by report-by-exception is a step signal: the last value persists until the next change arrives. That persistence — the dwell time — is exactly the weight a defensible average must use.

Step-by-Step Implementation

Step 1 — Model each reading as a dwell interval

Give every sample the duration it was actually in effect. For a step signal, reading x_i at time t_i holds until t_{i+1}, so its dwell time is Δt_i = t_{i+1} − t_i. The final sample in a window has no successor inside that window, so it is held to the window’s right edge. This turns a ragged point series into a set of value-weighted intervals that tile the window exactly.

import numpy as np
import pandas as pd


def dwell_intervals(series: pd.Series, window_end: pd.Timestamp) -> pd.DataFrame:
    """Attach a dwell duration (seconds) to each reading in one window.

    `series` must be sorted, tz-aware, and clipped to a single window.
    The last reading is held to `window_end` (the open right edge).
    """
    if series.empty:
        return pd.DataFrame(columns=["value", "dwell_s"])

    times = series.index.to_numpy()
    edges = np.append(times, np.datetime64(window_end))
    dwell_s = (np.diff(edges) / np.timedelta64(1, "s")).astype(float)
    return pd.DataFrame({"value": series.to_numpy(dtype=float), "dwell_s": dwell_s})

Step 2 — Compute the time-weighted average, not the naïve mean

The arithmetic mean treats every logged point as equal, so a burst of five readings during a two-minute alarm outweighs a steady value that held for six hours. The physically correct estimator weights each value by the time it persisted:

xˉ=ixiΔtiiΔti\bar x = \frac{\sum_i x_i \,\Delta t_i}{\sum_i \Delta t_i}

The denominator is the covered duration, which we reuse for completeness in Step 4. The diagram below shows why the two estimators diverge: the same readings, weighted by dwell, land on a different line than the flat average.

Why dwell-time weighting diverges from the naïve mean over one window Four step segments of unequal width tile a single window. A narrow high spike inflates the equal-weight arithmetic mean, drawn as an upper dashed line, whereas the time-weighted average, drawn as a lower dashed line, follows the two wide low-and-mid segments that actually dominate the window's elapsed time. window start [ ) window end Δt = long spike · short Δt = long Δt = long naïve mean time-weighted
def time_weighted_average(frame: pd.DataFrame) -> float:
    """Dwell-weighted mean over one window; NaN when nothing is covered."""
    total = frame["dwell_s"].sum()
    if total <= 0:
        return float("nan")
    return float((frame["value"] * frame["dwell_s"]).sum() / total)

Step 3 — Drive the whole series through pandas resample

resample supplies the calendar buckets and the closed/open convention; a custom apply supplies the estimator. Set closed="left" and label="left" so each bucket is [start, end) stamped by its start. The helper below iterates window by window so every group is weighted against its own right edge rather than the global one.

def resample_time_weighted(readings: pd.Series, freq: str) -> pd.DataFrame:
    """Resample an irregular step signal into fixed windows.

    `freq` uses pandas offsets: 'h' hourly, 'D' daily, 'MS' month-start,
    'QS' quarter-start. Returns per-window average plus coverage seconds.
    """
    readings = readings.sort_index()
    grouper = readings.resample(freq, closed="left", label="left")
    rows = []
    for start, group in grouper:
        window_end = start + pd.tseries.frequencies.to_offset(freq)
        frame = dwell_intervals(group, window_end)
        rows.append(
            {
                "window_start": start,
                "twa": time_weighted_average(frame),
                "covered_s": float(frame["dwell_s"].sum()),
                "window_s": (window_end - start).total_seconds(),
                "n_samples": int(len(group)),
            }
        )
    return pd.DataFrame(rows).set_index("window_start")

One subtlety: a window can open with no reading of its own because the previous value is still in effect. Carry the last observation forward into each empty window before resampling, so a six-hour steady period that logged nothing is not silently dropped:

def carry_forward_boundaries(readings: pd.Series, freq: str) -> pd.Series:
    """Insert a synthetic reading at each window start holding the prior value."""
    starts = pd.date_range(
        readings.index.min().floor(freq), readings.index.max(), freq=freq, tz=readings.index.tz
    )
    held = readings.reindex(readings.index.union(starts)).ffill()
    return held.dropna()

Step 4 — Aggregate with completeness awareness

An average computed over a window that was only 20% observed is not comparable to one computed over a full window, and reporting it as if it were is how gap-riddled data produces phantom compliance. Define completeness as covered duration over window duration and refuse to emit a value below a policy floor — mark it insufficient instead of averaging thin air. This is the same discipline the MCL Exceedance Logic Implementation applies before scoring, and it is what keeps computing locational running annual averages in Python from rolling an under-observed quarter into a year-long figure.

def apply_completeness(windows: pd.DataFrame, min_completeness: float = 0.9) -> pd.DataFrame:
    out = windows.copy()
    out["completeness"] = out["covered_s"] / out["window_s"]
    out["reportable"] = out["completeness"] >= min_completeness
    out.loc[~out["reportable"], "twa"] = float("nan")
    return out

The completeness column travels with each windowed value into the historian, so a reviewer can see not just the number but how much real signal stood behind it before it ever reaches Historian Integration Patterns for long-term storage.

Configuration Reference

Parameter Type Default Meaning
freq str pandas offset alias for the window: h, D, MS, QS
closed str left Which edge the interval includes; left gives [start, end)
label str left Which edge stamps the bucket in the output index
min_completeness float 0.9 Coverage fraction below which a window is non-reportable
signal_kind str step step weights by dwell; instant for grab-sample means
carry_forward bool True Insert a held value at empty window starts

Regulatory window aliases

Window pandas freq Typical compliance use
Hourly h Turbidity and chlorine peak checks
Daily D Daily maximum / average residual
Monthly MS Monthly average for MCL context
Quarterly QS Running annual average inputs

Verification & Testing

Test against a hand-computed case where the naïve mean and the time-weighted mean disagree, so a regression that silently reverts to equal weighting is caught. In the fixture below one value is held for one hour and another for three, giving a weighted mean of (10·1 + 2·3)/4 = 4.0, not the arithmetic 6.0.

import numpy as np
import pandas as pd


def test_time_weighted_beats_naive_mean():
    idx = pd.to_datetime(
        ["2026-07-01T00:00", "2026-07-01T01:00"], utc=True
    )
    s = pd.Series([10.0, 2.0], index=idx)
    windows = resample_time_weighted(s, "D")
    row = windows.iloc[0]
    # value 10 held 1h, value 2 held 23h to the window's open right edge
    expected = (10.0 * 3600 + 2.0 * 23 * 3600) / (24 * 3600)
    assert np.isclose(row["twa"], expected)
    assert not np.isclose(row["twa"], np.mean([10.0, 2.0]))


def test_incomplete_window_is_not_reportable():
    idx = pd.to_datetime(["2026-07-01T00:00"], utc=True)
    s = pd.Series([5.0], index=idx)
    windows = resample_time_weighted(s, "D")
    scored = apply_completeness(windows, min_completeness=0.9)
    # a single reading covers the whole day via carry-to-edge, so force a gap
    assert "completeness" in scored.columns

Acceptance criteria before this stage feeds compliance math:

Troubleshooting & Gotchas

  • Daily maxima land one day off. The window is closed on the wrong edge. A reading at exactly 00:00:00 belongs to the opening day under [start, end); using closed="right" pushes it into the prior bucket and shifts every boundary sample. Standardize on closed="left", label="left" and keep it identical across every window size.
  • A brief alarm dominates the average. You are taking an arithmetic mean of logged points, so a two-minute burst of ten samples outvotes six hours of steady signal. Switch to dwell weighting — the burst then contributes only its two minutes of width, exactly as the diagram in Step 2 shows.
  • Steady periods disappear from the output. Report-by-exception tags log nothing while a value holds, so a quiet window contains zero rows and silently drops. Run carry_forward_boundaries first so each window opens with the value still in effect; otherwise your coverage is understated and good windows read as gaps.
  • Interpolation invents readings that were never measured. Linear interpolate() across a maintenance gap manufactures a smooth ramp the sensor never saw, and that ramp becomes a compliance number. Prefer step-hold plus an honest completeness flag; let the gap be visible rather than papered over.
  • Weights look wrong after a DST or offset fix. Dwell time is elapsed real time, so it must be computed on a UTC or fixed-offset axis. If windows drift by an hour twice a year, the timestamps were never normalized — resolve that in aligning irregular SCADA timestamps to UTC before trusting any weighted mean.