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
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
#[derive(Clone, Debug)]
pub struct Yin {
threshold: f64,
tau_max: usize,
tau_min: usize,
sample_rate: usize,
}
impl Yin {
pub fn init(threshold: f64, freq_min: f64, freq_max: f64, sample_rate: usize) -> Yin {
let tau_max = sample_rate / freq_min as usize;
let tau_min = sample_rate / freq_max as usize;
Yin {
threshold,
tau_max,
tau_min,
sample_rate,
}
}
pub fn estimate_freq(&self, audio_sample: &[f64]) -> f64 {
let sample_frequency = compute_sample_frequency(
audio_sample,
self.tau_min,
self.tau_max,
self.sample_rate,
self.threshold,
);
sample_frequency
}
}
fn diff_function(audio_sample: &[f64], tau_max: usize) -> Vec<f64> {
let mut diff_function = vec![0.0; tau_max];
for tau in 1..tau_max {
for j in 0..(audio_sample.len() - tau_max) {
let tmp = audio_sample[j] - audio_sample[j + tau];
diff_function[tau] += tmp * tmp;
}
}
diff_function
}
fn cmndf(raw_diff: &[f64]) -> Vec<f64> {
let mut cmndf_diff: Vec<f64> = raw_diff[1..]
.iter()
.enumerate()
.scan(0.0, |state, x| {
*state = *state + x.1;
let result = *x.1 * (x.0 + 1) as f64 / *state;
Some(result)
})
.collect();
cmndf_diff.insert(0, 1.0);
cmndf_diff
}
fn compute_diff_min(diff_fn: &[f64], min_tau: usize, max_tau: usize, harm_threshold: f64) -> usize {
let mut tau = min_tau;
while tau < max_tau {
if diff_fn[tau] < harm_threshold {
while tau + 1 < max_tau && diff_fn[tau + 1] < diff_fn[tau] {
tau += 1;
}
return tau;
}
tau += 1;
}
0
}
fn convert_to_frequency(sample_period: usize, sample_rate: usize) -> f64 {
let value: f64 = sample_period as f64 / sample_rate as f64;
1.0 / value
}
pub fn compute_sample_frequency(
audio_sample: &[f64],
tau_min: usize,
tau_max: usize,
sample_rate: usize,
threshold: f64,
) -> f64 {
let diff_fn = diff_function(&audio_sample, tau_max);
let cmndf = cmndf(&diff_fn);
let sample_period = compute_diff_min(&cmndf, tau_min, tau_max, threshold);
convert_to_frequency(sample_period, sample_rate)
}
#[cfg(test)]
mod tests {
use dasp::{signal, Signal};
fn produce_sample(sample_rate: usize, frequency: f64, noise_ratio: f64) -> Vec<f64> {
use rand::prelude::*;
let mut rng = thread_rng();
let mut signal = signal::rate(sample_rate as f64).const_hz(frequency).sine();
let sample: Vec<f64> = (0..sample_rate)
.map(|_| signal.next() + noise_ratio * rng.gen::<f64>())
.collect();
sample
}
use super::*;
#[test]
fn sanity_basic_sine() {
let sample = produce_sample(12, 4.0, 0.0);
let yin = Yin::init(0.1, 2.0, 5.0, 12);
let computed_frequency = yin.estimate_freq(&sample);
assert_eq!(computed_frequency, 4.0);
}
#[test]
fn sanity_non_multiple() {
let sample = produce_sample(44100, 4000.0, 0.0);
let yin = Yin::init(0.1, 3000.0, 5000.0, 44100);
let computed_frequency = yin.estimate_freq(&sample);
let difference = computed_frequency - 4000.0;
assert!(difference.abs() < 50.0);
}
#[test]
fn sanity_full_sine() {
let sample = produce_sample(44100, 441.0, 0.0);
let yin = Yin::init(0.1, 300.0, 500.0, 44100);
let computed_frequency = yin.estimate_freq(&sample);
assert_eq!(computed_frequency, 441.0);
}
#[test]
fn readme_doctest() {
let estimator = Yin::init(0.1, 10.0, 30.0, 80);
let mut example = vec![];
let mut prev_value = -1.0;
for i in 0..80 {
if i % 2 != 0 {
example.push(0.0);
} else {
prev_value *= -1.0;
example.push(prev_value);
}
}
let freq = estimator.estimate_freq(&example);
assert_eq!(freq, 20.0);
}
}