激光与光电子学进展, 2021, 58 (8): 0810024, 网络出版: 2021-04-16  

任意步长两步相移法的鲁棒高精度相位解算方法 下载: 578次

Robust High-Precision Phase Solution Method Based on Two-Step Phase-Shifting Method with Arbitrary Step Length
作者单位
1 南京理工大学理学院, 江苏 南京 210094
2 东南大学土木工程学院, 江苏 南京 211189
3 山东科技大学计算机科学与工程学院, 山东 青岛 266590
摘要
两步相移法是均衡栅线投影技术的高速与高精度的重要方法。但是目前的求解算法的精度较低或算法复杂度较高。提出一种基于变量分组优化的两步相移求解方法。该方法将原始迭代变量分组为线性组变量和非线性组变量,针对确定的非线性组变量,通过最小二乘法获得线性组变量显式的最优解。通过对非线性组进行参数优化,获得全局的最优解。利用数值模拟与实验对本文方法进行了验证。结果表明本文方法有效地降低了多元非线性优化的相位值陷入局部最优的可能性,且降低了算法的复杂度。
Abstract
The two-step phase-shifting method is an important method to balance high speed and high precision in fringe projection profilometry. However, the current solution algorithm has low accuracy or high algorithm complexity. This paper presents a two-step phase-shifting solution method based on variable grouping optimization. In this method, the original iterative variables are divided into linear group variables and nonlinear group variables. For the determined nonlinear group variables, the explicit optimal solution of linear group variables can be obtained by the least square method. By optimizing the parameters of the nonlinear group, the global optimal solution is obtained. The method in this paper is verified by the numerical simulation and experiment. The results show that the method in this paper effectively reduces the possibility of falling into the local optimum of the phase value obtained by the multivariate nonlinear optimization and reduces the complexity of the algorithm.

1 引言

随着工业、医学等领域对于测量需求的不断增长,光学测量由于其高精度、非接触等优点得到飞速发展[1-4]。 其中栅线投影技术(FPP)具有测量精度高、测量视场大、测量速度快等优点[1-2,5],在三维形貌测量中得到了广泛的应用。

目前栅线投影中的相移解算方法主要包括傅里叶方法[6-7]和相移方法[5,8]。其中相移方法通过引入多张图像,简便地得到高精度的相位值。但是传统的相移方法至少需要3幅相移角已知的相移图像[8],为了强化其抗噪性,在实际测量中往往选取至少4步相移。过多的相移步数将会占据较大的时域空间,但理想的相移方法要求在单次测量时域内物体的形貌保持不变,故对于高速运动的物体,测量中会产生一定的相位求解误差[9]。为了进行高速动态测量,采用步数较少的相移方法是必要的,国内外学者对于两步相移法进行了一定研究。Zhang等[10]提出了一种基于两张相移图像和一张背景光图像的2+1相移算法;康新等[11]提出了一种基于局部窗口信息来消去背景光及折射项,并直接采用反正切函数进行相位求解的方法;李宝顺等[12]提出了一种基于三角波两步相移的方法;张晓璇等[13]提出了一种将傅里叶方法与相移方法相结合的相位解算方法;Mao等[14]提出了一种基于几何关系及迭代优化的相位解算方法;Yang等[15]提出了一种利用图像梯度信息进行变量消除,从而使得相位可解的方法;Zhang等[16]提出了一种基于希尔伯特变换的相位求解方法;Flores等[17]提出了一种基于Lissajous椭圆拟合的相位估计方法。

目前两步相移方法多以π[11-13]为相位角进行测量,但由于数字光处理(DLP)投影技术可以灵活地控制相移的步长,因此基于3π/2[14]、π/2[15]等不同步长的相移方法被提出。本文提出了一种适用于任意相移步长的相位求解方法[18-19]。通过优化参量的分组减少非线性优化的参量个数,使得求解过程更加稳定高效。相关模拟与实验证明了该方法的有效性。

2 基本理论

2.1 任意步长相移方法的优化目标正则化

对于任意步长的两步相移法,经过物体表面高度调制的栅线分布为

