CN115127683B - 一种动态干涉仪干涉图参数失配的相位提取的方法 - Google Patents

Abstract

本发明公开了一种动态干涉仪同步移相干涉图参数失配的相位提取方法，该方法为：首先采集参数失配的空间同步移相干涉图，分割为四幅子干涉图；然后使用条纹细化法计算各干涉图的对比度和载频误差，将对比度、载频误差和预设的移相间隔的初值代入线性方程组计算干涉图的相位、背景以及调制度；接着设定干涉图间的相位、背景、调制度正确，根据方程组计算的对比度、载频误差以及移相间隔的偏差值，更新初值。重复两个过程直至收敛判定函数优于阈值，停止迭代，得到最终的对比度、载频误差以及移相间隔，最后代入最小二乘方程组计算相位。本发明降低了动态干涉仪相位提取的约束性，干涉图参数失配相位提取的准确性和可靠性高，简单高效、适用范围广。

Description

一种动态干涉仪干涉图参数失配的相位提取的方法

技术领域

本发明涉及光干涉计量技术领域，特别是一种动态干涉仪同步移相干涉图参数失配的相位提取方法。

背景技术

空间移相干涉可以动态实时捕获多幅移相干涉图，提取相位空间分布的参数，用于测量动态变化的对象，在光学制造、天文、国防等领域具有广泛的应用。解算相位时，若干涉图之间的对比度、载频以及移相间隔不一致，会引入相位误差。为了保证相位的复原精度，需要对失配干涉图进行参数计算，然后再进行相位提取，从而保障干涉仪的测量精度以及准确度。

其中如AIA[Z.Wang and B.Han,“Advanced iterative algorithm for randomlyphase-shifted interferograms with intra-and inter-frame intensityvariations,”Opt.Lasers Eng.45(2),274–280(2007).]和QCA[J.Vargas andC.O.S.Sorzano,"Quadrature Component Analysis for interferometry,"Opt.LasersEng 51,637–641(2013)]法，可以快速准确地从包含随机平移相移且载频对比度一致的干涉图中提取相位。Deck[L.L.Deck,“Model-based phase shifting interferometry,”Appl.Opt.53(21),4628–4636(2014).]公开了一种基于多参数物理模型的MPSI方法，采用泰勒级数展开，通过线性回归得到相位，有效抑制了干涉图载频误差。该方法补偿了载频以及移相误差，将背景以及调制度作为时间无关量，但不能补偿对比度。Liu等人[Q.Liu,Y.Wang,J.He,and F.Ji,“Phase shift extraction and wavefront retrieval frominterferograms with background and contrast fluctuations,”J.Opt.17(2),025704(2015).]在干涉图参数失配提取相位中做出了大量工作，在MPSI基础上通过条件近似，将非线性最小二乘化简为线性形式计算，得到对比度修正后的相位提取结果，但是不能补偿条纹载频误差。后续他又提出通过三次近似线性的方程组，交替迭代补偿载频误差，收敛速度快，但只适用于高载频干涉图。除此之外，空间载频法[E.McDonnell and L.Deck,“Solutions for environmentally robust interferometric optical testing”inOptical Manufacturing and Testing XIII,R.Rascher and D.W.Kim,eds.(SPIE,2020),p.16.]可以实现单帧干涉图计算，克服了干涉图间干涉参数不一致的影响，但是参考光束与测试光束间的大量倾斜会引入非共光路的回程误差，降低相位提取精度。上述的方法都仅仅补偿三者中部分参数失配，对于干涉图完全失配的相位提取不适用。

发明内容

本发明的目的在于提供一种准确性和可靠性高、简单高效、适用范围广的动态干涉仪同步移相干涉图参数失配的相位提取方法。

实现本发明目的的技术解决方案为：一种动态干涉仪同步移相干涉图参数失配的相位提取方法，包括以下步骤：

步骤1、使用动态干涉仪采集空间同步移相干涉图，并分割为四个子干涉图；

步骤2、使用等效对比度计算公式计算出每个子干涉图中的等效对比度信息，对各个子干涉图进行条纹细化，并进行插值拟合得到载频量，根据载频量计算各个子干涉图间的载频误差；

