<think>我们已知DLT系数（即投影矩阵M）描述了理想针孔相机模型下的投影关系。但当相机存在透镜畸变时，实际成像点会偏离理想投影位置。因此，修正DLT系数的核心思路是：**先通过畸变模型校正图像坐标，再用校正后的坐标进行DLT计算**。以下是具体步骤：
---
### **1. 理解畸变模型**
透镜畸变主要包括**径向畸变**和**切向畸变**：
- **径向畸变**：由镜头形状引起，成像点沿径向偏移
$$ \begin{cases}
x_d = x(1 + k_1 r^2 + k_2 r^4 + k_3 r^6) \\
y_d = y(1 + k_1 r^2 + k_2 r^4 + k_3 r^6)
\end{cases} \quad \text{其中} \quad r^2 = x^2 + y^2 $$
- **切向畸变**：由镜头安装不平行引起
$$ \begin{cases}
x_d = x + [2p_1 xy + p_2(r^2 + 2x^2)] \\
y_d = y + [p_1(r^2 + 2y^2) + 2p_2 xy]
\end{cases} $$
综合畸变模型为：
$$ \begin{bmatrix} u_d \\ v_d \end{bmatrix} = \begin{bmatrix} u \\ v \end{bmatrix} + \underbrace{ \begin{bmatrix}
u(k_1 r^2 + k_2 r^4 + k_3 r^6) + (2p_1 uv + p_2(r^2 + 2u^2)) \\
v(k_1 r^2 + k_2 r^4 + k_3 r^6) + (p_1(r^2 + 2v^2) + 2p_2 uv)
\end{bmatrix} }_{\text{畸变偏移}} \tag{1} $$
其中 $(u, v)$ 是理想坐标，$(u_d, v_d)$ 是实际观测坐标[^2][^4]。
---
### **2. 修正DLT系数的流程**
#### **步骤1：获取控制点数据**
- 物方点：$P_i = (X_i, Y_i, Z_i)$
- 实际图像点：$p_i = (u_{d_i}, v_{d_i})$（含畸变）
#### **步骤2：估计初始DLT系数（忽略畸变）**
用含畸变的图像点直接求解初始投影矩阵 $M_0$：
1. 对每个点构建方程：
$$ \begin{bmatrix}
X_i & Y_i & Z_i & 1 & 0 & 0 & 0 & 0 & -u_{d_i}X_i & -u_{d_i}Y_i & -u_{d_i}Z_i \\
0 & 0 & 0 & 0 & X_i & Y_i & Z_i & 1 & -v_{d_i}X_i & -v_{d_i}Y_i & -v_{d_i}Z_i
\end{bmatrix} \begin{bmatrix} \mathbf{m} \end{bmatrix} = \begin{bmatrix} u_{d_i} \\ v_{d_i} \end{bmatrix} $$
（其中 $\mathbf{m}$ 是 $M_0$ 的向量形式）
2. 用SVD求解超定方程组，得到 $M_0$
#### **步骤3：估计畸变系数**
固定 $M_0$，联合优化畸变参数 $\mathbf{k} = (k_1, k_2, p_1, p_2)$：
1. 用 $M_0$ 计算理想坐标 $(u_i, v_i)$：
$$ \begin{bmatrix} u_i \\ v_i \\ 1 \end{bmatrix} \simeq M_0 \begin{bmatrix} X_i \\ Y_i \\ Z_i \\ 1 \end{bmatrix} $$
2. 建立畸变校正方程（反向映射）：
$$ \begin{cases}
u_i = u_{d_i} + \Delta u(u_{d_i}, v_{d_i}, \mathbf{k}) \\
v_i = v_{d_i} + \Delta v(u_{d_i}, v_{d_i}, \mathbf{k})
\end{cases} \tag{2} $$
3. 将 $(2)$ 式线性化（以 $\Delta u$ 为例）：
$$ \Delta u \approx \frac{\partial \Delta u}{\partial k_1} \delta k_1 + \frac{\partial \Delta u}{\partial k_2} \delta k_2 + \frac{\partial \Delta u}{\partial p_1} \delta p_1 + \frac{\partial \Delta u}{\partial p_2} \delta p_2 $$
其中雅可比矩阵为：
$$ J = \begin{bmatrix}
\frac{\partial \Delta u}{\partial k_1} & \frac{\partial \Delta u}{\partial k_2} & \frac{\partial \Delta u}{\partial p_1} & \frac{\partial \Delta u}{\partial p_2} \\
\frac{\partial \Delta v}{\partial k_1} & \frac{\partial \Delta v}{\partial k_2} & \frac{\partial \Delta v}{\partial p_1} & \frac{\partial \Delta v}{\partial p_2}
\end{bmatrix} = \begin{bmatrix}
u_d r^2 & u_d r^4 & 2u_d v_d & r^2 + 2u_d^2 \\
v_d r^2 & v_d r^4 & r^2 + 2v_d^2 & 2u_d v_d
\end{bmatrix} $$
4. 构建最小二乘问题：
$$ \min_{\delta \mathbf{k}} \sum_i \left\| \begin{bmatrix} u_i - u_{d_i} \\ v_i - v_{d_i} \end{bmatrix} - J \delta \mathbf{k} \right\|^2 $$
解得 $\delta \mathbf{k}$ 并更新 $\mathbf{k} \leftarrow \mathbf{k} + \delta \mathbf{k}$
#### **步骤4：联合优化DLT系数和畸变参数**
将 $M$ 和 $\mathbf{k}$ 同时作为优化变量：
1. 定义重投影误差：
$$ E = \sum_i \left\| \begin{bmatrix} u_{d_i} \\ v_{d_i} \end{bmatrix} - \text{distort} \left( M \begin{bmatrix} X_i \\ Y_i \\ Z_i \\ 1 \end{bmatrix}, \mathbf{k} \right) \right\|^2 $$
（其中 $\text{distort}(\cdot)$ 表示应用畸变模型）
2. 使用Levenberg-Marquardt算法优化 $E$
---
### **3. 代码实现框架（Python）**
```python
import numpy as np
from scipy.optimize import least_squares
# 步骤1: 加载控制点 (world_pts, image_pts)
world_pts = np.loadtxt('world_points.txt') # Nx3
image_pts = np.loadtxt('image_points.txt') # Nx2
# 步骤2: 初始DLT求解 (忽略畸变)
def dlt_calib(world_pts, image_pts):
A = []
for i in range(len(world_pts)):
X, Y, Z = world_pts[i]
u, v = image_pts[i]
A.append([X, Y, Z, 1, 0, 0, 0, 0, -u*X, -u*Y, -u*Z, -u])
A.append([0, 0, 0, 0, X, Y, Z, 1, -v*X, -v*Y, -v*Z, -v])
U, S, Vt = np.linalg.svd(np.array(A))
M = Vt[-1].reshape(3,4)
return M / M[2,3] # 归一化
M_init = dlt_calib(world_pts, image_pts)
# 步骤3: 估计畸变系数 (k1,k2,p1,p2)
def distort_correct(params, M, world_pts, image_pts):
k1, k2, p1, p2 = params
errors = []
for i in range(len(world_pts)):
# 计算理想投影
P = np.append(world_pts[i], 1)
uv = M @ P
u_ideal, v_ideal = uv[0]/uv[2], uv[1]/uv[2]
# 计算畸变偏移
u_obs, v_obs = image_pts[i]
r2 = u_obs**2 + v_obs**2
delta_u = u_obs * (k1*r2 + k2*r2**2) + 2*p1*u_obs*v_obs + p2*(r2 + 2*u_obs**2)
delta_v = v_obs * (k1*r2 + k2*r2**2) + p1*(r2 + 2*v_obs**2) + 2*p2*u_obs*v_obs
# 计算误差 (理想点 = 观测点 + 反向偏移)
errors.append(u_ideal - (u_obs + delta_u))
errors.append(v_ideal - (v_obs + delta_v))
return errors
k_init = [0,0,0,0]
res = least_squares(distort_correct, k_init, args=(M_init, world_pts, image_pts))
k_opt = res.x
# 步骤4: 联合优化
def full_reprojection_error(params, world_pts, image_pts):
# 解包参数: M(12维) + k(4维)
M = params[:12].reshape(3,4)
k = params[12:]
errors = []
for i in range(len(world_pts)):
# 计算无畸变投影
P = np.append(world_pts[i], 1)
uv = M @ P
u_ideal, v_ideal = uv[0]/uv[2], uv[1]/uv[2]
# 应用畸变模型
u_obs, v_obs = image_pts[i]
r2 = u_obs**2 + v_obs**2
u_dist = u_ideal + u_ideal*(k[0]*r2 + k[1]*r2**2) + 2*k[2]*u_ideal*v_ideal + k[3]*(r2 + 2*u_ideal**2)
v_dist = v_ideal + v_ideal*(k[0]*r2 + k[1]*r2**2) + k[2]*(r2 + 2*v_ideal**2) + 2*k[3]*u_ideal*v_ideal
# 计算重投影误差
errors.append(u_obs - u_dist)
errors.append(v_obs - v_dist)
return errors
params_init = np.hstack([M_init.flatten(), k_opt])
res_joint = least_squares(full_reprojection_error, params_init, args=(world_pts, image_pts))
params_opt = res_joint.x
M_opt = params_opt[:12].reshape(3,4)
k_opt_joint = params_opt[12:]
print("修正后的DLT矩阵:\n", M_opt)
print("畸变系数(k1,k2,p1,p2):", k_opt_joint)
```
---
### **关键注意事项**
1. **畸变中心假设**：通常假设畸变中心与主点重合，若需优化中心坐标，需扩展参数[^4]。
2. **数值稳定性**：高次径向畸变 ($k_2, k_3$) 易导致数值不稳定，建议用正则化约束。
3. **初始值依赖**：联合优化需良好初始值（步骤2和3的结果）。
4. **畸变方向**：代码中实现的是**正向畸变模型**（从理想点到畸变点），但实际常用**反向映射**（畸变校正）[^4]。
---
### **相关问题**
1. 如何从单应矩阵分解相机内参？
2. 鱼眼相机模型与普通畸变模型有何区别？
3. 在已知畸变参数后，如何实时校正摄像头图像？
4. 张正友标定法如何集成畸变校正？
[^1]: 相机标定之DLT算法学习
[^2]: 相机标定(原理、过程、代码实现、畸变矫正)
[^3]: 3D视觉相机模型与标定详解
[^4]: 相机标定中的坐标变换，内外参求解，畸变校正</think>### 相机存在透镜畸变时修正DLT系数的方法
当相机存在透镜畸变时，DLT系数（投影矩阵 $M$）计算需增加畸变校正环节。核心思路是：**先通过畸变模型校正图像坐标，再用校正后的坐标求解DLT系数**。以下是具体步骤：
---
#### **1. 畸变模型概述**
透镜畸变主要包括两类[^2][^4]：
1. **径向畸变**（Radial Distortion）
$$ \begin{cases}
x_d = x(1 + k_1 r^2 + k_2 r^4 + k_3 r^6) \\
y_d = y(1 + k_1 r^2 + k_2 r^4 + k_3 r^6)
\end{cases} \quad \text{其中} \quad r^2 = x^2 + y^2 $$
- $k_1, k_2, k_3$：径向畸变系数
- $(x,y)$：理想坐标（无畸变）
- $(x_d,y_d)$：实际观测坐标
2. **切向畸变**（Tangential Distortion）
$$ \begin{cases}
x_d = x + [2p_1 xy + p_2(r^2 + 2x^2)] \\
y_d = y + [p_1(r^2 + 2y^2) + 2p_2 xy]
\end{cases} $$
- $p_1, p_2$：切向畸变系数
---
#### **2. 修正DLT系数的流程**
##### **步骤1：求解初始DLT系数（忽略畸变）**
1. 使用含畸变的图像点直接计算初始投影矩阵 $M_0$
$$ \min_{M} \sum_i \left\| \begin{bmatrix} u_{d_i} \\ v_{d_i} \end{bmatrix} - \text{proj}(M, P_i) \right\|^2 $$
- $P_i$：物方点 $(X_i, Y_i, Z_i)$
- $(u_{d_i}, v_{d_i})$：实际图像点（含畸变）
##### **步骤2：估计畸变系数**
1. 用 $M_0$ 计算理想坐标 $(u_i, v_i)$：
$$ \begin{bmatrix} u_i \\ v_i \\ 1 \end{bmatrix} \simeq M_0 \begin{bmatrix} X_i \\ Y_i \\ Z_i \\ 1 \end{bmatrix} $$
2. 建立畸变校正方程：
$$ \begin{bmatrix} u_{d_i} - u_i \\ v_{d_i} - v_i \end{bmatrix} = \begin{bmatrix}
u_i r_i^2 & u_i r_i^4 & 2u_i v_i & r_i^2 + 2u_i^2 \\
v_i r_i^2 & v_i r_i^4 & r_i^2 + 2v_i^2 & 2u_i v_i
\end{bmatrix} \begin{bmatrix} k_1 \\ k_2 \\ p_1 \\ p_2 \end{bmatrix} $$
- $r_i^2 = u_i^2 + v_i^2$
3. 用最小二乘法求解畸变系数 $\mathbf{k} = [k_1, k_2, p_1, p_2]^T$
##### **步骤3：校正图像坐标**
1. 应用反向畸变模型校正坐标：
$$ \begin{bmatrix} u_i' \\ v_i' \end{bmatrix} = \mathcal{D}^{-1} \left( \begin{bmatrix} u_{d_i} \\ v_{d_i} \end{bmatrix}, \mathbf{k} \right) $$
- $\mathcal{D}^{-1}$：畸变校正函数（需迭代求解[^4]）
##### **步骤4：重新计算DLT系数**
1. 用校正后的坐标 $(u_i', v_i')$ 求解精确投影矩阵 $M$：
$$ \min_{M} \sum_i \left\| \begin{bmatrix} u_i' \\ v_i' \end{bmatrix} - \text{proj}(M, P_i) \right\|^2 $$
---
#### **3. 联合优化方法（推荐）**
将DLT系数和畸变参数同时优化：
```python
import numpy as np
from scipy.optimize import least_squares
def reprojection_error(params, world_pts, image_pts):
# 解包参数: M(12维) + 畸变系数(4维)
M = params[:12].reshape(3,4)
k = params[12:16]
errors = []
for i in range(len(world_pts)):
# 物方点转齐次坐标
P = np.append(world_pts[i], 1)
# 投影到理想图像平面
uv = M @ P
u_ideal, v_ideal = uv[0]/uv[2], uv[1]/uv[2]
# 应用畸变模型
r2 = u_ideal**2 + v_ideal**2
u_dist = u_ideal*(1 + k[0]*r2 + k[1]*r2**2) + 2*k[2]*u_ideal*v_ideal + k[3]*(r2 + 2*u_ideal**2)
v_dist = v_ideal*(1 + k[0]*r2 + k[1]*r2**2) + k[2]*(r2 + 2*v_ideal**2) + 2*k[3]*u_ideal*v_ideal
# 计算与观测点的误差
errors.append(image_pts[i,0] - u_dist)
errors.append(image_pts[i,1] - v_dist)
return np.array(errors)
# 初始化参数 (M0 + 零畸变)
params_init = np.hstack([M_init.flatten(), [0,0,0,0]])
# 非线性优化
res = least_squares(reprojection_error, params_init, args=(world_pts, image_pts))
params_opt = res.x
# 提取最终结果
M_opt = params_opt[:12].reshape(3,4)
k_opt = params_opt[12:16]
```
---
#### **4. 关键注意事项**
1. **畸变中心**：通常假设主点 $(c_x, c_y)$ 是畸变中心，需提前标定[^4]
2. **数值稳定性**：
- 径向畸变系数 $k_1, k_2$ 需正则化约束
- 切向畸变 $p_1, p_2$ 量级通常较小 ($10^{-4}$)
3. **初始值依赖**：
- 步骤1的 $M_0$ 应使用RANSAC抵抗异常点
- 畸变系数初始值设为0
4. **高次项处理**：
- 普通镜头：保留 $k_1, k_2$ 即可
- 鱼眼镜头：需 $k_3$ 甚至更高阶项
> **实验建议**：实际标定时优先使用张正友法（棋盘格标定），其内置了畸变校正的鲁棒实现[^4]。
---
### **相关问题**
1. 如何从单幅图像估计透镜畸变参数？
2. 鱼眼相机模型与普通畸变模型有何区别？
3. 在线标定（实时畸变校正）如何实现？
4. 深度学习能否替代传统畸变校正方法？
[^1]: 相机标定之DLT算法学习
[^2]: 相机标定(原理、过程、代码实现、畸变矫正)
[^3]: 3D视觉相机模型与标定详解
[^4]: 相机标定中的坐标变换，内外参求解，畸变校正
本文介绍了透镜的径向和切向畸变原理，以及如何用OpenCV进行摄像机标定，包括内参矩阵和畸变向量的计算，强调了通过多角度观察结构获取相机参数的重要性。
扫一扫
3658
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