In=1(x,y)=R(x,y){A(x,y)+B(x,y)cos[2πfx+φ(x,y)]}In=2(x,y)=R(x,y){A(x,y)+B(x,y)cos[2πfx+φ(x,y)+ζ]},(1)

式中:In=1(x,y)、In=2(x,y)分别为第1、2张图像的像素点的灰度值;n为相移步数;(x,y)为像素坐标;φ(x,y)为相位值;B(x,y)为栅线幅度;R(x,y)为物体表面的折射率;A(x,y)为投影仪的基础光强与环境的背景光强之和; B(x,y)A(x,y)为栅线条纹的对比度;f为栅线的频率;ζ为相移的步长。

首先借助传统方法的思想[12],对物体表面的反射率项进行消除,即

In=1(x,y)In=2(x,y)=1+B(x,y)A(x,y)cos[2πfx+φ(x,y)]1+B(x,y)A(x,y)cos[2πfx+φ(x,y)+ζ](2)

为了简化书写,省略坐标y并且将对比度参数标识为新的形式,即

In=1(x)In=2(x)=1+BA(x)cos[2πfx+φ(x)]1+BA(x)cos[2πfx+φ(x)+ζ],(3)

式中:BA(x)为栅线对比度。

进一步引入两点假设,即当选取的图像窗口较小时,假设:1)窗口内像素对应的相位值为线性分布;2)窗口内像素对应的对比度为恒定值。

进一步省略x坐标,此时方程组共有3个参数:对比度、相位值、相位变化线性参数,满足

In=1(i)In=2(i)=1+BA(i)cos[φn+iΔφ]1+BA(i)cos[φn+ζ+iΔφ],(4)

式中:Δφ为与栅线频率f有关的线性变化参数,在周期性的栅线分布中可以认为Δφ=2πf;φn为窗口中心点处的相位值;i为窗口内以窗口中心为原点的相对坐标。

一般选择关于该像素点对称的窗口,故i=-W,…,-1,0,1,…,W,W为窗口半径。此时由于有3个参数,则W必须大于等于1才能进行求解。换言之,至少需引入三个像素点才能实现相位求解。此时相位求解可转化为非线性最小二乘(LS)极值问题,即

minBA,φn,Δφi=-WW{In=1(i)[1+BA(i)cos(φn+ζ+iΔφ)]-In=2(i)1+BA(i)cos(φn+iΔφ)}2(5)

注意到这个最优化目标函数是含有多个位置参量的形式,从数值角度来说其在优化过程中并不稳定,具有较多的局部最优解,这进一步加大了陷入错误相位解算值的可能性。故本文对上述目标函数进行以下推导:

In=1(i)1+BAcos(φn+ζ+iΔφ)-In=2(i)1+BAcos(φn+iΔφ)=In=1(i)1+BAcosφn+ζ2+iΔφ+ζ2-In=2(i)1+BAcos(φn+ζ2+iΔφ-ζ2)=In=1(i)-In=2(i)1+BAcosφn+ζ2+iΔφcosζ2-In=1(i)+In=2(i)BAsinφn+ζ2+iΔφsinζ2=I-(i)cosζ2BAcosφn+ζ2+iΔφ-I+(i)sinζ2BAsinφn+ζ2+iΔφ+I-(i)=I-(i)-I-(i)cosζ22+I+(i)sinζ22BAcosφn+iΔφ+ζ2-arctanI-(i)cosζ2I+(i)sinζ2=I-(i)-I*(i)BAcosφn+iΔφ+φ*(i),(6)

其中

I-(i)=In=1(i)-In=2(i)I+(i)=In=1(i)+In=2(i)I*(i)=I-(i)cosζ22+I+(i)sinζ22φ*(i)=ζ2-arctanI-(i)cosζ2I+(i)sinζ2(7)

由于ζIn=1(i)与In=2(i)全部为已知量,故I*(i)与φ*(i)也为已知量。综上所述,非线性最小二乘的目标函数可以表示为

minBA,φn,Δφi=-WW{I-(i)-I*(i)BAcosφn+iΔφ+φ*(i)}2(8)