步骤3、将计算的对比度、载频误差以及预设的移相间隔代入干涉图帧之间相关的方程组，计算干涉图表达式中的背景、调制度以及相位分布；

步骤4、设定计算背景、调制度以及相位分布为实际值，代入干涉图帧内相关的方程组，计算对比度、载频误差以及移相间隔小偏差值，对对比度、载频误差以及移相间隔进行更新；

步骤5、根据计算的背景、调制度、相位分布、对比度、载频误差以及移相间隔，计算收敛判定函数值，进行判定，如果不满足收敛条件则重复步骤3～步骤5，直至收敛判定函数满足收敛条件；

步骤6、根据收敛条件下得到的对比度、载频误差以及移相间隔，代入干涉图模型进行最小二乘法计算，确定计算相位。

本发明与现有技术相比，其显著优点在于：(1)使用条纹细化法拟合得到载频误差，更方便简单；(2)通过对对比度、载频误差、移相间隔进行解耦得到干涉图所有的失配参数，准确性和可靠性高；(3)方法简单高效，不需要任何附加的辅助硬件，适用于绝大多数的动态干涉仪。

附图说明

图1是本发明动态干涉仪同步移相干涉图参数失配的相位提取方法的流程图。

图2是实施例中动态干涉仪采集的空间同步移相干涉图。

图3是实施例中四幅子干涉图参数失配干涉图的相位解算结果图。

图4是实施例中点光源异位动态斐索干涉仪采集得到空间同步移相干涉图。

图5是实施例中参数失配干涉图的相位解算结果图。

具体实施方式

本发明动态干涉仪同步移相干涉图参数失配的相位提取方法，包括以下步骤：

步骤1、使用动态干涉仪采集空间同步移相干涉图，并分割为四个子干涉图；

步骤2、使用等效对比度计算公式计算出每个子干涉图中的等效对比度信息，对各个子干涉图进行条纹细化，并进行插值拟合得到载频量，根据载频量计算各个干涉图间的载频误差；

步骤3、将计算的对比度、载频误差以及预设的移相间隔代入干涉图帧之间相关的方程组，计算干涉图表达式中的背景、调制度以及相位分布；

步骤4、设定计算背景、调制度以及相位分布为实际值，代入干涉图帧内相关的方程组，计算对比度、载频误差以及移相间隔小偏差值，对对比度、载频误差以及移相间隔进行更新；

步骤5、根据计算的背景、调制度、相位分布、对比度、载频误差以及移相间隔，计算收敛判定函数值，进行判定，如果不满足收敛条件则重复步骤3、4、5，直至收敛判定函数满足收敛条件；

步骤6、根据收敛条件下得到的对比度、载频误差以及移相间隔，代入干涉图模型进行最小二乘法计算，确定计算相位。

进一步地，步骤1所述的采集空间同步移相干涉图，并分割为四个子干涉图，具体如下：

在动态干涉仪中采集得到空间同步移相干涉图，对其进行分割，得到四个子干涉图I₀、I₁、I₂和I₃。

进一步地，步骤2所述的使用等效对比度计算公式计算出每个子干涉图中的等效对比度信息，对各个子干涉图进行条纹细化，并进行插值拟合得到载频量，根据载频量计算各个干涉图间的载频误差，具体如下：

步骤2.1、分别对四个子干涉图I₀、I₁、I₂和I₃，根据等效对比度计算公式计算对应的等效对比度γ₀、γ₁、γ₂、γ₃；其中干涉图的等效对比度计算公式为：

式中γ表示干涉图的对比度，I(r)表示干涉图在空间位置矢量r下的灰度，表示干涉图的灰度均值，M×N表示参与对比度计算的干涉图尺寸大小，即像素数，Ω为干涉图有效区域；

步骤2.2、对各个子干涉图进行条纹细化处理，并进行插值拟合相位，根据倾斜估计载频量，根据载频量计算各个子干涉图间的载频误差。

进一步地，步骤3所述的将计算的对比度、载频误差以及预设的移相间隔代入干涉图帧之间相关的方程组，计算干涉图表达式中的背景、调制度以及相位分布，具体如下：

步骤3.1、预设移相间隔为(0°，90°，180°，270°)；

