概念
有关ICP问题中求解旋转矩阵和平移向量的知识在论文“Least-Squares Fitting of Two 3-D Point Sets”中有详细介绍。
ICP算法的思路就是:找到两组点云集合中距离最近的点对,根据估计的变换关系(
R
,
t
R,t
R,t)来计算距离最近点对经过变换之后的误差,进过不断的迭代直至误差小于某一阈值或者达到迭代次数来确定最终的变换关系。
- ICP算法能够使两个不同坐标系下的点集匹配到一个坐标系中,这个过程就是配准,配准的操作就是找到从坐标系1变换到坐标系2的刚性变换;
- ICP的本质就是配准,但有不同的配准方案,ICP算法本质是基于最小二乘的最优配准方法。该方法重复进行选择对应关系对,计算最优刚体变换,直到满足正确配准的收敛精度要求;
- ICP算法的目的就是找到待匹配点云数据与参考点云数据之间的旋转参数 R R R 和平移参数 t t t,使得两点数据之间满足某种度量准则下的最优匹配;
- 每次操作的都是目标点集,使目标点集不断靠近参考点集,因此求 R R R 和 t t t 也是每次考虑目标点集中每个点在参考点集中的最近点。
功能
通过给定的俩组空间点云坐标信息来解算俩个坐标系之间的位姿变换信息( R , t R,t R,t)。
步骤
两个点集 P 1 , P 2 P_1,P_2 P1,P2,每一步迭代,都朝着距离最小的目标进行。
- 筛选点对:由 P 1 P_1 P1中的点,在 P 2 P_2 P2中搜索出其最近的点,组成一个点对;找出两个点集中所有的点对。点对集合相当于进行有效计算的两个新点集。
- 跟据点集对,即两个新点集,计算两个质心。
- 由新点集,计算出下一步计算的旋转矩阵 R R R,和平移矩阵 t t t(其实来源于质心的差异)。
- 得到旋转矩阵和平移矩阵 R , t R,t R,t,就可以计算点集 P 1 P_1 P1进行刚体变换之后的新点集 P 2 P_2 P2,由计算 P 1 P_1 P1到 P 2 P_2 P2的距离平方和,以连续两次距离平方和之差绝对值,作为是否收敛的依据。若小于阈值,就收敛,停止迭代。
- 重复1-5,直到收敛或达到既定的迭代次数。
其中,计算旋转矩阵 R R R和平移向量 t t t时,需要矩阵方面的运算。
下面主要针对已知精确的空间点云集匹配关系的条件下,求解旋转矩阵
R
R
R和评议向量
t
t
t的过程。
已有精准匹配的点云集
p
i
,
p
i
′
(
i
=
1
⋯
n
)
p_i,p'_i(i=1\cdots n)
pi,pi′(i=1⋯n),俩者之间存在关系等式
p
i
′
=
R
p
i
+
t
+
N
i
p'_i=Rp_i+t+N_i
pi′=Rpi+t+Ni其中
R
R
R是旋转矩阵,
t
t
t是平移向量,
N
i
N_i
Ni是对应的噪声(为了简洁我们默认没有噪声)。
为了求解精确的旋转矩阵
R
R
R和平移向量
t
t
t,最小化投影误差平方和
∑
i
=
1
n
∣
∣
p
i
′
−
(
R
p
i
+
t
)
∣
∣
2
\sum^n_{i=1}||p'_i-(Rp_i+t)||^2
i=1∑n∣∣pi′−(Rpi+t)∣∣2获取来个点集的质心
p
=
1
n
∑
i
=
1
n
p
i
p
′
=
1
n
∑
i
=
1
n
p
i
′
p=\frac1 n\sum^n_{i=1}p_i\\p'=\frac1 n\sum^n_{i=1}p'_i
p=n1i=1∑npip′=n1i=1∑npi′进行去质心处理
q
i
=
p
i
−
p
q
i
′
=
p
i
′
−
p
′
q_i=p_i-p\\q'_i=p'_i-p'
qi=pi−pqi′=pi′−p′于是最小化投影误差信息
∑
i
=
1
n
∣
∣
p
i
′
−
(
R
p
i
+
t
)
∣
∣
2
=
∑
i
=
1
n
∣
∣
p
i
′
−
R
p
i
−
t
−
p
+
R
p
′
+
p
−
R
p
′
∣
∣
2
=
∑
i
=
1
n
∣
∣
(
p
i
′
−
p
′
−
R
(
p
i
−
p
)
)
+
(
p
′
−
R
p
−
t
)
∣
∣
2
=
∑
i
=
1
n
∣
∣
p
i
′
−
p
′
−
R
(
p
i
−
p
)
∣
∣
2
+
∣
∣
p
′
−
R
p
−
t
∣
∣
2
+
2
(
p
i
′
−
p
′
−
R
(
p
i
−
p
)
)
T
(
p
′
−
R
p
−
t
)
\sum^n_{i=1}||p'_i-(Rp_i+t)||^2\\=\sum^n_{i=1}||p'_i-Rp_i-t-p+Rp'+p-Rp'||^2\\=\sum^n_{i=1}||(p'_i-p'-R(p_i-p))+(p'-Rp-t)||^2\\=\sum^n_{i=1}||p'_i-p'-R(p_i-p)||^2+||p'-Rp-t||^2+2(p'_i-p'-R(p_i-p))^T(p'-Rp-t)
i=1∑n∣∣pi′−(Rpi+t)∣∣2=i=1∑n∣∣pi′−Rpi−t−p+Rp′+p−Rp′∣∣2=i=1∑n∣∣(pi′−p′−R(pi−p))+(p′−Rp−t)∣∣2=i=1∑n∣∣pi′−p′−R(pi−p)∣∣2+∣∣p′−Rp−t∣∣2+2(pi′−p′−R(pi−p))T(p′−Rp−t)值得注意的是
(
p
i
′
−
p
′
−
R
(
p
i
−
p
)
)
(p'_i-p'-R(p_i-p))
(pi′−p′−R(pi−p))部分在求和之后为0,因此上述等式为
∑
i
=
1
n
∣
∣
p
i
′
−
p
′
−
R
(
p
i
−
p
)
∣
∣
2
+
∣
∣
p
′
−
R
p
−
t
∣
∣
2
\sum^n_{i=1}||p'_i-p'-R(p_i-p)||^2+||p'-Rp-t||^2
i=1∑n∣∣pi′−p′−R(pi−p)∣∣2+∣∣p′−Rp−t∣∣2最小化该等式右侧项可以用来获取平移向量
t
=
p
′
−
R
p
t=p'-Rp
t=p′−Rp通过最小化
∑
i
=
1
n
∣
∣
p
i
′
−
p
′
−
R
(
p
i
−
p
)
∣
∣
2
=
∑
i
=
1
n
∣
∣
q
′
−
R
q
∣
∣
2
\sum^n_{i=1}||p'_i-p'-R(p_i-p)||^2\\=\sum^n_{i=1}||q'-Rq||^2
i=1∑n∣∣pi′−p′−R(pi−p)∣∣2=i=1∑n∣∣q′−Rq∣∣2来求解旋转矩阵信息
R
R
R
∑
i
=
1
n
∣
∣
q
′
−
R
q
∣
∣
2
=
∑
i
=
1
n
(
q
i
′
T
q
i
′
+
q
i
T
R
T
R
q
i
−
2
q
′
T
R
q
i
)
\sum^n_{i=1}||q'-Rq||^2\\=\sum^n_{i=1}(q'^T_iq'_i+q^T_iR^TRq_i-2q'^TRq_i)
i=1∑n∣∣q′−Rq∣∣2=i=1∑n(qi′Tqi′+qiTRTRqi−2q′TRqi)使
∑
i
=
1
n
q
i
′
T
R
q
i
\sum^n_{i=1}q'^T_iRq_i
i=1∑nqi′TRqi最大。
假设
W
=
∑
i
=
1
n
q
i
′
q
i
T
⟹
W
=
U
Σ
V
T
W=\sum^n_{i=1}q'_iq_i^T\\\implies W=U\Sigma V^T
W=i=1∑nqi′qiT⟹W=UΣVT为了能使得
∑
i
=
1
n
q
i
′
T
R
q
i
\sum^n_{i=1}q'^T_iRq_i
∑i=1nqi′TRqi最大,
R
=
U
V
T
R=UV^T
R=UVT使得
M
A
X
∑
i
=
1
n
q
i
′
T
R
q
i
=
U
Σ
U
T
MAX\sum^n_{i=1}q'^T_iRq_i=U\Sigma U^T
MAXi=1∑nqi′TRqi=UΣUT有关这一部分的证明在论文“Least-Squares Fitting of Two 3-D Point Sets”中
(
11
)
−
(
16
)
(11)-(16)
(11)−(16)有详细说明。
至此,我们获得了旋转矩阵
R
=
U
V
T
R=UV^T
R=UVT的求解公式以及平移向量
t
=
p
′
−
R
p
t=p'-Rp
t=p′−Rp的求解公式。