这一形式比前文的目标函数更加简洁,更加有利于优化过程的建立。

2.2 变量局部线性化

在相移法中,可以将非线性方程进行线性化,从而利用线性的最小二乘法进行求解。注意到,(8)式所示的目标函数与相移方法的目标函数的形式类似,借鉴相移方法的思想,若给定Δφ的值,则有

I-(i)-I*(i)BAcosφn+iΔφ+φ*(i)=I-(i)-I*(i)BAcosiΔφ+φ*(i)cosφn+I*(i)BAsiniΔφ+φ*(i)sinφn=I-(i)-I*(i)cosiΔφ+φ*(i)Bc-I*(i)siniΔφ+φ*(i)Bs,(9)

式中:Bc=BAcos φn;Bs=-BAsin φn。由于I*(i)BAsin iΔφ+φ*(i)I*(i)cos[iΔφ+φ*(i)]均为已知量,则BAφn这两个非线性参量被成功地转化为线性的BsBc。那么对其的优化可以表示为

I*(i)BAsiniΔφ+φ*(-W)I*(i)BAcosiΔφ+φ*(-W)I*(i)BAsiniΔφ+φ*(-W+1)I*(i)BAcosiΔφ+φ*(-W+1)I*(i)BAsiniΔφ+φ*(W-1)I*(i)BAcosiΔφ+φ*(W+1)I*(i)BAsiniΔφ+φ*(W)I*(i)BAcosiΔφ+φ*(W)BsBc=I-(-W)I-(-W+1)I-(W-1)I-(W)(10)

将(10)式中三个矩阵依次记作NxI,则(10)式可以简写为

Nx=I(11)

根据线性最小二乘理论,其最优解为

x=(NTN)-1NTI(12)

最终该窗口中心点的相位值为

φn=arctan-BsBc(13)

2.3 参量分组优化

优化过程易陷入局部最优解的问题仍未解决。上述分别对优化方程进行正则化并对部分参量进行线性化。那么参与迭代的变量可以从3个简化为1个,即只需要对Δφ进行迭代。对于每一个给定Δφ,通过显式表达即可获得最优解,从而获得整个最优化系统中当前值与目标值的误差。优化流程如图1所示。

图 1. 变量分组优化基本流程

Fig. 1. Flow chart of variable grouping optimization

下载图片 查看所有图片

从局部最优分析可以看出,本文方法对于初值的选取并不敏感,选择贪心类算法对频率f进行非线性优化。由于T=1/f,故直接对周期T进行优化即可,其中初始步长L为1。按照贪心的原则,选择两侧最小的值,若两侧值的误差均比原始值的误差大,则步长缩减为原来的一半,即

{Tnew;Lnew}=Told-Lold;Lnew   Err(Told-Lold)<Err(Told+Lold)<Err(Told)Told+Lold;Lnew  Err(Told+Lold)<Err(Told-Lold)<Err(Told)Told;Lold/2     Err(Told)>Err(Told±Lold),(14)

式中:Tnew为更新后的周期值;Told为当前周期值;Lnew为更新后的步长;Err(·)为当前优化误差。

初值的选取可以采用傅里叶方法、包络方法等一系列算法,由于初值的不敏感性,本文对图像进行傅里叶变换,选取其一阶主频对应的周期作为初始值。

3 模拟与实验

首先通过数值方法验证本文迭代方法不会陷入过大的局部最优解。通过生成周期T=60的标准正弦曲线,为其添加均值为0、标准差为1的正态分布噪声。相移幅度过小时图像的差异性会被噪声所掩盖,且在求解过程中会引入奇异性问题,故本文生成了相移角度在 π/6,11π/6区间的不同相移图像。通过100次模拟求平均值,在不同相移步长下代入周期值得到的相位平均绝对值误差(MAE)如图2所示。可以明显看出,当代入周期值大幅度波动时,其仅有T=60这个最优解,并没有其他额外的局部最优解出现。由于优化函数为凸函数,这意味着本文所提出的方法将不具有初值敏感性,对于初值的选取将不再重要,所提方法可以通过已有方法得到一个大致的周期初值,甚至通过肉眼估计一个初始值以实现稳定收敛。