步骤3.2、将计算的对比度、载频误差以及预设的移相间隔代入四幅子干涉图模型方程，得到空间同步移相干涉图为：

式中A(r)为干涉图背景，γ_n为干涉图对比度，下标n为干涉图序号，表示待测相位，r为空间位置矢量，x空间位置坐标中的横坐标，y为纵坐标，a_n和b_n为相移项中两个正交方向上的倾斜误差因子，δ_n为干涉图间的移相量；

步骤3.3、定义：Δ_n(r,n)＝a_n·x+b_n·y+δ_n，式(2)可以转换为：

I(r,n)＝A(r)+C(r)γ_ncos(Δ_n(r,n))+S(r)γ_nsin(Δ_n(r,n)) (3)

由最小二乘准则建立方程组：

[A(r) C(r) S(r)]^T＝Q^-1B (4)

式中：

其中I_mes(r,θ_n)为采集干涉图光强。使用最小二乘法计算未知的背景、调制度A(r)；

步骤3.4、使用C(r)、S(r)计算相位

进一步地，步骤4所的述设定计算背景、调制度以及相位分布为实际值，代入干涉图帧内相关的方程组，计算对比度、载频误差以及移相间隔小偏差值，对对比度、载频误差以及移相间隔进行更新，具体如下：

步骤4.1从已知相位干涉图计算对比度、载频误差以及移相量：

将式(2)经过一阶泰勒展开为小偏差值的线性化形式：

I_k(r,n)＝I_n(r,n)+γ′_nA(r)cos(Δ(r,n))-(a′_nx+b′_ny+δ′_n)γ_nA(r)sin(Δ(r,n)) (6)

式中I_k为k次迭代后结果，γ′_n、a′_n、b′_n、δ′_n分别表示对比度、载频误差、移相间隔对应参数的小偏差值，式中

根据式(6)使最小二乘准则建立方程组，计算干涉图参数的小偏差值γ′(θ_n)、δ′(θ_n)、a′(θ_n)、b′(θ_n)：

由式(7)可以计算干涉图参数的小偏差值γ′(θ_n)、δ′(θ_n)、a′(θ_n)、b′(θ_n)；

步骤4.2、更新对比度、载频误差以及移相间隔的初值新估计：

式中(k)表示第k次迭代结果，μ_n为步长系数,n＝1、2、3、4，μ₁＝μ₂＝μ₃＝μ₄＝1。

进一步地，步骤5所述的根据计算的背景、调制度、相位分布、对比度、载频误差以及移相间隔计算收敛判定函数值，根据计算的背景、调制度、相位分布、对比度、载频误差以及移相间隔，计算收敛判定函数值，进行判定，如果不满足收敛条件则重复步骤3、4、5，直至收敛判定函数满足收敛条件，具体如下：

步骤5.1、由于式(6)是小偏差值的线性系统，因此定义评价函数为：

其中I_mes(x_i,n)为采集干涉图光强，M×N表示参与对比度计算的干涉图尺寸大小，即像素数；

对评价函数进行修改，得到该线性系统的收敛判定函数：

max{χ^(k-1)(n)-χ^(k)(n)}≤ε

式中max{}为取最大值函数，ε为收敛临界阈值，取值选取10^-6；

步骤5.2、每次迭代计算出新的对比度、载频误差以及移相间隔后，将计算数据代入收敛判定函数中，若不满足收敛临界阈值，则返回步骤3重新迭代，直至收敛判定函数满足收敛临界阈值。

进一步地，步骤6所述的根据收敛条件下得到的对比度、载频误差以及移相间隔，代入干涉图模型进行最小二乘法计算，确定计算相位，具体如下：

将计算出的四个子干涉图中的对比度、载频误差、移相间隔代入式(4)计算出C(r)、S(r)，然后根据式(5)得到最终同步移相干涉图的相位分布。

下面结合附图和具体实施例对本发明作进一步的详细说明。

实施例1

结合图1，本实施例提供了一种动态干涉仪同步移相干涉图参数失配的相位提取方法，包括以下步骤：

步骤1、使用动态干涉仪采集空间同步移相干涉图，并分割为四个子干涉图，具体如下：

