CN110532721B - 基于优化量子粒子群算法准确提取未知相移量的方法 - Google Patents

Abstract

本发明涉及一种基于优化量子粒子群算法准确提取未知相移量的方法，其首先基于衍射场的统计性质近似计算四步相移中各相移量的菲涅尔衍射场统计平均值；其次，根据相位随机条件求取中间值，用其表示近似相移量；最后，设置平均误差函数作为适应值函数，并改进量子粒子群算法以进一步逼近相移量的真实值，其中，量子粒子群算法中添加高斯变异操作，使算法符合粒子的初始分布规律。本发明构思合理，理论完善，不仅可以实现四步相移中未知相移量的近似预估计，提高计算方向的准确性，而且可以提高计算未知相移量的真实值的准确度和速度，有效克服传统迭代方法中出现在真实值附近的“锯齿解现象”，具有极快的速度和优秀的求解质量等优点。

Description

基于优化量子粒子群算法准确提取未知相移量的方法

技术领域

本发明涉及一种基于优化量子粒子群算法准确提取未知相移量的方法，属于相移干涉测量以及相位提取等领域，特别是适用于通过迭代算法来提取未知相移量的领域。

背景技术

光学三维测量技术是一种非接触、高精度投影结构光或自然光的表面微观形貌测量技术，这种技术已经广泛应用在电子、汽车、机械加工、纺织等现代工业中。相移干涉测量(Phase-shifting interferometry，PSI)是光学三维测量技术中的一种主动测量方式，其具备光学三维测量技术的优点，已广泛应用于从相位测量到显微镜系统的各种应用中。传统的PSI需要特殊且恒定的相移量(2π/N,N为大于等于3的整数)，例如多步定步长方法和多步等步长方法，这些方法要求每次的相移量都是一个固定值，而在实验过程中，由于受到相移器件标定误差、机械振动以及空气扰动等因素的影响，这个要求很难在实际中完全满足，导致波前重建及测量精度会受到影响。因此，为了克服这一不足，研究能够提取未知相移信息的广义相移干涉术(GPSI)对提高PSI的实际应用价值具有重要意义。

广义相移干涉术(GPSI)不需要严格控制每一步的相移量，而是从干涉图中提取任意未知相移值，再由相应的波前重建方法得到物光波，由于其相移值可以是任意未知的，所以GPSI可以减少对相移器件精度和环境稳定性的依赖，提高对系统误差的免疫能力。目前已有的相移提取算法可以分为迭代算法和非迭代算法两种类别，迭代算法根据迭代原理的不同可以分为交替式迭代算法和渐进式迭代算法，两种迭代方式都需要物光波相位分布或者误差函数和相移值之间的多次迭代，以逐渐逼近准确的相移值和物光波相位分布，其实验精度受到迭代步数的影响，而且运算时间较长，很难做到实时应用。非迭代算法根据提取相移值的原理的不同可以分为基于衍射场相位分布统计特性的相移值求解法、干涉图极值法、代数拟合法等，其中，基于衍射场相位分布统计特性的相移值求解法的基本依据是利用菲涅尔衍射场的振幅分布和相位分布相互独立且衍射场的相位分布近似随机变化这一特性，该方法会计算出相移量的近似值，实际上，由于各种实际因素，例如全息图大小及采样点个数的有限等因素的影响，求取的近似值存在较大的误差，需要进一步使用以干涉图的匹配为依据的迭代方法逐步逼近真实值，但是，相移量迭代步长的选择对计算速度有较大的影响：步长太小，需要增多迭代次数达到真实值；步长太大，又容易错过真实值，还需要减小步长反向寻找，会出现真实值附近的“锯齿解现象”。因此需要探究能够更好的契合基于衍射场相位分布统计特性的相移值求解法的新迭代方法。

发明内容

