CN103983212A - 基于gpu并行计算小波变换的光学条纹图相位提取方法 - Google Patents

Abstract

一种用于结构光投影的三维面形测量的基于GPU并行计算小波变换的光学条纹图相位提取方法。所述方法使用以CPU作为主机端、GPU作为设备端的计算平台，主要包括以下计算步骤：CPU主机端接收光学条纹图像，根据图形尺寸生成所需小波变换的小波矩阵，并传输至GPU设备端；GPU设备端接收图形和小波矩阵，计算小波变换系数矩阵，并传回CPU主机端；CPU主机端接收小波变换矩阵，求解小波脊的包裹相位，解包裹后获取光学条纹图的连续相位。本方法在提高了数据处理速度的同时，保持了小波变换光学条纹图相位提取方法的精度。

Description

基于GPU并行计算小波变换的光学条纹图相位提取方法

技术领域

本发明涉及结构光投影三维面形测量，光学图像处理方法，尤其涉及一种用于结构光投影的三维面形测量的基于GPU(图形处理单元)并行计算小波变换的光学条纹图像相位提取方法。

背景技术

在面结构光投影三维面形测量技术中，CCD记录结构光被待测物体三维面形调制后的光学条纹图，对光学条纹图进行高精度相位提取，是实现实时三维面形测量的关键。光学条纹图的相位提取方法包括：常用的相移法(在先技术1：Shaoyan Gai,Feipeng Da.A novel phase shifting method based on strip marker[J].Opt.&Laser inEngng.2010.48(2):205-211)和傅立叶变换法(在先技术2：卢世江，李勇，王辉等，含孤立物体场景的高速、高密度三维面形采集[J].光电工程，2012,39(1):34-39)已经在动态物体的实时测量中得到了广泛的应用。小波变换光学条纹图相位提取方法(在先技术3：李思坤,苏显渝,陈文静.一种新的小波变换空间载频条纹相位重建方法[J].中国激光,2010,37(12):3060-3065)相对相移法而言，只需要一幅变形光学条纹图就能得到物体三维面形信息，且不需要精确的相移以及CCD同步，速度更快；相对傅立叶变换法而言，小波变换具备局部化、多分辨率分析能力和更强的噪声抑制能力。但是小波变换的计算涉及多个尺度因子下的行列循环和卷积运算，数据处理过程依循行列、多尺度的顺序逐次循环计算。多次循环计算过程中，每一列像素点相位的计算相互之间独立，各个尺度之间的计算也各不相关，体现出了很高的计算密度。计算速度无法满足实时、动态测量的需求。

发明内容

本发明的目的在于提供一种用于结构光投影的三维面形测量的基于GPU并行计算小波变换的光学条纹图相位提取方法，该方法保持了原有小波变换光学条纹图相位提取方法的精度的同时，提高了数据处理速度。

本发明的技术解决方案如下：

一种用于结构光投影的三维面形测量的基于GPU并行计算小波变换的光学条纹图相位提取方法，特征在于该方法使用CPU作为主机端和GPU作为设备端的计算平台对CCD获取光学条纹图进行处理，包括以下步骤：

1)在待测物体表面的斜上方设置投影仪，在待测物体表面的斜上方的同一水平面设置所述的投影仪和CCD相机，所述的CCD相机与CPU的输入端相连，该CPU与GPU相连，所述的CCD相机与待测物体所在的参考平面的距离为L₀，所述的投影系统与所述的CCD相机之间的距离为d，投影仪投影正弦结构条纹到被测物体表面，CCD获取光学条纹图，CPU读入CCD获取的光学条纹图，并根据光学条纹图分辨率以及小波变换中用到的尺度因子大小和数目n，生成各个尺度对应的小波矩阵，并将n个小波矩阵与图像传输至设备端GPU的内存；

2)所述的设备端GPU根据小波变换轮廓术逐个计算每个尺度因子下光学条纹图的小波变换系数矩阵，并将所有小波变换系数矩阵传回主机端CPU，计算过程如下：

