qualia_core_db/solvers/transforms/
ztransform.rs1use super::fourier::Cplx;
6
7#[inline]
8fn cmul(a: Cplx, b: Cplx) -> Cplx {
9 (a.0 * b.0 - a.1 * b.1, a.0 * b.1 + a.1 * b.0)
10}
11fn cinv(z: Cplx) -> Option<Cplx> {
13 let d = z.0 * z.0 + z.1 * z.1;
14 if d == 0.0 {
15 return None;
16 }
17 Some((z.0 / d, -z.1 / d))
18}
19fn csub(a: Cplx, b: Cplx) -> Cplx {
20 (a.0 - b.0, a.1 - b.1)
21}
22
23pub fn z_transform_finite(x: &[f64], z: Cplx) -> Option<Cplx> {
25 let zinv = cinv(z)?;
26 let mut acc = (0.0, 0.0);
27 let mut zpow = (1.0, 0.0); for &xn in x {
29 acc = (acc.0 + xn * zpow.0, acc.1 + xn * zpow.1);
30 zpow = cmul(zpow, zinv);
31 }
32 Some(acc)
33}
34
35pub fn unit_step_z(z: Cplx) -> Option<Cplx> {
38 let zinv = cinv(z)?;
39 let denom = csub((1.0, 0.0), zinv); let id = cinv(denom)?;
41 Some(id)
42}
43
44pub fn geometric_z(a: f64, z: Cplx) -> Option<Cplx> {
46 let zinv = cinv(z)?;
47 let azinv = (a * zinv.0, a * zinv.1);
48 let denom = csub((1.0, 0.0), azinv);
49 cinv(denom)
50}
51
52#[cfg(test)]
53mod tests {
54 use super::*;
55 const EPS: f64 = 1e-9;
56
57 #[test]
58 fn finite_sequence_evaluates() {
59 let v = z_transform_finite(&[1.0, 2.0, 3.0], (2.0, 0.0)).unwrap();
61 assert!((v.0 - 2.75).abs() < EPS && v.1.abs() < EPS);
62 let d = z_transform_finite(&[1.0], (3.0, -1.0)).unwrap();
64 assert!((d.0 - 1.0).abs() < EPS && d.1.abs() < EPS);
65 assert!(z_transform_finite(&[1.0], (0.0, 0.0)).is_none());
66 }
67
68 #[test]
69 fn closed_forms_match_truncated_sums() {
70 let z = (2.0, 0.0);
72 let closed = geometric_z(0.5, z).unwrap();
73 let seq: Vec<f64> = (0..60).map(|n| 0.5_f64.powi(n)).collect();
74 let approx = z_transform_finite(&seq, z).unwrap();
75 assert!((closed.0 - approx.0).abs() < 1e-6 && (closed.1 - approx.1).abs() < 1e-6);
76 let us = unit_step_z(z).unwrap();
78 let ones = vec![1.0; 60];
79 let ua = z_transform_finite(&ones, z).unwrap();
80 assert!((us.0 - ua.0).abs() < 1e-6);
81 }
82}