网站优化哪家最专业,网站建设搭建专业网站平台公司,长春网站建设制作,热门话题推荐文章目录 1. 前置知识1.1. 非齐次线性方程组求解1.1.1. 传统求解方法1.1.2. 奇异值分解法1.1.3. 牛顿法或者梯度下降法 1.2. 齐次线性方程组的最小二乘解1.3. 非线性方程组的最小二乘解 2. 相机标定2.1. 相机内参数求解2.1.1. 求解 u 0 u_0 u0 和 v 0 v_0 v02.1.2. 求解 … 文章目录 1. 前置知识1.1. 非齐次线性方程组求解1.1.1. 传统求解方法1.1.2. 奇异值分解法1.1.3. 牛顿法或者梯度下降法 1.2. 齐次线性方程组的最小二乘解1.3. 非线性方程组的最小二乘解 2. 相机标定2.1. 相机内参数求解2.1.1. 求解 u 0 u_0 u0 和 v 0 v_0 v02.1.2. 求解 θ \theta θ2.1.3. 求解 α \alpha α 和 β \beta β 2.2. 相机外参数求解2.2.1. 求解旋转矩阵 R \mathbf{R} R2.2.2. 求解平移向量 T \mathbf{T} T 2.3. 径向畸变情况下的相机标定 3. Faugeras 定理 课程视频链接
计算机视觉之三维重建深入浅出SfM与SLAM核心算法——2.摄像机标定。 1. 前置知识
1.1. 非齐次线性方程组求解 在相机标定任务中需要求解线性方程组以确定相机内外参数 x 1 , x 2 , . . . , x q x_1, x_2, ..., x_q x1,x2,...,xq。通过采样一系列点可以列出该方程组如图所示此时通常有方程个数 p p p 参数个数 q q q构成超定方程组 A x y Ax y Axy。 当矩阵 A A A 列满秩时求解线性方程组 A x y Ax y Axy A ∈ R m × n , x ∈ R n , y ∈ R m \mathbf{A} \in \mathbb{R}^{m × n}, x \in \mathbb{R}^{n}, y \in \mathbb{R}^{m} A∈Rm×n,x∈Rn,y∈Rm的方法有如下三种
1.1.1. 传统求解方法
先证引理若 A \mathbf{A} A 是列满秩则 A T A \mathbf{A}^T\mathbf{A} ATA 可逆。 证假设存在向量 x ∈ R n x \in \mathbb{R}^n x∈Rn 满足 A T A x 0 \mathbf{A}^T \mathbf{A} x 0 ATAx0两边同时左乘 x T x^T xT x T A T A x 0 ⇒ ∣ ∣ A x ∣ ∣ 2 2 0 x^T\mathbf{A}^T \mathbf{A} x 0 \Rightarrow ||\mathbf{A} x||^2_2 0 xTATAx0⇒∣∣Ax∣∣220由范数的非负性可得 A x 0 \mathbf{A} x 0 Ax0由于 A \mathbf{A} A 是列满秩列向量线性无关方程 A x 0 \mathbf{A} x 0 Ax0 仅有零解 x 0 x 0 x0所以 A T A \mathbf{A}^T\mathbf{A} ATA 可逆。 根据引理 A x y ⇒ A T A x A T y ⇒ x ( A T A ) − 1 A T y \mathbf{A} x y \Rightarrow \mathbf{A}^T \mathbf{A} x \mathbf{A}^T y \Rightarrow x (\mathbf{A}^T \mathbf{A})^{-1} A^T y Axy⇒ATAxATy⇒x(ATA)−1ATy。
1.1.2. 奇异值分解法
设 A ∈ R m × n \mathbf{A} \in \mathbb{R}^{m × n} A∈Rm×n其奇异值分解(SVD)为 A U Σ V T \mathbf{A} \mathbf{U} \mathbf{\Sigma} \mathbf{V}^T AUΣVT其中 U ∈ R m × m \mathbf{U} \in \mathbb{R}^{m × m} U∈Rm×m 是左奇异矩阵正交阵 U U T I \mathbf{U} \mathbf{U}^T \mathbf{I} UUTI V ∈ R n × n \mathbf{V} \in \mathbb{R}^{n × n} V∈Rn×n 是右奇异矩阵正交阵 V V T I \mathbf{V} \mathbf{V}^T \mathbf{I} VVTI Σ ∈ R m × n \mathbf{\Sigma} \in \mathbb{R}^{m × n} Σ∈Rm×n 是对角阵对角线元素 σ 1 ≥ σ 2 ≥ ⋅ ⋅ ⋅ ≥ σ r \sigma_1 ≥ \sigma_2 ≥ \cdot\cdot\cdot ≥ \sigma_r σ1≥σ2≥⋅⋅⋅≥σr 为奇异值 r r a n k ( A ) r rank(\mathbf{A}) rrank(A)其余元素为 0 0 0。
则有 A x U Σ V T x 0 \mathbf{A} x \mathbf{U} \mathbf{\Sigma} \mathbf{V}^Tx 0 AxUΣVTx0左乘 V T \mathbf{V}^T VT U T U Σ V T 0 ⇒ Σ V T x 0 \mathbf{U}^T \mathbf{U} \mathbf{\Sigma} \mathbf{V}^T 0 \Rightarrow \mathbf{\Sigma} \mathbf{V}^T x 0 UTUΣVT0⇒ΣVTx0令 y V T x y \mathbf{V}^T x yVTx则 Σ y 0 \mathbf{\Sigma} y 0 Σy0。不妨设 Σ ( Σ r 0 0 0 ) \mathbf{\Sigma} \begin{pmatrix} \Sigma_r 0 \\ 0 0 \end{pmatrix} Σ(Σr000) 且 U T b c ∈ R m \mathbf{U}^T b c \in \mathbb{R}^m UTbc∈Rm则有 Σ r d i a g ( σ 1 , σ 2 , . . . , σ r ) \mathbf{\Sigma_r} diag(\sigma_{1}, \sigma_2, ..., \sigma_r) Σrdiag(σ1,σ2,...,σr)同时将 y y y 分块为 y ( y r y n − r ) y \begin{pmatrix} y_r \\ y_{n - r} \end{pmatrix} y(yryn−r)其中 y r ∈ R r y_r \in \mathbb{R}^r yr∈Rr y n − r ∈ R n − r y_{n - r} \in \mathbb{R}^{n - r} yn−r∈Rn−r则有 ( Σ r 0 0 0 ) ( y r y n − r ) ( c r c n − r ) ⇒ y ( y r y n − r ) ( Σ r − 1 c r 0 ) ⇒ x V ( Σ r − 1 c r 0 ) \begin{pmatrix} \mathbf{\Sigma_r} 0 \\ 0 0 \end{pmatrix} \begin{pmatrix} y_r \\ y_{n - r} \end{pmatrix} \begin{pmatrix} c_r \\ c_{n - r} \end{pmatrix} \Rightarrow y \begin{pmatrix} y_r \\ y_{n - r} \end{pmatrix} \begin{pmatrix} \mathbf{\Sigma_{r}^{-1}} c_r \\ 0 \end{pmatrix} \Rightarrow x \mathbf{V} \begin{pmatrix} \mathbf{\Sigma_{r}^{-1}} c_r \\ 0 \end{pmatrix} (Σr000)(yryn−r)(crcn−r)⇒y(yryn−r)(Σr−1cr0)⇒xV(Σr−1cr0)其中 c r ∈ R r c_r \in \mathbb{R}^r cr∈Rr c m − r ∈ R n − r c_{m - r} \in \mathbb{R}^{n - r} cm−r∈Rn−r。
1.1.3. 牛顿法或者梯度下降法
令 f ( x ) ∥ A x − b ∥ 2 f(\mathbf{x}) \| \mathbf{A} x - \mathbf{b} \|^2 f(x)∥Ax−b∥2即可利用牛顿法或者梯度下降法求解 x x x。
1.2. 齐次线性方程组的最小二乘解 值得注意的是求解齐次线性方程组的时候令 ∥ x ∥ 1 \|x\| 1 ∥x∥1实际上 λ x \lambda x λx 都是齐次方程组的解。
1.3. 非线性方程组的最小二乘解 2. 相机标定
相机的内、外参数矩阵描述了三维世界到二维像素的映射关系相机标定的任务就是求解摄像机内、外参数矩阵。 如下图所示相机标定问题抽象为已知世界坐标系中 P 1 , P 2 , . . . , P n P_1, P_2, ..., P_n P1,P2,...,Pn 的位置及这些点在像素坐标系上投影的位置 p 1 , p 2 , . . . , p n p_1, p_2, ..., p_n p1,p2,...,pn计算出相机的内、外参数矩阵。这里用小写 p p p 表示像素坐标系大写 P P P 表示世界坐标系。 值得注意的是选取的点 P i P_i Pi 不能位于同一平面上。 假设世界坐标系中的点 P i P_i Pi 投影到像素坐标系的点 p i p_i pi则有 p i ( u i v i ) ( m 1 P i m 3 P i m 2 P i m 3 P i ) p_i \begin{pmatrix} u_i \\ v_i \end{pmatrix} \begin{pmatrix} \dfrac{m_1 P_i}{m_3 P_i} \\ \dfrac{m_2 P_i}{m_3 P_i} \end{pmatrix} pi(uivi) m3Pim1Pim3Pim2Pi 即有 { u i m 1 P i m 3 P i v i m 2 P i m 3 P i ⇒ { u i ( m 3 P i ) m 1 P i v i ( m 3 P i ) m 2 P i ⇒ { m 1 P i − u i ( m 3 P i ) 0 m 2 P i − v i ( m 3 P i ) 0 \begin{cases} u_i \dfrac{m_1P_i}{m_3 P_i} \\ v_i \dfrac{m_2 P_i}{m_3 P_i} \end{cases} \Rightarrow \begin{cases} u_i(m_3 P_i) m_1 P_i \\ v_i(m_3 P_i) m_2 P_i \end{cases} \Rightarrow \begin{cases} m_1 P_i - u_i(m_3 P_i) 0 \\ m_2 P_i - v_i(m_3 P_i) 0 \end{cases} ⎩ ⎨ ⎧uim3Pim1Pivim3Pim2Pi⇒{ui(m3Pi)m1Pivi(m3Pi)m2Pi⇒{m1Pi−ui(m3Pi)0m2Pi−vi(m3Pi)0由于 M \mathbf{M} M 中有 11 个未知量所以至少需要 6 对点才能求解出 M \mathbf{M} M 中的所有未知量。在实验中使用多于 6 对点以获得更加鲁棒的结果。 假设有 n n n 对点 ( p i , P i ) (p_i, P_i) (pi,Pi)则有方程组 { − u 1 ( m 3 P 1 ) m 1 P 1 0 − v 1 ( m 3 P 1 ) m 2 P 1 0 ⋮ − u n ( m 3 P n ) m 1 P n 0 − v n ( m 3 P n ) m 2 P n 0 ⇒ P m → 0 \begin{cases} -u_{1} (m_{3} P_{1}) m_{1} P_{1} 0 \\ -v_{1} (m_{3} P_{1}) m_{2} P_{1} 0 \\ \quad \quad \quad \quad \quad \vdots \\ -u_{n} (m_{3} P_{n}) m_{1} P_{n} 0 \\ -v_{n} (m_{3} P_{n}) m_{2} P_{n} 0 \end{cases} \Rightarrow \mathbf{P} \overrightarrow{\boldsymbol{m}} 0 ⎩ ⎨ ⎧−u1(m3P1)m1P10−v1(m3P1)m2P10⋮−un(m3Pn)m1Pn0−vn(m3Pn)m2Pn0⇒Pm 0其中 P def ( P 1 T 0 T − u 1 P 1 T 0 T P 1 T − v 1 P 1 T ⋮ P n T 0 T − u n P n T 0 T P n T − v n P n T ) 2 n × 12 \mathbf{P} \overset{\text{def}}{} \begin{pmatrix} P_1^T 0^T -u_1 P_1^T \\ 0^T P_1^T -v_1 P_1^T \\ \vdots \\ P_n^T 0^T -u_n P_n^T \\ 0^T P_n^T -v_n P_n^T \end{pmatrix}_{2n \times 12} Pdef P1T0TPnT0T0TP1T⋮0TPnT−u1P1T−v1P1T−unPnT−vnPnT 2n×12 且 m → def ( m 1 T m 2 T m 3 T ) 12 × 1 \overrightarrow{\boldsymbol{m}} \stackrel{\text{def}}{} \left( \begin{array}{c} \mathrm{m}_{1}^{T} \\ \mathrm{~m}_{2}^{T} \\ \mathrm{~m}_{3}^{T} \end{array} \right)_{12 \times 1} m def m1T m2T m3T 12×1。 使用 SVD 求解上述超定齐次线性方程组 P m → 0 \mathbf{P} \overrightarrow{\boldsymbol{m}} 0 Pm 0如下图所示 值得注意的是实际的投影矩阵为 ρ M \rho \mathbf{M} ρM如下图所示 假设旋转矩阵 R ( r 1 T r 2 T r 3 T ) \mathbf{R} \begin{pmatrix} r_1^T \\ r_2^T \\ r_3^T \end{pmatrix} R r1Tr2Tr3T 平移向量 T ( t x t y t z ) \mathbf{T} \begin{pmatrix} t_x \\ t_y \\ t_z \end{pmatrix} T txtytz 则有 ρ M K [ R , T ] ( α − α cot θ u 0 0 β sin θ v 0 0 0 1 ) ( r 1 T t x r 2 T t y r 3 T t z ) ( α r 1 T − α cot θ r 2 T u 0 r 3 T α t x − α cot θ t y u 0 t z β sin θ r 2 T v 0 r 3 T β sin θ t y v 0 t z r 3 T t z ) \begin{align*} \rho \mathbf{M} \mathbf{K}[\mathbf{R}, T] \begin{pmatrix} \alpha -\alpha \cot \theta u_0 \\ 0 \dfrac{\beta}{\sin \theta} v_0 \\ 0 0 1 \end{pmatrix} \begin{pmatrix} r_1^T t_x \\ r_2^T t_y \\ r_3^T t_z \end{pmatrix} \\ \begin{pmatrix} \alpha r_1^T - \alpha \cot \theta\, r_2^T u_0 r_3^T \alpha t_x - \alpha \cot \theta\, t_y u_0 t_z \\ \dfrac{\beta}{\sin \theta} r_2^T v_0 r_3^T \dfrac{\beta}{\sin \theta} t_y v_0 t_z \\ r_3^T t_z \end{pmatrix} \end{align*} ρMK[R,T] α00−αcotθsinθβ0u0v01 r1Tr2Tr3Ttxtytz αr1T−αcotθr2Tu0r3Tsinθβr2Tv0r3Tr3Tαtx−αcotθtyu0tzsinθβtyv0tztz
2.1. 相机内参数求解
2.1.1. 求解 u 0 u_0 u0 和 v 0 v_0 v0
不妨设 M [ A , b ] \mathbf{M} [\mathbf{A}, b] M[A,b]则根据 ρ [ A , b ] [ ρ A , ρ b ] K [ R , T ] [ K R , K T ] \rho [\mathbf{A}, \mathbf{b}] [\rho \mathbf{A}, \rho \mathbf{b}] \mathbf{K} [\mathbf{R}, \mathbf{T}] [\mathbf{K} \mathbf{R}, \mathbf{K}\mathbf{T}] ρ[A,b][ρA,ρb]K[R,T][KR,KT] 有 ρ A K R \mathbf{\rho \mathbf{A}} \mathbf{K} \mathbf{R} ρAKR即有 ρ A ρ ( a 1 T a 2 T a 3 T ) ( α r 1 T − α cot θ r 2 T u 0 r 3 T β sin θ r 2 T v 0 r 3 T r 3 T ) K R \rho A \rho \begin{pmatrix} a_1^T \\ a_2^T \\ a_3^T \end{pmatrix} \begin{pmatrix} \alpha r_1^T - \alpha \cot \theta\, r_2^T u_0 r_3^T \\ \dfrac{\beta}{\sin \theta} r_2^T v_0 r_3^T \\ r_3^T \end{pmatrix} KR ρAρ a1Ta2Ta3T αr1T−αcotθr2Tu0r3Tsinθβr2Tv0r3Tr3T KR注意这里的 a 1 T , a 2 T , a 3 T a_1^T, a_2^T, a_3^T a1T,a2T,a3T 都是已知的。 由于旋转矩阵 R \mathbf{R} R 是正交阵所以有 r i T ⋅ r j T 0 r_i^T \cdot r_j^T 0 riT⋅rjT0 ( i ≠ j i ≠ j ij)这里的 ⋅ \cdot ⋅ 表示点积 且 ∥ r i ∥ 1 \|r_i\|1 ∥ri∥1。 一方面因为 ρ a 3 T r 3 T ⇒ ∥ ρ a 3 T ∥ ∥ r 3 T ∥ 1 \rho a_3^T r_3^T \Rightarrow \| \rho a_3^T\| \|r_3^T\|1 ρa3Tr3T⇒∥ρa3T∥∥r3T∥1所以 ρ ± 1 ∥ a 3 ∥ (1) \rho ±\dfrac{1}{\|a_3\|} \tag{1} ρ±∥a3∥1(1)这里还无法确定正负号。 另一方面 ρ a 1 T ⋅ ρ a 3 T ρ 2 a 1 T ⋅ a 3 T ( α r 1 T − α cot θ r 2 T u 0 r 3 T ) ⋅ r 3 T α r 1 T ⋅ r 3 T − α cot θ r 2 T ⋅ r 3 T u 0 r 3 T ⋅ r 3 T u 0 \begin{align*} \rho a_1^T \cdot \rho a_3^T \rho^2a_1^T \cdot a_3^T (\alpha r_1^T - \alpha \cot \theta\, r_2^T u_0 r_3^T) \cdot r_3^T \\ \alpha r_1^T \cdot r_3^T - \alpha \cot \theta\, r_2^T \cdot r_3^T u_0 r_3^T \cdot r_3^T \\ u_0 \end{align*} ρa1T⋅ρa3Tρ2a1T⋅a3T(αr1T−αcotθr2Tu0r3T)⋅r3Tαr1T⋅r3T−αcotθr2T⋅r3Tu0r3T⋅r3Tu0 ρ a 2 T ⋅ ρ a 3 T ρ 2 a 2 T ⋅ a 3 T ( β sin θ r 2 T v 0 r 3 T ) ⋅ r 3 T β sin θ r 2 T ⋅ r 3 T v 0 r 3 T ⋅ r 3 T v 0 \begin{align*} \rho a_2^T \cdot \rho a_3^T \rho^2a_2^T \cdot a_3^T \left( \frac{\beta}{\sin \theta} r_2^T v_0 r_3^T \right) \cdot r_3^T \\ \frac{\beta}{\sin \theta} r_2^T \cdot r_3^T v_0 r_3^T \cdot r_3^T \\ v_0 \end{align*} ρa2T⋅ρa3Tρ2a2T⋅a3T(sinθβr2Tv0r3T)⋅r3Tsinθβr2T⋅r3Tv0r3T⋅r3Tv0所以有 { u 0 ρ 2 a 1 ⋅ a 3 v 0 ρ 2 a 2 ⋅ a 3 (2) \begin{cases} u_0 \rho^2 a_1 \cdot a_3 \\ v_0 \rho^2 a_2 \cdot a_3 \end{cases} \tag{2} {u0ρ2a1⋅a3v0ρ2a2⋅a3(2)
2.1.2. 求解 θ \theta θ
由于旋转矩阵 R \mathbf{R} R 是正交阵所以有 r 1 T × r 2 T r 3 T , r 1 T × r 3 T r 2 T , r 3 T × r 1 T r 2 T r_1^T × r_2^T r_3^T, r_1^T × r_3^T r_2^T, r_3^T × r_1^T r_2^T r1T×r2Tr3T,r1T×r3Tr2T,r3T×r1Tr2T这里的 × × × 表示叉积。 则有 ρ a 1 T × ρ a 3 T ρ 2 ( a 1 × a 3 ) α ( r 1 T × r 3 T ) − α c o t θ ( r 2 T × r 3 T ) u 0 ( r 3 T × r 3 T ) − α r 2 − α c o t θ r 1 \begin{align*} \rho a_1^T × \rho a_3^T \rho^2(a_1 × a_3) \alpha (r_1^T × r_3^T) - \alpha cot\theta (r_2^T × r_3^T) u_0 (r_3^T × r_3^T) \\ -\alpha r_2 - \alpha cot\theta r_1 \end{align*} ρa1T×ρa3Tρ2(a1×a3)α(r1T×r3T)−αcotθ(r2T×r3T)u0(r3T×r3T)−αr2−αcotθr1两边同时取模则有 ρ 2 ∥ a 1 × a 3 ∥ ∥ α r 2 α c o t θ r 1 ∥ \rho^2 \|a_1 × a_3\| \|\alpha r_2 \alpha cot\theta r_1\| \\ ρ2∥a1×a3∥∥αr2αcotθr1∥又因为 ∥ α r 2 α c o t θ r 1 ∥ 2 ( α r 2 α c o t θ r 1 ) T ⋅ ( α r 2 α c o t θ r 1 ) α 2 α 2 c o t 2 θ α 2 ( 1 c o s 2 θ s i n 2 θ ) α 2 s i n 2 θ \begin{align*} \|\alpha r_2 \alpha cot\theta r_1\|^2 (\alpha r_2 \alpha cot\theta r_1)^T \cdot (\alpha r_2 \alpha cot\theta r_1) \alpha^2 \alpha^2cot^2\theta \\ \alpha^2(1 \dfrac{cos^2\theta}{sin^2\theta}) \dfrac{\alpha^2}{sin^2\theta} \end{align*} ∥αr2αcotθr1∥2(αr2αcotθr1)T⋅(αr2αcotθr1)α2α2cot2θα2(1sin2θcos2θ)sin2θα2即有 ρ 2 ∥ a 1 × a 3 ∥ α s i n θ (3) \rho^2 \|a_1 × a_3\| \dfrac{\alpha}{sin\theta} \tag{3} ρ2∥a1×a3∥sinθα(3)注意这里的 α f x d x 0 \alpha \dfrac{f_x}{dx} 0 αdxfx0 且 0 θ ≤ π 2 0 \theta ≤ \dfrac{\pi}{2} 0θ≤2π。 又因为 ρ a 2 T × ρ a 3 T ρ 2 ( a 2 × a 3 ) β s i n θ r 2 T × r 3 T v 0 r 3 T × r 3 T β s i n θ r 1 \begin{align*} \rho a_2^T × \rho a_3^T \rho^2(a_2 × a_3) \dfrac{\beta}{sin\theta} r_2^T × r_3^T v_0 r_3^T × r_3^T \\ \dfrac{\beta}{sin\theta} r_1 \end{align*} ρa2T×ρa3Tρ2(a2×a3)sinθβr2T×r3Tv0r3T×r3Tsinθβr1所以 ρ 2 ∥ a 2 × a 3 ∥ β s i n θ (4) \rho^2 \|a_2 × a_3\| \dfrac{\beta}{sin\theta} \tag{4} ρ2∥a2×a3∥sinθβ(4)从而 { ρ 2 ( a 1 × a 3 ) ⋅ ρ 2 ( a 2 × a 3 ) ( α r 2 − α c o t θ r 1 ) ⋅ β s i n θ − α β c o t θ s i n θ − α β c o s θ s i n 2 θ ρ 2 ∥ a 1 × a 3 ∥ ⋅ ρ 2 ∥ a 2 × a 3 ∥ α β s i n 2 θ \begin{cases} \rho^2(a_1 × a_3) \cdot \rho^2(a_2 × a_3) (\alpha r_2 - \alpha cot\theta r_1) \cdot \dfrac{\beta}{sin\theta}-\dfrac{\alpha \beta cot\theta}{sin\theta} -\dfrac{\alpha \beta cos\theta}{sin^2\theta} \\ \rho^2 \|a_1 × a_3\| \cdot \rho^2\|a_2×a_3\| \dfrac{\alpha \beta}{sin^2\theta} \end{cases} ⎩ ⎨ ⎧ρ2(a1×a3)⋅ρ2(a2×a3)(αr2−αcotθr1)⋅sinθβ−sinθαβcotθ−sin2θαβcosθρ2∥a1×a3∥⋅ρ2∥a2×a3∥sin2θαβ两式相除有 c o s θ − ( a 1 × a 3 ) ⋅ ( a 2 × a 3 ) ∥ a 1 × a 3 ∥ ⋅ ∥ a 2 × a 3 ∥ (5) cos\theta-\dfrac{(a_1 × a_3) \cdot(a_2 × a_3)}{\|a_1 × a_3\| \cdot \|a_2 × a_3\|} \tag{5} cosθ−∥a1×a3∥⋅∥a2×a3∥(a1×a3)⋅(a2×a3)(5)
2.1.3. 求解 α \alpha α 和 β \beta β
根据公式 ( 3 ) (3) (3) 和 ( 4 ) (4) (4)可得 { α ρ 2 ∥ a 1 × a 3 ∥ s i n θ β ρ 2 ∥ a 2 × a 3 ∥ s i n θ (6) \begin{cases} \alpha \rho^2 \| a_1 × a_3 \| sin\theta \\ \beta \rho^2 \|a_2 × a_3 \| sin \theta \end{cases} \tag{6} {αρ2∥a1×a3∥sinθβρ2∥a2×a3∥sinθ(6) 至此 5 个相机内参数已经全部求出分别对应公式 ( 1 ) , ( 2 ) , ( 5 ) , ( 6 ) (1), (2), (5), (6) (1),(2),(5),(6)。
2.2. 相机外参数求解
2.2.1. 求解旋转矩阵 R \mathbf{R} R
因为 ρ a 3 T r 3 T \rho a_3^T r_3^T ρa3Tr3T 且 ρ ± 1 ∥ a 3 ∥ \rho ±\dfrac{1}{\|a_3\|} ρ±∥a3∥1所以 r 3 ± a 3 ∥ a 3 ∥ r_3 ± \dfrac{a_3}{\|a_3\|} r3±∥a3∥a3又因为 ρ 2 ( a 2 × a 3 ) β s i n θ r 1 \rho^2(a_2 × a_3) \dfrac{\beta}{sin\theta} r_1 ρ2(a2×a3)sinθβr1即 a 2 × a 3 a_2 × a_3 a2×a3 与 r 1 r_1 r1 同方向所以 r 1 a 2 × a 3 ∥ a 2 × a 3 ∥ r_1 \dfrac{a_2 × a_3}{\|a_2 × a_3\|} r1∥a2×a3∥a2×a3最后有 r 2 r 3 × r 1 r_2 r_3 × r_1 r2r3×r1综上有 { r 3 ± a 3 ∥ a 3 ∥ r 1 a 2 × a 3 ∥ a 2 × a 3 ∥ r 2 r 3 × r 1 (7) \begin{cases} r_3 ± \dfrac{a_3}{\|a_3\|} \\ r_1 \dfrac{a_2 × a_3}{\|a_2 × a_3\|} \\ r_2 r_3 × r_1 \end{cases} \tag{7} ⎩ ⎨ ⎧r3±∥a3∥a3r1∥a2×a3∥a2×a3r2r3×r1(7)
2.2.2. 求解平移向量 T \mathbf{T} T
因为 ρ b K T \rho b \mathbf{K} \mathbf{T} ρbKT所以 T ρ K − 1 b \mathbf{T} \rho \mathbf{K}^{-1} b TρK−1b由于 K ( α − α cot θ u 0 0 β sin θ v 0 0 0 1 ) \mathbf{K} \begin{pmatrix} \alpha -\alpha \cot \theta u_0 \\ 0 \dfrac{\beta}{\sin \theta} v_0 \\ 0 0 1 \end{pmatrix} K α00−αcotθsinθβ0u0v01 显然 K \mathbf{K} K 是满秩矩阵所以 K \mathbf{K} K 是可逆的。 至此相机的内外参数都求解出来了如下图所示
2.3. 径向畸变情况下的相机标定
前面讨论的内容没有考虑径向畸变的情况实际上 M P i \mathbf{M} P_{i} MPi 计算的是理想情况下无径向畸变的像素坐标 ( u ′ v ′ ) \begin{pmatrix} u \\ v \end{pmatrix} (u′v′)畸变前的像素坐标与畸变后的像素坐标 ( u v ) \begin{pmatrix} u \\ v \end{pmatrix} (uv) 满足如下关系 { u u ′ ( 1 k 1 r 2 k 2 r 4 k 3 r 6 ) v v ′ ( 1 k 1 r 2 k 2 r 4 k 3 r 6 ) \begin{cases} u u (1 k_1 r^2 k_2 r^4 k_3 r^6) \\ v v (1 k_1 r^2 k_2 r^4 k_3 r^6) \end{cases} {uu′(1k1r2k2r4k3r6)vv′(1k1r2k2r4k3r6)其中 k 1 , k 2 , k 3 k_1, k_2, k_3 k1,k2,k3 是畸变系数 r 2 ( u ′ ) 2 ( v ′ ) 2 r^2 (u)^2 (v)^2 r2(u′)2(v′)2。 不妨令 λ ( 1 k 1 r 2 k 2 r 4 k 3 r 6 ) \lambda (1 k_1 r^2 k_2 r^4 k_3 r^6) λ(1k1r2k2r4k3r6)则有 ( λ 0 0 0 λ 0 0 0 1 ) M P i ( u i v i ) p i \begin{pmatrix} \lambda 0 0 \\ 0 \lambda 0 \\ 0 0 1 \end{pmatrix} \mathbf{M} P_i \begin{pmatrix} u_i \\ v_i \end{pmatrix} p_i λ000λ0001 MPi(uivi)pi记 Q ( λ 0 0 0 λ 0 0 0 1 ) M ( q 1 q 2 q 3 ) \mathbf{Q} \begin{pmatrix} \lambda 0 0 \\ 0 \lambda 0 \\ 0 0 1 \end{pmatrix} \mathbf{M} \begin{pmatrix} q_1 \\ q_2 \\ q_3 \end{pmatrix} Q λ000λ0001 M q1q2q3 则有 p i ( u i v i ) ( q 1 P i q 3 P i q 2 P i q 3 P i ) ⇒ { u i q 3 P i q 1 P i v i q 3 P i q 2 P i ⇒ { u i q 3 P i − q 1 P i 0 v i q 3 P i − q 2 P i 0 p_i \begin{pmatrix} u_i \\ v_i \end{pmatrix} \begin{pmatrix} \dfrac{q_1 P_i}{q_3 P_i} \\ \dfrac{q_2 P_i}{q_3 P_i} \end{pmatrix} \Rightarrow \begin{cases} u_i q_3 P_i q_1 P_i \\ v_i q_3 P_i q_2 P_i \end{cases} \Rightarrow \begin{cases} u_i q_3 P_i - q_1 P_i 0 \\ v_i q_3 P_i - q_2 P_i 0 \end{cases} pi(uivi) q3Piq1Piq3Piq2Pi ⇒{uiq3Piq1Piviq3Piq2Pi⇒{uiq3Pi−q1Pi0viq3Pi−q2Pi0我们无法用之前的方法求解出 k 1 , k 2 , k 3 , m 1 , m 2 , m 3 k_1, k_2, k_3, m_1, m_2, m_3 k1,k2,k3,m1,m2,m3因为方程组是非线性的。因此采用非线性求解方法进行求解如下图所示 从初始解开始迭代若初始解与实际相距较远可能会很慢。为了加快算法求解的速度可以先求解出 m 1 m_1 m1 和 m 2 m_2 m2这样后续需要求解的变量就从 k 1 , k 2 , k 3 , m 1 , m 2 , m 3 k_1, k_2, k_3, m_1, m_2, m_3 k1,k2,k3,m1,m2,m3 减少到 k 1 , k 2 , k 3 , m 3 k_1, k_2, k_3, m_3 k1,k2,k3,m3求解 m 1 m_1 m1 和 m 2 m_2 m2 的方法如下图所示
3. Faugeras 定理 第一点证明如下 先证必要性如果 M \mathbf{M} M 是透视投影矩阵则相机内参矩阵 K \mathbf{K} K 满足 d e t ( K ) ≠ 0 det(\mathbf{K}) ≠ 0 det(K)0旋转矩阵 R \mathbf{R} R 满足 d e t ( R ) ≠ 0 det(\mathbf{R}) ≠ 0 det(R)0所以 d e t ( A ) d e t ( K R ) d e t ( K ) d e t ( R ) ≠ 0 det(\mathbf{A}) det(\mathbf{K} \mathbf{R}) det(\mathbf{K}) det(\mathbf{R}) ≠ 0 det(A)det(KR)det(K)det(R)0。 再证充分性由于 d e t ( A ) ≠ 0 det(\mathbf{A}) ≠ 0 det(A)0所以 A \mathbf{A} A 可逆将 A \mathbf{A} A 进行 QR 分解有 A K R \mathbf{A} \mathbf{K} \mathbf{R} AKR其中 K \mathbf{K} K 是可逆的上三角矩阵满足内参矩阵结构 R \mathbf{R} R 是正交矩阵满足旋转矩阵性质。进而可以还原完整的相机投影矩阵 M [ A , b ] K [ R , K − 1 b ] K [ R , T ] \mathbf{M} [\mathbf{A}, b] \mathbf{K} [\mathbf{R}, \mathbf{K}^{-1}b] \mathbf{K}[\mathbf{R}, T] M[A,b]K[R,K−1b]K[R,T]证毕。
第二点证明 由公式 ( 5 ) (5) (5) 可知 c o s θ − ( a 1 × a 3 ) ⋅ ( a 2 × a 3 ) ∥ a 1 × a 3 ∥ ⋅ ∥ a 2 × a 3 ∥ cos\theta-\dfrac{(a_1 × a_3) \cdot(a_2 × a_3)}{\|a_1 × a_3\| \cdot \|a_2 × a_3\|} cosθ−∥a1×a3∥⋅∥a2×a3∥(a1×a3)⋅(a2×a3)。因此 M \mathbf{M} M 为零倾斜透视投影矩阵 ⇔ \Leftrightarrow ⇔ θ π 2 \theta \dfrac{\pi}{2} θ2π ⇔ \Leftrightarrow ⇔ ( a 1 × a 3 ) ⋅ ( a 2 × a 3 ) 0 (a_1 × a_3) \cdot(a_2 × a_3) 0 (a1×a3)⋅(a2×a3)0。 证毕。
第三点证明 根据公式 ( 6 ) (6) (6) 可知 { α ρ 2 ∥ a 1 × a 3 ∥ s i n θ β ρ 2 ∥ a 2 × a 3 ∥ s i n θ \begin{cases} \alpha \rho^2 \| a_1 × a_3 \| sin\theta \\ \beta \rho^2 \|a_2 × a_3 \| sin \theta \end{cases} {αρ2∥a1×a3∥sinθβρ2∥a2×a3∥sinθ则宽高比为 1 ⇔ \Leftrightarrow ⇔ ∥ a 1 × a 3 ∥ ∥ a 2 × a 3 ∥ \| a_1 × a_3 \| \|a_2 × a_3 \| ∥a1×a3∥∥a2×a3∥ ⇔ \Leftrightarrow ⇔ ( a 1 × a 3 ) ⋅ ( a 1 × a 3 ) ( a 2 × a 3 ) ⋅ ( a 2 × a 3 ) (a_1 × a_3) \cdot (a_1 × a_3) (a_2 × a_3) \cdot (a_2 × a_3) (a1×a3)⋅(a1×a3)(a2×a3)⋅(a2×a3)。 结合第二点的证明则可以完成第三点的证明。 证毕。