[education/kstars] kstars: PID Auto-Tune: fix 2x gain error, RA periodic-error contamination, spurious warning

Jasem Mutlaq <[email protected]> Tue, 4 Aug 2026 15:39:18 +0000 (UTC)
Newsgroups gmane.comp.kde.cvs
Message-ID <[email protected]>
Git commit aff7639baef6701fc90b0b8be369efccacceb875 by Jasem Mutlaq.
Committed on 04/08/2026 at 15:38.
Pushed by mutlaqja into branch 'master'.

PID Auto-Tune: fix 2x gain error, RA periodic-error contamination, spurious warning

Four issues were found in computeAndApplyAxisGain() by reviewing a real
protocol run against the mount's independently-measured calibration
(ms/arcsec), documented in PID_AUTOTUNE_ISSUES.pdf:

1. K was 2x the true process gain. extractPulseCurve() baseline-subtracts
   each curve, so paired EAST/WEST (or NORTH/SOUTH) plateaus end up
   opposite-signed; differencing them therefore added their magnitudes, but
   the result was divided by a single pulse's magnitude. Confirmed against
   the guide's own calibration (133.8 ms/arcsec vs 68.0 ms/arcsec as coded,
   136.0 ms/arcsec at half that value).

2. RA's process gain was contaminated by periodic error. Pairing was by
   array order, not by firing time, so a paired EAST/WEST pulse could land
   ~50s apart against a ~120s PE period, roughly opposite points of the
   cycle; PE then adds a systematic bias instead of cancelling, and repeats
   at nearly the same phase across repetitions (355s repeat interval is
   close to 3 PE cycles), so the result looked repeatable while still being
   wrong. Confirmed by lack of magnitude-independence: RA's fitted gain
   differed ~1.8x between 500ms and 1000ms pulses (DEC did not).

3. The C++ gain-lock had no equivalent of the offline trainer's confidence
   gate, so it applied RA's contaminated fit anyway.

4. applyPIDAutoTuneGainLock()'s own gain write tripped refreshFingerprint()'s
   "guide settings changed during the protocol" warning, reading as a user
   error rather than the protocol's own expected action.

Fixes all four at the root: computeAndApplyAxisGain() no longer pairs and
differences opposite-direction curves at all. Each pulse-response curve is
fit independently as a step plus a local linear trend (drift and periodic
error both show up as that trend and are removed per-curve), and the step
is read off by extrapolating the trend back to the pulse instant. This
removes both the differencing-induced 2x error and the need for temporally
adjacent pairing. Fits are grouped by pulse magnitude, and since the process
gain must be magnitude-independent for a linear system, the gain is refused
if per-magnitude medians disagree by more than 25% (RA's contaminated data
would fail this check); the largest-magnitude group is used when they
agree, since it has the best signal-to-contamination ratio. The fingerprint
warning is fixed by re-baselining directly after a gain-lock write instead
of going through refreshFingerprint(), which preserves the one-time warning
for a genuine later user edit in the same run rather than consuming it on
the protocol's own change.

offline_trainer/pid_autotune.py's advisory calculation is updated to mirror
the same per-curve/magnitude-grouped approach, so the live gain-lock and
the offline cross-check agree; opsaiconfig.ui/kstars.kcfg's option
descriptions are corrected to say the gain is genuinely applied and locked
(they previously and inaccurately said "advisory only, never applied
automatically").

Patch by Pavan ([email protected]), from
https://invent.kde.org/pavansg/kstars/-/commit/5015a8bf95767858825819bde8db4230b7226d2f,
applied as-is after review and a clean rebuild.

Co-Authored-By: Pavan <[email protected]>

M  +113  -76   kstars/ekos/guide/aiguideprotocol.cpp
M  +96   -35   kstars/ekos/guide/offlinetrainer/pid_autotune.py
M  +1    -1    kstars/ekos/guide/opsaiconfig.ui
M  +1    -1    kstars/kstars.kcfg

https://invent.kde.org/education/kstars/-/commit/aff7639baef6701fc90b0b8be369efccacceb875

diff --git a/kstars/ekos/guide/aiguideprotocol.cpp b/kstars/ekos/guide/aiguideprotocol.cpp
index fb1015ed88..4a1a225a63 100644
--- a/kstars/ekos/guide/aiguideprotocol.cpp
+++ b/kstars/ekos/guide/aiguideprotocol.cpp
@@ -22,6 +22,7 @@
 #include <QSet>
 #include <QVector>
 #include <algorithm>
+#include <limits>
 #include <cmath>
 
 namespace Ekos
@@ -1289,26 +1290,6 @@ PulseCurve extractPulseCurve(const QJsonObject &session, const QString &axisKey)
     return curve;
 }
 