本发明的目的是克服传统相移法中每步相移是固定值的限制，结合智能算法——量子粒子群算法，设计一种与基于衍射场相位分布统计特性的相移值求解法结合更好的迭代方法。本发明提出基于优化量子粒子群算法准确提取未知相移量的方法，以基于衍射场相位分布统计特性的相移值求解法为算法基础，在计算出相移量的近似值后，使用优化量子粒子群算法进行迭代计算以求取相移量的真实值，其中，在量子粒子群算法中增加一次高斯变异操作，以增加初始粒子群的多样性，提高求解速度。使用从(3)式推导出的相移量的表达式(7)式恢复出记录平面P_H上的物波场，从而得到A_o和

的准确值。本发明实验操作简单方便，算法流程完善，求取相移量精度高，能更好的满足实际应用。

本发明提供的基于优化量子粒子群算法准确提取未知相移量的方法，具体步骤如下：

第1、依据实验光路图搭建光路，控制偏振片来引入任意相移量a1,a2,a3，用CMOS在时间序列上采集四帧全息图，并记录下物光和参考光的强度；假设物光波前经菲涅耳衍射后在记录平面P_H上的复振幅分布为

假设与记录平面P_H垂直的轴向平面波为参考波，其复振幅分布为

其中，A_o(x,y)和A_r(x,y)分别表示物光和参考光的振幅分布，可分别取物光光强I_o(x,y)和参考光光强I_r(x,y)的平方根，

和

分别表示物光和参考光的位相分布；记录平面P_H上的总光场为(为表达简便，省略(x,y))：

U＝O+R (1)

强度分布为：

则可以将四帧全息图分别表示为：

第2、对采集到的四帧全息图之间进行减法操作，并对得到的差值的绝对值进行整幅图取平均运算，从而形成三个中间值：

其中，<|...|＞是四帧全息图之间作差后取绝对值，然后作取平均运算。当系统满足物与记录器件的横向尺寸均显著小于物平面P_O与记录平面P_H之间的距离，即满足菲涅尔近似条件

时，物光波经菲涅尔衍射后，在记录平面P_H上，待测物体的相位

是随机分布的，或者说在[0，2π]内是均匀取值的，称该条件为相位随机条件。若满足该条件，不管相移量α取何值，均有近似条件

成立，即(4)式中的

因此(4)式中等号右边的<|...|＞项近似相等，记为中间变量c，即

则有：

p＝csin(α₁/2),q＝csin(α₂/2),r＝csin(α₃/2) (5)

当物光和参考光的强度被记录时，则c可以直接计算得出近似值，再经过简单的变换计算可得：

α₁＝2arcsin(p/c),α₂＝2arcsin(q/c),α₃＝2arcsin(r/c) (6)

即求得四步GPSI相邻步长之间的相移量的近似值α₁,α₂,α₃；

第3、将所求得的α₁,α₂,α₃作为初始值，使用优化量子粒子群算法进一步逼近它们的真实值：

首先，将(6)式计算得到的初始值α₁,α₂,α₃和测量得到的参考光强度I_r＝A_r ²，以及记录的四帧全息图，带入到(7)式中(由(3)式推导出)：

计算可得振幅A_o和相位

的初始近似值；

其次，构造一个使振幅A_o、相位

和相移量α₁,α₂,α₃相关联的平均误差函数ε作为量子粒子群算法的适应值函数：

其中，“∑”表示对整帧图像求和；

再次，根据p_ed(m)＝p_ed(m)±N(0,ξ₁·Δn_range)对α₁,α₂,α₃进行第一次高斯变异操作，优化初始粒子群的分布；

最后，按照优化的量子粒子群算法计算得出使平均误差函数ε最小的一组α₁,α₂,α₃，该组α₁,α₂,α₃则是相移过程中的真实值或最接近真实值的最优值；

第4、将第3中求得的α₁,α₂,α₃的最优值，测量得到的参考光强度I_r＝A_r ²，以及记录的四帧全息图再次带入(7)式(由(3)式推导出)，可以恢复出记录平面P_H上的物波场，从而得到振幅A_o和相位

的准确值。

第5、本发明中应用的优化量子粒子群算法原理如下：

