fix: scale sinc kernel by bandwidth

Author: roderickvd

dasp_interpolate/src/sinc/mod.rs

@@ -93,9 +93,11 @@ where
 
         (0..max_depth).fold(Self::Frame::EQUILIBRIUM, |mut v, n| {
             v = {
-                let a = PI * bandwidth * (phil + n as f64);
-                let first = if a == 0.0 { 1.0 } else { sin(a) / a };
-                let second = 0.5 + 0.5 * cos(a / (depth as f64 * bandwidth));
+                let t = phil + n as f64;
+                let a = PI * bandwidth * t;
+                let b = PI * t / depth as f64;
+                let first = if a.abs() < f64::EPSILON { bandwidth } else { bandwidth * sin(a) / a };
+                let second = 0.5 + 0.5 * cos(b);
                 v.zip_map(self.frames[nl - n], |vs, r_lag| {
                     vs.add_amp(
                         (first * second * r_lag.to_sample::<f64>())
@@ -105,9 +107,11 @@ where
                 })
             };
 
-            let a = PI * bandwidth * (phir + n as f64);
-            let first = if a == 0.0 { 1.0 } else { sin(a) / a };
-            let second = 0.5 + 0.5 * cos(a / (depth as f64 * bandwidth));
+            let t = phir + n as f64;
+            let a = PI * bandwidth * t;
+            let b = PI * t / depth as f64;
+            let first = if a.abs() < f64::EPSILON { bandwidth } else { bandwidth * sin(a) / a };
+            let second = 0.5 + 0.5 * cos(b);
             v.zip_map(self.frames[nr + n], |vs, r_lag| {
                 vs.add_amp(
                     (first * second * r_lag.to_sample::<f64>())

examples/aliasing.rs

@@ -0,0 +1,53 @@
+use dasp::{interpolate::sinc::Sinc, ring_buffer, signal, Signal};
+use hound::{WavSpec, WavWriter};
+
+fn main() {
+    // Scaled version of issue #174:
+    // Original: 900 Hz @ 2000 Hz → 1500 Hz (Nyquist 750 Hz, aliases to 600 Hz)
+    // Scaled 8x: 7200 Hz @ 16000 Hz → 12000 Hz (Nyquist 6000 Hz, aliases to 4800 Hz)
+    let my_signal = signal::rate(16000.0).const_hz(7200.0).sine();
+
+    let spec = WavSpec {
+        channels: 1,
+        sample_rate: 12000,
+        bits_per_sample: 16,
+        sample_format: hound::SampleFormat::Int,
+    };
+
+    let ring_buffer = ring_buffer::Fixed::from([0.0; 100]);
+    let sinc = Sinc::new(ring_buffer);
+    let mut my_signal = my_signal.from_hz_to_hz(sinc, 16000.0, 12000.0);
+
+    let buffer_size = 100;
+    let source_hz = 16000.0;
+    let target_hz = 12000.0;
+    let samples_to_skip = (buffer_size as f64 / (source_hz / target_hz)).ceil() as usize;
+
+    for _ in 0..samples_to_skip {
+        my_signal.next();
+    }
+
+    // Check first 1000 samples to verify signal
+    let check_samples: Vec<f64> = (&mut my_signal).take(1000).collect();
+    let check_rms = (check_samples.iter().map(|s| s * s).sum::<f64>() / check_samples.len() as f64).sqrt();
+    eprintln!("RMS of resampled signal: {:.4}", check_rms);
+    eprintln!("Max amplitude: {:.4}", check_samples.iter().map(|s| s.abs()).fold(0.0f64, f64::max));
+
+    let mut writer = WavWriter::create("my_output.wav", spec).unwrap();
+
+    // Write the check samples first
+    for &frame in &check_samples {
+        let clamped = frame.max(-1.0).min(1.0);
+        writer
+            .write_sample((clamped * i16::MAX as f64) as i16)
+            .unwrap();
+    }
+
+    // Then write the rest
+    for frame in my_signal.take(12000 * 4 - 1000) {
+        let clamped = frame.max(-1.0).min(1.0);
+        writer
+            .write_sample((clamped * i16::MAX as f64) as i16)
+            .unwrap();
+    }
+}