图 2. 相位误差与代入周期值的关系图

Fig. 2. Phase error varying with used period

下载图片 查看所有图片

表 2. 各相位求解算法的计算精度

Table 2. Calculation accuracy of different phase solution methods

PhasestepMethodPhase MAE /radMAE of our method /rad
π/2Method in Ref. [15]0.63060.0569
πMethod in Ref. [11]0.13690.0445
/2Method in Ref. [14]0.09410.0554

查看所有表

下面模拟存在物体表面反射率差异、背景光差异的真实测量环境。其中为了模拟反射率差异的随机性,采用Lena测试图对其进行模拟。光圈和曝光时间对于图像的增益也是线性的,这与折射率对于图像灰度值的增益效果一致,故总可以假设折射率的分布为[0,1]。注意此时的反射率是一种当然反射率,由材料决定,但并不反映材料的真实反射率。对于一般测量物体而言,加工精度或材料的不均质性导致的反射率差异往往是较小的(往往不会超过0.1),为了验证本文算法的强鲁棒性,将Lena图的灰度映射至[0.8,1.0]。同时,为了刻画环境光的影响,假设一个位于图像左下角处的点光源的光强增益等价分布为[-0.1,0.1]。注意该分布只是一种为了模拟光线分布差异性的等价形式,真实的环境光必然是大于0的。其中物体的表面高度采用二元高斯概率密度函数来计算。为了验证本文方法对于噪声的敏感性,添加了均值为0、方差为2个灰度等级的随机噪声。图像中每个像素点的相位可以表示为

表 1. 各相位求解算法的表现

Table 1. Performance of different phase solution methods

MethodPhase MAE /radConsuming timeProportion of LOS /%
(a)0.059961<10
(b)0.038584<10
(c)0.027161880.10244.7
(d)0.02213662.19800

查看所有表

φ(x,y)=2πfpx+z·dL-z,(15)

式中:fp为投影栅线的频率;z为被测物体表面的高度;L为投影仪光心至被测物体坐标系的距离;d为投影仪光心至相机光心的距离。在本文中,取fp=1/20,L=80,d=30。图3为两步相移法的模拟结果。

图 3. 两步相移法的模拟结果。(a)环境光强;(b)表面反射率;(c)高度分布;(d)栅线图像

Fig. 3. Simulation results of two-step phase-shifting method. (a) Environmental light intensity distribution; (b) surface reflectance; (c) height distribution; (d) fringe line image

下载图片 查看所有图片

下面以相移步长为3π/2为例验证本文所提出的算法的有效性。对比仅考虑理想光源方法、局部信息图解法、直接迭代优化方法[14]的求解结果,第256行的求解误差如图4所示,其全图整体精度以及求解时间等信息如表1所示。其中直接迭代优化方法对初值具有敏感性,本文首先通过文献[ 12]的方法获取初值,以降低其陷入过于偏激的局部最优解的可能性、加快优化,其中优化过程采用Matlab优化工具箱。可以明显看出,本文算法的精度远高于前两种方法,且略高于第三种算法,同时本文算法的计算耗时远少于第三种方法。此外,若收敛相位的绝对值相位误差大于1 rad,则认为算法在优化过程中陷入了局部最优解(LOS)。如表1所示,本文方法相比于其他方法可以有效避免此类现象的出现。

图 4. 不同相位求解方法的误差。(a)仅考虑理想光源的方法;(b)图解法;(c)直接迭代优化法;(d)本文算法

Fig. 4. Error of different phase solution methods. (a) Method only considering ideal light source; (b) graphical method; (c) method of direct iterative optimization; (d) method of this paper

下载图片 查看所有图片