粒子群算法(PSO)的基本思想来源于对鸟类群集行为的模拟，根据粒子群中单个粒子对目标环境的适应值函数的大小将单个粒子移动到其适应值较小的位置。PSO算法在优化过程中不会对单个粒子使用进化算子操作，而是将其中的每个粒子看作是D维可行性搜索空间中一个已经忽略了体积的粒子点，并且能够以一定的速度在可行性解空间中进行有效的搜索。PSO算法虽然具有很多优点，但是由于粒子的收敛是以轨道的形式存在的，并且粒子的收敛速度总是在有限的区域内。而在本发明应用的QPSO算法中，可以用量子空间中的一个粒子来描述粒子群中的单个粒子，从而能够体现粒子具有的不确定性。粒子不存在速度的限制，只有位置变化，粒子位置的更替由以下三个方程决定：

x_ed(m+1)＝p_ed(m)±β(m)|m_best(m)-x_ed(m)|ln[1/u_ed(m)] (11)

式中，u_ed，

介于(0,1)之间的随机数；m_best：所有个体当前最优位置的中心点；m：迭代次数；M：粒子个数；p_ed：粒子群的吸引点，介于p_ed和p_gd之间的随机位置；β：压缩扩张因子。

粒子的当前最优位置p_ed和粒子群全局最优位置p_ed的更新表达式为

在QPSO算法的迭代过程中，压缩-扩张因子β会影响到单个粒子的收敛行为，从而会对整个算法的收敛性产生有意义的影响。因此，基于单个粒子的层面上，一种自适应的参数控制方法被引入。为了实现这种压缩-扩张因子β的自适应赋值，评价函数F(m)被引入来评价单个粒子与粒子群最好位置的接近程度，表达式为

F(m)＝[f_e(m)-f_gbest(m)]/min{abs[f_e(m)],abs[f_gbest(m)]} (14)

式中，min{abs[f_e(m)],abs[f_gbest(m)]}：abs[f_e(m)]和abs[f_gbest(m)]中较小的值；f_gbest(m)：当前最好粒子的适应值；f_e(m)：第e个粒子的适应值；m：迭代次数。

在公式(14)中，F(m)的函数值越小，则表明该粒子越接近当前最优位置，同时该粒子在迭代过程中的搜索范围越小，容易出现算法的早熟现象；F(m)的函数值越大，则该粒子与当前最优位置距离越远，从而影响该粒子的收敛速度。因此，F(m)的函数值小，则给β设置较大的值，使算法搜索到更大的空间；F(m)的函数值大，则给β设置较小的值，使粒子具有一定程度的聚集性。由此，β可以看成是logF,e,m的函数，即为β(logF,e,m)。β的自适应参数调节方法为：对于当前迭代次数的粒子群，当logF的值从小到大排列时，β(logF,e,m)可以设置成从0.3-0.001的分段减小常数。

为了提高算法的收敛速率，将(6)式得到的α₁,α₂,α₃作为粒子群的parent-M，对parent-M进行第一次高斯变异操作，以引入其他M-1个粒子。

p_ed(m)＝p_ed(m)±N(0,ξi·Δn_range) (15)

式中，N(0,ξ_i·Δn_range)：均值为0，方差为ξ_i·Δn_range的正态分布随机数；Δn_range：初始α₁,α₂,α₃的调制范围[-0.005,0.005]；ξ_i：变异常数，取(0，1)的随机常数。第6、本发明中应用的优化量子粒子群算法具体实施步骤如下：

Step1.设置式(6)求得的α₁,α₂,α₃为初始值，根据式(15)对α₁,α₂,α₃进行高斯变异操作，求取M-1个粒子；

Step2.通过适应值函数式(8)计算每个粒子的适应值，并比较适应值的大小，粒子的初始位置即为它的最好初始位置p_ed(0)，p_ed(0)中适应值较小的粒子的位置即为粒子群的最好位置的初始值p_gd(0)；

Step3.按照式(9)、(10)分别计算粒子群的初始吸引点p_ed(0)和初始平均最优位置m_best(0)；

Step4.按照式(14)及β与logF,e,m之间的关系选择压缩-扩张因子β；

