1/**
2 * The theoretical bits/sample: the entropy rate of the stationary Gaussian
3 * process y = h * x with x ~ N(0, σ²) i.i.d., quantized at unit step. In the
4 * fine-quantization (high-resolution) regime the discrete entropy rate
5 * approaches the differential entropy rate (Kolmogorov):
6 *
7 * R = ½ log₂(2πe) + ∫₀^{1/2} log₂ S(f) df, S(f) = σ² |H(f)|²
8 *
9 * with f in cycles/sample. With no filter this reduces to ½ log₂(2πe σ²).
10 * The formula ignores dither and degrades where S(f) falls to the order of
11 * the quantization step or below — exactly the regime the app probes.
12 */
14const INTEGRATION_POINTS = 8192
16export function entropyRateBits(kernel: Float64Array, sigma: number): number {
17 const L = kernel.length
18 let integral = 0
19 for (let k = 0; k < INTEGRATION_POINTS; k++) {
20 // Midpoint rule keeps f = 0 (where a bandpass H vanishes) off the grid;
21 // the log singularity at isolated zeros is integrable.
22 const f = (0.5 * (k + 0.5)) / INTEGRATION_POINTS
23 let re = 0
24 let im = 0
25 for (let i = 0; i < L; i++) {
26 re += kernel[i] * Math.cos(2 * Math.PI * f * i)
27 im -= kernel[i] * Math.sin(2 * Math.PI * f * i)
28 }
29 const S = sigma * sigma * (re * re + im * im)
30 integral += Math.log2(Math.max(S, 1e-300))
31 }
32 integral *= 0.5 / INTEGRATION_POINTS
33 return 0.5 * Math.log2(2 * Math.PI * Math.E) + integral
34}