scan_core/transform/
perspective.rs

1use crate::Point;
2
3/// Port of perspective-transform (js)
4/// Solves for coefficients mapping src points to dst points.
5pub struct PerspectiveTransform {
6    coeffs: [f32; 9],
7    coeffs_inv: [f32; 9],
8}
9
10impl PerspectiveTransform {
11    pub fn new(src: &[Point], dst: &[Point]) -> Option<Self> {
12        if src.len() < 4 || dst.len() != src.len() {
13            return None;
14        }
15        
16        let coeffs = get_normalization_coefficients(src, dst, false)?;
17        let coeffs_inv = get_normalization_coefficients(src, dst, true)?;
18        
19        Some(Self {
20            coeffs,
21            coeffs_inv,
22        })
23    }
24    
25    pub fn transform(&self, x: f32, y: f32) -> Point {
26        let c = self.coeffs;
27        let den = c[6] * x + c[7] * y + 1.0;
28        Point {
29            x: (c[0] * x + c[1] * y + c[2]) / den,
30            y: (c[3] * x + c[4] * y + c[5]) / den,
31        }
32    }
33    
34    pub fn transform_inverse(&self, x: f32, y: f32) -> Point {
35        let c = self.coeffs_inv;
36        let den = c[6] * x + c[7] * y + 1.0;
37        Point {
38            x: (c[0] * x + c[1] * y + c[2]) / den,
39            y: (c[3] * x + c[4] * y + c[5]) / den,
40        }
41    }
42}
43
44fn get_normalization_coefficients(src: &[Point], dst: &[Point], is_inverse: bool) -> Option<[f32; 9]> {
45    let (src, dst) = if is_inverse { (dst, src) } else { (src, dst) };
46    let n = src.len();
47    
48    // We want to solve Mh = B for h = [h0..h7] (assuming h8=1).
49    // M is (2*n x 8), B is (2*n x 1).
50    // To solve via Least Squares for n >= 4:
51    // Solve (M^T * M) * h = (M^T * B).
52    // Let A = M^T * M (8x8) and Y = M^T * B (8x1).
53    
54    let mut a = [[0.0f32; 8]; 8];
55    let mut y = [0.0f32; 8];
56    
57    for i in 0..n {
58        let sx = src[i].x;
59        let sy = src[i].y;
60        let dx = dst[i].x;
61        let dy = dst[i].y;
62        
63        // Rows of M:
64        // [sx, sy, 1, 0,  0,  0, -dx*sx, -dx*sy]  = dx
65        // [0,  0,  0, sx, sy, 1, -dy*sx, -dy*sy]  = dy
66        
67        let r0 = [sx, sy, 1.0, 0.0, 0.0, 0.0, -dx * sx, -dx * sy];
68        let r1 = [0.0, 0.0, 0.0, sx, sy, 1.0, -dy * sx, -dy * sy];
69        
70        // Update A = Sum(r_i^T * r_i) and Y = Sum(r_i^T * b_i)
71        for row in 0..8 {
72            for col in 0..8 {
73                a[row][col] += r0[row] * r0[col] + r1[row] * r1[col];
74            }
75            y[row] += r0[row] * dx + r1[row] * dy;
76        }
77    }
78    
79    // Now solve 8x8 system A * h = Y using Gaussian Elimination
80    let mut matrix = [[0.0f32; 9]; 8];
81    for i in 0..8 {
82        for j in 0..8 {
83            matrix[i][j] = a[i][j];
84        }
85        matrix[i][8] = y[i];
86    }
87    
88    // Gaussian elimination (same as before)
89    let size = 8;
90    for i in 0..size {
91        let mut pivot_row = i;
92        for j in i + 1..size {
93            if matrix[j][i].abs() > matrix[pivot_row][i].abs() {
94                pivot_row = j;
95            }
96        }
97        matrix.swap(i, pivot_row);
98        let pivot = matrix[i][i];
99        if pivot.abs() < 1e-9 { return None; }
100        
101        for j in i..=size { matrix[i][j] /= pivot; }
102        for k in 0..size {
103            if k != i {
104                let factor = matrix[k][i];
105                for j in i..=size { matrix[k][j] -= factor * matrix[i][j]; }
106            }
107        }
108    }
109    
110    let mut res = [0.0; 9];
111    for i in 0..8 {
112        res[i] = matrix[i][8];
113    }
114    res[8] = 1.0;
115    
116    Some(res)
117}
scan_core/transform/perspective.rs

scan_core/transform/
perspective.rs