Step5.根据式(11)更新粒子的位置x_ed(m+1)；

Step6.计算每个粒子的适应值，按照式(12)、(13)分别更新粒子当前的最好位置为p_ed(m)，并且把p_ed(m)中最好的分配给粒子群当前的最好位置p_gd(m)；

Step7.按照式(15)对粒子的当前最优位置p_ed(m)进行高斯变异操作，并且计算变异前后粒子位置的适应值。如果变异操作以后这个粒子位置具有较小的适应值，那么更新它在当前迭代中的最优位置和适应值。假设

f[p_ed(m)+N(0,ξ₂·Δn_range)]＜f[p_ed(m)]，那么

p_ed(m)＝p_ed(m)±N(0,ξ₂·Δn_range)和f[p_ed(m)]＝f[p_ed(m)+N(0,ξ₂·Δn_range)]；

Step8.按照式(9)、(10)分别计算当前迭代粒子群的吸引点p_ed(m)和平均最优位置m_best(m)，按照式(11)计算粒子的当前位置x_ed(m+1)；

Step9.按照适应值函数公式(8)计算每个粒子当前位置的适应值，并与前一次迭代获得的目前粒子全局最优位置进行比较，如果粒子的当前位置较好，则更新粒子的当前最优位置及其适应值。即为，如果适应值函数f[x_ed(m+1)]＜f[p_ed(m)]，那么粒子当前迭代的最好位置p_ed(m+1)＝x_ed(m+1)，并且f[p_ed(m+1)]＝f[x_ed(m+1)]；

Step10.计算粒子目前迭代的最好位置p_ed(m)中适应值最好的粒子的位置分给目前全局最优位置p_gd(m+1)，即

并且f[p_ed(m+1)]＝f[p_ed(m)]；

Step11.比较当前迭代中粒子群全局最优的位置与其历史的全局最优位置的适应值大小，若当前粒子群的全局最优位置较好，则更新粒子群的全局最优位置及其适应值；

Step12.结束条件为迭代次数达到预设值M_max，或适应值ε≤0.0001，判断该算法是否满足结束条件，若条件已满足，则结束，否则返回Step4继续迭代。

本发明的优点及有益效果如下：

本发明将非迭代算法和智能算法有机结合，提出一种基于优化量子粒子群算法准确提取未知相移量的方法。本方法以基于衍射场相位分布统计特性的相移值求解法为算法基础，结合优化量子粒子群算法进行迭代计算以求取相移量的真实值，并且在量子粒子群算法中增加一次高斯变异操作，以增加初始粒子群的多样性，提高求解速度。本发明克服了传统相移法中每步相移都是固定值的限制，解决了传统的基于衍射场相位分布统计特性的相移值求解法在真实值附近出现的“锯齿解现象”。实验操作简单方便，算法流程完善，求取相移量精度高，能更好的满足实际应用，是一次将智能算法引入传统算法体系的很好尝试。

附图说明

图1为本发明模拟实验中的装置示意图，其中P_H表示记录平面，P_O表示物平面，BS表示分光棱镜；

图2为本发明模拟实验中人为设计的带有“光”字的振幅型物体图；

图3为本发明模拟实验中参考光的强度图；

图4为本发明模拟实验中四帧带有随机相移量的全息图；

图5为本发明模拟实验中以最优相移量重建的物光波前图。

具体实施方式

本部分通过MATLAB仿真例子对本发明进一步进行说明。

实施例：MATLAB仿真例子由于许多条件的限制，在实际应用中一个具体物体表面的准确表达式难以获得，故计算机模拟已成为验证一个算法的有效性的标准手段。本发明模拟实验中的装置示意图如图1所示。

