采用整体最小二乘法的条纹图配准方法

2016-12-12穆治亚樊孝贺何丁龙韩冬松于国栋

中国光学 2016年6期

关键词：刻划十字条纹

穆治亚，艾华，樊孝贺，何昕，何丁龙，韩冬松，于国栋

(1.中国科学院长春光学精密机械与物理研究所，吉林长春 130033；2.中国人民解放军防化研究院，北京 102205;3.白城兵器试验中心，吉林白城 137001)

采用整体最小二乘法的条纹图配准方法

穆治亚1*，艾华1，樊孝贺1，何昕1，何丁龙1，韩冬松2，于国栋3

(1.中国科学院长春光学精密机械与物理研究所，吉林长春 130033；2.中国人民解放军防化研究院，北京 102205;3.白城兵器试验中心，吉林白城 137001)

针对共光路菲索型动态干涉仪采集到的4幅条纹图的空间一致性问题，本文提出一种新的条纹图配准思路，将条纹图与图像配准在关联度上人为分开，搭建了专用的图像配准装置，避免将条纹与十字丝标志点混在一幅图像内而影响干涉仪的测量精度。首先通过本文搭建的条纹图配准装置对4台CMOS相机进行物理配准，然后利用整体最小二乘法对采集到的同一十字丝刻划板图像进行十字丝提取、交点计算以及旋转量计算，实现共光路菲索型动态干涉仪条纹图像的点点对应。最后通过试验对比验证，证明了本文算法的配准精度优于模板重心法的配准结果，互相关度达到96%以上。

共光路动态干涉仪；条纹图配准；整体最小二乘法

1 引言

共光路斐索型动态干涉仪[1]光机结构装调完毕后，由于在安装时不能保证空间相对位置绝对的精确，存在微量的旋转、上下平移、左右平移等问题；并且由于相机的摆放位置不准确，还有可能出现采集到的条纹图存在镜像、倒像等情况；另外在成像过程中不同的CMOS感光特性也有细微的差别；共光路动态干涉仪在进行面形重构时，相位信息需要从同步采集的4幅条纹图中获得，4幅条纹图中相同位置的像素点所对应的必须是光学面形上同一处细小区域反射光的干涉光强[2]。解决以上问题的经典方法是进行图像配准，通过对图像进行旋转、平移等操作，使4幅不同相位调制量的条纹图中，相同位置的像素点对应同一空间位置。由于条纹图中的特征信息较少，单纯对采集到的条纹图进行配准，难度较大，并且配准精度很低。因此需要提出一种合适的配准方法来实现4幅条纹图的像素点点对应。

目前能够解决空间一致性问题的方法主要有两种：一种是基于条纹图的图像配准技术，通过对条纹图的特征点进行提取，然后通过配准条纹图像来实现4幅条纹图的点点位置精确对应，此种方法的缺点是，由于采集到的条纹图像几乎没有相同的特征点，特征点提取精度不高导致图像的配准精度很低；第二种是通过CMOS相机的标定技术，通过特殊的配准装置，使4台CMOS相机同时采集含有十字丝的刻划板，通过肉眼检测的方式来确定十字丝的位置重合标定，此方法属于人工检测装调技术，具有很大的操作误差，配准精度依然不高。

南京理工大学左芬[3]、徐晨[4]等人提出两种位置匹配方法：(1)采用在条纹图中加入一个标志点来对4幅条纹进行配准，此方法的缺点是在条纹中人为加入标志点，会更改条纹像素值，最后会影响面型计算结果，并且将条纹与标志点混在一起也加大了提取标志点的难度。(2)利用计算条纹图区域重心法来对4幅条纹图像进行位置一致性确定，此方法的缺点是有效条纹区域的确定本身就有误差，在有误差的区域内再次进行重心法计算，会继续增加误差，使得中心点计算结果更加不精确，并且此方法无法标定条纹图像的旋转量。内蒙古工业大学田枫等人[5]提出一种基于统计分析方法的圆形域同步移相干涉图位置配准技术，其精度优势在于条纹图轮廓的精确划分，而通常条纹图对比度不高，影响条纹图轮廓提取，最终影响精度计算。

本文将综合以上方法的优缺点，提出一种新的配准思路，将条纹图与位置配准在关联度上人为分开，搭建了基于十字丝的专用图像配准装置，避免将条纹与标志点混在一幅图像内，然后利用最小二乘法对十字丝的两条直线进行拟合，配准得到CMOS相机的统一坐标系。此种方法的配准精度高，可以实现4幅条纹图的精确配准。

