diff --git a/CMakeLists.txt b/CMakeLists.txt index 65434a6..c83f6f5 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -1112,13 +1112,14 @@ add_test(NAME master_gain_tests COMMAND master_gain_tests) # filter: four targets along the module's own seams, so each asserts one domain. # filter_params_tests — the rate-free control mappings (cutoff/Q/drive) and their inverses. # filter_morph_tests — the pure morph-weight algebra under both morph laws; no DSP is run. -# filter_state_tests — numerical stability, the denormal flush, and the state lifecycle. +# filter_state_tests — numerical stability, the denormal flush, bounded-output/self-oscillation +# under full drive, and the state lifecycle — none of it needs the measurement harness below. # filter_tests — the frequency response: pins the SVF coefficients against an independent # derivation, holds the morph endpoints to the analytic 2-pole targets, and measures the # HP-BP-LP corner flatness, the HP-notch-LP null, and rate/level invariance and drive # stability by driving real sines. The seams above were chosen so this file alone owns the # analytic reference and the steady-state gain measurement — a forked copy of a measurement -# reference is a worse defect than a long file, which is why it sits over the ~600-line bar. +# reference is a worse defect than a long file. # NEITHER SDK. add_executable(filter_params_tests tests/test_filter_params.cpp) target_link_libraries(filter_params_tests PRIVATE filter) diff --git a/src/core/instrument/engine/filter/CLAUDE.md b/src/core/instrument/engine/filter/CLAUDE.md index 3c016e3..d4577eb 100644 --- a/src/core/instrument/engine/filter/CLAUDE.md +++ b/src/core/instrument/engine/filter/CLAUDE.md @@ -90,10 +90,15 @@ persisted field lands on it rather than on the SEM leg. and LP **together** across the whole sweep, `bp == 0` throughout. The notch is not tuned in: HP and LP sit at exactly +90° and −90° at the corner, so equal weights cancel there by construction. Here the corner magnitude deliberately goes to **zero** at the centre — - measured worst case −88 dB across every rate/cutoff/Q, typically −110 to −145 dB. The fold - makes that structural rather than a runtime near-miss: `m2 = lp - hp` is **exactly** `0.0f` - at the centre, because `cos` and `sin` of π/4 differ by about an ulp of *double*, nine - orders below float's spacing there, so they narrow to one float. + measured worst case −88 dB on the shipped `{250, 1000, 4000}` Hz cutoff grid, typically −110 to + −145 dB. Over the full control range (20 Hz – 20 kHz, Q 0.1 – 10) the worst residual is + shallower — −69.8 dB at 192 kHz / 30 Hz / Q=10 — from float conditioning in the folded + `x − k·v1` term as `fc/sr → 1e-4` at high Q; it is Q-dependent (Q=0.1 holds −110 dB everywhere) + and still an excellent notch, not a broadband defect. `test_filter.cpp`'s null test covers this + full range with a Q-scaled threshold rather than the flat −74 dB the shipped grid alone would + justify. The fold makes the centre's cancellation structural rather than a runtime near-miss: + `m2 = lp - hp` is **exactly** `0.0f` at the centre, because `cos` and `sin` of π/4 differ by + about an ulp of *double*, nine orders below float's spacing there, so they narrow to one float. SEM's zero is at the **notch frequency**, not a broadband level sag — off the corner the pair is still equal-power, so neither law's legs dip. Measuring that requires dividing by each diff --git a/src/core/instrument/engine/filter/filter_saturate.h b/src/core/instrument/engine/filter/filter_saturate.h index b78fc95..5c2425d 100644 --- a/src/core/instrument/engine/filter/filter_saturate.h +++ b/src/core/instrument/engine/filter/filter_saturate.h @@ -11,7 +11,9 @@ namespace reasampler::instrument::engine::filter { // // Three properties are load-bearing and none of them are tuning: // - depth == 0 makes this ALGEBRAICALLY the identity (x / sqrt(1) == x, exact in IEEE), so -// drive = 0 is bit-exact linear with no branch and no special case on the hot path. +// drive = 0 is bit-exact linear whether or not the caller special-cases it. (voice_filter.h +// gates the call on drive != 0 anyway, but as a perf optimization, not because correctness +// needs it.) // - |softLimit(x, d)| <= |x| for every d, so dropping it into the resonance state update can // only ever shrink the state. The filter therefore cannot gain energy from the drive stage: // stability at any Q and any cutoff is structural, not a tuned margin, and it can never diff --git a/tests/test_filter.cpp b/tests/test_filter.cpp index 8a4bbf4..226884e 100644 --- a/tests/test_filter.cpp +++ b/tests/test_filter.cpp @@ -222,21 +222,31 @@ static void testCornerMagnitudeIsFlatAtQAcrossTheHighBandLowSweep() { // The SEM's centre is a genuine null, not merely a dip: the corner magnitude falls to the float // noise floor because HP and LP sit at exactly +90 and -90 degrees there, so equal weights cancel -// by construction. Measured worst case across this whole grid is 3.8e-05 (-88 dB); the typical -// figure is -110 to -145 dB. The settle window has to clear the resonator's ring-down before the -// residual means anything — at 0.15 s and Q=10 the leftover transient alone reads as -52 dB and -// would be mistaken for the floor. +// by construction. Grid spans the full control range (20 Hz - 20 kHz), not just three interior +// cutoffs: the residual is worse near the low-cutoff/high-rate corner (float conditioning in the +// folded x - k*v1 term as fc/sr -> 1e-4 at high Q) and is Q-dependent, so the threshold scales +// with Q rather than repeating a flat bound sized off the shallow grid. Measured worst case on +// this wider grid: 2.6e-06 (-111.7 dB) at Q=0.1, 7.0e-05 (-83.1 dB) at Q=sqrt(2), 3.2e-04 +// (-69.8 dB) at Q=10, all at 192 kHz / 30 Hz — still an excellent notch, not a broadband defect. +// The settle window has to clear the resonator's ring-down before the residual means anything — +// at 0.15 s and Q=10 the leftover transient alone reads as -52 dB and would be mistaken for the +// floor. static void testHighNotchLowCentreIsATrueNullAtTheCorner() { for (int r = 0; r < kRateCount; ++r) { - for (double fc : {250.0, 1000.0, 4000.0}) { + for (double fc : {20.0, 30.0, 50.0, 250.0, 1000.0, 4000.0, 16000.0, 20000.0}) { + if (fc > kRates[r] * 0.48) continue; for (float res : {0.0f, 0.5f, 1.0f}) { + const double q = filterQFromNorm(res); + // Sized against measurement (margins 6.6x/1.55x/2.2x at Q=0.1/sqrt(2)/10 on this + // grid), not copied from the corner figure alone. + const double threshold = 1e-5 + 7e-5 * q; const double got = measuredGain(at(fc, res, kCentre, 0.0f, MorphLaw::HighNotchLow), kRates[r], fc, 0.25, 2.0, 0.5); - if (!(got < 2e-4)) { + if (!(got < threshold)) { std::printf("FAIL line %d: SEM notch at sr %.0f fc %.0f res %.1f is %.3e " - "(%.1f dB) — not a null\n", + "(%.1f dB) — not a null (threshold %.3e)\n", __LINE__, kRates[r], fc, res, got, - 20.0 * std::log10(got + 1e-300)); + 20.0 * std::log10(got + 1e-300), threshold); ++g_fail; } } @@ -451,82 +461,31 @@ static void testDriveZeroIsBitIdenticalToTheLinearKernel() { // The complaint the rewrite answers: resonance must not track how hard the sample hits the // filter unless the user asked for it. At drive 0 the response is identical over a 1000:1 level -// range; the tap this replaced moved by 14% over the same span. +// range; the tap this replaced moved by 14% over the same span. Runs under both laws; the centre +// is skipped under HighNotchLow because analyticMag has no notch formula to compare against there +// — level invariance at drive 0 is structural for any linear combination of the SVF's taps, so +// skipping one morph position on one law loses no real coverage. static void testDriveZeroResponseIsLevelInvariant() { const double sr = 48000.0, fc = 1000.0; - for (float morph : {kHighPass, kBandPass, kLowPass}) { - const double q = filterQFromNorm(1.0f); - const double want = analyticMag(morph, fc, fc, q, sr); - for (double amp : {0.001, 0.01, 0.1, 1.0}) { - const double got = measuredGain(at(fc, 1.0f, morph), sr, fc, amp); - if (!(std::fabs(got / want - 1.0) <= kAgreement)) { - std::printf("FAIL line %d: morph %.1f amp %g gain %.6f vs analytic %.6f " - "(%.3f%%)\n", - __LINE__, morph, amp, got, want, (got / want - 1.0) * 100.0); - ++g_fail; - } - } - } -} - -// Drive is bounded by construction, not by tuning: softLimit is a contraction, so the state -// update can only ever shrink the state and the filter cannot gain energy from it. This sweeps -// the corners that would expose a tuned margin instead. -static void testFullDriveStaysBoundedAtEveryCutoffResonanceAndRate() { - unsigned rng = 0x2468aceu; - auto noise = [&rng]() { - rng = rng * 1664525u + 1013904223u; - return static_cast(static_cast(rng >> 9) - (1 << 22)) / - static_cast(1 << 22); - }; - - for (int r = 0; r < kRateCount; ++r) { - const double sr = kRates[r]; - for (MorphLaw law : kBothLaws) { - for (int ci = 0; ci <= 8; ++ci) { - for (int mi = 0; mi <= 4; ++mi) { - for (float res : {0.0f, 0.5f, 1.0f}) { - VoiceFilter f; - f.prepare({ci / 8.0f, res, mi / 4.0f, 1.0f, law}, sr); - f.reset(); - for (int i = 0; i < 4000; ++i) { - const float y = f.process(0, noise()); - if (!std::isfinite(y) || std::fabs(y) > 8.0f) { - std::printf("FAIL line %d: %s sr=%.0f cutoff=%.2f morph=%.2f " - "res=%.1f full drive produced %g\n", - __LINE__, lawName(law), sr, ci / 8.0, mi / 4.0, res, y); - ++g_fail; - return; - } - } - } + for (MorphLaw law : kBothLaws) { + for (float morph : {kHighPass, kBandPass, kLowPass}) { + if (morph == kBandPass && law != MorphLaw::HighBandLow) continue; + const double q = filterQFromNorm(1.0f); + const double want = analyticMag(morph, fc, fc, q, sr); + for (double amp : {0.001, 0.01, 0.1, 1.0}) { + const double got = measuredGain(at(fc, 1.0f, morph, 0.0f, law), sr, fc, amp); + if (!(std::fabs(got / want - 1.0) <= kAgreement)) { + std::printf("FAIL line %d: %s morph %.1f amp %g gain %.6f vs analytic %.6f " + "(%.3f%%)\n", + __LINE__, lawName(law), morph, amp, got, want, + (got / want - 1.0) * 100.0); + ++g_fail; } } } } } -// Full drive at full resonance with no input must still go quiet. A nonlinearity in the loop is -// exactly where a self-oscillator would hide, and softLimit's sub-unit slope is what forbids it. -static void testFullDriveDoesNotSelfOscillate() { - for (int r = 0; r < kRateCount; ++r) { - const double sr = kRates[r]; - for (MorphLaw law : kBothLaws) { - for (float morph : {kHighPass, kBandPass, kLowPass}) { - VoiceFilter f; - f.prepare(at(1000.0, 1.0f, morph, 1.0f, law), sr); - f.reset(); - const int excite = static_cast(sr * 0.01); - for (int i = 0; i < excite; ++i) { - f.process(0, static_cast(std::sin(2.0 * kPi * 1000.0 * i / sr))); - } - for (int i = 0; i < static_cast(sr * 0.5); ++i) f.process(0, 0.0f); - CHECK(f.isSilent()); - } - } - } -} - // Drive has to actually do something at the top of its travel, and do it monotonically — the // brief's "extreme, not politely warm". Measured at the corner, where the resonance state is // what the limiter sees. @@ -548,56 +507,33 @@ static void testDriveCompressesTheResonantPeakMonotonically() { } } -static void testSoftLimitIsOddMonotoneBoundedAndExactAtZeroDepth() { - for (double x : {-3.0, -0.5, 0.0, 1e-9, 0.25, 7.0}) { - // Depth 0 is the identity by algebra, so drive 0 needs no special case on the hot path. - CHECK(softLimit(static_cast(x), 0.0f) == static_cast(x)); - } - CHECK_NEAR(softLimit(1.5f, 2.0f), -softLimit(-1.5f, 2.0f), 1e-9); - - for (float depth : {0.5f, 4.0f, 64.0f}) { - // The two properties the stability argument rests on, over the whole excursion range a - // resonating state can reach. Monotonicity is NOT asserted here: far past the knee the - // curve is asymptotically flat, so the true increment between adjacent samples falls - // below float epsilon and rounding can walk it backwards by an ulp. - for (int i = -400; i <= 400; ++i) { - const float x = static_cast(i) * 0.05f; - const float y = softLimit(x, depth); - CHECK(std::fabs(y) <= std::fabs(x)); // a contraction — the stability argument - CHECK(std::fabs(y) < 1.0f / depth + 1e-6f); // bounded by the knee - } - // Strictly increasing across the knee, which is where the shaping actually happens. - const float knee = 1.0f / depth; - float prev = -1e30f; - for (int i = -20; i <= 20; ++i) { - const float y = softLimit(static_cast(i) * 0.1f * knee, depth); - CHECK(y > prev); - prev = y; - } - } -} - // --------------------------------------------------------------------------- // Sample-rate invariance // --------------------------------------------------------------------------- // The rate must enter only through g = tan(pi*fc/sr), so the response at a given cutoff and Q is // the same filter at every rate. The retired feedback tap made this false: it closed the loop -// once per SAMPLE, so emphasis ran 5.02 at 48k against 8.52 at 192k. +// once per SAMPLE, so emphasis ran 5.02 at 48k against 8.52 at 192k. Runs under both laws; the +// centre is skipped under HighNotchLow because analyticMag has no notch formula to compare +// against there — SEM centre behavior across rates is covered by +// testHighNotchLowCentreIsATrueNullAtTheCorner instead. static void testResponseIsRateInvariantAtEveryMorph() { - for (float morph : {kHighPass, kBandPass, kLowPass}) { - for (float res : {0.2f, 0.5f, 1.0f}) { - const double q = filterQFromNorm(res); - for (double fc : {250.0, 1000.0, 4000.0}) { - for (int r = 0; r < kRateCount; ++r) { - const double got = measuredGain(at(fc, res, morph), kRates[r], fc); - const double want = analyticMag(morph, fc, fc, q, kRates[r]); - if (!(std::fabs(got / want - 1.0) <= kAgreement)) { - std::printf("FAIL line %d: morph %.1f res %.1f fc %.0f at %.0f Hz: %.6f " - "vs analytic %.6f (%.3f%%)\n", - __LINE__, morph, res, fc, kRates[r], got, want, - (got / want - 1.0) * 100.0); - ++g_fail; + for (MorphLaw law : kBothLaws) { + for (float morph : {kHighPass, kBandPass, kLowPass}) { + if (morph == kBandPass && law != MorphLaw::HighBandLow) continue; + for (float res : {0.2f, 0.5f, 1.0f}) { + const double q = filterQFromNorm(res); + for (double fc : {250.0, 1000.0, 4000.0}) { + for (int r = 0; r < kRateCount; ++r) { + const double got = measuredGain(at(fc, res, morph, 0.0f, law), kRates[r], fc); + const double want = analyticMag(morph, fc, fc, q, kRates[r]); + if (!(std::fabs(got / want - 1.0) <= kAgreement)) { + std::printf("FAIL line %d: %s morph %.1f res %.1f fc %.0f at %.0f Hz: " + "%.6f vs analytic %.6f (%.3f%%)\n", + __LINE__, lawName(law), morph, res, fc, kRates[r], got, want, + (got / want - 1.0) * 100.0); + ++g_fail; + } } } } @@ -646,10 +582,7 @@ int main() { testDriveZeroIsBitIdenticalToTheLinearKernel(); testDriveZeroResponseIsLevelInvariant(); - testFullDriveStaysBoundedAtEveryCutoffResonanceAndRate(); - testFullDriveDoesNotSelfOscillate(); testDriveCompressesTheResonantPeakMonotonically(); - testSoftLimitIsOddMonotoneBoundedAndExactAtZeroDepth(); testResponseIsRateInvariantAtEveryMorph(); testLowCutoffHighRateCornerHoldsTheAnalyticPeak(); diff --git a/tests/test_filter_morph.cpp b/tests/test_filter_morph.cpp index aec6601..ab2a012 100644 --- a/tests/test_filter_morph.cpp +++ b/tests/test_filter_morph.cpp @@ -1,7 +1,9 @@ // Standalone tests for the pure morph domain: normalized position -> tap weights under both // morph laws, and the fold of those weights into the kernel's three multipliers. Algebra only — -// no filter is run here. Interior expectations are derived from the intended law (in radicals, -// so they share not even a trig call with the implementation) rather than read back out of it. +// no filter is run here. Most interior expectations are derived from the intended law in radicals, +// sharing not even a trig call with the implementation; one check evaluates std::cos/std::sin +// directly at the same argument the implementation does, but a radical-derived check of the same +// leg sits right beside it, so no coverage rests solely on the shared call. // The MEASURED consequences of each law — HP-BP-LP's flat corner, HP-notch-LP's null — live in // test_filter.cpp, where a filter is actually driven. diff --git a/tests/test_filter_state.cpp b/tests/test_filter_state.cpp index 91d7d7b..d690599 100644 --- a/tests/test_filter_state.cpp +++ b/tests/test_filter_state.cpp @@ -1,12 +1,14 @@ // Standalone tests for the running filter's NUMERICAL behaviour and state lifecycle — bounded -// output under a live parameter sweep, the denormal flush, DC handling, the impulse response -// against the coefficients, and reset/prepare/per-channel state rules. Split from test_filter.cpp -// along the one seam that costs nothing: none of these need the frequency-response measurement -// harness, so the analytic reference lives in exactly one file and cannot fork. +// output under a live parameter sweep, full-drive stability and self-oscillation, the softLimit +// shaper's own properties, the denormal flush, DC handling, the impulse response against the +// coefficients, and reset/prepare/per-channel state rules. Split from test_filter.cpp along the +// one seam that costs nothing: none of these need the frequency-response measurement harness, so +// the analytic reference lives in exactly one file and cannot fork. #include "../src/core/instrument/engine/filter/filter_coeffs.h" #include "../src/core/instrument/engine/filter/filter_morph.h" #include "../src/core/instrument/engine/filter/filter_params.h" +#include "../src/core/instrument/engine/filter/filter_saturate.h" #include "../src/core/instrument/engine/filter/voice_filter.h" #include @@ -32,6 +34,9 @@ static constexpr float kCentre = 0.5f; static constexpr float kLowPass = 1.0f; static const MorphLaw kBothLaws[] = {MorphLaw::HighBandLow, MorphLaw::HighNotchLow}; +static const char* lawName(MorphLaw law) { + return law == MorphLaw::HighBandLow ? "HP-BP-LP" : "HP-notch-LP"; +} // The rates the invariance claims are made over. static const double kRates[] = {44100.0, 48000.0, 88200.0, 96000.0, 192000.0}; @@ -77,11 +82,102 @@ static void testFullRangeCutoffSweepAtAudioRateStaysBounded() { } } -// The flush tests the ENVELOPE — both integrators — not one sample. ic1 and ic2 are in -// quadrature, so a resonator swings each through zero twice a cycle; flushing on a single one -// injects a step in phase with the resonance, which the resonance amplifies, and the filter -// limit-cycles at the floor forever instead of going quiet. Re-verified for TPT rather than -// assumed to carry over from the retired Direct Form I state. +// Drive is bounded by construction, not by tuning: softLimit is a contraction, so the state +// update can only ever shrink the state and the filter cannot gain energy from it. This sweeps +// the corners that would expose a tuned margin instead. +static void testFullDriveStaysBoundedAtEveryCutoffResonanceAndRate() { + unsigned rng = 0x2468aceu; + auto noise = [&rng]() { + rng = rng * 1664525u + 1013904223u; + return static_cast(static_cast(rng >> 9) - (1 << 22)) / + static_cast(1 << 22); + }; + + for (int r = 0; r < kRateCount; ++r) { + const double sr = kRates[r]; + for (MorphLaw law : kBothLaws) { + for (int ci = 0; ci <= 8; ++ci) { + for (int mi = 0; mi <= 4; ++mi) { + for (float res : {0.0f, 0.5f, 1.0f}) { + VoiceFilter f; + f.prepare({ci / 8.0f, res, mi / 4.0f, 1.0f, law}, sr); + f.reset(); + for (int i = 0; i < 4000; ++i) { + const float y = f.process(0, noise()); + if (!std::isfinite(y) || std::fabs(y) > 8.0f) { + std::printf("FAIL line %d: %s sr=%.0f cutoff=%.2f morph=%.2f " + "res=%.1f full drive produced %g\n", + __LINE__, lawName(law), sr, ci / 8.0, mi / 4.0, res, y); + ++g_fail; + return; + } + } + } + } + } + } + } +} + +// Full drive at full resonance with no input must still go quiet. A nonlinearity in the loop is +// exactly where a self-oscillator would hide, and softLimit's sub-unit slope is what forbids it. +static void testFullDriveDoesNotSelfOscillate() { + for (int r = 0; r < kRateCount; ++r) { + const double sr = kRates[r]; + for (MorphLaw law : kBothLaws) { + for (float morph : {kHighPass, kBandPass, kLowPass}) { + VoiceFilter f; + f.prepare(at(1000.0, 1.0f, morph, 1.0f, law), sr); + f.reset(); + const int excite = static_cast(sr * 0.01); + for (int i = 0; i < excite; ++i) { + f.process(0, static_cast(std::sin(2.0 * kPi * 1000.0 * i / sr))); + } + for (int i = 0; i < static_cast(sr * 0.5); ++i) f.process(0, 0.0f); + CHECK(f.isSilent()); + } + } + } +} + +// The softLimit shaper's own properties, independent of any running filter: bit-exact at depth 0, +// odd, a contraction over the whole excursion range, bounded by the knee, and increasing where the +// shaping actually happens. +static void testSoftLimitIsOddMonotoneBoundedAndExactAtZeroDepth() { + for (double x : {-3.0, -0.5, 0.0, 1e-9, 0.25, 7.0}) { + // Depth 0 is the identity by algebra, so correctness doesn't require special-casing it — + // voice_filter.h gates the call anyway, but as a perf optimization (see its comment). + CHECK(softLimit(static_cast(x), 0.0f) == static_cast(x)); + } + CHECK_NEAR(softLimit(1.5f, 2.0f), -softLimit(-1.5f, 2.0f), 1e-9); + + for (float depth : {0.5f, 4.0f, 64.0f}) { + // The two properties the stability argument rests on, over the whole excursion range a + // resonating state can reach. Monotonicity is NOT asserted here: far past the knee the + // curve is asymptotically flat, so the true increment between adjacent samples falls + // below float epsilon and rounding can walk it backwards by an ulp. + for (int i = -400; i <= 400; ++i) { + const float x = static_cast(i) * 0.05f; + const float y = softLimit(x, depth); + CHECK(std::fabs(y) <= std::fabs(x)); // a contraction — the stability argument + CHECK(std::fabs(y) < 1.0f / depth + 1e-6f); // bounded by the knee + } + // Strictly increasing across the knee, which is where the shaping actually happens. + const float knee = 1.0f / depth; + float prev = -1e30f; + for (int i = -20; i <= 20; ++i) { + const float y = softLimit(static_cast(i) * 0.1f * knee, depth); + CHECK(y > prev); + prev = y; + } + } +} + +// The flush tests the ENVELOPE — both integrators — not one sample: isSilent() means "both are +// exactly zero," so both have to reach zero for that check to mean anything, and the conjunctive +// test is the cheapest guarantee of that (see voice_filter.h's flush comment). The stronger +// limit-cycle rationale belongs to the retired Direct Form I state, where the flushed variables +// were the actual filter OUTPUT rather than integrator state — it does not reproduce here. static void checkFlushGoesSilent(double sr, float morph, float drive, MorphLaw law) { // The decay to the floor is a fixed WALL-CLOCK time, so the budget scales with the rate. const int budget = static_cast(sr * 0.5); @@ -234,6 +330,9 @@ static void testChannelStateIsIndependent() { int main() { testFullRangeCutoffSweepAtAudioRateStaysBounded(); + testFullDriveStaysBoundedAtEveryCutoffResonanceAndRate(); + testFullDriveDoesNotSelfOscillate(); + testSoftLimitIsOddMonotoneBoundedAndExactAtZeroDepth(); testStateFlushesToZeroWithoutStallingInDenormals(); testHighPassSustainedDCDoesNotReRing(); testImpulseResponseMatchesTheKernel();