第一、模拟实验中使用如图2所示的人为设计的带有“光”字的振幅型物体作为测试物体，波长为λ＝632.8nm、波数为k＝2π/λ的平面波正入射到测量物体上，携带测量物体的信息衍射至CMOS靶面。根据抽样定理，设置记录平面P_H到物平面P_O的距离为Z＝20cm，衍射面的尺寸为5mm，参考光是正入射的呈余弦变化的平面波，如图3所示，其表达式为R＝exp[j*k(x_ocosα+y_ocosβ)]，其中，x_o和y_o分别表示对衍射面进行x方向和y方向的均分，设置α＝π/2和β＝π/2.02则分别表示参考光与x轴和y轴的夹角。

第二、由物平面P_O出射的光经菲涅尔衍射至记录平面P_H，它与被BS反射后的到达记录平面P_H上的参考光干涉叠加形成干涉图。记录过程中，物光和参考光在记录平面P_H上进行干涉，得到全息图I₀，然后改变参考光的相位，进行三次随机相移并得到另外三幅全息图I₁，I₂，I₃(四帧带有随机相移量的全息图如图4所示)，模拟实验中相移量为作为参照数据，将三次随机相移量分别预设为α₁＝0.4472rad,α₂＝0.4106rad,α₃＝0.4879rad。

第三、四帧全息图之间进行减法操作，并对得到的差值的绝对值进行取平均运算，从而形成三个中间值：分别记为p、q、r。本实施例中的模拟实验的设计满足菲涅尔衍射近似条件

即满足相位随机条件，经过公式(5)、(6)对p、q、r进行近似变换操作，即求得四步GPSI相邻步长之间的相移量的近似值α₁＝0.4494rad,α₂＝0.4111rad,α₃＝0.4893rad；

第四、将第三求得的α₁,α₂,α₃的近似值作为初始值代入优化量子粒子群算法，经过step1-12的迭代计算，求解出相移量α₁,α₂,α₃的真实值或最接近真实值的最优值。

表1中列出了本实施例的模拟实验中四步GPSI的迭代结果，其中标注为(Preset)的第一行是预设的三个随机相移量α₁,α₂,α₃，标注为(Initial)的第二行是由基于衍射场相位分布统计特性的相移值求解法简化后的统计方法计算得到的α₁,α₂,α₃的近似值，使用这一近似值作为初始值代入优化量子粒子群算法求解的迭代过程则由第三行及下面的各行表示。从表1中可以得出α₁,α₂,α₃分别在第12次、第7次、第10次迭代时得到最优值，相移量在迭代过程中变化量很小，仅为0.0001×N，是因为本实施例的模拟实验中为消除传统算法在准确值附近存在的“锯齿解现象”，设置的迭代步长的调控范围仅为[-0.005,0.005]。

表1基于优化量子粒子群算法准确提取未知相移量的计算实例

第五、将求得的α₁,α₂,α₃的最优解，测量得到的参考光强度I_r＝A_r ²，以及四帧全息图平均光强带入公式(7)，可以恢复出记录平面P_H上的物光波前，如图5所示。由图5所示的物光波前图可以看出，物光的再现图与原图基本保持一致，算法可以实现准确提取未知相移量，由于本发明中量子粒子群算法这一智能算法与传统非迭代算法的结合，使得算法的可靠性得到较大提升。实施例中的模拟实验充分验证了算法的可行性和合理性。

本发明的技术方案的主要贡献具体如下：

(1)提出将智能算法(量子粒子群算法)与传统算法(基于衍射场相位分布统计特性的相移值求解法)相结合，准确提取未知相移量的算法新思路；

(2)简化基于衍射场相位分布统计特性的相移值求解法的算法步骤，使得算法省去了计算较为复杂的三角变换以及二次求取误差函数的过程；

(3)优化量子粒子群算法，在量子粒子群算法中添加一次高斯变异步骤，使得初始粒子群分布更符合未知相移量的真实分布特征；

(4)在量子粒子群算法中设置[-0.005,0.005]的迭代步长调控范围，消除了传统算法在准确值附近存在的“锯齿解现象”，算法的连贯性和可靠性大幅度提升；

(5)设计了基于优化量子粒子群算法准确提取未知相移量的方法的详细步骤和完整算法体系，并且通过模拟实验验证了算法的可行性和合理性。

