sva_samples/
grid.rs

1// Concern: the sample grid: a sample's instant, the sample an instant or edge lands on, spans of samples | Non-concern: what a sample holds | IO: (n or t) -> t or n
2
3/// No pull, stored chunk or machine block crosses a multiple.
4pub const BLOCK: i64 = 1 << 12;
5
6pub fn block_end(n: i64) -> i64 {
7    (n.div_euclid(BLOCK) + 1).saturating_mul(BLOCK)
8}
9
10/// Samples `[start, end)` of the grid, whose sample 0 is t = 0.
11#[derive(Clone, Copy, Debug, PartialEq, Eq, Hash)]
12pub struct Extent {
13    pub start: i64,
14    pub end: i64,
15}
16
17impl Extent {
18    /// `i64::MIN` and `i64::MAX` stand for no edge.
19    pub const EVERYWHERE: Extent = Extent {
20        start: i64::MIN,
21        end: i64::MAX,
22    };
23
24    pub const NOWHERE: Extent = Extent { start: 0, end: 0 };
25
26    pub fn new(start: i64, end: i64) -> Extent {
27        assert!(start <= end, "an extent [{start}, {end}) runs backwards");
28        Extent { start, end }
29    }
30
31    pub fn from(start: i64) -> Extent {
32        Extent::new(start, i64::MAX)
33    }
34
35    pub fn is_bounded(&self) -> bool {
36        self.start != i64::MIN && self.end != i64::MAX
37    }
38
39    pub fn contains(&self, n: i64) -> bool {
40        self.start <= n && n < self.end
41    }
42
43    pub fn intersect(self, other: Extent) -> Extent {
44        let start = self.start.max(other.start);
45        let end = self.end.min(other.end);
46        match start < end {
47            true => Extent { start, end },
48            false => Extent::NOWHERE,
49        }
50    }
51
52    pub fn hull(self, other: Extent) -> Extent {
53        match (self.is_empty(), other.is_empty()) {
54            (true, _) => other,
55            (_, true) => self,
56            _ => Extent {
57                start: self.start.min(other.start),
58                end: self.end.max(other.end),
59            },
60        }
61    }
62
63    /// Every sample moved `by` later; an unbounded edge stays unbounded.
64    pub fn shifted(self, by: i64) -> Extent {
65        if self.is_empty() {
66            return self;
67        }
68        let edge = |n: i64| match n {
69            i64::MIN | i64::MAX => n,
70            n => n.saturating_add(by).clamp(i64::MIN + 1, i64::MAX - 1),
71        };
72        Extent {
73            start: edge(self.start),
74            end: edge(self.end),
75        }
76    }
77
78    pub fn secs(rate: u32, start_secs: f64, end_secs: f64) -> Extent {
79        let at = |secs: f64| (secs * f64::from(rate)).round() as i64;
80        Extent::new(at(start_secs), at(end_secs))
81    }
82
83    pub fn len(&self) -> usize {
84        debug_assert!(self.is_bounded(), "an unbounded extent has no length");
85        (self.end - self.start) as usize
86    }
87
88    pub fn is_empty(&self) -> bool {
89        self.end == self.start
90    }
91
92    pub fn start_secs(&self, rate: u32) -> f64 {
93        self.start as f64 / f64::from(rate)
94    }
95
96    pub fn span_secs(&self, rate: u32) -> f64 {
97        self.len() as f64 / f64::from(rate)
98    }
99}
100
101/// Sample `n` stands at `a*n/d` samples of `rate`, in lowest terms, `a, d > 0`: every grid
102/// starts at t = 0, so a sample index is the one clock every read and key counts in.
103#[derive(Clone, Copy, Debug, PartialEq, Eq, Hash, PartialOrd, Ord)]
104pub struct Grid {
105    pub rate: u32,
106    pub a: i128,
107    pub d: i128,
108}
109
110impl Grid {
111    pub const fn of(rate: u32) -> Grid {
112        Grid { rate, a: 1, d: 1 }
113    }
114
115    pub fn finer(rate: u32, scale: usize) -> Grid {
116        Grid {
117            rate,
118            a: 1,
119            d: scale as i128,
120        }
121    }
122
123    pub fn is_rate(&self) -> bool {
124        (self.a, self.d) == (1, 1)
125    }
126
127    pub(crate) fn exact(&self, n: i64) -> Option<(i128, i128)> {
128        let num = self.a.checked_mul(i128::from(n))?;
129        Some((num, self.d.checked_mul(i128::from(self.rate))?))
130    }
131
132    /// One quotient: correctly rounded while both integers are under 2^53; past that each
133    /// integer rounds once converting and the quotient once more.
134    pub fn instant(&self, n: i64) -> f64 {
135        if self.is_rate() {
136            return n as f64 / f64::from(self.rate);
137        }
138        let num = self.a.saturating_mul(i128::from(n));
139        num as f64 / self.d.saturating_mul(i128::from(self.rate)) as f64
140    }
141
142    /// The first sample whose exact instant is at or past `edge` read as the decimal it prints
143    /// as: the one rule every crop and indicator edge meets the grid by.
144    pub fn first_at(&self, edge: f64) -> i64 {
145        let beyond = if edge < 0.0 { i64::MIN } else { i64::MAX };
146        if !edge.is_finite() {
147            return beyond;
148        }
149        match self.decimal_steps(edge) {
150            Some(n) => i64::try_from(n).unwrap_or(beyond),
151            None => self.step_at(edge, Round::Ceil).unwrap_or(beyond),
152        }
153    }
154
155    /// An edge `first_at` meets at sample `n` exactly, beside `n`'s instant.
156    pub fn edge_at(&self, n: i64) -> Option<f64> {
157        let t = self.instant(n);
158        [t, t.next_down(), t.next_up()]
159            .into_iter()
160            .find(|edge| self.first_at(*edge) == n)
161    }
162
163    /// `ceil(edge * d * rate / a)`, `edge` as the shortest decimal that prints it.
164    fn decimal_steps(&self, edge: f64) -> Option<i128> {
165        const LIMIT: i128 = 1 << 120;
166        let text = format!("{edge:e}");
167        let (mantissa, exp) = text.split_once('e')?;
168        let (whole, frac) = mantissa.split_once('.').unwrap_or((mantissa, ""));
169        let digits: i128 = format!("{whole}{frac}").parse().ok()?;
170        let shift = exp.parse::<i32>().ok()? - frac.len() as i32;
171        let ten = 10i128
172            .checked_pow(shift.unsigned_abs())
173            .filter(|p| *p < LIMIT)?;
174        let (num, den) = match shift >= 0 {
175            true => (digits.checked_mul(ten)?, 1),
176            false => (digits, ten),
177        };
178        let top = num.checked_mul(self.d.checked_mul(i128::from(self.rate))?)?;
179        let bottom = den.checked_mul(self.a)?;
180        Some(top.div_euclid(bottom) + i128::from(top.rem_euclid(bottom) != 0))
181    }
182
183    /// Whether sample `n` lies in `[l, r)` by `first_at`, which only a tie of its instant needs.
184    pub(crate) fn inside(&self, n: i64, l: f64, r: f64) -> bool {
185        let t = self.instant(n);
186        let reached = |edge: f64| match t == edge {
187            true => n >= self.first_at(edge),
188            false => t > edge,
189        };
190        reached(l) && !reached(r)
191    }
192
193    pub fn position(&self, n: i64) -> f64 {
194        self.a.saturating_mul(i128::from(n)) as f64 / self.d as f64
195    }
196
197    pub fn sr(&self) -> f64 {
198        f64::from(self.rate) * self.d as f64 / self.a as f64
199    }
200
201    pub fn count(&self, t: f64) -> f64 {
202        t * f64::from(self.rate) * self.d as f64 / self.a as f64
203    }
204
205    /// The step `t` falls nearest, rounded exactly from `t`'s own binary value; `None` where
206    /// that is past what the integers hold.
207    pub fn step_at(&self, t: f64, round: Round) -> Option<i64> {
208        match self.is_rate() {
209            true => rated(t, self.rate, round).or_else(|| self.wide(t, round)),
210            false => self.wide(t, round),
211        }
212    }
213
214    fn wide(&self, t: f64, round: Round) -> Option<i64> {
215        #[cfg(test)]
216        WIDE.with(|w| w.set(w.get() + 1));
217        if !t.is_finite() {
218            return None;
219        }
220        let bits = t.abs().to_bits();
221        let (exp, frac) = ((bits >> 52) as i32, (bits & ((1 << 52) - 1)) as i128);
222        let (mantissa, shift) = match exp {
223            0 => (frac, 1074),
224            e => (frac | (1 << 52), 1075 - e),
225        };
226        let zeros = mantissa.trailing_zeros().min(127) as i32;
227        let (mantissa, shift) = match mantissa {
228            0 => (0, 0),
229            m => (m >> zeros, shift - zeros),
230        };
231        let mantissa = if t < 0.0 { -mantissa } else { mantissa };
232        let (num, den) = match shift {
233            s if s <= 0 => (mantissa.checked_mul(1i128.checked_shl((-s) as u32)?)?, 1),
234            s if s < 127 => (mantissa, 1i128 << s),
235            _ => return None,
236        };
237        let top = num.checked_mul(self.d.checked_mul(i128::from(self.rate))?)?;
238        let bottom = self.a.checked_mul(den)?;
239        let (floor, rem) = (top.div_euclid(bottom), top.rem_euclid(bottom));
240        let k = match round {
241            Round::Floor => floor,
242            Round::Ceil => floor + i128::from(rem != 0),
243            Round::Even => match (2 * rem).cmp(&bottom) {
244                std::cmp::Ordering::Less => floor,
245                std::cmp::Ordering::Greater => floor + 1,
246                std::cmp::Ordering::Equal => floor + floor.rem_euclid(2),
247            },
248        };
249        i64::try_from(k).ok()
250    }
251}
252
253/// `t*rate` rounded from the double nearest it and its exact error: halves of `t` times `rate`
254/// exactly, summed exactly; the error decides only a sum on a whole or half step. `None` where
255/// `wide` may refuse or a product could round.
256fn rated(t: f64, rate: u32, round: Round) -> Option<i64> {
257    let held = t == 0.0 || (2f64.powi(-74)..2f64.powi(51)).contains(&t.abs());
258    if !held || rate >= 1 << 26 {
259        return None;
260    }
261    let r = f64::from(rate);
262    let split = 134_217_729.0 * t;
263    let hi = split - (split - t);
264    let (x, y) = (hi * r, (t - hi) * r);
265    let p = x + y;
266    let back = p - x;
267    let e = (x - (p - back)) + (y - back);
268    if p.abs() >= 2f64.powi(51) {
269        return None;
270    }
271    let f = p.floor();
272    let k = match round {
273        Round::Floor if p == f && e < 0.0 => f - 1.0,
274        Round::Floor => f,
275        Round::Ceil if p == p.ceil() && e > 0.0 => p + 1.0,
276        Round::Ceil => p.ceil(),
277        Round::Even => match (p - f).total_cmp(&0.5) {
278            std::cmp::Ordering::Less => f,
279            std::cmp::Ordering::Greater => f + 1.0,
280            std::cmp::Ordering::Equal if e > 0.0 => f + 1.0,
281            std::cmp::Ordering::Equal if e < 0.0 => f,
282            std::cmp::Ordering::Equal => f + f.rem_euclid(2.0),
283        },
284    };
285    Some(k as i64)
286}
287
288#[cfg(test)]
289thread_local! {
290    pub(crate) static WIDE: std::cell::Cell<u64> = const { std::cell::Cell::new(0) };
291}
292
293#[derive(Clone, Copy, Debug, PartialEq, Eq, Hash, PartialOrd, Ord)]
294pub enum Round {
295    Even,
296    Floor,
297    Ceil,
298}
299
300#[cfg(test)]
301mod tests {
302    use super::{Grid, Round, rated};
303
304    /// 0.8333333333333334 s lies past 5/6 s, sample 40000's at 48 kHz, though the doubles tie.
305    #[test]
306    fn an_edge_in_seconds_tying_an_instant_is_decided_by_its_decimal() {
307        let cd = Grid::of(44_100);
308        assert_eq!(cd.first_at(0.1), 4410);
309        assert!(cd.inside(4410, 0.1, 1.0));
310        let grid = Grid::of(48_000);
311        let edge = 60.0 / 72.0;
312        assert_eq!(grid.instant(40_000), edge);
313        assert_eq!(grid.first_at(edge), 40_001);
314        assert!(!grid.inside(40_000, edge, 2.0));
315        assert!(grid.inside(40_000, 0.0, edge));
316    }
317
318    fn instants(rate: u32) -> Vec<f64> {
319        let r = f64::from(rate);
320        let mut out = vec![0.0, -0.0, 1e-300, -1e-300, 1e300, f64::MIN_POSITIVE];
321        let mut seed = 0x9e37_79b9_7f4a_7c15u64;
322        for k in -2_000i64..2_000 {
323            for half in [0.0, 0.5] {
324                let at = (k as f64 + half) / r;
325                let mut near = at;
326                for _ in 0..3 {
327                    near = near.next_up();
328                    out.push(near);
329                }
330                near = at;
331                for _ in 0..3 {
332                    near = near.next_down();
333                    out.push(near);
334                }
335                out.push(at);
336                out.push(at * 1e9);
337            }
338            seed = seed.wrapping_mul(6_364_136_223_846_793_005).wrapping_add(1);
339            out.push(f64::from_bits(seed >> 2) * if seed & 1 == 0 { 1.0 } else { -1.0 });
340        }
341        out
342    }
343
344    #[test]
345    fn a_step_rounds_alike_in_doubles_and_in_wide_integers() {
346        for rate in [1, 8_000, 44_100, 48_000, 88_200, 96_000, 192_000] {
347            let grid = Grid::of(rate);
348            for t in instants(rate) {
349                for round in [Round::Floor, Round::Ceil, Round::Even] {
350                    if let Some(k) = rated(t, rate, round) {
351                        assert_eq!(Some(k), grid.wide(t, round), "{t:e} at {rate} {round:?}");
352                    }
353                }
354            }
355        }
356    }
357}
sva_samples/grid.rs

sva_samples/
grid.rs