-// Linear interpolation of (t, v) at queryT; clamps to the curve's endpoints outside its range.
-double interpAt(const QVector<double> &t, const QVector<double> &v, double queryT)
-{
-    if (t.isEmpty())
-        return 0.0;
-    if (queryT <= t.first())
-        return v.first();
-    if (queryT >= t.last())
-        return v.last();
-    for (int i = 1; i < t.size(); ++i)
-    {
-        if (t[i] >= queryT)
-        {
-            const double frac = (queryT - t[i - 1]) / std::max(t[i] - t[i - 1], 1e-9);
-            return v[i - 1] + frac * (v[i] - v[i - 1]);
-        }
-    }
-    return v.last();
-}
-
 double medianOf(QVector<double> values)
 {
     if (values.isEmpty())
@@ -1326,16 +1307,20 @@ constexpr int MIN_FITS_TO_APPLY = 4;
 // Same default as offline_trainer/pid_autotune.py's SIMC_LAMBDA_L_FACTOR.
 constexpr double SIMC_LAMBDA_L_FACTOR = 3.0;
 constexpr double INTEGRAL_GAIN_CONSERVATIVE_FRACTION = 0.25;
+// K is arcsec of motion per ms of pulse, so it must not depend on the pulse magnitude used
+// to measure it. When the per-magnitude estimates disagree by more than this, the response
+// is contaminated (drift/PE that did not cancel in the pairing) and the fit is not trusted.
+constexpr double K_MAGNITUDE_CONSISTENCY_TOLERANCE = 0.25;
 
 } // namespace
 