2 基于整体最小二乘法的十字丝配准方法

2.1 基于刻划板十字丝图像的物理配准

由于相移分光产生的4个不同相位的干涉条纹由4个CMOS相机分别接收，因此在图像处理时要确保处理的图像像素要在严格的相同位置上，这就需要对CMOS相机进行位置配准，为此本文搭建了基于刻划板十字丝图像的物理配准装置。

本文搭建的干涉条纹配准结构图如图1所示。调整过程如下：在镜头的目视位置安装十字刻划板，在镜头前端安装法线与镜头光轴平行的平面反射镜，通过均匀光源照亮刻划板后并通过镜头放大成像在4个CMOS相机上，使其中一个像在接近靶面中心位置后固定，调整其他3个CMOS相机上的像与第一个像位置重合，即保证了4个CMOS相机靶面位置空间位置的配准。在配准完成的CMOS相机进行空间位置固定之后，CMOS相机依然会出现不同影响的位置偏移及旋转，此时需要通过图像配准技术对十字丝图像进行精确配准。

图1 干涉条纹配准结构三维图 Fig.1 3D picture of inference fringe registration structure

利用4个CMOS相机同时采集十字刻划板图像，采集到的原始图像如图2所示，然后通过算法计算十字丝交点坐标的偏移量以及十字丝的旋转量来实现四幅干涉条纹图的位置配准。

图2 调整前刻划板采集图 Fig.2 Collecting picture of depiction board before regulating depiction

2.2 基于sobel算子的图像预处理

由于原始图像十字丝特征不明显，首先我们需要对图2中的原始图像进行边缘处理，提高十字丝的对比度，常见的边缘提取有Prewitt算子[6]、拉普拉斯算子[7]、candy算子[8]等。由于这些算法局限性较高，不能保证提取目标的完整边缘，本文将利用Sobel算子对采集到的刻划板图像进行预处理，然后在其基础上利用整体最小二乘拟合对十字丝进行直线提取操作，最后联立直线方程，实现求出两直线的交点坐标。

Sobel算子本质上是一种梯度幅度检测算子[6]，其原理是：以像素点为中心的3×3窗口内进行灰度计算。然后根据该点的极值状态来进行边缘提取检测。令f(x,y)为像素点的灰度值，fx为水平方向上的梯度，fy为垂直方向上的梯度。

Sobel边缘检测算子定义为：

式中，fx=[f(i-1,x-1)+2f(i-1,j)+

f(i-1, j+1)]-[f(i+1,j-1)+

2f(i+1,j)+f(i+1, j+1)]

fy=[f(i-1,y-1)+2f(i,j-1)+

f(i+1, j-1)]-[f(i-1, j+1)+

2f(i, j+1)+f(i+1, j+1)].

Sobel边缘检测算子综合了图像每个像素点的上下左右邻点灰度的和，接近模板中心的权值。适当选取阈值T就可以实现图像边缘检测。

对0°相位的十字丝图像预处理结果如图3所示。由结果可以看出，预处理后的十字丝图像对比度高，十字丝更加明显，对于之后的直线拟合起到关键性作用。

图3 预处理结果及十字丝放大显示图 Fig.3 Result of pretreatment and the enlarged image of cross wire

2.3 基于整体最小二乘法的直线拟合

在图像预处理之后，需要对图3中的十字丝进行直线拟合，求取两直线交点的坐标，为此我们将求交点坐标转换为直线拟合问题。直线拟合问题可以描述为：假设给定n个坐标点，直线拟合的最终目的是通过计算得到一条最佳拟合直线，此直线能够尽量多的通过或者接近这些坐标点。其本质是求直线斜率和截距的最佳估计，最经典最常用的拟合方法是最小二乘法。普通最小二乘法[9]的缺陷在于自变量和因变量的选择对直线拟合的结果影响很大，造成拟合精度差异较大。经典最小二乘法要求误差仅来源于模型中的因变量，解释变量无误差或相比因变量误差而言其误差大小可以忽略不计，而在很多工程应用中，模型中解释变量有误差且是不可忽略的，因此20世纪80年代提出了整体最小二乘法[10](Total Least Squares)，其思想是建立在变量含误差EIV(Errors-In-Variables)模型上，能够很好地解决所有数据被随机误差影响的问题，其精度比经典最小二乘法更高。本文在最小二乘直线拟合方法的基础上，提出了基于整体最小二乘直线拟合方法，完成了十字丝两直线交点坐标的计算。