在动态干涉仪中采集得到空间同步移相干涉图，对其进行分割，得到四个子干涉图I₀、I₁、I₂和I₃，如图2所示。

步骤2、使用等效对比度计算公式计算出每个子干涉图中的等效对比度信息，对各个子干涉图进行条纹细化，并进行插值拟合得到载频量，根据载频量计算各个子干涉图间的载频误差，具体如下：

步骤2.1、分别对四个子干涉图I₀、I₁、I₂和I₃，根据等效对比度计算公式计算对应的等效对比度γ₀、γ₁、γ₂、γ₃；其中干涉图的等效对比度计算公式为：

式中γ表示干涉图的对比度，I(r)表示干涉图在空间位置矢量r下的灰度，表示干涉图的灰度均值，M×N表示参与对比度计算的干涉图尺寸大小，即像素数，Ω为干涉图有效区域；

步骤2.2、对各个子干涉图进行条纹细化处理，并进行插值拟合相位，根据倾斜估计载频量，根据载频量计算各个子干涉图间的载频误差。

步骤3、将计算的对比度、载频误差以及预设的移相间隔代入干涉图帧之间相关的方程组，计算干涉图表达式中的背景、调制度以及相位分布，具体如下：

步骤3.1、空间移相干涉图采集的干涉图一般在90°移相位置附近，因此预设的移相间隔为(0°，90°，180°，270°)；

步骤3.2、将计算的对比度、载频误差以及预设的移相间隔代入四幅子干涉图模型方程，得到空间同步移相干涉图为：

式中A(r)为干涉图背景，γ_n为干涉图对比度，下标n为干涉图序号，表示待测相位，r为空间位置矢量，x空间位置坐标中的横坐标，y为纵坐标，a_n和b_n为相移项中两个正交方向上的倾斜误差因子，δ_n为干涉图间的移相量；

步骤3.3、在已知干涉图对比度、载频误差、移相步长后，可以进一步计算未知相位。在此采用经典的最小二乘方法来实现相位计算，即用适定和超定方程组求解未知相位。定义：Δ_n(r,n)＝a_n·x+b_n·y+δ_n，则式(2)可以写为：

I(r,n)＝A(r)+C(r)γ_ncos(Δ_n(r,n))+S(r)γ_nsin(Δ_n(r,n)) (3)

由最小二乘准则建立方程组：

[A(r) C(r) S(r)]^T＝Q^-1B (4)

式中:

其中I_mes(r,θ_n)为采集干涉图光强。使用最小二乘法计算未知的背景、调制度A(r)；

步骤3.4、由C(r)、S(r)计算相位分布

步骤4、设定计算背景、调制度以及相位分布为实际值，代入干涉图帧内相关的方程组，计算对比度、载频误差以及移相间隔小偏差值，对对比度、载频误差以及移相间隔进行更新，具体如下：

步骤4.1从已知相位干涉图计算对比度、载频误差以及移相量：

将式(2)经过一阶泰勒展开为小偏差值的线性化形式：

I_k(r,n)＝I_n(r,n)+γ′_nA(r)cos(Δ(r,n))-(a′_nx+b′_ny+δ′_n)γ_nA(r)sin(Δ(r,n)) (6)

式中I_k为k次迭代后结果，γ′_n、a′_n、b′_n、δ′_n分别表示对比度、载频误差、移相间隔对应参数的小偏差值，式中

根据式(6)使用最小二乘准则可以建立方程组，计算干涉图参数的小偏差值γ′(θ_n)、δ′(θ_n)、a′(θ_n)、b′(θ_n)。

由式(7)计算干涉图参数的小偏差值γ′(θ_n)、δ′(θ_n)、a′(θ_n)、b′(θ_n)；

步骤4.2、对比度、载频误差以及移相间隔的偏差值计算出来后，可以进一步更新对比度、载频误差以及移相间隔的初值新估计：

式中(k)表示第k次迭代结果，μ_n为步长系数,n＝1、2、3、4，μ₁＝μ₂＝μ₃＝μ₄＝1。