-// Live C++ port of offline_trainer/pid_autotune.py's per-axis SIMC-style calculation,
-// simplified to avoid a nonlinear curve fit: the plateau amplitude is estimated as the
-// mean of the last third of a paired pulse's response curve rather than fit, since on
-// every rig tested so far the step response is effectively flat well before the
-// response-frame window ends and the fitted tau has never resolved below the sampling
-// floor anyway (pid_autotune_plan.md ยง9.1) -- this is a reasonable same-night stand-in,
-// not a replacement for the offline trainer's fuller fit.
+// Live C++ counterpart of offline_trainer/pid_autotune.py's per-axis SIMC-style
+// calculation; both estimate the step the same way so the applied gain and the trainer's
+// advisory agree. Each response is fitted on its own as a step plus a local linear trend
+// and the step is read off the fitted line, which removes drift, periodic error and the
+// residual tracking-rate offset an RA pulse leaves behind. Note tau and the dead time L
+// are not resolvable at guide cadence -- the mount reaches its plateau within one frame --
+// so the SIMC arithmetic reduces to a conservative constant over the process gain.
 bool AIGuideProtocol::computeAndApplyAxisGain(const QString &axis, double msPerArcsec)
 {
     if (msPerArcsec <= 0.0)
@@ -1350,85 +1335,97 @@ bool AIGuideProtocol::computeAndApplyAxisGain(const QString &axis, double msPerA
     const double pixelScale = m_SysIdData.value("equipment").toObject()
                               .value("pixel_scale_arcsec_per_px").toDouble(1.0);
 
-    QMap<double, QList<QJsonObject>> posByMag, negByMag;
+    QMap<double, QList<QJsonObject>> curvesByMag;
     const QJsonArray sessions = m_SysIdData.value("sessions").toArray();
     for (const auto &s : sessions)
     {
         const QJsonObject so = s.toObject();
         if (so.value("type").toString() != "pulse_response" || so.value("pulse_axis").toString() != axis)
             continue;
-        const double mag = so.value("pulse_magnitude_ms").toDouble();
         const QString dir = so.value("pulse_direction").toString();
-        if (dir == posDir)
-            posByMag[mag].append(so);
-        else if (dir == negDir)
-            negByMag[mag].append(so);
+        if (dir == posDir || dir == negDir)
+            curvesByMag[so.value("pulse_magnitude_ms").toDouble()].append(so);
     }
 
-    QVector<double> kSamples, lSamples, tFirstSamples;
+    QMap<double, QVector<double>> kByMag, lByMag, tFirstByMag;
     QSet<int> signSet;
 
-    for (auto it = posByMag.constBegin(); it != posByMag.constEnd(); ++it)
+    // Each response is fitted on its own as a step plus a local linear trend, and the step
+    // is read off at the pulse instant. Detrending removes anything moving at a steady rate
+    // during the window -- drift, periodic error, and the residual tracking-rate offset an
+    // RA pulse leaves behind -- so the estimate no longer needs the opposite-direction pulse
+    // to have been fired close enough in time for those to cancel by subtraction.
+    for (auto it = curvesByMag.constBegin(); it != curvesByMag.constEnd(); ++it)
     {
         const double mag = it.key();
-        const QList<QJsonObject> &posList = it.value();
-        const QList<QJsonObject> &negList = negByMag.value(mag);
-        const int pairs = std::min(posList.size(), negList.size());
-        for (int i = 0; i < pairs; ++i)
+        for (const QJsonObject &so : it.value())
         {
-            const PulseCurve cp = extractPulseCurve(posList.at(i), axisKey);
-            const PulseCurve cn = extractPulseCurve(negList.at(i), axisKey);
-            if (cp.t.size() < 5 || cn.t.size() < 5)
+            const PulseCurve curve = extractPulseCurve(so, axisKey);
+            if (curve.t.size() < 6)
                 continue;
 
-            QVector<double> tArr, diff;
-            for (int j = 0; j < cp.t.size(); ++j)
+            // Frames whose exposure overlapped the pulse only integrate part of the motion
+            // and do not sit on the post-step trend line. At least one frame is always
+            // dropped; at a fast cadence more than one can fall inside the pulse.
+            int first = 1;
+            while (first < curve.t.size() - 4 && curve.t[first] < mag / 1000.0)
+                ++first;
+
+            double sumT = 0.0, sumV = 0.0, sumTT = 0.0, sumTV = 0.0;
+            const int n = curve.t.size() - first;
+            for (int j = first; j < curve.t.size(); ++j)
             {
-                if (cp.t[j] < cn.t.first() || cp.t[j] > cn.t.last())
-                    continue;
-                tArr.append(cp.t[j]);
-                diff.append(cp.pos[j] - interpAt(cn.t, cn.pos, cp.t[j]));
+                sumT  += curve.t[j];
+                sumV  += curve.pos[j];
+                sumTT += curve.t[j] * curve.t[j];
+                sumTV += curve.t[j] * curve.pos[j];
             }
-            if (tArr.size() < 5)
+            const double denom = n * sumTT - sumT * sumT;
+            if (std::abs(denom) < 1e-12)
                 continue;
-
-            // Plateau amplitude/noise estimate from the last third of samples.
-            const int tailCount = std::max(1, static_cast<int>(tArr.size()) / 3);
-            double tailSum = 0.0;
-            for (int j = tArr.size() - tailCount; j < tArr.size(); ++j)
-                tailSum += diff[j];
-            const double pFit = tailSum / tailCount;
+            const double slope     = (n * sumTV - sumT * sumV) / denom;
+            const double intercept = (sumV - slope * sumT) / n;
 
             double varSum = 0.0;
-            for (int j = tArr.size() - tailCount; j < tArr.size(); ++j)
-                varSum += (diff[j] - pFit) * (diff[j] - pFit);
-            const double residualStd = std::sqrt(varSum / tailCount);
+            for (int j = first; j < curve.t.size(); ++j)
+            {
+                const double resid = curve.pos[j] - (intercept + slope * curve.t[j]);
+                varSum += resid * resid;
+            }
+            // Two parameters were fitted, so n-2 is the unbiased residual estimate.
+            const double residualStd = std::sqrt(varSum / std::max(n - 2, 1));
 
-            if (std::abs(pFit) < 2.0 * std::max(residualStd, 1e-6))
-                continue; // noise-dominated, skip -- same gate as pulse_response_fit.py
+            if (std::abs(intercept) < 2.0 * std::max(residualStd, 1e-6))
+                continue; // noise-dominated, skip
 
-            double lSample = tArr.last();
+            double lSample = curve.t.last();
             const double threshold = std::max(2.5 * residualStd, 1e-6);
-            for (int j = 0; j < tArr.size(); ++j)
+            for (int j = 0; j < curve.t.size(); ++j)
             {
-                if (std::abs(diff[j]) > threshold)
+                if (std::abs(curve.pos[j]) > threshold)
                 {
-                    lSample = tArr[j];
+                    lSample = curve.t[j];
                     break;
                 }
             }
 
-            kSamples.append(std::abs(pFit) * pixelScale / mag);
-            lSamples.append(lSample);
-            tFirstSamples.append(tArr.first());
-            signSet.insert(pFit > 0 ? 1 : -1);
+            kByMag[mag].append(std::abs(intercept) * pixelScale / mag);
+            lByMag[mag].append(lSample);
+            tFirstByMag[mag].append(curve.t.first());
+            // The two directions must move the star opposite ways; anything else is noise.
+            const int dirSign = (so.value("pulse_direction").toString() == posDir) ? 1 : -1;
+            signSet.insert(dirSign * (intercept > 0 ? 1 : -1));
         }
     }
 
-    if (kSamples.size() < MIN_FITS_TO_APPLY)
+    int totalFits = 0;
+    for (auto it = kByMag.constBegin(); it != kByMag.constEnd(); ++it)
+        totalFits += it.value().size();
+
+    if (totalFits < MIN_FITS_TO_APPLY)
     {
         emit protocolLog(QString("PID Auto-Tune [%1]: only %2 usable pulse-response fit(s) (need >= %3) -- keeping current gain.")
-                         .arg(axis).arg(kSamples.size()).arg(MIN_FITS_TO_APPLY));
+                         .arg(axis).arg(totalFits).arg(MIN_FITS_TO_APPLY));
         return false;
     }
     if (signSet.size() > 1)
@@ -1438,9 +1435,45 @@ bool AIGuideProtocol::computeAndApplyAxisGain(const QString &axis, double msPerA
         return false;
     }
 
-    const double K   = medianOf(kSamples);
-    const double L   = medianOf(lSamples);
-    const double tau = std::max(medianOf(tFirstSamples), L);
+    // K must be magnitude-independent; disagreement means the pairing failed to cancel
+    // drift/PE, and pooling the samples would produce a median between two populations
+    // that no single measurement supports.
+    double kMin = std::numeric_limits<double>::max(), kMax = 0.0;
+    QStringList perMag;
+    for (auto it = kByMag.constBegin(); it != kByMag.constEnd(); ++it)
+    {
+        const double km = medianOf(it.value());
+        kMin = std::min(kMin, km);
+        kMax = std::max(kMax, km);
+        perMag << QString("%1ms: %2\"/ms (n=%3)").arg(it.key(), 0, 'f', 0).arg(km, 0, 'f', 5).arg(it.value().size());
+    }
+    if (kByMag.size() < 2)
+    {
+        emit protocolLog(QString("PID Auto-Tune [%1]: only one usable pulse magnitude (%2) -- the "
+                                 "cross-magnitude check cannot run, keeping current gain.")
+                         .arg(axis, perMag.join(", ")));
+        return false;
+    }
+    if (kMin <= 0.0 || (kMax - kMin) / kMin > K_MAGNITUDE_CONSISTENCY_TOLERANCE)
+    {
+        emit protocolLog(QString("PID Auto-Tune [%1]: process gain disagrees across pulse magnitudes "
+                                 "(%2) -- responses are contaminated, keeping current gain.")
+                         .arg(axis, perMag.join(", ")));
+        return false;
+    }
+
+    // The largest pulse has the best signal-to-contamination ratio, so prefer its estimate.
+    const double bestMag = kByMag.lastKey();
+    if (kByMag.value(bestMag).size() < MIN_FITS_TO_APPLY)
+    {
+        emit protocolLog(QString("PID Auto-Tune [%1]: only %2 usable fit(s) at %3ms (need >= %4) -- "
+                                 "keeping current gain.")
+                         .arg(axis).arg(kByMag.value(bestMag).size()).arg(bestMag, 0, 'f', 0).arg(MIN_FITS_TO_APPLY));
+        return false;
+    }
+    const double K   = medianOf(kByMag.value(bestMag));
+    const double L   = medianOf(lByMag.value(bestMag));
+    const double tau = std::max(medianOf(tFirstByMag.value(bestMag)), L);
 
     if (K <= 0.0)
     {
@@ -1468,7 +1501,7 @@ bool AIGuideProtocol::computeAndApplyAxisGain(const QString &axis, double msPerA
     emit protocolLog(QString("PID Auto-Tune [%1]: K=%2\"/ms  L=%3s  tau=%4s  (n=%5 fits) "
                              "-- gain %6 -> %7 (locked for the rest of this session)")
                      .arg(axis).arg(K, 0, 'f', 5).arg(L, 0, 'f', 2).arg(tau, 0, 'f', 2)
-                     .arg(kSamples.size()).arg(oldGain, 0, 'f', 3).arg(proportionalGain, 0, 'f', 3));
+                     .arg(kByMag.value(bestMag).size()).arg(oldGain, 0, 'f', 3).arg(proportionalGain, 0, 'f', 3));
     return true;
 }
 
@@ -1490,8 +1523,12 @@ void AIGuideProtocol::applyPIDAutoTuneGainLock()
     const bool raApplied  = computeAndApplyAxisGain("RA",  cal.raPulseMillisecondsPerArcsecond());
     const bool decApplied = computeAndApplyAxisGain("DEC", cal.decPulseMillisecondsPerArcsecond());
 
+    // The gains just written are the protocol's own doing, so re-baseline the fingerprint
+    // directly rather than through refreshFingerprint(), which would report them as a user
+    // changing settings mid-run -- and would leave that warning suppressed for the rest of
+    // the protocol, hiding a genuine change made later.
     if (raApplied || decApplied)
-        refreshFingerprint();
+        m_SysIdData["model_fingerprint"] = buildFingerprint();
 }
 
 }
\ No newline at end of file
diff --git a/kstars/ekos/guide/offlinetrainer/pid_autotune.py b/kstars/ekos/guide/offlinetrainer/pid_autotune.py
index 12c66d7e79..c9b72d558d 100644
--- a/kstars/ekos/guide/offlinetrainer/pid_autotune.py
+++ b/kstars/ekos/guide/offlinetrainer/pid_autotune.py
@@ -26,7 +26,6 @@ SPDX-License-Identifier: GPL-2.0-or-later
 
 import numpy as np
 
-from pulse_response_fit import fit_pulse_response
 
 # Default SIMC design parameter: lambda = max(tau, SIMC_LAMBDA_L_FACTOR * dead_time).
 # Larger -> slower/more robust closed loop; smaller -> faster/less margin. Kept as
@@ -93,6 +92,66 @@ def _calibration_ms_per_arcsec(sysid: dict, axis: str) -> float:
     return float(np.median(values)) if values else 0.0
 
 
+
+def _step_samples_by_magnitude(sysid: dict, axis: str) -> dict:
+    """
+    Step amplitude of each pulse response, in pixels, keyed by pulse magnitude.
+
+    Each response is fitted on its own as a step plus a local linear trend and the
+    step is read off at the pulse instant. Detrending per curve removes anything
+    moving at a steady rate during the window -- drift, periodic error, and the
+    residual tracking-rate offset an RA pulse leaves behind -- so the estimate does
+    not depend on the opposite-direction pulse having been fired close enough in
+    time for those to cancel by subtraction. This mirrors
+    AIGuideProtocol::computeAndApplyAxisGain() so the advisory and the gain the
+    protocol applies live are computed the same way.
+    """
+    axis_key = "ra_raw_px" if axis.upper() == "RA" else "dec_raw_px"
+    pos_dir, neg_dir = ("EAST", "WEST") if axis.upper() == "RA" else ("NORTH", "SOUTH")
+
+    by_mag, dead_times, first_times, signs = {}, {}, {}, set()
+    for s in sysid["sessions"]:
+        if s.get("type") != "pulse_response" or s.get("pulse_axis", "").upper() != axis.upper():
+            continue
+        direction = s.get("pulse_direction", "")
+        if direction not in (pos_dir, neg_dir):
+            continue
+        mag = float(s.get("pulse_magnitude_ms", 0.0))
+        frames = s.get("response_frames", [])
+        # The pulse lands mid-exposure, so the first frame integrates only part of
+        # the motion and does not sit on the post-step trend line.
+        if mag <= 0.0 or len(frames) < 6:
+            continue
+        base_frames = s.get("baseline_frames", [])
+        baseline = float(np.mean([f.get(axis_key, 0.0) for f in base_frames])) if base_frames else 0.0
+        t = np.array([f.get("t", 0.0) for f in frames], dtype=float)
+        pos = np.array([f.get(axis_key, 0.0) for f in frames], dtype=float) - baseline
+
+        # Frames whose exposure overlapped the pulse only integrate part of the motion.
+        # At least one is always dropped; at a fast cadence more than one can fall inside.
+        first = 1
+        while first < len(t) - 4 and t[first] < mag / 1000.0:
+            first += 1
+        slope, intercept = np.polyfit(t[first:], pos[first:], 1)
+        resid = pos[first:] - (intercept + slope * t[first:])
+        # Two parameters were fitted, so n-2 is the unbiased residual estimate.
+        residual_std = float(np.sqrt(np.sum(resid ** 2) / max(len(resid) - 2, 1)))
+        if abs(intercept) < 2.0 * max(residual_std, 1e-6):
+            continue  # noise-dominated
+
+        by_mag.setdefault(mag, []).append(abs(intercept))
+        dead_times.setdefault(mag, []).append(
+            _estimate_dead_time_s(t, pos, max(residual_std, 1e-6)))
+        first_times.setdefault(mag, []).append(float(t[0]))
+        signs.add((1 if direction == pos_dir else -1) * (1 if intercept > 0 else -1))
+
+    return {"by_mag": by_mag, "dead_times": dead_times,
+            "first_times": first_times, "signs": signs}
+
+# Mirrors K_MAGNITUDE_CONSISTENCY_TOLERANCE in aiguideprotocol.cpp.
+K_MAGNITUDE_CONSISTENCY_TOLERANCE = 0.25
+
+
 def _recommend_axis_pid_gain(sysid: dict, axis: str, guide_exp: float,
                              lambda_l_factor: float, verbose: bool) -> dict:
     """
@@ -105,15 +164,6 @@ def _recommend_axis_pid_gain(sysid: dict, axis: str, guide_exp: float,
     data yet (Options::aIPIDAutoTune() was off, or this mount class's
     protocol doesn't collect it).
     """