为了验证本文算法对任意步长相移在实际测量中的有效性,采用步长为π/2,π,3π/2的相移方法进行验证,分别对比文献[ 11]、文献[ 15]、文献[ 14]中所提出的方法。选择人脸面具中的眼鼻区域作为被测物体。本文采用12步相移,选取12幅图像中的第1&7、1&4、1&10幅图像作为以上两步相移的测量图像,从而通过计算获得相位值。同时利用完整的12步相移法的测量结果作为每个像素点处的精确值。各相位求解算法的计算结果如表2图5所示。可以明显看出,文献[ 11]的方法含数值求导过程,其对于噪声是敏感的,噪声较大会导致严重的误差。且文献[ 11]与文献[ 15]中的方法均需要较大的局部窗口信息,物体高度变化较为明显时会降低测量分辨率。文献[ 14]中的方法虽然精度有所提升,但其时间复杂度极高,不利于求解。而本文方法在任何相移步长下均可以在相对较低的时间复杂度下进行高精度且鲁棒的相位求解。

图 5. 验证实验结果。(a)包裹相位的精确值分布;(b)相移步长为π/2时文献[ 15]中方法的求解值;(c)相移步长为π时文献[ 11]中方法的求解值;(d)相移步长为3π/2时文献[ 14]中方法的求解值;(e)~(g)相移步长分别为π/2,π,3π/2时本文方法的求解值;(h)图像第一行中相位精确值分布;(i)~(k)图像第一行相移步长分别为π/2,π,3π/2时的求解误差

Fig. 5. Results of verification experiment. (a) Accurate value distribution of wrapping phase; (b) solution value of method in Ref. [15] when phase shift step is π/2; (c) solution value of method in Ref. [11] when phase shift step is π; (d) solution value of method in Ref. [14] when phase shift step is 3π/2; (e)-(g) solution values of proposed method when phase shift steps are π/2,π, and 3π/2, respectively; (h) distribution of phase precise value in first row of image; (i)-(k) solution errors in first row of image when phase shift steps are π/2, π, and 3π/2, respectively

下载图片 查看所有图片

4 结论

提出了一种变量分组优化的两步相移相位求解方法,并将待优化参量分组为线性组与非线性组。当给定非线性值时,对于线性组的参量可以直接通过最小二乘法的显式获得,即获得系统的最小误差。故所提方法仅需要对非线性组的参量进行迭代优化即可,这极大减少了参与迭代的参数数量,从而提升了计算速度,且降低了迭代陷入局部最优解的可能性。实验结果表明,所提方法的精度高、计算速度快,进一步提升了两步相移法的可应用性。

参考文献

[1] Liu C, Chen L J, He X Y, et al. Coaxial projection profilometry based on speckle and fringe projection[J]. Optics Communications, 2015, 341: 228-236.

[2] Zuo C, Chen Q, Gu G H, et al. High-speed three-dimensional profilometry for multiple objects with complex shapes[J]. Optics Express, 2012, 20(17): 19493-19510.

[3] 邵新星, 陈振宁, 戴云彤, 等. 数字图像相关方法若干关键问题研究进展[J]. 实验力学, 2017, 32(3): 305-325.

    Shao X X, Chen Z N, Dai Y T, et al. Research progress of several key problems in digital image correlation method[J]. Journal of Experimental Mechanics, 2017, 32(3): 305-325.

[4] 吴文杰, 刘聪, 徐志洪, 等. 二值化数字散斑功率谱理论研究[J]. 光学学报, 2020, 40(3): 0312002.

    Wu W J, Liu C, Xu Z H, et al. Theoretical study on binary digital speckle power spectrum[J]. Acta Optica Sinica, 2020, 40(3): 0312002.

[5] Zhang S. Recent progresses on real-time 3D shape measurement using digital fringe projection techniques[J]. Optics and Lasers in Engineering, 2010, 48(2): 149-158.

[6] Chen W J, Su X Y, Cao Y, et al. Fourier transform profilometry based on a fringe pattern with two frequency components[J]. Optik - International Journal for Light and Electron Optics, 2008, 119(2): 57-62.

[7] Su X Y, Chen W J. Fourier transform profilometry: review[J]. Optics and Lasers in Engineering, 2001, 35(5): 263-284.