整体最小二乘法将自变量和因变量的误差加入到直线拟合方程中，其表达式为：

相应的误差方程可以按照EIV模型描述，如式(3)所示:

式中，设计矩阵A与观测向量l的误差分别由EA、El表示。

在设计矩阵A中，可以看到有一列元素的值恒为1，这是混合最小二乘法的求解问题[11]，因此，应按混合最小二乘求解。

令：A=[A1A2]，其中，A1=[1 1 … 1]T，A2=[x1x2…xm]T，此时构造一个增广矩阵C，并对其进行QR三角分解：

式中，Q为正交矩阵，R为上三角矩阵，则:

R=QTC=QT=

式中，R11、R12、R22、R1l、R2l为标量，可将上式写成如下的形式：

R22δa=R2l.

在对上述参数进行求解时，先利用整体最小二乘法对方程(7)进行求解以获取参数δa，把δa带入方程(6)，然后再按照普通的最小二乘法进行求解参数δb。

构造增广矩阵CR=[R22R2l]，并进行奇异值分解：

式中，U=[U1U2]，V=[N1N2]，Σ=diag(σ1,σ2)，σ1>σ2。

可得参数δa的整体最小二乘解如下式所示：

得到了参数δa，δb的解，需要对其残差和精度进行理论上的评价。

把增广矩阵CR写成下述形式：

R的改正量即为：

设计矩阵和观测向量的改正数为：

此时可得出残差：

在正交方向上的距离残差为：

距离残差在几何意义上就是最小的正交距离的残差平方和，所以，样本精度的评定应该按照正交距离的残差的平方和来计算:

由于整体最小二乘法同时考虑了两个坐标轴的误差，并且在两个坐标轴方向上的正交距离残差都比普通最小二乘法小，其单位权重误差比普通的最小二乘法更小。

3 实验结果与分析

本文搭建了共光路动态干涉仪的CMOS相机图像配准试验平台，平台选择KEYENCE公司的VH-Z50L长距离高性能变焦镜头，其观察距离在85 mm时能够达到500的放大率，该装置能够辅助配准干涉条纹。在缩短CMOS相机像面和PBS距离后，能够选取放大率更大的镜头。利用4个CMOS相机同时拍摄刻划板实物，刻划板上包含一个十字丝标志。表1为利用整体最小二乘法得到的配准结果，包含平移量和旋转量。图4为经过整体最小二乘配准后的十字丝刻划板图，图5为配准完成后利用配准得到的旋转平移公式变换后的4幅干涉条纹图。

图4 整体最小二乘法配准后刻划板图 Fig.4 Depiction board picture after registration using total least square

图5 按配准得到的公式变换后的条纹图 Fig.5 Fringe picture of transform using the registration result

图像序列十字丝中心点平移量直线斜率旋转量/(°)0相位(1030,1045)(6,21)0.01360.7792180相位(1016,1021)(-8,-3)0.00130.074590相位(1023,1025)(-1,1)0.00080.0458270相位(1019,1027)(-5,3)0.00110.0630

由图2可以看出，初始采集到的4幅十字丝图像具有明显的位置差异，利用本文提出的整体最小二乘法得到4个十字丝中心点位置坐标，并将四幅配准图的中心点平移到(1 024 ,1 024)坐标处，然后利用整体最小二乘法得到四幅配准图的旋转量，以0°相位配准图为标准，将其他3幅配准图转换到0°相位配准图的位置，表2计算出配准后的4幅条纹图两两互相关度结果表明，经过整体最小二乘方法进行配准后，两两条纹图之间的互相关度结果达到了0.97左右。结果表明，本算法达到了预期的处理目的。

表2 互相关度计算结果

图6 手动提取模板 Fig.6 Cyclostyle by manpower

图7 有效区域提取结果 Fig.7 Results collection of valid region

为了更好的验证本算法的优越性，利用重心法配准结果[12-15]与本文提出的配准算法进行比较。图6为重心法条纹图手动提取模板，图7为利用提取模板得到的有效条纹区域，图8为配准结果，由表2可以看出，重心法由于只考虑了条纹中心配准结果，未考虑图像的旋转量，其互相关度结果只在0.92以下，较本文算法有一定的差距。

图8 重心法配准结果 Fig.8 Registration result of barycenter method

4 结论

