cathar 0.5.4

Audio toolkit in pure Rust — denoise, de-hum, de-click, de-clip, de-reverb, normalise, and more.
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71

//! Anti-aliased sample-rate conversion (Kaiser-windowed sinc).

/// Modified Bessel function of the first kind, order 0 — for the Kaiser window.
fn bessel_i0(x: f64) -> f64 {
    let mut sum = 1.0f64;
    let mut term = 1.0f64;
    let half_sq = (x / 2.0) * (x / 2.0);
    for k in 1..=30 {
        term *= half_sq / (k as f64 * k as f64);
        sum += term;
        if term < 1e-13 * sum {
            break;
        }
    }
    sum
}

/// Resample one channel from `from_rate` to `to_rate` with a Kaiser-windowed
/// sinc (arbitrary ratio). The cutoff tracks the lower of the two Nyquist
/// limits, so downsampling is anti-aliased and upsampling adds no imaging;
/// the filter support widens at low cutoffs to keep the stopband sharp. Returns
/// the input unchanged when the rates already match.
pub fn resample(signal: &[f32], from_rate: u32, to_rate: u32) -> Vec<f32> {
    if from_rate == to_rate || from_rate == 0 || to_rate == 0 || signal.is_empty() {
        return signal.to_vec();
    }
    let ratio = to_rate as f64 / from_rate as f64;
    let out_len = (signal.len() as f64 * ratio).round().max(1.0) as usize;

    // Cutoff as a fraction of the source rate (≤ 0.5); downsampling lowers it to
    // the destination Nyquist.
    let cutoff = 0.5 * ratio.min(1.0);
    // Fixed number of sinc lobes per side; support (in source samples) grows as
    // the cutoff falls so the filter always spans the same zero crossings.
    const LOBES: f64 = 16.0;
    let half_width = LOBES / (2.0 * cutoff);
    let beta = 9.0;
    let i0_beta = bessel_i0(beta);

    let n = signal.len() as isize;
    let mut out = vec![0.0f32; out_len];
    for (i, o) in out.iter_mut().enumerate() {
        let center = i as f64 / ratio; // position in source samples
        let first = (center - half_width).ceil() as isize;
        let last = (center + half_width).floor() as isize;
        let mut acc = 0.0f64;
        let mut norm = 0.0f64;
        for idx in first..=last {
            if idx < 0 || idx >= n {
                continue;
            }
            let dx = center - idx as f64;
            // Low-pass sinc 2·fc·sinc(2·fc·dx) = sin(π·t)/(π·dx), t = 2·fc·dx.
            let t = 2.0 * cutoff * dx;
            let sinc = if dx.abs() < 1e-9 {
                2.0 * cutoff
            } else {
                (std::f64::consts::PI * t).sin() / (std::f64::consts::PI * dx)
            };
            // Kaiser window over |dx| ≤ half_width.
            let r = dx / half_width;
            let w =
                if r.abs() < 1.0 { bessel_i0(beta * (1.0 - r * r).sqrt()) / i0_beta } else { 0.0 };
            let k = sinc * w;
            acc += signal[idx as usize] as f64 * k;
            norm += k;
        }
        *o = (acc / norm.max(1e-12)) as f32;
    }
    out
}