debug: use stateful decimating FIR in pre-demod path

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
+}

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))
+	}
+}

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