CN110500968B - 基于稀疏傅里叶变换的数字莫尔干涉相位实时测量方法 - Google Patents
基于稀疏傅里叶变换的数字莫尔干涉相位实时测量方法 Download PDFInfo
- Publication number
- CN110500968B CN110500968B CN201910782763.9A CN201910782763A CN110500968B CN 110500968 B CN110500968 B CN 110500968B CN 201910782763 A CN201910782763 A CN 201910782763A CN 110500968 B CN110500968 B CN 110500968B
- Authority
- CN
- China
- Prior art keywords
- phase
- spectrum
- moire
- dimensional
- frequency spectrum
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 60
- 238000001228 spectrum Methods 0.000 claims abstract description 106
- 230000008707 rearrangement Effects 0.000 claims abstract description 38
- 238000005259 measurement Methods 0.000 claims abstract description 27
- 230000009467 reduction Effects 0.000 claims abstract description 21
- 230000008569 process Effects 0.000 claims abstract description 20
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 15
- 238000005070 sampling Methods 0.000 claims abstract description 15
- 238000012545 processing Methods 0.000 claims abstract description 8
- 230000015572 biosynthetic process Effects 0.000 claims abstract description 7
- 230000010363 phase shift Effects 0.000 claims abstract description 7
- 238000003786 synthesis reaction Methods 0.000 claims abstract description 7
- 238000002156 mixing Methods 0.000 claims abstract description 3
- 230000003595 spectral effect Effects 0.000 claims description 36
- 238000001914 filtration Methods 0.000 claims description 11
- 230000008859 change Effects 0.000 claims description 9
- 238000000691 measurement method Methods 0.000 claims description 5
- 238000006073 displacement reaction Methods 0.000 claims description 3
- 238000002360 preparation method Methods 0.000 claims description 3
- 230000002441 reversible effect Effects 0.000 claims description 3
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 claims 1
- 230000003287 optical effect Effects 0.000 abstract description 10
- 238000004364 calculation method Methods 0.000 abstract description 7
- 238000005516 engineering process Methods 0.000 abstract description 4
- 230000009286 beneficial effect Effects 0.000 abstract description 3
- 230000000694 effects Effects 0.000 description 3
- 239000007788 liquid Substances 0.000 description 3
- 238000004088 simulation Methods 0.000 description 3
- 238000010586 diagram Methods 0.000 description 2
- 230000006872 improvement Effects 0.000 description 2
- 238000003754 machining Methods 0.000 description 2
- 238000011084 recovery Methods 0.000 description 2
- 230000000717 retained effect Effects 0.000 description 2
- 230000009466 transformation Effects 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 230000010365 information processing Effects 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000004441 surface measurement Methods 0.000 description 1
- 230000001052 transient effect Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
- G01B11/00—Measuring arrangements characterised by the use of optical techniques
- G01B11/24—Measuring arrangements characterised by the use of optical techniques for measuring contours or curvatures
- G01B11/2441—Measuring arrangements characterised by the use of optical techniques for measuring contours or curvatures using interferometry
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/10—Image enhancement or restoration using non-spatial domain filtering
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/20—Image enhancement or restoration using local operators
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20024—Filtering details
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20048—Transform domain processing
- G06T2207/20056—Discrete and fast Fourier transform, [DFT, FFT]
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Length Measuring Devices By Optical Means (AREA)
- Investigating Or Analysing Materials By Optical Means (AREA)
Abstract
本发明涉及基于稀疏傅里叶变换的数字莫尔移相干涉相位实时测量技术，属于光学测量领域。本发明首先构造虚拟干涉图，进行莫尔合成。以一维稀疏傅里叶变换(SFFT)为基础，设计二维的莫尔条纹图重排规律，通过空域的重排将频谱中少量的频谱幅值大值点稀疏地分散在频谱中。其次，设置二维的SFFT的窗函数，进行二维信号的空域混叠和降采样FFT，将原本的N×N大小的图像降至B×B大小，减少了后续计算量。再经过位置坐标还原和估值还原完成N×N大小的频谱还原。最后将还原频谱用于数字莫尔移相方法的解相过程。本发明可以降低解相过程的计算量，提高数字莫尔移相解相方法的算法速度，有助于实现精密加工、工程实践、生物医学领域相位的实时测量。
Description
技术领域
本发明涉及一种基于稀疏傅里叶变换的数字莫尔移相干涉相位实时测量方法，属于光学测量领域。
背景技术
许多物理量的测量可以转化为相位的测量。比如精密加工光学元件表面面形测量中，利用干涉法可以把被测元件面形变化转化为反射激光相位变化进行测量；气液流场、生物细胞形态测量中，透明气液及组织的厚度也可通过微分干涉等转化为透射光相位变化进行测量。在上述测量过程中，由于被测物理量往往随时间发生变化，有时变化速度较高，因此需要一种可以进行相位实时测量的技术，方能解决物理量瞬态分布的测量问题，为精密加工、工程实践、生物医学研究提供有效的测量手段。
干涉法的基本原理是利用干涉这种物理现象将相位空间分布转化为干涉图的强度空间分布，继而对记录的干涉图像进行处理，恢复相位。利用多幅移相干涉图进行解相的方法具有精度高、分辨率高的优点，但是系统包含移相机构，相对复杂，且需要采集多幅干涉图，难以适应实时测量的需求。单幅干涉图解相方法只需采集一幅干涉图，具有抗振性好、速度快的优势，能够满足日渐迫切的实时测量需求。公告号为CN 1587950的专利“一种用部分补偿透镜实现非球面面形的干涉测量方法”公开了一种用于非球面测量的数字莫尔移相干涉测量方法，这种方法利用数字莫尔移相取代实际的移相过程，是一种单幅干涉图相位测量方法。
数字莫尔移相方法用于测量一幅实际干涉图的相位，为了实现这一测量过程，首先构造k幅虚拟移相干涉图。具体来说，去除直流项的实际干涉图分布iR(x,y)和构造的虚拟干涉图分布iV(x,y)如式(1)所示。
其中bR(x,y)是实际干涉图调制量，pR(x,y)为实际干涉图待测相位，f是x方向添加的空间载波；bV(x,y)是虚拟干涉图调制量，pV(x,y)为构造的虚拟相位，δk为附加的移相相位值，δk＝1，2,3,4＝0,π/2,π,3π/2，对应四步移相算法。实际干涉图与虚拟干涉图进行莫尔合成，得到莫尔条纹分布iM(x,y)：
经过快速傅里叶变换得到莫尔条纹的频谱，在频谱中通过低通滤波，滤除bR(x,y)bV(x,y)cos[pR(x,y)+pV(x,y)+4πfx-δk]部分，只保留低频部分：
iLk(x,y)＝bR(x,y)bV(x,y)cos[pR(x,y)-pV(x,y)+δk] (3)
设Δ(x,y)＝pR(x,y)-pV(x,y)，为实际干涉图相位与虚拟干涉图相位之差。使用传统四步移相算法可以求得Δ(x,y)的分布：
通过pR(x,y)＝Δ(x,y)+pV(x,y)即可求得实际相位分布。
上述测量方法的相位求解过程需要经过四次傅里叶变换，而由于傅里叶变换的计算较为复杂，其运算速度往往是影响相位求解速度的重要因素。目前传统且常用的傅里叶变换算法是快速傅里叶变换(FFT)，这种傅里叶变换方法在信息处理、图像处理和频域分析中有广泛的应用。对于某些信号，信号的频谱仅在少量的频率处存在较大的幅值，其余频率处幅值近似为零。这种信号可以被称为频域稀疏的信号。但是FFT算法中，所有数据点都具有等价且等量的运算量，没有注意或利用到频谱本身可能存在的稀疏性。2012年，麻省理工的研究人员提出稀疏快速傅里叶变换(SFFT)，利用频谱的稀疏特性简化计算量，相较FFT算法提高了运算效率。对于数字莫尔移相解相过程，二维的稀疏傅里叶变换是一个有价值的加速干涉图相位求解速度的途径。
发明内容
本发明的目的是为了解决现有技术的干涉图相位解算时间长，无法实现实时测量的问题，提供一种基于稀疏傅里叶变换的数字莫尔移相干涉相位实时测量方法，该方法通过二维稀疏傅里叶变换，降低解相过程中傅里叶变换的计算量，提升运算效率，为实现相位实时测量提供技术支持。
本发明的目的是通过以下技术方案实现的。
基于稀疏傅里叶变换的数字莫尔移相干涉相位实时测量方法，以一维SFFT变换方法为基础，设计二维的空域干涉图重排规律，通过空域的重排将频谱中少量的频谱幅值大值点稀疏地分散在频谱中。设置二维的SFFT的窗函数，作用是在频域将被分散后的频谱大值点独立选取出来。经过窗函数处理后，进行二维信号的空域混叠降采样和FFT，将原本的N×N大小的图像降至B×B大小，N>B，降采样后的FFT过程减少了后续计算量。经过位置坐标还原和估值还原即可完成N×N大小的频谱重构，最后将重构频谱用于数字莫尔移相方法的解相过程。
基于稀疏傅里叶变换的数字莫尔移相干涉相位实时测量方法，包括如下步骤：
步骤1：使用数字莫尔移相解相方法完成解相准备
通过数字莫尔移相解相方法，按照公式(1)和公式(2)对实际干涉图与构造的虚拟干涉图进行莫尔合成，并得到4幅二维的莫尔条纹图iM(x,y)，每幅莫尔条纹图的大小为N×N。以下步骤2到步骤5对每一幅莫尔条纹图分别进行，不再赘述。
步骤2：设计二维SFFT的莫尔条纹图空域重新排列规律
取σx和σy分别为x,y方向的重排系数，对二维SFFT空域莫尔条纹图iM(x,y)进行重新排列得到重排信号i'(x',y')，x',y'为重排信号的坐标。
i′(x′,y′)＝iM(σxx,σyy) (5)
经过空域重新排列，原始莫尔条纹图像的K个频谱大值点被分散至整个频谱的各个位置。
步骤3：设置二维窗函数
理想的窗函数为一个空域的滤波函数，特点是在频域呈方形分布，方形内为通带，值为1，方形外为禁带，值为0。
以二维sinc函数Sc(x′,y′)为基础构建窗函数。对于sinc函数的频谱来说，除了中心的主峰外，其余位置处也存在较多的次级峰。为了消去次级峰，仅保留主峰作为窗函数的通带，需要在sinc函数上叠加一个高斯函数g(x′,y′)。将二维sinc函数和二维高斯函数配合使用得到二维的空域窗函数wf(x′,y′)。
其中b和c为系数，分别用于控制二维sinc函数的幅值和频谱宽度；d和μ为系数，分别用于控制高斯函数的幅值和频谱宽度。
使用空域窗函数wf(x′,y′)对重排信号i′(x′,y′)进行滤波处理得到加窗的重排图像信号j(x′,y′)。
j(x′,y′)＝wf(x′,y′)·i′(x′,y′) (7)
窗函数wf(x′,y′)的作用是将分散的频谱大值点独立选择出来。
步骤4：进行二维信号的空域混叠和降采样FFT
将加窗的重排图像信号j(x′,y′)进行空域混叠，即将原大小为N×N的j(x′,y′)处理为大小为B×B的降采样混叠图像信号z(u,v)。
对z(u,v)进行FFT得到B×B点的频谱Z(α,β)，Z(α,β)中包含有需要恢复频谱的大值点。
步骤5：频谱还原
对于降采样到大小为B×B的频谱来说，不需要将所有频谱的值还原，仅需还原K个频谱大值点到大小为N×N的频谱，将其余点置零即可，此处K≤B。
这个步骤分两步完成，第一是位置还原，第二是估值还原。
5.1位置还原
位置还原包括两步，第一步是频谱上采样，第二步是频谱逆重排。
频谱上采样是指将大小为B×B的频谱还原为大小为N×N的频谱。由于频域降采样过程为等间隔抽取频谱点，那么由Z(α,β)的频谱坐标(α,β)还原为重排后降采样前频谱J(λ′,θ′)的频谱坐标(λ′,θ′)为：
λ′＝(N/B)·α，θ′＝(N/B)·β (9)
频谱逆重排是指还原由于步骤2信号重排带来的频谱坐标变化。由频谱重排后坐标(λ′,θ′)计算频谱重排前坐标(λ,θ)的方法为：
λ＝rem(σm·λ′,N),θ＝rem(σn·θ′,N) (10)
得到频谱大值点在原始频谱图中的坐标，即频谱重排前坐标；其中rem表示取余操作；
5.2估值还原
二维状态下的大值点还原公式为(11)所示：
上式中Ir(λ,θ)为还原的二维频谱，Z(α,β)为经降采样FFT后的得到的B×B的频谱，WF为窗函数wf的频谱，oσ(λ)为λ方向偏移量，oσ(θ)为θ方向偏移量，定义如公式(12)所示。
通过估值还原频谱大值点的幅值。
重复步骤2-5对每一幅莫尔条纹图分别进行解算；
步骤6：将还原频谱用于数字莫尔移相解相过程
通过低通滤波对还原的莫尔条纹频谱Ir(λ,θ)进行处理，得到低频部分频谱分布。经过逆傅里叶变换和公式(4)即可得到实际干涉图与虚拟干涉图的相位差Δ(x,y)，结合已知的pV(x,y)，可以得到实际干涉图相位pR(x,y)＝Δ(x,y)+pV(x,y)。
有益效果
1、本发明的基于稀疏傅里叶变换的数字莫尔移相干涉相位实时测量方法，将其扩展为适用于数字莫尔移相解相方法的二维SFFT算法。二维SFFT算法的时间复杂度为O(N2)，而FFT算法的时间复杂度为O(2N2·log2N)，可以看到SFFT的时间复杂度更小，简化了傅里叶变换的计算量。
2、将本发明应用于数字莫尔移相解相方法的傅里叶变换过程，可以降低解相过程的计算量，提高数字莫尔移相解相方法的算法速度，有助于实现相位的实时测量。
附图说明
图1是基于系数傅里叶变换的数字莫尔移相解相方法整体流程图。
图2是仿真中的实际干涉图。
图3是仿真中的虚拟干涉图。
图4是实际干涉图的原始相位分布。
图5是莫尔合成图原始频谱分布。
图6是重排后的莫尔合成图频谱分布。
图7是sinc函数的频谱分布。
图8是高斯函数的频谱分布。
图9是窗函数的频谱分布
图10降采样至B×B大小的频谱分布。
图11重构频谱分布。
图12频谱的重构误差分布。
图13求解得到的实际干涉图相位分布。
图14相位求解误差分布。
具体实施方式
下面结合附图与实施例对本发明作进一步说明。
实施例1
本实例是基于MATLAB仿真以验证方法的可行性。
本实施例以精密加工光学表面的面形测量为例，说明如何利用稀疏傅里叶变换实现面形的实时测量。
首先搭建斐索或泰曼格林型干涉仪，将正在加工的光学表面作为被测面。从干涉仪发出的准直激光束入射光学表面并被其反射，反射光以相位空间分布的形式携带光学表面的面形信息。反射光与标准面反射光发生干涉，将相位空间分布转化为光强空间分布，被面阵探测器记录，形成干涉图。在加工过程中，光学表面的面形持续变化，最终接近设计面形。利用面阵探测器持续采集干涉图，则这一系列干涉图的相位变化将反应光学表面的面形变化。以下说明如何进行单幅干涉图的实时解相。
如图1所示，本实施例公开的基于稀疏傅里叶变换的数字莫尔移相干涉相位实时测量方法，具体实现步骤如下：
步骤1：使用数字莫尔移相解相方法完成解相准备
通过数字莫尔移相解相方法，按照公式(1)和公式(2)对实际干涉图与构造的虚拟干涉图进行莫尔合成，并得到4幅二维的莫尔条纹图iM(x,y)，每幅莫尔条纹图的大小为N×N。通过步骤2到步骤5分别对每一幅莫尔条纹图进行解算。
在数字莫尔移相方法中实际干涉图的分布如图2所示，其大小为1024×1024。这是一个平面光楔反射激光形成的干涉图。针对其面形分布特点，构造虚拟干涉图分布如图3所示，实际干涉图的原始相位分布如图4所示。为了相对清晰的展示莫尔合成图原始频谱，仅展示在中心位置处的部分频谱如图5所示。
步骤2：设计二维SFFT的莫尔条纹图空域重新排列规律
取σx和σy分别为x,y方向的重排系数，对二维SFFT空域莫尔条纹图iM(x,y)进行重新排列得到重排信号i'(x',y')，x',y'为重排信号的坐标。
i′(x′,y′)＝iM(σxx,σyy) (5)
经过空域重新排列，原始莫尔条纹图像的K个频谱大值点被分散至整个频谱的各个位置。此处K＝32。
为了相对清晰的展示经过重排后的频谱分布，仅展示频谱中心位置处的频谱分布如图6所示，可以看到原本距离较近的大值点被分散开。
步骤3：设置二维窗函数
理想的窗函数为一个空域的滤波函数，特点是在频域呈方形分布，方形内为通带，值为1，方形外为禁带，值为0。
以二维sinc函数Sc(x′,y′)为基础构建窗函数，sinc函数的频谱分布放大展示如图7所示。对于sinc函数的频谱来说，除了中心的主峰外，其余位置处也存在较多的次级峰。为了消去次级峰，仅保留主峰作为窗函数的通带，需要在sinc函数上叠加一个高斯函数g(x′,y′)，如图8所示。将二维sinc函数和二维高斯函数配合使用得到二维的空域窗函数wf(x′,y′)，放大展示如图9所示。
其中b和c为系数，分别用于控制二维sinc函数的幅值和频谱宽度；d和μ为系数，分别用于控制高斯函数的幅值和频谱宽度。
使用空域窗函数wf(x′,y′)对重排信号i′(x′,y′)进行滤波处理得到加窗的重排图像信号j(x′,y′)。
j(x′,y′)＝wf(x′,y′)·i′(x′,y′) (7)
窗函数wf(x′,y′)的作用是将分散的频谱大值点独立选择出来。
步骤4：进行二维信号的空域混叠降采样和FFT
将加窗的重排图像信号j(x′,y′)进行空域混叠，即将原大小为N×N的j(x′,y′)处理为大小为B×B的混叠图像信号z(u,v)。
对z(u,v)进行FFT得到B×B点的频谱Z(α,β)，此处B＝256，如图10所示，Z(α,β)中包含有需要恢复频谱的大值点。
步骤5：频谱还原
对于降采样到大小为B×B的频谱来说，不需要将所有频谱的值还原，仅需还原K个频谱大值点到大小为N×N的频谱，将其余点置零即可，此处K＝32<B。
本步骤分两步完成，第一是位置还原，第二是估值还原。
5.1位置还原
位置还原包括两步，第一步是频谱上采样，第二步是频谱逆重排。
频谱上采样是指将大小为B×B的频谱还原为大小为N×N的频谱。由于频域降采样过程为等间隔抽取频谱点，那么由Z(α,β)的频谱坐标(α,β)还原为重排后降采样前频谱J(λ′,θ′)的频谱坐标(λ′,θ′)为：
λ′＝(N/B)·α，θ′＝(N/B)·β (9)
频谱逆重排是指还原由于步骤2信号重排带来的频谱坐标变化。由频谱重排后坐标(λ′,θ′)计算频谱重排前坐标(λ,θ)的方法为：
λ＝rem(σm·λ′,N),θ＝rem(σn·θ′,N) (10)
得到频谱大值点在原始频谱图中的坐标，即频谱重排前坐标；其中rem表示取余操作；
5.2估值还原
二维状态下的大值点还原公式为(11)所示：
上式中Ir(λ,θ)为还原的二维频谱，Z(α,β)为经降采样FFT后的得到的B×B的频谱，WF为窗函数wf的频谱，oσ(λ)为λ方向偏移量，oσ(θ)为θ方向偏移量，定义如公式(12)所示。
通过估值还原频谱大值点的幅值。
得到的还原频谱如图11所示，与原始频谱相比，频谱还原误差如图12所示。可以看到，存在一定的误差，但误差相较原始频谱均小1至3个数量级，说明还原频谱基本准确。
重复步骤2-5对每一幅莫尔条纹图分别进行解算；
步骤六：将还原频谱用于数字莫尔移相解相过程
通过低通滤波对还原的莫尔条纹频谱Ir(λ,θ)进行处理，得到低频部分频谱分布。经过逆傅里叶变换和公式(4)即可得到实际干涉图与虚拟干涉图的相位差Δ(x,y)，结合已知的pV(x,y)，可以得到实际干涉图相位pR(x,y)＝Δ(x,y)+pV(x,y)，如图13所示。图13与图4点对点相减求解误差如图14所示，可以看到恢复的相位分布基本正确，也存在一定误差，但这个误差除了频谱重构误差带来的影响，也有数字莫尔移相方法本身的误差影响。因此整体来说，将稀疏傅里叶变换算法用于数字莫尔移相解相方法是可行的。
上述基于SFFT的数字莫尔移相干涉相位实时测量技术在除了FFT以外的步骤与原有方法相同，计算复杂度相同。但FFT算法的时间复杂度为O(2N2·log2N)，二维SFFT算法的时间复杂度为O(N2)，当N＝1024时，本实施例与原有的基于FFT的数字莫尔移相相位测量方法相比，两者计算时间差别约为20倍。考虑到目前FFT是数字莫尔移相干涉算法中计算量最大、占时最久的步骤，这一计算时间的提升将极大有利于实时测量的实现。
除了实施例1列举的精密加工光学表面的面形测量领域以外，本发明公开的技术对于气液流场、生物细胞形态测量过程中产生的相位实时测量问题同样适用。
以上所述的具体描述，对发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限定本发明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。
Claims (2)
1.基于稀疏傅里叶变换的数字莫尔移相干涉相位实时测量方法，其特征在于：包括如下步骤：
步骤1、使用数字莫尔移相解相方法完成解相准备
通过数字莫尔移相解相方法，对实际干涉图与构造的虚拟干涉图进行莫尔合成，并得到4幅二维的莫尔条纹图iM(x,y)，每幅莫尔条纹图的大小为N×N；通过步骤2到步骤5对每一幅莫尔条纹图分别进行解算，再将还原频谱用于数字莫尔移相解相过程，即能够得到实际干涉图相位，实现实时测量；
步骤2、设计二维SFFT的莫尔条纹图空域重新排列规律
σx和σy分别为x,y方向的重排系数，对二维SFFT空域的莫尔条纹分布iM(x,y)进行重新排列得到重排信号i'(x',y')，x',y'为重排信号的坐标；
i′(x′,y′)＝iM(σxx,σyy) (1)
经过空域重新排列，原始莫尔条纹图像的K个频谱大值点被分散至整个频谱的各个位置；
步骤3、设置二维窗函数
理想的窗函数为一个空域的滤波函数，特点是在频域呈方形分布，方形内为通带，值为1，方形外为禁带，值为0；
以二维sinc函数Sc(x′,y′)为基础构建窗函数；对于sinc函数的频谱来说，除了中心的主峰外，其余位置处也存在较多的次级峰；为了消去次级峰，仅保留主峰作为窗函数的通带，需要在sinc函数上叠加一个高斯函数g(x′,y′)；将二维sinc函数和二维高斯函数配合使用得到二维的空域窗函数wf(x′,y′)；
其中b和c为系数，分别用于控制二维sinc函数的幅值和频谱宽度；d和μ为系数，分别用于控制高斯函数的幅值和频谱宽度；
使用空域窗函数wf(x′,y′)对重排信号i′(x′,y′)进行滤波处理得到加窗的重排图像信号j(x′,y′)；
j(x′,y′)＝wf(x′,y′)·i′(x′,y′) (3)
窗函数wf(x′,y′)的作用是将分散的频谱大值点独立选择出来；
步骤4、进行二维信号的空域混叠和降采样FFT
将加窗的重排图像信号j(x′,y′)进行空域混叠，即将原大小为N×N的j(x′,y′)处理为大小为B×B的降采样混叠图像信号z(u,v)；
对z(u,v)进行FFT得到B×B点的频谱Z(α,β)，Z(α,β)中包含有需要恢复频谱的大值点；
步骤5、频谱还原
对于降采样到大小为B×B的频谱来说，不需要将所有频谱的值还原，仅需还原K个频谱大值点到大小为N×N的频谱，将其余点置零即可，此处K≤B；
5.1位置还原
位置还原包括两步，第一步是频谱上采样，第二步是频谱逆重排；
频谱上采样是指将大小为B×B的频谱还原为大小为N×N的频谱；由于频域降采样过程为等间隔抽取频谱点，那么由Z(α,β)的频谱坐标(α,β)还原为重排后降采样前频谱J(λ′,θ′)的频谱坐标(λ′,θ′)为：
λ′＝(N/B)·α，θ′＝(N/B)·β (5)
频谱逆重排是指还原由于步骤2信号重排带来的频谱坐标变化；由频谱重排后的频谱坐标(λ′,θ′)计算频谱重排前坐标(λ,θ)的方法为：
λ＝rem(σx·λ′,N),θ＝rem(σy·θ′,N) (6)
得到频谱大值点在原始频谱图中的坐标，即频谱重排前坐标；其中rem表示取余操作；
5.2估值还原
二维状态下的大值点还原公式为(7)所示：
上式中Ir(λ,θ)为还原的二维频谱，Z(α,β)为经降采样FFT后的得到的B×B的频谱，WF为窗函数wf的频谱，oσ(λ)为λ方向偏移量，oσ(θ)为θ方向偏移量，定义如公式(8)所示；
通过估值还原频谱大值点的幅值；
重复步骤2-5对每一幅莫尔条纹图分别进行解算；
步骤6、将还原频谱用于数字莫尔移相解相过程
通过低通滤波对还原的莫尔条纹频谱Ir(λ,θ)进行处理，得到低频部分频谱分布；经过逆傅里叶变换和公式(9)即可得到实际干涉图与虚拟干涉图的相位差Δ(x,y)，结合已知的构造虚拟相位pV(x,y)，得到实际干涉图相位pR(x,y)＝Δ(x,y)+pV(x,y)；
2.如权利要求1所述的基于稀疏傅里叶变换的数字莫尔移相干涉相位实时测量方法，其特征在于：实现步骤1所述对实际干涉图与构造的虚拟干涉图进行莫尔合成的方法为：
去除直流项的实际干涉图分布iR(x,y)和构造的虚拟干涉图分布iV(x,y)如式(10)所示；
其中bR(x,y)是实际干涉图调制量，pR(x,y)为实际干涉图待测相位，f是x方向添加的空间载波；bV(x,y)是虚拟干涉图调制量，pV(x,y)为构造的虚拟相位，δk为附加的移相相位值，δk＝1，2,3,4＝0,π/2,π,3π/2，对应四步移相算法； 实际干涉图与虚拟干涉图进行莫尔合成，得到莫尔条纹分布iM(x,y)：
Applications Claiming Priority (2)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN2019106244977 | 2019-07-11 | ||
| CN201910624497 | 2019-07-11 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110500968A CN110500968A (zh) | 2019-11-26 |
| CN110500968B true CN110500968B (zh) | 2021-04-20 |
Family
ID=68589035
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910782763.9A Active CN110500968B (zh) | 2019-07-11 | 2019-08-23 | 基于稀疏傅里叶变换的数字莫尔干涉相位实时测量方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110500968B (zh) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111238397B (zh) * | 2020-02-10 | 2021-03-02 | 北京理工大学 | 光学元件面形的瞬态数字莫尔移相干涉测量装置和方法 |
Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101876538A (zh) * | 2010-05-07 | 2010-11-03 | 中国科学院光电技术研究所 | 一种接近式纳米光刻中的间隙测量方法 |
| CN102410819A (zh) * | 2011-08-08 | 2012-04-11 | 苏州大学 | 一种测量膜基反射镜三维面形的方法 |
| CN104813648A (zh) * | 2012-12-11 | 2015-07-29 | 富士胶片株式会社 | 图像处理装置、摄像装置、图像处理方法及图像处理程序 |
| CN105572649A (zh) * | 2015-12-11 | 2016-05-11 | 中北大学 | 基于稀疏傅里叶变换的雷达目标检测方法 |
| US9390468B1 (en) * | 2009-06-22 | 2016-07-12 | Lockheed Martin Corporation | Imaging systems and methods for generating image data |
| CN107388963A (zh) * | 2017-07-13 | 2017-11-24 | 北京理工大学 | 将小波分析和低通滤波结合的数字莫尔条纹相位提取方法 |
| CN109101462A (zh) * | 2018-09-13 | 2018-12-28 | 南京邮电大学 | 可降低时间损耗的傅里叶变换算法 |
Family Cites Families (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US7043082B2 (en) * | 2000-01-06 | 2006-05-09 | Canon Kabushiki Kaisha | Demodulation and phase estimation of two-dimensional patterns |
| US8345242B2 (en) * | 2008-10-28 | 2013-01-01 | Molecular Imprints, Inc. | Optical system for use in stage control |
| US8983029B2 (en) * | 2013-02-22 | 2015-03-17 | Shimadzu Corporation | Radiographic apparatus and method for the same |
| US20150071105A1 (en) * | 2013-09-12 | 2015-03-12 | The Provost, Fellows, Foundation Scholars, & the Other Members of Board, of The College of the Holy | Method and System for Compensating for Interference Due to Carrier Frequency Offset in an OFDM Communication System |
-
2019
- 2019-08-23 CN CN201910782763.9A patent/CN110500968B/zh active Active
Patent Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US9390468B1 (en) * | 2009-06-22 | 2016-07-12 | Lockheed Martin Corporation | Imaging systems and methods for generating image data |
| CN101876538A (zh) * | 2010-05-07 | 2010-11-03 | 中国科学院光电技术研究所 | 一种接近式纳米光刻中的间隙测量方法 |
| CN102410819A (zh) * | 2011-08-08 | 2012-04-11 | 苏州大学 | 一种测量膜基反射镜三维面形的方法 |
| CN104813648A (zh) * | 2012-12-11 | 2015-07-29 | 富士胶片株式会社 | 图像处理装置、摄像装置、图像处理方法及图像处理程序 |
| CN105572649A (zh) * | 2015-12-11 | 2016-05-11 | 中北大学 | 基于稀疏傅里叶变换的雷达目标检测方法 |
| CN107388963A (zh) * | 2017-07-13 | 2017-11-24 | 北京理工大学 | 将小波分析和低通滤波结合的数字莫尔条纹相位提取方法 |
| CN109101462A (zh) * | 2018-09-13 | 2018-12-28 | 南京邮电大学 | 可降低时间损耗的傅里叶变换算法 |
Non-Patent Citations (6)
| Title |
|---|
| Fast two-step layer-based method for computer generated hologram using subsparse 2D fast Fourier transform;JIA JIA等;《OPTICS EXPRESS》;20180625;第26卷(第13期);第17487-17497页 * |
| SFFT based phase demodulation for faster interference fringes analysis;Chen-Yu Lee等;《Proc. of SPIE》;20160316;第9754卷;第975417-1~ 975417-7页 * |
| Simple and Practical Algorithm for Sparse Fourier Transform;Haitham Hassanieh等;《Proc.SODA’12》;20121231;正文第1-12页 * |
| Study on influential factors of measurement precision using Fourier transform profilometry;Wei Shang等;《Optik》;20180430;第159卷;第133-142页 * |
| 基于一维空间数字相关的实时相移标定方法;罗印龙等;《中国激光》;20080804;第32卷(第8期);第1117-1122页 * |
| 大梯度相位单幅干涉图的解相方法;叶明哲等;《中国激光》;20190531;第46卷(第5期);第0504002-1~0504002-3页第1-2节 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN110500968A (zh) | 2019-11-26 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| US4791584A (en) | Sub-nyquist interferometry | |
| Surrel | Fringe analysis | |
| Fillard | Sub-pixel accuracy location estimation from digital signals | |
| AU659271B2 (en) | Spatial wavefront evaluation by intensity relationship | |
| CN108872153B (zh) | 基于非均匀傅里叶变换的平行平板光学均匀性的测量方法 | |
| CN107917676B (zh) | 一种基于条纹图像频谱分析的干涉测量方法 | |
| Fercher | Computer-generated holograms for testing optical elements: error analysis and error compensation | |
| CN112525070B (zh) | 基于非均匀采样校正的抗振动白光干涉测量方法 | |
| CN112097678B (zh) | 基于频率盲估计的多表面面形测量方法 | |
| Lassahn et al. | Multiphase fringe analysis with unknown phase shifts | |
| CN111521112B (zh) | 一种傅里叶及窗口傅里叶变换的联合相位重构算法 | |
| Lehmann et al. | Analysis of interference microscopy in the spatial frequency domain | |
| CN110500968B (zh) | 基于稀疏傅里叶变换的数字莫尔干涉相位实时测量方法 | |
| CN114088010B (zh) | 一种干涉三维形貌解算方法 | |
| CN108170636B (zh) | 一种从线性载频移相干涉图中恢复出真实相位的方法 | |
| Colonna de Lega | Processing of non-stationary interference patterns: adapted phase-shifting algorithms and wavelet analysis. Application to dynamic deformation measurements by holographic and speckle interferometry | |
| CN117109450B (zh) | 一种自适应参数谱估计的高分辨光谱共焦干涉薄膜测厚法 | |
| US6693715B2 (en) | Fringe analysis method using fourier transform | |
| CN111562088A (zh) | 基于采样函数的平行平板光学参数的测量方法 | |
| Takeda | Spatial carrier heterodyne techniques for precision interferometry and profilometry: an overview | |
| Taylor et al. | An error analysis of holographic strains determined by cubic splines: The errors are analyzed which arise when computer-simulated holographic data of a loaded cantilever beam are interpolated and differentiated via cubic splines to yield strains and bending moments | |
| JPH0652162B2 (ja) | 干渉測定方法 | |
| Loheide | Innovative evaluation method for shearing interferograms | |
| Kara-Mohammed et al. | Imaging of particles with 3D full parallax mode with two-color digital off-axis holography | |
| Larkin | Topics in multi-dimensional signal demodulation |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |