feat: a de-noise math node rejects impulses with a median and a one-pole

One op, not two: a running median followed by a one-pole at `fc`. The median kills an impulse
outright -- a peak shorter than half its window cannot move it at all -- and the pole then
takes the broadband hash the median leaves behind. Chaining the existing `LowPass` after a
median would let the user wire them the wrong way round, and a pole *before* a median feeds
the impulse into the filter's state for its whole settling time, which is the thing the node
exists to avoid.

`DeNoise` is **appended** at `math_op_at(12)` (with `Div` given its explicit `case 11`, so the
`default` no longer swallows it): a session stores an op as its index through that accessor,
so inserting it beside the filters in the dropdown would reinterpret every node in every
existing session file. It sits last in the dropdown instead, which is the right side of that
trade.

The median window is `MathNode::med_n`, odd and 3..9, enforced by `median_n()` rather than by
the dialog -- a file holding anything else still loads -- and the ring is *pooled* in
`MathRuntime::median`, sized in `math_prepare()`, because nine doubles as a member would push
`MathNodeRt` across two cache lines for every node, filter or not. The pole reuses `b0`/`z1`,
which `math_prepare()` already carries across a re-prepare, so a decimation change does not
make its trace jump.

Warm-up takes the median of what the window holds rather than NaN: in this application NaN
means *a hole in the data*, and spending that meaning on a start-up transient would make a gap
ambiguous. A non-finite input stays a hole and leaves both the window and the pole alone.

src/math_channel.cpp