-    _, _, fit_info = fit_pulse_response(sysid, axis, guide_exp, False, return_fits=True)
-    fits = fit_info["fits"]
-
-    if fit_info["sign_consistent"] is False:
-        return {"confidence": "unavailable",
-                "reason": "pulse-response signs inconsistent across directions -- fits are noise, not mechanics"}
-    if not fits:
-        return {"confidence": "unavailable", "reason": "no usable pulse-response step-response fits"}
-
     cal_ms_per_arcsec = _calibration_ms_per_arcsec(sysid, axis)
     if cal_ms_per_arcsec <= 0.0:
         return {"confidence": "unavailable",
@@ -121,28 +171,39 @@ def _recommend_axis_pid_gain(sysid: dict, axis: str, guide_exp: float,
 
     pixel_scale = _effective_pixel_scale(sysid)
 
-    K_samples, tau_samples, L_samples, t_first_samples = [], [], [], []
-    for f in fits:
-        if f["pulse_magnitude_ms"] <= 0.0:
-            continue
-        K_samples.append(abs(f["P_fit_px"]) * pixel_scale / f["pulse_magnitude_ms"])
-        t_first_samples.append(f["t_first_s"])
-        L_samples.append(_estimate_dead_time_s(f["t_arr"], f["pos_arr"], f["residual_std_px"]))
-        if f["tau_fit_s"] <= 9.8:  # same "pinned at bound, degenerate with drift" guard as the spring fit
-            tau_samples.append(f["tau_fit_s"])
-
-    if not K_samples:
-        return {"confidence": "unavailable", "reason": "no fit had a usable pulse magnitude"}
-
-    K = float(np.median(K_samples))  # arcsec of steady-state response per ms of pulse
-    L = float(np.median(L_samples))
-    # tau can't be resolved below the sampling floor either: a fit whose tau_fit
-    # landed below t_first ("spring already released") is only known to be
-    # <= t_first, not physically ~0 -- floor it the same way L is floored, so a
-    # spuriously tiny fitted tau can't produce a dangerously aggressive Kc.
-    tau_raw = float(np.median(tau_samples)) if tau_samples else float(np.median(t_first_samples))
-    tau = max(tau_raw, L)
-    resolution_limited = bool(np.isclose(L, float(np.median(t_first_samples)), rtol=0.05))
+    steps = _step_samples_by_magnitude(sysid, axis)
+    by_mag = steps["by_mag"]
+    if not by_mag:
+        return {"confidence": "unavailable",
+                "reason": "no usable pulse-response step-response fits"}
+    if len(steps["signs"]) > 1:
+        return {"confidence": "unavailable",
+                "reason": "pulse-response signs inconsistent across directions -- fits are noise, not mechanics"}
+
+    n_fits = sum(len(v) for v in by_mag.values())
+    # K is arcsec per ms of pulse, so it must not depend on the magnitude used to
+    # measure it. Disagreement means the responses are contaminated and no median
+    # over the pooled samples would be meaningful.
+    k_by_mag = {mag: float(np.median(v)) * pixel_scale / mag for mag, v in by_mag.items()}
+    k_lo, k_hi = min(k_by_mag.values()), max(k_by_mag.values())
+    detail = ", ".join(f"{m:.0f}ms: {k:.5f}" for m, k in sorted(k_by_mag.items()))
+    if len(k_by_mag) < 2:
+        return {"confidence": "unavailable",
+                "reason": f"only one usable pulse magnitude ({detail}) -- the cross-magnitude check cannot run"}
+    if k_lo <= 0.0 or (k_hi - k_lo) / k_lo > K_MAGNITUDE_CONSISTENCY_TOLERANCE:
+        return {"confidence": "unavailable",
+                "reason": f"process gain disagrees across pulse magnitudes ({detail}) -- responses are contaminated"}
+
+    # The largest pulse has the best signal-to-contamination ratio.
+    best_mag = max(by_mag)
+    if len(by_mag[best_mag]) < MIN_FITS_FOR_MEDIUM_CONFIDENCE:
+        return {"confidence": "unavailable",
+                "reason": f"only {len(by_mag[best_mag])} usable fit(s) at {best_mag:.0f}ms"}
+    K = k_by_mag[best_mag]                      # arcsec of response per ms of pulse
+    L = float(np.median(steps["dead_times"][best_mag]))
+    t_first = float(np.median(steps["first_times"][best_mag]))
+    tau = max(t_first, L)
+    resolution_limited = bool(np.isclose(L, t_first, rtol=0.05))
 
     lam = max(tau, lambda_l_factor * L)
     Kc = (1.0 / K) * tau / (lam + L)               # ms of pulse per arcsec of error
@@ -150,11 +211,11 @@ def _recommend_axis_pid_gain(sysid: dict, axis: str, guide_exp: float,
 
     proportional_gain = Kc / cal_ms_per_arcsec
     integral_gain = INTEGRAL_GAIN_CONSERVATIVE_FRACTION * proportional_gain
-    confidence = "low" if (resolution_limited or len(fits) < MIN_FITS_FOR_MEDIUM_CONFIDENCE) else "medium"
+    confidence = "low" if (resolution_limited or n_fits < MIN_FITS_FOR_MEDIUM_CONFIDENCE) else "medium"
 
     if verbose:
         print(f"  [{axis} PID] K={K:.5f} arcsec/ms  tau={tau:.2f}s  L={L:.2f}s  lambda={lam:.2f}s "
-              f"(n={len(fits)} fits, cal={cal_ms_per_arcsec:.1f}ms/arcsec)")
+              f"(n={n_fits} fits @ {best_mag:.0f}ms, cal={cal_ms_per_arcsec:.1f}ms/arcsec)")
         print(f"  [{axis} PID] Recommended proportional_gain={proportional_gain:.3f}  "
               f"integral_gain={integral_gain:.3f}  confidence={confidence}")
 
@@ -168,7 +229,7 @@ def _recommend_axis_pid_gain(sysid: dict, axis: str, guide_exp: float,
         "lambda_s":                   lam,
         "tau_i_s":                    tau_i,
         "calibration_ms_per_arcsec":  cal_ms_per_arcsec,
-        "n_fits":                     len(fits),
+        "n_fits":                     n_fits,
         "resolution_limited":         resolution_limited,
     }
 
diff --git a/kstars/ekos/guide/opsaiconfig.ui b/kstars/ekos/guide/opsaiconfig.ui
index dfe1bee09f..1dafecb77d 100644
--- a/kstars/ekos/guide/opsaiconfig.ui
+++ b/kstars/ekos/guide/opsaiconfig.ui
@@ -79,7 +79,7 @@
          <string>PID Auto-Tune (Step-Response Gain Calibration)</string>
         </property>
         <property name="toolTip">
-         <string>Sends known pulses at the start of the protocol and measures the mount's step response to recommend a base RA/DEC guiding gain for any mount type (offline trainer only, never applied automatically). Adds roughly 10 minutes to the protocol; enabled by default.</string>
+         <string>Sends known pulses at the start of the protocol and measures the mount's step response to derive a base RA/DEC guiding gain for any mount type. The measured gain is applied and locked for the rest of the session; the current gain is kept if the responses do not pass the consistency checks. Adds roughly 20 minutes to the protocol; enabled by default.</string>
         </property>
        </widget>
       </item>
diff --git a/kstars/kstars.kcfg b/kstars/kstars.kcfg
index 5393127cbe..e0de84d097 100644
--- a/kstars/kstars.kcfg
+++ b/kstars/kstars.kcfg
@@ -3124,7 +3124,7 @@
       </entry>
       <entry name="AIPIDAutoTune" type="Bool">
          <label>Run PID Auto-Tune (step-response gain calibration) as part of the AI data collection protocol.</label>
-         <whatsthis>Sends known pulses at the start of the protocol and measures the mount's step response to recommend a base RA/DEC guiding gain for any mount type. Advisory only -- never applied automatically. Adds roughly 10 minutes to the protocol. Enabled by default.</whatsthis>
+         <whatsthis>Sends known pulses at the start of the protocol and measures the mount's step response to derive a base RA/DEC guiding gain for any mount type. The measured gain is applied and locked for the rest of the session, so the data that follows is collected under it; the current gain is kept unchanged if the responses do not pass the consistency checks. Adds roughly 20 minutes to the protocol. Enabled by default.</whatsthis>
          <default>true</default>
       </entry>
       <entry name="AIShadowMode" type="Bool">