CPD(Coherent Point Drift)相干点漂移算法 python 扩展库使用

一、CPD算法原理

1.点集配准

点集配准是将两个两个点集对其的空间变换过程。简单来说就是有两个点集，二维的三维的或者高维的点集。它们是对同一个物体或者广泛的说是对某种事物的呈现，但由于视角不同等原因，两个点集并不相同。就好像我们人站在不同的角度透视的观察一个物体，不同的角度看到的不一样，因为我们看到的特征之间的相对位置关系并不一样，但我们仍能认出这个物体来，是因为特征没有变。实际在应用角度来讲，就是有两个数据点集（二维的，三维的或者高维的），把其中一个作为参照，把另一个点集中的点的位置往对应的点的位置不断靠近。当然两个数据集的点的数量可能不一样。而CPD算法就是这样的一个过程，最终得到的输出是变换后的点集。

2.CPD算法

这个东东的原理我确实没看懂，概率论学得太差，不懂，不做解释。大家如果有需要的话可以去看一下论文（并不是很推荐，除非数学逆天，不然就是浪费时间）论文链接

二、CPD算法pyton软件包使用

1.pycpd安装

CPD算法本来是有matlab的软件包的。但确实也不好找，我之前找了半天也没找到。然后就想着看看python版本有没有，果然就被我找到了。这个给大家说一下如何找python扩展库，PyPI，直接搜索就可找到这个网站。里面有各种python的扩展库。这里放一下链接here.然后里面搜索CPD（cpd）就可以找到了。里面有一个 pycpd 的扩展库，直接pip install pycpd 就可以了。

2.pycpd使用

在PyPI网站里面给出了这个开源项目的github地址，里面有教程如何使用，以及全部的源代码。link

pycpd 里面有主要几个不同的类:
AffineRegistration(仿射）
DeformableRegistraion(可变性）
RigidRegistration(刚性变换）
分别对应不同类型的匹配。

然后主要有一个register函数，返回要得到的结果

def register(self,callback=lambda **kwargs: None):
	'''
	省略许多代码
	'''
	return self.TY, self.get_registration_parameters()

其中TY就是返回的变换点集的结果，是一个点集。
而get_registration_parameters()返回的则是完成这一变换的参数。不同的类有不同返回结果。
下面是RigidRegistration类get_registration_parameters（）源代码

'''
    R: numpy array (semi-positive definite)
        DxD rotation matrix. Any well behaved matrix will do,
        since the next estimate is a rotation matrix.

    t: numpy array
        1xD initial translation vector.

    s: float (positive)
        scaling parameter. 
'''
def get_registration_parameters(self):
        """
        Return the current estimate of the rigid transformation parameters.

        """
        return self.s, self.R, self.t

R是旋转矩阵，t是转换向量，s是缩放参数。
但我没理解算法，所以也不知道这些参数是干嘛的。

然后是一个简单的示例（github里面有数据，我就是用的他的数据）：

import matplotlib.pyplot as plt
import pycpd
import numpy as np

X = np.loadtxt('../data/bunny_source.txt')
Y = np.loadtxt('../data/bunny_target.txt')
fig = plt.figure()
ax = fig.gca(projection="3d")
ax.scatter(X[:,0],X[:,1],X[:,2])
ax.scatter(Y[:,0],Y[:,1],Y[:,2])
plt.show()

结果如下，实际上这个两个兔子，分别在坐标轴的不同位置。

下面分别绘制X，和Y

fig = plt.figure()
ax = fig.gca(projection="3d")
ax.scatter(x[:,0],x[:,1],x[:,2])

fig = plt.figure()
ax = fig.gca(projection="3d")
ax.scatter(Y[:,0],Y[:,1],Y[:,2])

可以看到X和Y的的形状一模一样，只是位置不同。

reg =pypd.RigidRegistration(**{'X': X, 'Y': Y})
data,(s,r,t) = reg.register()

以X点集作为基准，把Y点集往X点集变换，data为最终得到的结果。
下面绘制data点集

fig = plt.figure()
ax = fig.gca(projection="3d")
ax.scatter(data[:,0],data[:,1],data[:,2])

可以看到data点集描绘的形状也是兔子，而且在坐标轴上的位置几乎和X一模一样。
下面我们开看一下具体的数据

    print("X\n",X)
    print("data\n",data)

可以看到每个点都很接近。

    ax = fig.gca(projection="3d")
    ax.scatter(X[:,0],X[:,1],X[:,2])
    ax.scatter(data[:, 0],data[:, 1],data[:, 2])
    plt.show()

可以看到两个点集几乎重叠在一起了。

下面返回的几个参数（我也不懂）

data,(s,R,t) = reg.register(callback)
    print('s: ',s)
    print('R:\n',R)
    print('t:\n',t)

三、总结

大概这就是pycpd的基本使用方法，还有一点要说明的是，就是github的例子里面有一个可是化的函数，直观的显示了这个变化的过程，感觉很直观，大家可以去看一下。