idsp/
pll.rs

1use core::num::Wrapping as W;
2use dsp_fixedpoint::Q32;
3use dsp_process::{Process, SplitProcess};
4
5use crate::ClampWrap;
6
7/// Type-2, order-3 sampled phase, discrete time PLL
8///
9/// This PLL tracks the frequency and phase of an input signal with respect to the sampling clock.
10/// The open loop transfer function is type 2 (DC double integrator) from input phase to output phase.
11///
12/// The transfer functions (for phase and frequency) contain an additional zero at Nyquist.
13///
14/// The PLL locks to any frequency (i.e. it locks to the alias in the first Nyquist zone) and is
15/// stable for any numerically valid gain (in units of the sample rate: 7e-5 to 5e-2).
16/// It has a single parameter that determines the loop bandwidth.
17/// The gains can be changed freely between updates.
18///
19/// All math is naturally wrapping 32 bit integer. Phase and frequency are understood modulo that
20/// overflow in the first Nyquist zone. Expressing the IIR equations in other ways (e.g. single
21/// (T)-DF-{I,II} biquad/IIR) would break on overflow (i.e. every cycle).
22///
23/// There are no floating point rounding errors. The integer quantization/truncation
24/// error is fed back (first order noise shaping).
25///
26/// This PLL clamps phase wraps accumulated during (frequency) lock acquisition.
27///
28/// The phase detector is symmetric (additive): the loop filter has negative gain.
29/// The output will compensate the input phase: it will settle to the complement.
30/// The output phase increment (the loop filter output, the frequency) is the negative
31/// of the input increment.
32#[derive(Debug, Clone, Default)]
33pub struct PLL {
34    /// Lead lag coefficients
35    ///
36    /// `f0 += b0*y0 + b1*y1 + a1*f1`
37    pub ba: [Q32<32>; 3],
38}
39
40impl PLL {
41    /// Return Pll from zero/pole/gain
42    pub fn from_zpk(zero: f32, pole: f32, gain: f32) -> Self {
43        Self {
44            ba: [gain, -gain * zero, -(1.0 - pole)].map(Q32::from_f32),
45        }
46    }
47
48    /// Given a crossover, create a PLL
49    ///
50    /// About 1.5 dB peaking, 62 deg phase margin for split=4
51    pub fn from_bandwidth(bw: f32, split: f32) -> Self {
52        let a = bw * 2.0 * core::f32::consts::PI;
53        let z = 1.0 - a / split;
54        let p = 1.0 - a * split;
55        let k = -a * a * split;
56        Self::from_zpk(z, p, k)
57    }
58}
59
60/// PLL state
61#[derive(Debug, Clone, Default)]
62pub struct PLLState {
63    /// Input phase difference clamp
64    pub clamp: ClampWrap<W<i32>>,
65    /// Loop filter state: after clamp
66    pub z0: i32,
67    /// After nyquist zero
68    pub y0: i32,
69    /// After lead-lag
70    pub f0: i64,
71    /// After DC pole
72    pub f: W<i64>,
73    /// Current output phase
74    pub y: W<i32>,
75}
76
77impl PLLState {
78    /// Return the current phase estimate
79    pub fn phase(&self) -> W<i32> {
80        self.y
81    }
82
83    /// Return the current frequency estimate
84    pub fn frequency(&self) -> W<i32> {
85        W((self.f.0 >> 32) as _)
86    }
87}
88
89impl SplitProcess<W<i32>, W<i32>, PLLState> for PLL {
90    fn process(&self, state: &mut PLLState, x: W<i32>) -> W<i32> {
91        // advance output phase, oscillator DC pole
92        state.y += state.frequency();
93        // phase error
94        let z0 = state.clamp.process(x + state.y).0 >> 1;
95        // nyquist zero
96        let y0 = z0 + state.z0;
97        state.z0 = z0;
98        // lead lag, wide state
99        state.f0 +=
100            (self.ba[0] * y0 + self.ba[1] * state.y0 + self.ba[2] * (state.f0 >> 32) as i32).inner
101                + ((self.ba[2].inner as i64 * state.f0 as u32 as i64) >> 32);
102        state.y0 = y0;
103        // DC pole, frequency, wide state
104        state.f += W(state.f0);
105        state.y
106    }
107}
108
109#[cfg(test)]
110mod tests {
111    use crate::Accu;
112
113    use super::*;
114    use core::num::Wrapping as W;
115
116    #[test]
117    fn converge_pll() {
118        let p = PLL::from_bandwidth(5e-2, 4.0);
119        println!("{p:?}");
120        let mut s = PLLState::default();
121        let a = Accu::<W<i32>>::new(W(0x0), W(0x71f63049));
122        let n = 1 << 9;
123        for (i, x) in a.take(n).enumerate() {
124            let y = p.process(&mut s, x);
125            println!("x: {x:#010x} y+x: {:#010x}", y + x);
126            if i > n / 2 {
127                assert!((a.step + s.frequency()).0.abs() <= 1);
128                assert!((x + y).0.abs() <= 4);
129            }
130        }
131    }
132
133    #[test]
134    fn converge_narrow() {
135        let p = PLL::from_bandwidth(8e-5, 4.0);
136        println!("{p:?}");
137        let mut s = PLLState::default();
138        let a = Accu::<W<i32>>::new(W(0x0), W(0x140_1235));
139        let n = 1 << 18;
140        for (i, x) in a.take(n).enumerate() {
141            let y = p.process(&mut s, x);
142            println!("x: {x:#010x} y+x: {:#010x}", y + x);
143            if i > n / 2 {
144                assert!((a.step + s.frequency()).0.abs() <= 1 << 16);
145                assert!((x + y).0.abs() <= 1 << 16);
146            }
147        }
148    }
149}