本文首先介绍了条纹图位置关系对检测结果造成的影响，在此基础上提出采用硬件装调与图像处理相结合的方法完成共光路动态干涉仪的条纹图配准工作。以硬件装调为主，使得对4幅条纹图进行的数字信号处理操作尽可能少，以保证其配准的同时减少信息丢失。在目前国内的光机加工水平下，在安装时，难以保证空间相对位置绝对的精确，仍然会存在微量的旋转、俯仰、横摆、上下平移、左右平移等问题，本文采用整体最小二乘法来对十字丝进行配准，提高条纹图的配准精度。最后通过一系列的试验验证了本文算法的正确性，并与重心配准法进行了性能比较，验证了本文算法的优越性。

[1] 王永红，冯家亚，王鑫，等.基于狭缝光阑的剪切散斑干涉动态测量[J].光学精密工程，2015,23(3):645-651. WANG Y H,FENG J Y,WANG X,etal.. Shearing speckle interferometry based on slit aperture for dynamic measurement[J].Opt.PrecisionEng.,2015,23(3)：645-651.(in Chinese)

[2] 王永红，李骏睿，孙建飞,等.散斑干涉相位条纹图的频域滤波处理[J].中国光学，2013,7(3):389-395. WANG Y H,LI J R,Sun J F,etal.. Frequency domain filtering for phase fringe patterns of digital speckle pattern interfeinmetry[J].ChineOptics,2013,7(3):389-395.(in Chinese)

[3] 左芬.同步移相干涉测量的抗振技术研究[D].南京:南京理工大学，2008. ZUO F. Study on the simultaneous phase-shifting anti-vibrational interferometry[D]. Nanjing:Nanjing University of Science and Technology,2008.(in Chinese)

[4] 徐晨.动态干涉测量技术与应用研究[D].南京:南京理工大学，2009. XU C. Study of dynamic interferometry technology and application[D]. Nanjing:Nanjing University of Science and Technology,2009.(in Chinese)

[5] 田枫、白福忠、吴亚琴,等.基于统计分析方法的同步移相干涉图位置配准[J].光学学报，2014，6:0626001.1-7. TIAN F,BAI F ZH,WU Y Q,etal. Position registration for simultaneous phase-shifting interferograms based on statistical analysis method[J].ActaOpticaSinica,2014，6:0626001.1-7.(in Chinese)

[6] KAUFMANN G H,GALIZZI G E. Phase measurement in temporal speckle pattern interferometry:comparison between the phase-shifting and the Fourier transform methods[J].AppliedOptics,2002,41(34):7254-7263.

[7] WANG F,PENG X S,LIU S Y. The method to filter the noise in the speckle fringe pattern based on spin filtering[J].StrongLaserandIonBeam,2009,21(z):230-234.

[8] WANG M X,FANG Y,HU H P. The independent component analysis method based on two-dimensional wavelet transform and its application in image separation[J].ElectronicsandInformationJ.,2006，28(3):471-475.

[9] 武迎春,曹益平,肖焱山.任意相移最小二乘法迭代的在线三维检测[J].光学精密工程，2014,22(5):1347-1353. WU Y CH,CAO Y P,XIAO Y SH. On-line three-dimensional inspection using randomly phase-shifting fringe based on least-square iteration[J].Opt.PrecisionEng.,2014,22(5):1347-1353.(in Chinese)

[10] 丁克良，沈云中，欧吉坤.整体最小二乘法直线拟合[J].辽宁工程技术大学学报(自然科学版)，2010,19(1):44-47. DING K L,SHENG Y ZH,OU J K. Methods of line-fitting based on total least-squares[J].J.LiaoningTechnicalUniversity(Natural Science),2010,19(1):44-47.(in Chinese)

[11] 张海华，刘春.顾及粗差的混合最小二乘平差实验分析[J].现代测绘，2010，5:8-12. ZHANG H H,LIU CH. Mixed total least squares and least squares adjustment experimental anal ysis with gross error considered[J].ModernSurveyingandMapping,2010，5:8-12.(in Chinese)

[12] 郭交，李真芳，保铮.基于相位梯度的干涉条纹图配准方法[J].西安电子科技大学学报(自然科学版)，2011，1:110-116. GUO J,LI ZH F,BAO ZH. Coregistration of multiple-pass interferogram of a single-baseline InSAR system based on phase gradients[J].J.XidianUniversity,2011，1:110-116.(in Chinese)