步骤5、根据计算的背景、调制度、相位分布、对比度、载频误差以及移相间隔，计算收敛判定函数值，进行判定，如果不满足收敛条件则重复步骤3、4、5，直至收敛判定函数满足收敛条件，具体如下：

步骤5.1、计算收敛判定函数值，由于式(6)是小偏差值的线性系统，因此可以定义评价函数：

其中I_mes(x_i,n)为采集干涉图光强，M×N表示参与对比度计算的干涉图尺寸大小，即像素数。为了进一步提高迭代的单调收敛性，对评价函数进行修改得到该线性系统的收敛判定函数：

max{χ^(k-1)(n)-χ^(k)(n)}≤ε

式中max{}为取最大值函数，ε为收敛临界阈值，取值选取10^-6；

步骤5.2、每次迭代计算出新的对比度、载频误差以及移相间隔后，将计算数据代入收敛判定函数中，若不满足收敛临界阈值，则返回步骤3重新迭代，直至收敛判定函数满足收敛临界阈值。

步骤6、根据收敛条件下得到的对比度、载频误差以及移相间隔，代入干涉图模型进行最小二乘法计算，确定计算相位，具体如下：

将计算出的四个子干涉图中的对比度、载频误差、移相间隔代入式(4)后计算出C(r)、S(r)，然后根据式(5)得到最终同步移相干涉图的相位分布，如图3所示。

本实施例利用基于点光源异位空间移相的动态斐索干涉仪，实现了参数失配空间同步移相干涉图的相位提取，该干涉仪利用一个2×2的扩展光源阵列实现了空间移相，得到了参数失配的空间同步移相干涉图如图4所示，得到最终干涉图相位分布如图5所示。

由图4、图5可以看出，本发明动态干涉仪干涉图参数失配的相位提取的方法，利用干涉图模型对失配空间同步移相干涉仪进行迭代参数求解，实现对比度、载频误差、移相间隔的准确计算，最后利用最小二乘法进行相位的准确计算；使用等效对比度计算公式，使用条纹细化法估计载频误差，提高初值估计的速度以及准确性；该方法简单高效，不需要任何附加的辅助硬件，适用于绝大多数的动态干涉仪。

Claims (2)

1.一种动态干涉仪同步移相干涉图参数失配的相位提取方法，其特征在于，包括以下步骤：

步骤1、使用动态干涉仪采集空间同步移相干涉图，并分割为四个子干涉图；

步骤2、使用等效对比度计算公式计算出每个子干涉图中的等效对比度信息，对各个子干涉图进行条纹细化，并进行插值拟合得到载频量，根据载频量计算各个子干涉图间的载频误差；

步骤3、将计算的对比度、载频误差以及预设的移相间隔代入干涉图帧之间相关的方程组，计算干涉图表达式中的背景、调制度以及相位分布；

步骤4、设定计算背景、调制度以及相位分布为实际值，代入干涉图帧内相关的方程组，计算对比度、载频误差以及移相间隔小偏差值，对对比度、载频误差以及移相间隔进行更新；

步骤5、根据计算的背景、调制度、相位分布、对比度、载频误差以及移相间隔，计算收敛判定函数值，进行判定，如果不满足收敛条件则重复步骤3～步骤5，直至收敛判定函数满足收敛条件；

步骤6、根据收敛条件下得到的对比度、载频误差以及移相间隔，代入干涉图模型进行最小二乘法计算，确定计算相位；

步骤1所述的使用动态干涉仪采集空间同步移相干涉图，并分割为四个子干涉图，具体如下：

在动态干涉仪中采集得到空间同步移相干涉图，对其进行分割，得到四个子干涉图I₀、I₁、I₂和I₃；

步骤2所述的使用等效对比度计算公式计算出每个子干涉图中的等效对比度信息，对各个子干涉图进行条纹细化，并进行插值拟合得到载频量，根据载频量计算各个子干涉图间的载频误差，具体如下：

步骤2.1、分别对四个子干涉图I₀、I₁、I₂和I₃，根据等效对比度计算公式计算对应的等效对比度γ₀、γ₁、γ₂、γ₃；其中干涉图的等效对比度计算公式为：

式中γ表示干涉图的对比度，I(r)表示干涉图在空间位置矢量r下的灰度，表示干涉图的灰度均值，M×N表示参与对比度计算的干涉图尺寸大小，即像素数，Ω为干涉图有效区域；

