Skip to main content

max / audiofiles

DSP: fix export scaling/dither/MFCC correctness, add 4GB WAV guard, NaN hardening - Export integer scaling: scale by 2^(N-1) (not 2^(N-1)-1) so full-scale -1.0 maps to the most-negative code and +1.0 saturates to MAX; the clamp bounds are now both reachable and meaningful. Applies to WAV (encode.rs) and AIFF. - Dither: draw from xorshift64's high 32 bits and divide by 2^32 so next_f32 is half-open [0,1) and the TPDF dither is proper zero-mean triangular; add tests. - MFCC: take the log of the power spectrum (|X|^2, was magnitude), use an orthonormal DCT-II (was a flat *2.0 that left C0 dwarfing C1..C12), and give degenerate low mel bands their center bin instead of a constant dead value. Bump FEATURE_VERSION 1 -> 2 so existing vectors recompute via feature backfill. - WAV: reject >4GB exports before writing (mirrors the AIFF guard); the streaming path is exactly where this overflow is reachable. - Harden: clamp gain scale finite (avoid 0*inf=NaN on huge dB); cluster k-means uses total_cmp (no NaN-panic); clarify the fade "logarithmic" curve comment. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Author: Max Johnson <me@maxj.phd> · 2026-06-22 00:54 UTC
Commit: 87d3b7efc29b2e0ee2a13ff3849fd64b5b0bd4fe
Parent: 9e8e5a2
9 files changed, +132 insertions, -30 deletions
@@ -17,7 +17,12 @@ pub const NUM_FEATURES: usize = 35;
17 17 /// Version stamp for the feature-extraction layout/params. Bump whenever the
18 18 /// 35-feature vector's composition changes so persisted vectors (and shared
19 19 /// `.afcl` bundles) can be invalidated/recomputed rather than silently mixed.
20 - pub const FEATURE_VERSION: u32 = 1;
20 + ///
21 + /// v2: MFCCs switched to power-spectrum input + orthonormal DCT-II scaling, and
22 + /// degenerate low mel bands now use their center bin instead of a constant dead
23 + /// value (the 26 MFCC dims of the vector changed numerically). Existing libraries
24 + /// recompute these via the normal feature backfill on next launch.
25 + pub const FEATURE_VERSION: u32 = 2;
21 26
22 27 // ── SampleClass enum (unchanged) ──
23 28
@@ -101,11 +101,7 @@ fn init_centroids(points: &[Vec<f64>], k: usize) -> Vec<Vec<f64>> {
101 101 .collect();
102 102
103 103 let first = (0..points.len())
104 - .min_by(|&a, &b| {
105 - sq_dist(&points[a], &mean)
106 - .partial_cmp(&sq_dist(&points[b], &mean))
107 - .unwrap()
108 - })
104 + .min_by(|&a, &b| sq_dist(&points[a], &mean).total_cmp(&sq_dist(&points[b], &mean)))
109 105 .unwrap();
110 106
111 107 let mut centroids = vec![points[first].clone()];
@@ -113,7 +109,7 @@ fn init_centroids(points: &[Vec<f64>], k: usize) -> Vec<Vec<f64>> {
113 109
114 110 while centroids.len() < k {
115 111 let next = (0..points.len())
116 - .max_by(|&a, &b| min_d[a].partial_cmp(&min_d[b]).unwrap())
112 + .max_by(|&a, &b| min_d[a].total_cmp(&min_d[b]))
117 113 .unwrap();
118 114 centroids.push(points[next].clone());
119 115 for (i, p) in points.iter().enumerate() {
@@ -146,11 +142,7 @@ pub fn cluster_library(db: &Database, k: usize, max_iters: usize) -> Result<Clus
146 142 let mut changed = false;
147 143 for (i, p) in points.iter().enumerate() {
148 144 let best = (0..k)
149 - .min_by(|&a, &b| {
150 - sq_dist(p, &centroids[a])
151 - .partial_cmp(&sq_dist(p, &centroids[b]))
152 - .unwrap()
153 - })
145 + .min_by(|&a, &b| sq_dist(p, &centroids[a]).total_cmp(&sq_dist(p, &centroids[b])))
154 146 .unwrap();
155 147 if best != assign[i] {
156 148 assign[i] = best;
@@ -193,9 +185,7 @@ pub fn cluster_library(db: &Database, k: usize, max_iters: usize) -> Result<Clus
193 185 let medoid = *idxs
194 186 .iter()
195 187 .min_by(|&&a, &&b| {
196 - sq_dist(&points[a], &centroids[c])
197 - .partial_cmp(&sq_dist(&points[b], &centroids[c]))
198 - .unwrap()
188 + sq_dist(&points[a], &centroids[c]).total_cmp(&sq_dist(&points[b], &centroids[c]))
199 189 })
200 190 .unwrap();
201 191 clusters.push(Cluster {
@@ -82,9 +82,13 @@ fn build_mel_filterbank(fft_size: usize, sample_rate: u32, num_filters: usize) -
82 82 let end = bin_points[i + 2];
83 83
84 84 if start >= end || center <= start || center >= end {
85 + // Degenerate at this FFT resolution: low mel bands span less than one
86 + // FFT bin and collapse. Use the single center bin (weight 1.0) so the
87 + // filter still captures real energy, instead of emitting a constant
88 + // ln(1e-10) dead feature that carries no information.
85 89 filters.push(MelFilter {
86 - start_bin: 0,
87 - weights: Vec::new(),
90 + start_bin: center.min(spectrum_len - 1),
91 + weights: vec![1.0],
88 92 });
89 93 continue;
90 94 }
@@ -127,14 +131,16 @@ pub fn compute_mfccs(
127 131 let mut all_mfccs: Vec<[f64; NUM_MFCC]> = Vec::with_capacity(magnitude_frames.len());
128 132
129 133 for frame in magnitude_frames {
130 - // Apply mel filterbank
134 + // Apply mel filterbank over the power spectrum (|X|^2), the standard MFCC
135 + // input — `frame[bin]` is magnitude, so square it.
131 136 let mut mel_energies = [0.0f64; NUM_MEL_FILTERS];
132 137 for (i, filter) in filterbank.filters.iter().enumerate() {
133 138 let mut energy = 0.0;
134 139 for (j, &w) in filter.weights.iter().enumerate() {
135 140 let bin = filter.start_bin + j;
136 141 if bin < frame.len() {
137 - energy += w * frame[bin];
142 + let mag = frame[bin];
143 + energy += w * mag * mag;
138 144 }
139 145 }
140 146 mel_energies[i] = energy.max(1e-10);
@@ -144,7 +150,9 @@ pub fn compute_mfccs(
144 150 let log_energies: [f64; NUM_MEL_FILTERS] =
145 151 std::array::from_fn(|i| mel_energies[i].ln());
146 152
147 - // DCT-II: keep first NUM_MFCC coefficients
153 + // Orthonormal DCT-II: keep first NUM_MFCC coefficients. The per-coefficient
154 + // scaling (sqrt(1/N) for k=0, sqrt(2/N) otherwise) keeps C0 on the same
155 + // scale as C1..C12, instead of the flat *2.0 that left C0 dwarfing the rest.
148 156 let n = NUM_MEL_FILTERS as f64;
149 157 let mut mfccs = [0.0f64; NUM_MFCC];
150 158 for (k, coeff) in mfccs.iter_mut().enumerate() {
@@ -154,7 +162,8 @@ pub fn compute_mfccs(
154 162 * (std::f64::consts::PI * k as f64 * (2.0 * i as f64 + 1.0) / (2.0 * n))
155 163 .cos();
156 164 }
157 - *coeff = sum * 2.0;
165 + let ortho = if k == 0 { (1.0 / n).sqrt() } else { (2.0 / n).sqrt() };
166 + *coeff = sum * ortho;
158 167 }
159 168
160 169 all_mfccs.push(mfccs);
@@ -15,8 +15,9 @@ impl FadeCurve {
15 15 match self {
16 16 FadeCurve::Linear => t,
17 17 FadeCurve::Logarithmic => {
18 - // Attempt a logarithmic curve: fast rise then flattening
19 - // log2(1 + t) / log2(2) = log2(1 + t)
18 + // Concave (fast rise, then flattening): log2(1 + t) maps [0,1] to
19 + // [0,1] since log2(1)=0 and log2(2)=1. A logarithmic-shaped
20 + // amplitude ramp, not an equal-loudness/dB fade.
20 21 (1.0 + t).log2()
21 22 }
22 23 FadeCurve::SCurve => {
@@ -6,7 +6,11 @@ pub fn apply_gain(samples: &mut [f32], db: f64) {
6 6 return;
7 7 }
8 8
9 + // A huge finite dB overflows the dB→linear conversion to +inf; `0.0 * inf` is
10 + // NaN, which `clamp` propagates, silently poisoning silent samples. Clamp the
11 + // scale finite so silence stays 0.0 and loud samples saturate to ±1.0.
9 12 let scale = 10.0_f64.powf(db / 20.0) as f32;
13 + let scale = if scale.is_finite() { scale } else { f32::MAX };
10 14 for sample in samples.iter_mut() {
11 15 *sample = (*sample * scale).clamp(-1.0, 1.0);
12 16 }
@@ -16,10 +16,44 @@ impl SimpleRng {
16 16 self.state ^= self.state << 13;
17 17 self.state ^= self.state >> 7;
18 18 self.state ^= self.state << 17;
19 - (self.state & 0xFFFF_FFFF) as u32
19 + // Take the high 32 bits: xorshift64's low bits are the weakest, so the
20 + // top word is the better-quality half to draw the dither sample from.
21 + (self.state >> 32) as u32
20 22 }
21 23
24 + /// Uniform sample in the half-open range `[0, 1)`. Dividing by `2^32` (not
25 + /// `u32::MAX`, which would make the range closed) keeps it half-open, so
26 + /// `next_f32 + next_f32 - 1.0` is a proper zero-mean TPDF value in `[-1, 1)`.
22 27 pub(crate) fn next_f32(&mut self) -> f32 {
23 - self.next_u32() as f32 / u32::MAX as f32
28 + self.next_u32() as f32 / 4_294_967_296.0 // 2^32
29 + }
30 + }
31 +
32 + #[cfg(test)]
33 + mod tests {
34 + use super::*;
35 +
36 + #[test]
37 + fn next_f32_is_half_open_unit_interval() {
38 + let mut rng = SimpleRng::new(1);
39 + for _ in 0..100_000 {
40 + let x = rng.next_f32();
41 + assert!((0.0..1.0).contains(&x), "out of [0,1): {x}");
42 + }
43 + }
44 +
45 + #[test]
46 + fn tpdf_dither_is_zero_mean_and_bounded() {
47 + // (uniform + uniform - 1) is triangular on [-1, 1) with mean ~0.
48 + let mut rng = SimpleRng::new(0x1234_5678);
49 + let n = 200_000;
50 + let mut sum = 0.0f64;
51 + for _ in 0..n {
52 + let d = rng.next_f32() + rng.next_f32() - 1.0;
53 + assert!((-1.0..1.0).contains(&d), "dither out of range: {d}");
54 + sum += d as f64;
55 + }
56 + let mean = sum / n as f64;
57 + assert!(mean.abs() < 0.01, "TPDF mean not ~0: {mean}");
24 58 }
25 59 }
@@ -28,9 +28,13 @@ pub(crate) fn write_sample<W: Write + Seek>(
28 28 bit_depth: u16,
29 29 rng: &mut SimpleRng,
30 30 ) -> Result<(), CoreError> {
31 + // Scale by 2^(N-1) so full-scale -1.0 maps exactly to the most-negative code
32 + // and +1.0 maps to 2^(N-1) (one above MAX), which the clamp pins to MAX.
33 + // Scaling by 2^(N-1)-1 (the old code) left the most-negative code unreachable
34 + // and the clamp's lower bound dead — an internal inconsistency.
31 35 match bit_depth {
32 36 8 => {
33 - let scale = i8::MAX as f32;
37 + let scale = -(i8::MIN as f32); // 128 = 2^7
34 38 let dither = (rng.next_f32() + rng.next_f32() - 1.0) / scale;
35 39 let clamped = ((sample + dither) * scale)
36 40 .round()
@@ -38,7 +42,7 @@ pub(crate) fn write_sample<W: Write + Seek>(
38 42 writer.write_sample(clamped)
39 43 }
40 44 16 => {
41 - let scale = i16::MAX as f32;
45 + let scale = -(i16::MIN as f32); // 32768 = 2^15
42 46 let dither = (rng.next_f32() + rng.next_f32() - 1.0) / scale;
43 47 let clamped = ((sample + dither) * scale)
44 48 .round()
@@ -46,7 +50,7 @@ pub(crate) fn write_sample<W: Write + Seek>(
46 50 writer.write_sample(clamped)
47 51 }
48 52 24 => {
49 - let scale = 8_388_607.0f32; // 2^23 - 1
53 + let scale = 8_388_608.0f32; // 2^23
50 54 let scaled = (sample * scale).round().clamp(-8_388_608.0, 8_388_607.0) as i32;
51 55 writer.write_sample(scaled)
52 56 }
@@ -69,6 +73,23 @@ pub(crate) fn write_sample<W: Write + Seek>(
69 73 /// signal-independent noise; the only effect is that exports are reproducible.
70 74 pub(crate) const DITHER_SEED: u64 = 0x9E3779B97F4A7C15;
71 75
76 + /// Reject WAV exports whose data would overflow the RIFF format's 32-bit size
77 + /// fields (about 4 GB). `total_samples` is the total f32 sample count, i.e.
78 + /// frames times channels. This mirrors the AIFF guard in `encode_aiff`; without
79 + /// it the WAV and streaming paths trust hound's 32-bit RIFF fields and would
80 + /// silently wrap on an export larger than 4 GB.
81 + pub(crate) fn wav_size_check(total_samples: u64, bit_depth: u16) -> Result<(), CoreError> {
82 + let bytes_per_sample = (bit_depth / 8) as u64;
83 + let data_size = total_samples.saturating_mul(bytes_per_sample);
84 + // The RIFF chunk-size field (u32) covers the 36-byte header remainder + data.
85 + if data_size + 36 > u32::MAX as u64 {
86 + return Err(CoreError::Export(format!(
87 + "WAV: file too large ({total_samples} samples, {bit_depth}-bit = {data_size} bytes, exceeds 4 GB RIFF limit)"
88 + )));
89 + }
90 + Ok(())
91 + }
92 +
72 93 /// Encode audio to a WAV file at the given path.
73 94 ///
74 95 /// - 8-bit: unsigned PCM (the WAV spec stores 8-bit as offset-128 unsigned),
@@ -84,6 +105,7 @@ pub fn encode_wav(audio: &ConvertedAudio, bit_depth: u16, dest: &Path) -> Result
84 105 "unsupported bit depth: {bit_depth} (expected 8, 16, 24, or 32)"
85 106 )));
86 107 }
108 + wav_size_check(audio.samples.len() as u64, bit_depth)?;
87 109 let spec = WavSpec {
88 110 channels: audio.channels,
89 111 sample_rate: audio.sample_rate,
@@ -150,6 +172,33 @@ mod tests {
150 172 }
151 173
152 174 #[test]
175 + fn wav_16bit_full_scale_hits_min_and_max() {
176 + // Scaling by 2^15 maps -1.0 exactly to i16::MIN and clamps +1.0 (and
177 + // beyond) to i16::MAX. The old 2^15-1 scaling left MIN unreachable.
178 + let dir = tempfile::tempdir().unwrap();
179 + let path = dir.path().join("fs.wav");
180 + // 32-bit avoids dither perturbing the boundary check below; use 16-bit
181 + // with dither off by reading the raw codes at the extremes.
182 + let audio = make_audio(vec![-1.0, 1.0, 2.0, -2.0], 1, 44100);
183 + encode_wav(&audio, 16, &path).unwrap();
184 + let mut reader = hound::WavReader::open(&path).unwrap();
185 + let samples: Vec<i16> = reader.samples::<i16>().map(|s| s.unwrap()).collect();
186 + assert_eq!(samples[0], i16::MIN, "-1.0 should map to i16::MIN");
187 + assert_eq!(samples[1], i16::MAX, "+1.0 should clamp to i16::MAX");
188 + assert_eq!(samples[2], i16::MAX, "over-range + should clamp to MAX");
189 + assert_eq!(samples[3], i16::MIN, "over-range - should clamp to MIN");
190 + }
191 +
192 + #[test]
193 + fn wav_size_check_rejects_over_4gb() {
194 + // ~4 GB / 2 bytes = ~2.1e9 samples at 16-bit. Just over the limit errors;
195 + // a small count passes.
196 + assert!(wav_size_check(8, 16).is_ok());
197 + let over = (u32::MAX as u64) / 2 + 1; // *2 bytes > u32::MAX
198 + assert!(wav_size_check(over, 16).is_err());
199 + }
200 +
201 + #[test]
153 202 fn wav_24bit_roundtrip() {
154 203 let dir = tempfile::tempdir().unwrap();
155 204 let path = dir.path().join("test_24.wav");
@@ -94,7 +94,7 @@ pub fn encode_aiff(audio: &ConvertedAudio, bit_depth: u16, dest: &Path) -> Resul
94 94 // allocator address, making AIFF output nondeterministic and any
95 95 // byte-parity test flaky by construction.
96 96 let mut rng = SimpleRng::new(super::encode::DITHER_SEED);
97 - let scale = i16::MAX as f32;
97 + let scale = -(i16::MIN as f32); // 32768 = 2^15; see encode::write_sample
98 98 for &sample in &audio.samples {
99 99 let dither = (rng.next_f32() + rng.next_f32() - 1.0) / scale;
100 100 let dithered = (sample + dither) * scale;
@@ -103,7 +103,7 @@ pub fn encode_aiff(audio: &ConvertedAudio, bit_depth: u16, dest: &Path) -> Resul
103 103 }
104 104 }
105 105 24 => {
106 - let scale = 8_388_607.0f32; // 2^23 - 1
106 + let scale = 8_388_608.0f32; // 2^23
107 107 for &sample in &audio.samples {
108 108 let scaled = (sample * scale)
109 109 .round()
@@ -82,6 +82,16 @@ pub(crate) fn try_stream_wav_export(
82 82 let out_ch = out_channel_count(target_channels, src_channels) as usize;
83 83 let dst_rate = target_rate.unwrap_or(src_rate);
84 84
85 + // Reject a >4 GB output before writing a single byte (the streaming path
86 + // exists for long files, so this overflow is reachable here in a way the
87 + // buffered path rarely hits).
88 + let out_frames = if dst_rate == src_rate || src_rate == 0 {
89 + n_frames
90 + } else {
91 + ((n_frames as f64) * (dst_rate as f64 / src_rate as f64)).round() as usize
92 + };
93 + super::encode::wav_size_check((out_frames as u64).saturating_mul(out_ch as u64), bit_depth)?;
94 +
85 95 let spec = WavSpec {
86 96 channels: out_ch as u16,
87 97 sample_rate: dst_rate,