光学条纹图中的结构光方向为图形的行方向，设为x方向，列方向设为y方向光学条纹图第y行的强度信息如公式(1)所示：

I'(x)_y＝I_1y+I_2ycos(2πfx+Δφ(x)_y)           (1)

其中，I_1y为背景光强度，I_2y为条纹的调制度，f为投影条纹的基频，Δφ(x)_y为由待测物体高度信息引起的调制相位；

采用Morlet复小波对第y行强度信息进行连续小波变换：

其中，W(a,b)_y为第y行小波变换系数，为母波函数Ψ(x)的子波函数，a为尺度因子，b为平移因子，是ψ_a,b(x)的复共轭函数，为卷积核，表示傅立叶变换，表示傅立叶逆变换；

小波变换系数矩阵W(a,b)是光学条纹图每一行强度信息小波变换之后所求的

小波系数合集： $W(a,b)=\left[\begin{array}{c}{W(a,b)}_1\\ \·\\ \·\\ \·\\ W{(a,b)}_y\\ \·\\ \·\\ \·\\ W(a,b{)}_M\end{array}\right],$ 式中M为光学条纹图的行数；

3)主机端CPU按下列式(3)从所述的n个小波系数矩阵W(a,b)求取小波脊矩阵，按式(4)从小波脊矩阵中求得包裹相位，解包裹之后得到连续的相位分布，再根据测量光路的几何参数按(5)式恢复待测物体形貌：

取小波变换系数矩阵第b列在所有尺度因子a方向模的极大值点，作为小波变换脊；

ridge(b)_y＝max[|W(a,b)_y|]             (3)

各元素的幅角即为所需包裹相位：

$\mathrm{\φ}\left({\mathrm{ridge}\left(b\right)}_y\right)=\arctan \left\{\frac{\mathrm{Im}\left[\mathrm{ridge}{\left(b\right)}_y\right]}{\mathrm{Re}\left[\mathrm{ridge}{\left(b\right)}_y\right]}\right\}---\left(4\right)$

在远心光路条件下，L₀？h(x,y)，被测物体高度分布与调制相位的关系为：

$h(x,y)=-\frac{L_0}{2\mathrm{\πfd}}\mathrm{\Δ\φ}(x,y)---\left(5\right).$

所述的生成小波矩阵如下：

1)尺度因子a对小波母函数进行了伸缩变换，小波函数长度随着尺度因子的取值不同而产生变化，若尺度因子为a时的小波函数长度为L_a，光学条纹图的尺寸为Mpixel×Npixel，其中M为y方向的光学条纹图尺寸，N为x方向光学条纹图尺寸，小波函数先在结构方向通过左右对称地赋0值扩展为长度为L_a+N-1数组，然后在y方向复制成大小为M×(L_a+N-1)的小波矩阵；

2)根据公式(2)小波系数的求解过程，包含了两次傅里叶变换加一次逆傅里叶变换运算，其中对于小波矩阵的傅里叶变换过程与测得的光学条纹图强度数据无关，在读取光学条纹图前在主机端CPU执行，将小波矩阵做傅里叶变换之后的作为输入传至设备端GPU，以减少设备端GPU运算量，减少解调过程中的运算时间。

与在先技术相比，本发明具有以下优点：

本发明采用GPU并行计算平台及小波变换的快速傅立叶变换算法，基于SIMD架构的图形处理器，包含多网格、线程块、线程的并行资源，对高计算密集度的图形处理有很高的计算加速性能。

本发明所述方法利用GPU多个网格和线程模块同时计算同一尺度下光学条纹图多列信息的小波变换，循环计算每个尺度因子下光学条纹图的小波变换系数，进而完成相位提取，提高了小波变换光学条纹图相位提取方法的速度。

加快了小波变换光学条纹图相位提取方法在三维面形测量中的相位提取速度，并保持了原有方法的计算精度。

附图说明