@@ -59,6 +59,28 @@ void biquad(MathOp op, double fc, double q, double dt, double& b0, double& b1, d
 
 /// One biquad step, transposed direct form II -- the form whose state words are the two things
 /// that actually need carrying between samples.
+/// One pole's coefficient for a corner of `fc` Hz at `dt` seconds a sample: `a = 1 - exp(-2*pi*fc*dt)`,
+/// clamped just under Nyquist so `a` cannot exceed 1. `fc == 0` returns 1, which is a pass-through -- so
+/// "no pole" costs the step no branch.
+double pole_a(double fc, double dt) {
+  if (!(fc > 0.0) || !(dt > 0.0)) return 1.0;
+  const double nyquist = 0.5 / dt;
+  const double use = fc < nyquist * 0.99 ? fc : nyquist * 0.99;
+  return 1.0 - std::exp(-2.0 * kPi * use * dt);
+}
+
+/// The median window's length: odd, and between 3 and 9.
+///
+/// Odd because an even window's median is the average of two samples, which is a value the signal never
+/// took -- and the whole point of a median here is that its output is a sample that was really there.
+/// Rounded *up*, so a request for 4 gets 5 rather than 3. The ceiling of 9 is what `math_step`'s fixed
+/// `double sorted[9]` is sized for.
+int median_n(int want) {
+  int n = want < 3 ? 3 : (want > 9 ? 9 : want);
+  if ((n & 1) == 0) ++n;
+  return n;
+}
+
 double biquad_step(double x, double b0, double b1, double b2, double a1, double a2, double& z1,
                    double& z2) {
   // A reserved hole is NaN in the source, and a NaN integrated into the state words is
@@ -90,11 +112,12 @@ const char* math_op_name(MathOp op) {
     case MathOp::Sub: return "difference";
     case MathOp::Mul: return "product";
     case MathOp::Div: return "ratio";
+    case MathOp::DeNoise: return "de-noise";
   }
   return "?";
 }
 
-int math_op_count() { return 12; }
+int math_op_count() { return 13; }
 
 MathOp math_op_at(int index) {
   switch (index) {
@@ -109,7 +132,8 @@ MathOp math_op_at(int index) {
     case 8: return MathOp::Sum;
     case 9: return MathOp::Sub;
     case 10: return MathOp::Mul;
-    default: return MathOp::Div;
+    case 11: return MathOp::Div;
+    default: return MathOp::DeNoise;
   }
 }
 
@@ -136,6 +160,8 @@ bool math_op_uses_filter(MathOp op) {
   return op == MathOp::LowPass || op == MathOp::HighPass || op == MathOp::BandPass;
 }
 bool math_op_uses_clamp(MathOp op) { return op == MathOp::Clamp; }
+bool math_op_uses_fc(MathOp op) { return math_op_uses_filter(op) || op == MathOp::DeNoise; }
+bool math_op_uses_median(MathOp op) { return op == MathOp::DeNoise; }
 
 void math_prepare(const MathChannel& channel, MathRuntime& rt, double dt,
                   const std::function<int(const std::string&)>& index_of) {
@@ -152,6 +178,17 @@ void math_prepare(const MathChannel& channel, MathRuntime& rt, double dt,
   }
   rt.dt = dt;
   rt.ready = true;
+  // The median windows, pooled: one pass to size the pool, then the node loop hands out the slices. A
+  // window starts from rest after a re-prepare -- the same choice the biquads' `z2` makes -- while the
+  // pole's `z1` is carried, so a decimation change does not make the trace jump.
+  std::size_t med_total = 0;
+  for (std::size_t i = 0; i < n; ++i) {
+    if (channel.nodes[i].op == MathOp::DeNoise) {
+      med_total += static_cast<std::size_t>(median_n(channel.nodes[i].med_n));
+    }
+  }
+  rt.median.assign(med_total, 0.0);
+  std::size_t med_next = 0;
   for (std::size_t i = 0; i < n; ++i) {
     const MathNode& node = channel.nodes[i];
     MathNodeRt& out = rt.node[i];
@@ -159,6 +196,15 @@ void math_prepare(const MathChannel& channel, MathRuntime& rt, double dt,
       out.src = index_of ? index_of(node.source) : -1;
       continue;
     }
+    if (node.op == MathOp::DeNoise) {
+      out.med_n = median_n(node.med_n);
+      out.med_at = med_next;
+      out.med_pos = 0;
+      out.med_have = 0;
+      med_next += static_cast<std::size_t>(out.med_n);
+      out.b0 = pole_a(node.fc, dt);
+      continue;
+    }
     if (!math_op_uses_filter(node.op)) continue;
     biquad(node.op, node.fc, node.q, dt, out.b0, out.b1, out.b2, out.a1, out.a2);
   }
@@ -221,6 +267,46 @@ double math_step(const MathChannel& channel, MathRuntime& rt, const float* value
         x = std::fabs(d) > 0.0 ? at(node.in_a) / d : 0.0;
         break;
       }
+      case MathOp::DeNoise: {
+        const double in = at(node.in_a);
+        // A hole stays a hole, and the window is left alone: folding a NaN into the median would smear one
+        // gap across `med_n` samples, and the whole trace after it in the pole. `med_n <= 0` is a node the
+        // runtime has not prepared -- the window pointer below would be a read past the pool, so it is
+        // checked rather than assumed.
+        if (!std::isfinite(in) || self.med_n <= 0) {
+          x = kNaN;
+          break;
+        }
+        const std::size_t len = static_cast<std::size_t>(self.med_n);
+        double* win = rt.median.data() + self.med_at;
+        win[self.med_pos] = in;
+        self.med_pos = self.med_pos + 1 == len ? 0 : self.med_pos + 1;
+        if (self.med_have < len) ++self.med_have;
+        // An insertion sort of at most nine doubles into a local. The window is tiny and fixed, so this
+        // beats anything that allocates, and it leaves the ring itself in arrival order -- which is what
+        // lets the oldest sample be overwritten in place. Before the window has filled, the median is taken
+        // over what there is: a shorter window, deliberately, rather than NaN. A derived trace that reads
+        // "no data" for its first few samples after every rebuild looks like a broken graph.
+        double sorted[9];
+        const std::size_t m = self.med_have;
+        for (std::size_t j = 0; j < m; ++j) {
+          const double v = win[j];
+          std::size_t p = j;
+          while (p > 0 && sorted[p - 1] > v) {
+            sorted[p] = sorted[p - 1];
+            --p;
+          }
+          sorted[p] = v;
+        }
+        const double med = sorted[m / 2];
+        // Median *then* pole, and the order is the design. The median kills an impulse outright; the pole
+        // then takes the broadband hash it leaves behind. The other way round feeds the pole the impulse,
+        // which smears it over the pole's whole settling time and leaves the median a wide bump it cannot
+        // reject.
+        self.z1 += self.b0 * (med - self.z1);
+        x = self.z1;
+        break;
+      }
     }
     self.value = x;
   }

src/math_channel.hpp

@@ -25,6 +25,10 @@ enum class MathOp {
   Sub,       ///< a - b
   Mul,       ///< a * b
   Div,       ///< a / b, zero when b is 0 (a scope trace must not go to inf/NAN)
+  /// A running median of `med_n` samples, then a one-pole at `fc` Hz. For impulsive interference -- an EM
+  /// peak a few samples wide -- which a linear filter can only smear over its own settling time. `fc = 0`
+  /// is the median alone.
+  DeNoise,
 };
 
 const char* math_op_name(MathOp op);
@@ -41,6 +45,11 @@ int math_op_arity(MathOp op);
 bool math_op_uses_gain(MathOp op);
 bool math_op_uses_filter(MathOp op);
 bool math_op_uses_clamp(MathOp op);
+/// Whether an op reads `fc`. Wider than `math_op_uses_filter()`, which keeps its present meaning of *a
+/// biquad*: `DeNoise` reads `fc` for a one-pole and has no `q`.
+bool math_op_uses_fc(MathOp op);
+/// Whether an op reads `med_n`.
+bool math_op_uses_median(MathOp op);
 
 /// One node of a math channel's graph. Node indices are positions in `MathChannel::nodes`, and a
 /// node may only read ones before it: the evaluator walks the list in order, which is what makes a
@@ -56,6 +65,10 @@ struct MathNode {
   double q = 0.707;     ///< BP quality factor.
   double lo = -1.0;     ///< `Clamp` bounds.
   double hi = 1.0;
+  /// `DeNoise`: the median window, in samples. Odd and between 3 and 9 -- `median_n()` in the
+  /// implementation is what enforces that, so a session file holding anything else still loads. An impulse
+  /// up to `med_n / 2` samples wide is rejected outright.
+  int med_n = 5;
   double x = 0.0;       ///< the editor's box position, logical pixels.
   double y = 0.0;
 
@@ -67,7 +80,7 @@ struct MathNode {
   bool operator==(const MathNode& other) const {
     return op == other.op && source == other.source && in_a == other.in_a && in_b == other.in_b &&
            gain == other.gain && offset == other.offset && fc == other.fc && q == other.q &&
-           lo == other.lo && hi == other.hi;
+           lo == other.lo && hi == other.hi && med_n == other.med_n;
   }
 };
 
@@ -110,10 +123,25 @@ struct MathNodeRt {
   /// per source node, per sample, behind a type-erased call. Names cannot change without a
   /// rebuild, and every rebuild re-prepares.
   int src = -1;
+  /// `DeNoise`: where this node's median window starts in `MathRuntime::median`, how long it is, the ring's
+  /// write position, and how many samples it holds. The window is *pooled* rather than held here because a
+  /// node's state belongs on one cache line, and nine doubles would push this struct across two for every
+  /// node in the graph, filter or not.
+  ///
+  /// The one-pole reuses `b0` for its coefficient and `z1` for its value, which is not a shortcut:
+  /// `math_prepare()` deliberately carries `z1` across a re-prepare, so the pole keeps its history through
+  /// a decimation change exactly as the biquads keep theirs.
+  std::size_t med_at = 0;
+  int med_n = 0;
+  std::size_t med_pos = 0;
+  std::size_t med_have = 0;
 };
 
 struct MathRuntime {
   std::vector<MathNodeRt> node;
+  /// Every `DeNoise` node's median window, end to end: one allocation for the channel, sized by
+  /// `math_prepare()`, a node's slice starting at its own `med_at`.
+  std::vector<double> median;
   /// The rate the coefficients were built for; a change rebuilds them.
   double dt = 0.0;
   bool ready = false;