注意，上述仅为本发明的较佳实施例及所运用技术原理。本领域技术人员会理解，本发明不限于这里所述的特定实施例，对本领域技术人员来说能够进行各种明显的变化、重新调整和替代而不会脱离本发明的保护范围。因此，虽然通过以上实施例对本发明进行了较为详细的说明，但是本发明不仅仅限于以上实施例，在不脱离本发明构思的情况下，还可以包括更多其他等效实施例，而本发明的范围由所附的权利要求范围决定。

Claims (1)

1.基于优化量子粒子群算法精确提取未知相移量的方法，其特征在于，包括以下步骤：

第1、依据实验光路图搭建光路，控制偏振片来引入任意相移量α₁,α₂,α₃，用CMOS在时间序列上采集四帧全息图，并记录下物光和参考光的强度；定义与光束传播方向垂直的横截面上像素点坐标为(x,y)，假设物光波前经菲涅耳衍射后在记录平面P_H上的复振幅分布为

假设与记录平面P_H垂直的轴向平面波为参考波，其复振幅分布为

其中，A_o(x,y)和A_r(x,y)分别表示物光和参考光的振幅分布，可分别取物光光强I_o(x,y)和参考光光强I_r(x,y)的平方根，

和

分别表示物光和参考光的位相分布；为表达简便，省略(x,y)，记录平面P_H上的总光场为：

U＝O+R (1)

强度分布为：

则可以将四帧全息图分别表示为：

其中，

第2、对采集到的四帧全息图之间进行减法操作，并对得到的差值的绝对值进行取平均运算，从而形成三个中间值：

其中，<|...|>是四帧全息图之间作差后取绝对值，然后作取平均运算，根据相位随机条件，不管相移量α取何值，均有近似条件

成立，即(4)式中的

因此(4)式中等号右边的<|...|>项近似相等，记为中间变量c，即

则有：

p＝csin(α₁/2),q＝csin(α₂/2),r＝csin(α₃/2) (5)

当物光和参考光的强度被记录时，则c可以直接计算得出近似值，再经过简单的变换计算可得：

α₁＝2arcsin(p/c),α₂＝2arcsin(q/c),α₃＝2arcsin(r/c) (6)

即求得四步GPSI相邻步长之间的相移量的近似值α₁,α₂,α₃；

第3、将所求得的α₁,α₂,α₃作为初始值，使用优化量子粒子群算法进一步逼近它们的真实值：

首先，将(6)式计算得到的初始值α₁,α₂,α₃和测量得到的参考光强度I_r＝A_r ²，以及记录的I₀、I₁、I₂和I₃四帧全息图，带入到由(3)式推导出的(7)式中：

计算可得振幅A_o和相位

的初始近似值；

其次，构造一个使振幅A_o、相位

和相移量α₁,α₂,α₃相关联的平均误差函数ε作为量子粒子群算法的适应值函数：

其中，“∑”表示对整帧图像求和；

再次，为了优化初始粒子群的分布，根据公式(9)对α₁,α₂,α₃进行一次高斯变异操作：

p_ed(m)＝p_ed(m)±N(0,ξ₁·Δn_range) (9)

其中：p_ed是量子粒子群算法中粒子的当前最优位置，m是量子粒子群算法迭代次数，N(0,ξ₁·Δn_range)是一个均值为0，方差为ξ₁·Δn_range的正态分布随机数，Δn_range是初始α₁,α₂,α₃的调制范围，设置为[-0.005,0.005]，ξ₁是设置的随机变异常数，设置为(0,1)范围内的随机常数；

最后，按照优化的量子粒子群算法计算得出使平均误差函数ε最小的一组α₁,α₂,α₃，该组α₁,α₂,α₃则是相移过程中的真实值或最接近真实值的最优值；

第4、将第3中求得的α₁,α₂,α₃的最优值，测量得到的参考光强度I_r＝A_r ²，以及记录的I₀、I₁、I₂和I₃四帧全息图再次带入由(3)式推导出的(7)式，可以恢复出记录平面P_H上的物波场，从而得到振幅A_o和相位

的准确值。

Images