图1为基于GPU并行计算的小波变换面形测量光路及处理硬件原理图；

图2为基于GPU并行计算的小波变换光学条纹图相位提取方法流程图；

图3为本发明实施例中所用计算机生成光学条纹图；

图4为传统方法与本发明中基于GPU并行计算的小波变换光学条纹图相位提取结果。

具体实施方式

为使本发明的内容、实施过程和优点更加清楚，下面结合附图对本发明实施方式进一步详细描述，但不应以此限制本发明的保护范围。

实施例过程为，在包含CPU作为主机端和GPU作为设备端的计算平台上，对基于GPU并行计算小波变换光学条纹图相位提取方法的加速效果以及解调精度进行了测试。

实施例：

三维轮廓术中光学条纹图测量光路如图1所示，在待测物体表面的斜上方设置投影仪，在待测物体表面的斜上方的同一水平面设置所述的投影仪和CCD相机，所述的CCD相机与CPU的输入端相连，该CPU与GPU相连，所述的CCD相机与待测物体所在的参考平面的距离为L₀，所述的投影系统与所述的CCD相机之间的距离为d，投影仪投影正弦结构条纹到被测物体表面，CCD获取光学条纹图，CPU读入CCD获取的光学条纹图，并根据光学条纹图分辨率以及小波变换中用到的尺度因子大小和数目n，生成各个尺度对应的小波矩阵，并将n个小波矩阵与图像传输至设备端GPU的内存；投影系统投影正弦结构条纹到被测物体表面，CCD获取光学条纹图。主机端CPU(CPU型号为：Intel Xeon E5-2650)读取图像，如图3所示(所用光学条纹图采用计算机模拟的方式生成，变形光学条纹图为使用peak函数生成的相位对正弦函数进行调制后得到，图形尺寸为1024×1024，小波尺度因子共选取91个，小波母函数为Morlet复小波)、生成并传输小波矩阵、接受小波系数矩阵并完成相位求取和解包裹。小波变换对光学条纹图的计算过程包含了91个尺度的小波系数求取，每一个尺度包含了对1024行数据的计算过程。按照发明内容所述的GPU并行计算小波系数求解过程，基于SIMD(Single Instruction Multiple Data)架构的设备端GPU(所用GPU型号为：NVIDIA Tesla C2075)同时对1024行强度信息执行同一尺度小波系数的求取。

对于光学条纹图的第y(y＝1，2,3，…，1024)行长度为1024的强度信息I'(x)_y，当尺度因子a确定时，其小波系数求取过程中的卷积核固定不变。按照发明内容所述的小波矩阵 $G_a\left(b\right)=\left[\begin{array}{ccccc}0& ...& g_a{\left(b\right)}_1& ...& 0\\ 0& ...& g_a{\left(b\right)}_2& ...& 0\\ \·& \·& \·& \·& \·\\ \·& \·& \·& \·& \·\\ \·& \·& \·& \·& \·\\ 0& ...& g_a{\left(b\right)}_{1023}& ...& 0\\ 0& ...& g_a{\left(b\right)}_{1024}& ...& 0\end{array}\right],$

其中，矩阵宽度为L_a+1024-1，L_a为尺度因子为a时的子波函数长度。且其横向一维傅里叶变换可在CPU端读取光学条纹图之前优先计算。

将光学条纹图扩展为矩阵 $I_a\left(x\right)=\left[\begin{array}{ccccc}0& ...& I^{\′}{\left(x\right)}_1& ...& 0\\ 0& ...& I^{\′}{\left(x\right)}_2& ...& 0\\ \·& \·& \·& \·& \·\\ \·& \·& \·& \·& \·\\ \·& \·& \·& \·& \·\\ 0& ...& I^{\′}{\left(b\right)}_{1023}& ...& 0\\ 0& ...& I^{\′}(x{)}_{1024}& ...& 0\end{array}\right],$

其中，矩阵宽度为L_a+1024-1，L_a为尺度因子为a时的子波函数长度。