[13] 田枫.自参考干涉波前探测中光瞳形状问题研究[D].内蒙古工业大学，2014. TIAN F. Research on pupil shape problem in self-referencing interferometer wavefrom detection[D]. Inner Mongolia University of technology,2014.(in Chinese)

[14] WU S,HE X,YANG L. Enlarging the angle of view in Michelson interferometer based stereography by embedding a 4f system[J].Appl.Opt.,2011,50(21):1669-1705.

[15] 王新华，黄玮，欧阳继红.多探测器拼接成像系统实时图像配准[J].中国光学，2015,8(2)：211-219. WANG X H,HUANG W,OUYANG J H. Real-time image registration of the multi-detectors mosaic imaging system[J].ChineOptics,2015,8(2):211-219.(in Chinese)

《发光学报》

—EI核心期刊 (物理学类；无线电电子学、电信技术类)

《发光学报》是中国物理学会发光分会与中国科学院长春光学精密机械与物理研究所共同主办的中国物理学会发光分会的学术会刊。该刊是以发光学、凝聚态物质中的激发过程为专业方向的综合性学术刊物。

《发光学报》于1980年创刊，曾于1992年， 1996年， 2000年和2004年连续四次被《中文核心期刊要目总览》评为“物理学类核心期刊”，并于2000年同时被评为“无线电电子学、电信技术类核心期刊”。2000年获中国科学院优秀期刊二等奖。现已被《中国学术期刊(光盘版)》、《中国期刊网》和“万方数据资源系统”等列为源期刊。英国《科学文摘》(SA)自1999年；美国《化学文摘》(CA)和俄罗斯《文摘杂志》(AJ)自2000年；美国《剑桥科学文摘社网站》自2002年；日本《科技文献速报》(CBST， JICST)自2003年已定期收录检索该刊论文； 2008年被荷兰“Elsevier Bibliographic Databases”确定为源期刊; 2010年被美国“EI”确定为源期刊。2001年在国家科技部组织的“中国期刊方阵”的评定中，《发光学报》被评为“双效期刊”。2002年获中国科学院2001～2002年度科学出版基金“择重”资助。2004年被选入《中国知识资源总库·中国科技精品库》。本刊内容丰富、信息量大，主要反映本学科专业领域的科研和技术成就，及时报道国内外的学术动态，开展学术讨论和交流，为提高我国该学科的学术水平服务。

《发光学报》自2011年改为月刊， A4开本， 144页，国内外公开发行。国内定价： 40元，全年480元，全国各地邮局均可订阅。《发光学报》欢迎广大作者、读者广为利用，踊跃投稿。

地址：长春市东南湖大路3888号

《发光学报》编辑部

邮编： 130033

电话： (0431)86176862， 84613407

E-mail: fgxbt@126.com

国内统一刊号： CN 22-1116/O4

国际标准刊号： ISSN 1000-7032

国内邮发代号： 12-312

国外发行代号： 4863BM

http://www.fgxb.org

Inference fringe image registration using total least square method

MU Zhi-ya1*, AI Hua1, FAN Xiao-he1, HE Xin1, HE Ding-long1, HAN Dong-song2, YU Guo-dong3

(1.ChangchunInstituteofOptics,FineMechanicsandPhysics,ChineseAcademyofSciences,Changchun130033,China;2.ResearchInstituteofChemicalDefenseofPLA,Beijing102205,China;3.BaichengOrdnanceTestCenter,Baicheng147002,China)

To solve the space coherence of four inference fringe images by dynamic interferometer with common light paths, we propose a new notion of inference fringe image registration. We detach the conjunction between the inference fringe image and the image registration, and erect the equipment for inference fringe image registration, avoiding the confusion between the inference fringe and the cross wire which can influence the measure precision of dynamic interferometer. First, we realize the physical registration of four CMOS cameras using the equipment for inference fringe image registration. Then, we go along the image registration for the four cross wire images of the depict board by cross wire extraction, intersection point caculation and rotution amount caculation using total least square method, and realize the parallelism between the pix and pix of the dynamic interferometer with common light paths. Last, experimental results show that the proposed registration algorithm can improve the accuracy of registration, which is superior to that of the method of bary center. The method proposed in this paper can achieve the cross-correlation value of over 96%.

dynamic interferometer with common paths;inference fringe image registration;total least square method

2016-06-01；

2016-08-15

吉林省自然科学基金资助项目(No.201115124) Supported by Jilin Provincial Natural Science Foundation of China(No.201115124)

2095-1531(2016)06-0625-08

TP391.4

10.3788/CO.20160906.0625