步骤2.2、对各个子干涉图进行条纹细化处理，并进行插值拟合相位，根据倾斜估计载频量，根据载频量计算各个子干涉图间的载频误差；

步骤3所述的将计算的对比度、载频误差以及预设的移相间隔代入干涉图帧之间相关的方程组，计算干涉图表达式中的背景、调制度以及相位分布，具体如下：

步骤3.1、预设移相间隔为(0°，90°，180°，270°)；

步骤3.2、将计算的对比度、载频误差以及预设的移相间隔代入四幅子干涉图模型方程，得到空间同步移相干涉图为：

式中A(r)为干涉图背景，γ_n为干涉图对比度，下标n为干涉图序号，表示待测相位，r为空间位置矢量，x空间位置坐标中的横坐标，y为纵坐标，a_n和b_n为相移项中两个正交方向上的倾斜误差因子，δ_n为干涉图间的移相量；

步骤3.3、定义Δ_n(r,n)＝a_n·x+b_n·y+δ_n，则式(2)转换为：

I(r,n)＝A(r)+C(r)γ_ncos(Δ_n(r,n))+S(r)γ_nsin(Δ_n(r,n)) (3)

由最小二乘准则建立方程组：

[A(r) C(r) S(r)]^T＝Q^-1B (4)

式中：

其中I_mes(r,θ_n)为采集干涉图光强；使用最小二乘法计算未知的背景、调制度；

步骤3.4、使用C(r)、S(r)计算相位分布

步骤4所述的设定计算背景、调制度以及相位分布为实际值，代入干涉图帧内相关的方程组，计算对比度、载频误差以及移相间隔小偏差值，对对比度、载频误差以及移相间隔进行更新，具体如下：

步骤4.1从已知相位干涉图计算对比度、载频误差以及移相量：

将式(2)经过一阶泰勒展开为小偏差值的线性化形式：

I_k(r,n)＝I_n(r,n)+γ′_nA(r)cos(Δ(r,n))-(a′_nx+b′_ny+δ′_n)γ_nA(r)sin(Δ(r,n)) (6)

式中I_k为k次迭代后结果，γ′_n、a′_n、b′_n、δ′_n分别表示对比度、载频误差、移相间隔对应参数的小偏差值，式中

根据式(6)使用最小二乘准则建立方程组，计算干涉图参数的小偏差值γ′(θ_n)、δ′(θ_n)、a′(θ_n)、b′(θ_n)：

式中m为干涉图中有效像素数；由式(7)计算干涉图参数的小偏差值γ′(θ_n)、δ′(θ_n)、a′(θ_n)、b′(θ_n)；

步骤4.2、更新对比度、载频误差以及移相间隔的初值新估计：

式中上标(k)表示第k次迭代结果，μ_n为步长系数，n＝1、2、3、4，μ₁＝μ₂＝μ₃＝μ₄＝1；

所述步骤5具体如下：

步骤5.1、由于式(6)是小偏差值的线性系统，因此定义评价函数为：

其中I_mes(x_i,n)为采集干涉图光强，M×N表示参与对比度计算的干涉图尺寸大小，即像素数；

对评价函数进行修改，得到该线性系统的收敛判定函数：

max{χ^(k-1)(n)-χ^(k)(n)}≤ε

式中max{}为取最大值函数，ε为收敛临界阈值，取值选取10^-6；

步骤5.2、每次迭代计算出新的对比度、载频误差以及移相间隔后，将计算数据代入收敛判定函数中，若不满足收敛临界阈值，则返回步骤3重新迭代，直至收敛判定函数满足收敛临界阈值。

2.根据权利要求1中所述的动态干涉仪同步移相干涉图参数失配的相位提取方法，其特征在于，步骤6所述的根据收敛条件下得到的对比度、载频误差以及移相间隔，代入干涉图模型进行最小二乘法计算，确定计算相位，具体如下：

将计算出四个子干涉图中的对比度、载频误差、移相间隔代入式(4)计算出C(r)、S(r)，然后根据式(5)得到最终同步移相干涉图的相位分布。