矩阵G'_a(b)、I_a(x)传输至设备端GPU后，GPU按照公式(2)所述计算尺度因子为a时的小波系数所得矩阵列数为L_a+1024-1，行数为1024，取每一行中间的1024项组成小波系数矩阵W'_a(b)，传回主机端CPU。主机端接收该尺度因子对应的小波系数矩阵后，将下一尺度因子对应的小波矩阵和光学条纹图扩展矩阵传给GPU，如此循环计算直到所有尺度因子上的小波系数矩阵计算完成。

当所有尺度因子上的小波系数计算完成并传输至CPU端后，主机端CPU接收了大小为91×1024×1024三维小波系数矩阵W(a,b)_y，根据公式(3)、(4)得到变形光学条纹图的小波脊及包裹相位；再对原始未调制的光学条纹图做上述操作，得到参考图包裹相位。两者求差值后得到相位差值，解包裹得到待测相位信息。

以上阐述了基于GPU并行计算的小波变换光学条纹图相位求解过程，为更加清晰的表明本发明所提出的计算方法在实例中所取得的加速效果以及对精度的影响，下面将对比传统仅使用CPU计算平台的串行计算方法进行说明。

一、GPU并行计算小波变换光学条纹图相位提取方法提速效果：

实施例中所用光学条纹图为计算机生成图如图3所示，图形像素尺寸为1024pixel×1024pixel，分别采用本发明基于GPU并行计算的小波变换光学条纹图相位提取方法和传统的小波变换光学条纹图相位提取方法对相位进行提取，当使用的尺度因子数目为91时；传统方法对变形图和原始图的计算时长约为100s～110s，本发明计算时长约为24s～28s，取得了近4倍的提速效果。

二、GPU并行计算小波变换光学条纹图相位提取方法的相位提取精度：

图4(b)给出了实施例中本发明方法所提取的相位，图4(a)为传统小波变换光学条纹图相位提取方法所提取的相位，附图4(c)为原始物高对应的相位，物高相位的极大值为16.2124rad。图4(d)为传统小波变换光学条纹图相位提取方法的相位提取误差，其误差极大值为0.4621rad；图4(e)为基于GPU并行计算小波变换光学条纹图相位提取方法的相位提取误差，其误差极大值为0.4539rad；图4(f)为传统小波变换光学条纹图相位提取方法提取的相位与基于GPU并行计算小波变换光学条纹图相位提取方法提取的相位之间的差值，其极大值为0.1869rad；两种方法提取的相位处在同一精度量级。

三、图形尺寸及尺度因子数目对加速效果的影响

表1给出了本发明方法对于不同尺寸的光学条纹图在多尺度因子下相位提取时间，并给出了相同光学条纹图使用传统仅使用CPU的小波变换光学条纹图相位提取方法的相位提取时间，表格最后一列给出了相应尺度因子数目下，各尺寸光学条纹图相位提取所取得的加速效果。可以看到，当尺度因子数目较小时，GPU有很高的加速效果(512pixel*512pixel图形的单尺度因子下加速效果达到31.26倍，相位提取时间仅需0.076s)。加速效果受图形尺寸影响较小，然而随着小波尺度因子数目的增加，加速效果逐渐降低并趋于稳定。

以上所述仅为本发明实施的实例，并不用于限制本发明的具体实施，凡在本发明的精神和原则之内所做的任何修改、等同替换、改进等均在本发明保护之列。

表1不同大小光学条纹图多尺度因子计算时间

Claims (2)

1.一种用于结构光投影的三维面形测量的基于GPU并行计算小波变换的光学条纹图相位提取方法，特征在于该方法使用CPU作为主机端和GPU作为设备端的计算平台对CCD获取光学条纹图进行处理，包括以下步骤：