[8] Stoilov G, Dragostinov T. Phase-stepping interferometry: five-frame algorithm with an arbitrary step[J]. Optics and Lasers in Engineering, 1997, 28(1): 61-69.

[9] Feng S J, Chen Q, Gu G H, et al. Fringe pattern analysis using deep learning[J]. Advanced Photonics, 2019, 1(2): 025001.

[10] Zhang S, Yau S T. High-speed three-dimensional shape measurement system using a modified two-plus-one phase-shifting algorithm[J]. Optical Engineering, 2007, 46(11): 113603.

[11] 康新, 何小元. 两步相移实现投影栅相位测量轮廓术[J]. 光学学报, 2003, 23(1): 75-79.

    Kang X, He X Y. Two-step phase-shifting technique for phase measurement profilometry by grating projection[J]. Acta Optica Sinica, 2003, 23(1): 75-79.

[12] 李宝顺, 蔡青青, 包亚萍, 等. 采用相位半角的两步相移法[J]. 光子学报, 2014, 43(11): 1110002.

    Li B S, Cai Q Q, Bao Y P, et al. Two-step phase-shifting algorithm by the use of half angle of phase[J]. Acta Photonica Sinica, 2014, 43(11): 1110002.

[13] 张晓璇, 王月敏, 黄淑君, 等. 一种两步相移相位解算方法[J]. 光子学报, 2017, 46(3): 0311005.

    Zhang X X, Wang Y M, Huang S J, et al. A two-step phase-shifting algorithm for phase calculation[J]. Acta Photonica Sinica, 2017, 46(3): 0311005.

[14] Mao J Q, Yin Y K, Meng X F, et al. Two-step phase shifting in fringe projection: modeling and analysis[J]. Proceedings of SPIE, 2018, 1067: 106780V.

[15] Yang F J, He X Y. Two-step phase-shifting fringe projection profilometry: intensity derivative approach[J]. Applied Optics, 2007, 46(29): 7172-7178.

[16] Zhang H, Zhao H, Zha J, et al. Two-shot fringe pattern phase demodulation using the extreme value of interference with Hilbert-Huang per-filtering[J]. Proceedings of SPIE, 2019, 11056: 1105646.

[17] Flores V H, Rivera M. Robust two-step phase estimation using the Simplified Lissajous Ellipse Fitting method with Gabor Filters Bank preprocessing[J]. Optics Communications, 2020, 461: 125286.

[18] 牛文虎, 钟丽云, 孙鹏, 等. 一种改进的施密特正交化两步相移算法[J]. 中国激光, 2015, 42(6): 0608002.

    Niu W H, Zhong L Y, Sun P, et al. An improved two-step phase-shifting algorithm based on gram-Schmidt orthonormalization[J]. Chinese Journal of Lasers, 2015, 42(6): 0608002.

[19] 徐建亮, 周明安, 方晓汾. 一种两步相移任意步距相位轮廓测量技术研究[J]. 计量学报, 2016, 37(5): 472-475.

    Xu J L, Zhou M A, Fang X F. Study on an two phase-shifting profilometry with an arbitrary steps algorithm[J]. Acta Metrologica Sinica, 2016, 37(5): 472-475.

尹卓异, 刘聪, 赖立钊, 何小元, 刘晓鹏, 徐志洪. 任意步长两步相移法的鲁棒高精度相位解算方法[J]. 激光与光电子学进展, 2021, 58(8): 0810024. Zhuoyi Yin, Cong Liu, Lizhao Lai, Xiaoyuan He, Xiaopeng Liu, Zhihong Xu. Robust High-Precision Phase Solution Method Based on Two-Step Phase-Shifting Method with Arbitrary Step Length[J]. Laser & Optoelectronics Progress, 2021, 58(8): 0810024.

引用该论文: TXT   |   EndNote

相关论文

加载中...

关于本站 Cookie 的使用提示

中国光学期刊网使用基于 cookie 的技术来更好地为您提供各项服务,点击此处了解我们的隐私策略。 如您需继续使用本网站,请您授权我们使用本地 cookie 来保存部分信息。
全站搜索
您最值得信赖的光电行业旗舰网络服务平台!