From 96cf4d9a06c9afaacbd91059c60fc68a1a3dd8eb Mon Sep 17 00:00:00 2001 From: Jan Svabenik Date: Tue, 24 Mar 2026 18:06:42 +0100 Subject: [PATCH] debug: use stateful decimating FIR in pre-demod path --- internal/dsp/decimating_fir.go | 81 +++++++++++++++++++++++++++++ internal/dsp/decimating_fir_test.go | 57 ++++++++++++++++++++ internal/recorder/streamer.go | 20 ++++--- 3 files changed, 151 insertions(+), 7 deletions(-) create mode 100644 internal/dsp/decimating_fir.go create mode 100644 internal/dsp/decimating_fir_test.go diff --git a/internal/dsp/decimating_fir.go b/internal/dsp/decimating_fir.go new file mode 100644 index 0000000..9b0358d --- /dev/null +++ b/internal/dsp/decimating_fir.go @@ -0,0 +1,81 @@ +package dsp + +// StatefulDecimatingFIRComplex combines FIR filtering and decimation into a +// single stateful stage. This avoids exposing FIR settling/transient output as +// ordinary block-leading samples before decimation. +type StatefulDecimatingFIRComplex struct { + taps []float64 + delayR []float64 + delayI []float64 + factor int + phase int // number of input samples until next output sample (0 => emit now) +} + +func NewStatefulDecimatingFIRComplex(taps []float64, factor int) *StatefulDecimatingFIRComplex { + if factor < 1 { + factor = 1 + } + t := make([]float64, len(taps)) + copy(t, taps) + return &StatefulDecimatingFIRComplex{ + taps: t, + delayR: make([]float64, len(taps)), + delayI: make([]float64, len(taps)), + factor: factor, + phase: 0, + } +} + +func (f *StatefulDecimatingFIRComplex) Reset() { + for i := range f.delayR { + f.delayR[i] = 0 + f.delayI[i] = 0 + } + f.phase = 0 +} + +func (f *StatefulDecimatingFIRComplex) Process(iq []complex64) []complex64 { + if len(iq) == 0 || len(f.taps) == 0 { + return nil + } + if f.factor <= 1 { + out := make([]complex64, len(iq)) + for i := 0; i < len(iq); i++ { + copy(f.delayR[1:], f.delayR[:len(f.taps)-1]) + copy(f.delayI[1:], f.delayI[:len(f.taps)-1]) + f.delayR[0] = float64(real(iq[i])) + f.delayI[0] = float64(imag(iq[i])) + var accR, accI float64 + for k := 0; k < len(f.taps); k++ { + w := f.taps[k] + accR += f.delayR[k] * w + accI += f.delayI[k] * w + } + out[i] = complex(float32(accR), float32(accI)) + } + return out + } + + out := make([]complex64, 0, len(iq)/f.factor+1) + n := len(f.taps) + for i := 0; i < len(iq); i++ { + copy(f.delayR[1:], f.delayR[:n-1]) + copy(f.delayI[1:], f.delayI[:n-1]) + f.delayR[0] = float64(real(iq[i])) + f.delayI[0] = float64(imag(iq[i])) + + if f.phase == 0 { + var accR, accI float64 + for k := 0; k < n; k++ { + w := f.taps[k] + accR += f.delayR[k] * w + accI += f.delayI[k] * w + } + out = append(out, complex(float32(accR), float32(accI))) + f.phase = f.factor - 1 + } else { + f.phase-- + } + } + return out +} diff --git a/internal/dsp/decimating_fir_test.go b/internal/dsp/decimating_fir_test.go new file mode 100644 index 0000000..821cb09 --- /dev/null +++ b/internal/dsp/decimating_fir_test.go @@ -0,0 +1,57 @@ +package dsp + +import ( + "math/cmplx" + "testing" +) + +func TestStatefulDecimatingFIRComplexStreamContinuity(t *testing.T) { + taps := LowpassFIR(90000, 512000, 101) + factor := 2 + + input := make([]complex64, 8192) + for i := range input { + input[i] = complex(float32((i%17)-8)/8.0, float32((i%11)-5)/8.0) + } + + one := NewStatefulDecimatingFIRComplex(taps, factor) + whole := one.Process(input) + + chunkedProc := NewStatefulDecimatingFIRComplex(taps, factor) + var chunked []complex64 + for i := 0; i < len(input); i += 733 { + end := i + 733 + if end > len(input) { + end = len(input) + } + chunked = append(chunked, chunkedProc.Process(input[i:end])...) + } + + if len(whole) != len(chunked) { + t.Fatalf("length mismatch whole=%d chunked=%d", len(whole), len(chunked)) + } + for i := range whole { + if cmplx.Abs(complex128(whole[i]-chunked[i])) > 1e-5 { + t.Fatalf("sample %d mismatch whole=%v chunked=%v", i, whole[i], chunked[i]) + } + } +} + +func TestStatefulDecimatingFIRComplexMatchesBlockPipelineLength(t *testing.T) { + taps := LowpassFIR(90000, 512000, 101) + factor := 2 + input := make([]complex64, 48640) + for i := range input { + input[i] = complex(float32((i%13)-6)/8.0, float32((i%7)-3)/8.0) + } + + stateful := NewStatefulDecimatingFIRComplex(taps, factor) + out := stateful.Process(input) + + filtered := ApplyFIR(input, taps) + dec := Decimate(filtered, factor) + + if len(out) != len(dec) { + t.Fatalf("unexpected output len got=%d want=%d", len(out), len(dec)) + } +} diff --git a/internal/recorder/streamer.go b/internal/recorder/streamer.go index d5147cc..be7fb70 100644 --- a/internal/recorder/streamer.go +++ b/internal/recorder/streamer.go @@ -116,11 +116,12 @@ type streamSession struct { // Stateful pre-demod anti-alias FIR (eliminates cold-start transients // and avoids per-frame FIR recomputation) - preDemodFIR *dsp.StatefulFIRComplex - preDemodDecim int // cached decimation factor - preDemodRate int // cached snipRate this FIR was built for - preDemodCutoff float64 // cached cutoff - preDemodDecimPhase int // stateful decimation phase (index offset into next frame) + preDemodFIR *dsp.StatefulFIRComplex + preDemodDecimator *dsp.StatefulDecimatingFIRComplex + preDemodDecim int // cached decimation factor + preDemodRate int // cached snipRate this FIR was built for + preDemodCutoff float64 // cached cutoff + preDemodDecimPhase int // retained for backward compatibility in snapshots/debug // AQ-2: De-emphasis config (µs, 0 = disabled) deemphasisUs float64 @@ -890,9 +891,10 @@ func (sess *streamSession) processSnippet(snippet []complex64, snipRate int) ([] } // Lazy-init or reinit stateful FIR if parameters changed - if sess.preDemodFIR == nil || sess.preDemodRate != snipRate || sess.preDemodCutoff != cutoff { + if sess.preDemodDecimator == nil || sess.preDemodRate != snipRate || sess.preDemodCutoff != cutoff || sess.preDemodDecim != decim1 { taps := dsp.LowpassFIR(cutoff, snipRate, 101) sess.preDemodFIR = dsp.NewStatefulFIRComplex(taps) + sess.preDemodDecimator = dsp.NewStatefulDecimatingFIRComplex(taps, decim1) sess.preDemodRate = snipRate sess.preDemodCutoff = cutoff sess.preDemodDecim = decim1 @@ -901,7 +903,7 @@ func (sess *streamSession) processSnippet(snippet []complex64, snipRate int) ([] decimPhaseBefore := sess.preDemodDecimPhase filtered := sess.preDemodFIR.ProcessInto(fullSnip, sess.growIQ(len(fullSnip))) - dec = dsp.DecimateStateful(filtered, decim1, &sess.preDemodDecimPhase) + dec = sess.preDemodDecimator.Process(fullSnip) logging.Debug("boundary", "snippet_path", "signal", sess.signalID, "overlap_applied", overlapApplied, "snip_len", len(snippet), "full_len", len(fullSnip), "filtered_len", len(filtered), "dec_len", len(dec), "decim1", decim1, "phase_before", decimPhaseBefore, "phase_after", sess.preDemodDecimPhase) } else { logging.Debug("boundary", "snippet_path", "signal", sess.signalID, "overlap_applied", overlapApplied, "snip_len", len(snippet), "full_len", len(fullSnip), "filtered_len", len(fullSnip), "dec_len", len(fullSnip), "decim1", decim1, "phase_before", 0, "phase_after", 0) @@ -1322,6 +1324,7 @@ type dspStateSnapshot struct { pilotLPFHi *dsp.StatefulFIRReal pilotLPFLo *dsp.StatefulFIRReal preDemodFIR *dsp.StatefulFIRComplex + preDemodDecimator *dsp.StatefulDecimatingFIRComplex preDemodDecim int preDemodRate int preDemodCutoff float64 @@ -1354,6 +1357,7 @@ func (sess *streamSession) captureDSPState() dspStateSnapshot { pilotLPFHi: sess.pilotLPFHi, pilotLPFLo: sess.pilotLPFLo, preDemodFIR: sess.preDemodFIR, + preDemodDecimator: sess.preDemodDecimator, preDemodDecim: sess.preDemodDecim, preDemodRate: sess.preDemodRate, preDemodCutoff: sess.preDemodCutoff, @@ -1386,6 +1390,7 @@ func (sess *streamSession) restoreDSPState(s dspStateSnapshot) { sess.pilotLPFHi = s.pilotLPFHi sess.pilotLPFLo = s.pilotLPFLo sess.preDemodFIR = s.preDemodFIR + sess.preDemodDecimator = s.preDemodDecimator sess.preDemodDecim = s.preDemodDecim sess.preDemodRate = s.preDemodRate sess.preDemodCutoff = s.preDemodCutoff @@ -1767,6 +1772,7 @@ func (st *Streamer) ResetStreams() { defer st.mu.Unlock() for _, sess := range st.sessions { sess.preDemodFIR = nil + sess.preDemodDecimator = nil sess.preDemodDecimPhase = 0 sess.stereoResampler = nil sess.monoResampler = nil