1)在待测物体表面的斜上方设置投影仪，在待测物体表面的斜上方的同一水平面设置所述的投影仪和CCD相机，所述的CCD相机与CPU的输入端相连，该CPU与GPU相连，所述的CCD相机与待测物体所在的参考平面的距离为L₀，所述的投影系统与所述的CCD相机之间的距离为d，投影仪投影正弦结构条纹到被测物体表面，CCD获取光学条纹图，CPU读入CCD获取的光学条纹图，并根据光学条纹图分辨率以及小波变换中用到的尺度因子大小和数目n，生成各个尺度对应的小波矩阵，并将n个小波矩阵与图像传输至设备端GPU的内存；

2)所述的设备端GPU根据小波变换轮廓术逐个计算每个尺度因子下光学条纹图的小波变换系数矩阵，并将所有小波变换系数矩阵传回主机端CPU，计算过程如下：

光学条纹图中的结构光方向为图形的行方向，设为x方向，列方向设为y方向光学条纹图第y行的强度信息如公式(1)所示：

I'(x)_y＝I_1y+I_2ycos(2πfx+Δφ(x)_y)         (1)其中，I_1y为背景光强度，I_2y为条纹的调制度，f为投影条纹的基频，Δφ(x)_y为由待测物体高度信息引起的调制相位；

采用Morlet复小波对第y行强度信息进行连续小波变换：

其中，W(a,b)_y为第y行小波变换系数，为母波函数Ψ(x)的子波函数，a为尺度因子，b为平移因子，是ψ_a,b(x)的复共轭函数，为卷积核，表示傅立叶变换，表示傅立叶逆变换；

小波变换系数矩阵W(a,b)是光学条纹图每一行强度信息小波变换之后所求的小波系数合集， $W(a,b)=\left[\begin{array}{c}{W(a,b)}_1\\ \·\\ \·\\ \·\\ W{(a,b)}_y\\ \·\\ \·\\ \·\\ W(a,b{)}_M\end{array}\right],$ 式中M为光学条纹图的行数；

3)主机端CPU按下列式(3)从所述的n个小波系数矩阵W(a,b)求取小波脊矩阵，按式(4)从小波脊矩阵中求得包裹相位，解包裹之后得到连续的相位分布，再根据测量光路的几何参数按(5)式恢复待测物体形貌：

取小波变换系数矩阵第b列在所有尺度因子a方向模的极大值点，作为小波变换脊；

ridge(b)_y＝max[|W(a,b)_y|]          (3)

各元素的幅角即为所需包裹相位：

$\mathrm{\φ}\left({\mathrm{ridge}\left(b\right)}_y\right)=\arctan \left\{\frac{\mathrm{Im}\left[\mathrm{ridge}{\left(b\right)}_y\right]}{\mathrm{Re}\left[\mathrm{ridge}{\left(b\right)}_y\right]}\right\}---\left(4\right)$

在远心光路条件下，L₀？h(x,y)，被测物体高度分布与调制相位的关系为：

$h(x,y)=-\frac{L_0}{2\mathrm{\πfd}}\mathrm{\Δ\φ}(x,y)---\left(5\right).$

2.根据权利要求1所述的基于GPU并行计算小波变换的光学条纹图相位提取方法，其特征在于，所述的生成小波矩阵如下：

(1)尺度因子a对小波母函数进行了伸缩变换，小波函数长度随着尺度因子的取值不同而产生变化，若尺度因子为a时的小波函数长度为L_a，光学条纹图的尺寸为Mpixel×Npixel，其中M为y方向的光学条纹图尺寸，N为x方向光学条纹图尺寸，小波函数先在结构方向通过左右对称地赋0值扩展为长度为L_a+N-1数组，然后在y方向复制成大小为M×(L_a+N-1)的小波矩阵；

(2)根据公式(2)小波系数的求解过程，包含了两次傅里叶变换加一次逆傅里叶变换运算，其中对于小波矩阵的傅里叶变换过程与测得的光学条纹图强度数据无关，在读取光学条纹图前在主机端CPU执行，将小波矩阵做傅里叶变换之后的作为输入传至设备端GPU，以减少设备端GPU运算量，减少解调过程中的运算时间。