Skip to main content

sva_samples/measure/harmony/
dissonance.rs

1// Concern: Sethares' (1993) Plomp-Levelt roughness per moment, and the curve's extremes | Non-concern: where the partials came from | IO: (moments, framing, floor) -> Dissonance or Crowded
2
3use super::{Framing, Moment, Partial};
4
5pub const MODEL: &str = "Sethares, W. A. (1993), Local consonance and the relationship \
6                         between timbre and scale, J. Acoust. Soc. Am. 94(3) 1218-1228: \
7                         eq. 1 d(x) = exp(-a x) - exp(-b x) with a = 3.5, b = 5.75; eq. 3 \
8                         s = d*/(s1 f_low + s2) with d* = 0.24 and the fit s1 = 0.021, s2 = 19 \
9                         as printed; eq. 4 the amplitude factor is the product v1 v2; eq. 5 \
10                         roughness sums v1 v2 d(s |f2 - f1|) over every unordered pair";
11pub const A: f64 = 3.5;
12pub const B: f64 = 5.75;
13pub const D_STAR: f64 = 0.24;
14pub const S1: f64 = 0.021;
15pub const S2: f64 = 19.0;
16
17fn d(x: f64) -> f64 {
18    (-A * x).exp() - (-B * x).exp()
19}
20
21fn s(f_low: f64) -> f64 {
22    D_STAR / (S1 * f_low + S2)
23}
24
25pub fn pair(f1: f64, f2: f64, v1: f64, v2: f64) -> f64 {
26    v1 * v2 * d(s(f1.min(f2)) * (f2 - f1).abs())
27}
28
29/// `pair` over every point of two partials' bounds: `s` falls as `f_low` rises, and `d` rises
30/// to its one maximum at ln(b/a)/(b - a) and falls after, so the box's ends bound it.
31fn pair_bounds(p: &Partial, q: &Partial) -> (f64, f64) {
32    let ((p1, p2), (q1, q2)) = (p.hz_bounds, q.hz_bounds);
33    let df = ((p1 - q2).max(q1 - p2).max(0.0), (p2 - q1).max(q2 - p1));
34    let x = (s(p2.min(q2)) * df.0, s(p1.min(q1)) * df.1);
35    let top = ((B / A).ln() / (B - A)).clamp(x.0, x.1);
36    let (va, vb) = (p.amplitude_bounds, q.amplitude_bounds);
37    (va.0 * vb.0 * d(x.0).min(d(x.1)), va.1 * vb.1 * d(top))
38}
39
40#[derive(Clone, Copy, Debug, PartialEq)]
41pub struct DissonanceFrame {
42    pub t_secs: f64,
43    pub roughness: f64,
44    /// Over the partials whose amplitude bounds reach the floor, each pair at its bounds.
45    pub roughness_bounds: (f64, f64),
46    pub partials: usize,
47}
48
49#[derive(Clone, Copy, Debug, PartialEq)]
50pub struct Extreme {
51    pub value: f64,
52    pub t_secs: f64,
53}
54
55#[derive(Clone, Debug, PartialEq)]
56pub struct Dissonance {
57    pub framing: Framing,
58    pub floor_db: f64,
59    pub frames: Vec<DissonanceFrame>,
60    pub mean: Option<f64>,
61    pub max: Option<Extreme>,
62    pub min: Option<Extreme>,
63}
64
65#[derive(Clone, Copy, Debug, PartialEq)]
66pub struct Crowded {
67    pub t_secs: f64,
68    pub partials: usize,
69    pub pairs: usize,
70    pub cap: usize,
71}
72
73pub fn measure(
74    moments: &[Moment],
75    framing: Framing,
76    floor_db: f64,
77    max_pairs: usize,
78) -> Result<Dissonance, Crowded> {
79    let floor = 10f64.powf(floor_db / 20.0);
80    let mut frames = Vec::with_capacity(moments.len());
81    for moment in moments {
82        let reach: Vec<&Partial> = moment
83            .partials
84            .iter()
85            .filter(|p| p.amplitude_bounds.1 >= floor)
86            .collect();
87        let pairs = reach.len() * reach.len().saturating_sub(1) / 2;
88        if pairs > max_pairs {
89            return Err(Crowded {
90                t_secs: moment.t_secs,
91                partials: reach.len(),
92                pairs,
93                cap: max_pairs,
94            });
95        }
96        let (mut roughness, mut low, mut high) = (0.0, 0.0, 0.0);
97        for (i, p) in reach.iter().enumerate() {
98            for q in &reach[i + 1..] {
99                let (lo, hi) = pair_bounds(p, q);
100                high += hi;
101                if p.amplitude_bounds.0 >= floor && q.amplitude_bounds.0 >= floor {
102                    low += lo;
103                }
104                if p.amplitude >= floor && q.amplitude >= floor {
105                    roughness += pair(p.hz, q.hz, p.amplitude, q.amplitude);
106                }
107            }
108        }
109        frames.push(DissonanceFrame {
110            t_secs: moment.t_secs,
111            roughness,
112            roughness_bounds: (low, high),
113            partials: reach.iter().filter(|p| p.amplitude >= floor).count(),
114        });
115    }
116    let extreme = |better: fn(f64, f64) -> bool| {
117        frames
118            .iter()
119            .fold(None, |held: Option<Extreme>, f| match held {
120                Some(e) if !better(f.roughness, e.value) => Some(e),
121                _ => Some(Extreme {
122                    value: f.roughness,
123                    t_secs: f.t_secs,
124                }),
125            })
126    };
127    Ok(Dissonance {
128        framing,
129        floor_db,
130        mean: (!frames.is_empty())
131            .then(|| frames.iter().map(|f| f.roughness).sum::<f64>() / frames.len() as f64),
132        max: extreme(|a, b| a > b),
133        min: extreme(|a, b| a < b),
134        frames,
135    })
136}