CN101968383A - 一种抗扰动的时频域波前检测方法 - Google Patents
一种抗扰动的时频域波前检测方法 Download PDFInfo
- Publication number
- CN101968383A CN101968383A CN 201010270213 CN201010270213A CN101968383A CN 101968383 A CN101968383 A CN 101968383A CN 201010270213 CN201010270213 CN 201010270213 CN 201010270213 A CN201010270213 A CN 201010270213A CN 101968383 A CN101968383 A CN 101968383A
- Authority
- CN
- China
- Prior art keywords
- sequence
- phase
- light intensity
- frequency
- matrix
- 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.)
- Granted
Links
- 238000001514 detection method Methods 0.000 title claims abstract description 9
- 238000001228 spectrum Methods 0.000 claims abstract description 25
- 238000000034 method Methods 0.000 claims abstract description 22
- 239000011159 matrix material Substances 0.000 claims abstract description 19
- 238000001914 filtration Methods 0.000 claims abstract description 16
- 230000007613 environmental effect Effects 0.000 claims abstract description 7
- 230000010363 phase shift Effects 0.000 claims description 15
- 238000004364 calculation method Methods 0.000 claims description 10
- 238000005070 sampling Methods 0.000 claims description 7
- 238000012545 processing Methods 0.000 claims description 6
- 238000000605 extraction Methods 0.000 claims description 4
- 238000005259 measurement Methods 0.000 abstract description 9
- 230000000694 effects Effects 0.000 abstract description 5
- 230000003287 optical effect Effects 0.000 abstract description 4
- 238000010586 diagram Methods 0.000 description 5
- 238000003672 processing method Methods 0.000 description 4
- 238000005516 engineering process Methods 0.000 description 3
- 238000005305 interferometry Methods 0.000 description 3
- 230000003595 spectral effect Effects 0.000 description 3
- 238000000819 phase cycle Methods 0.000 description 2
- 230000001629 suppression Effects 0.000 description 2
- 230000004075 alteration Effects 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 238000013500 data storage Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000036039 immunity Effects 0.000 description 1
- 238000013178 mathematical model Methods 0.000 description 1
- 238000010183 spectrum analysis Methods 0.000 description 1
- 230000002123 temporal effect Effects 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 238000013519 translation Methods 0.000 description 1
Images
Landscapes
- Spectrometry And Color Measurement (AREA)
Abstract
本发明涉及一种抗扰动的时频域波前检测方法，属于光学测量领域。本发明首先采集一系列干涉图，将每张干涉图中任意一个相同位置的像素点的光强序列减去该序列的均值后，作快速傅里叶变换，再进行宽带滤波移频，然后对得到频谱序列作快速傅里叶逆变换，并对结果求对数，提取结果的虚部并解包裹后得到一个时序的相位统计序列；再然后对得到的相位统计序列作线性拟合，将拟合所得序列的第一个值存入矩阵中，作为该像素点初始相位的计算值；最后，重复上述操作，将其他像素点的初始相位计算值均存入矩阵后，再经空域解包裹重建波前原始相位的分布。本发明的方法简便易行，免受环境振动的影响，具有很好的抗扰动效果。
Description
技术领域
本发明涉及一种抗扰动的时频域波前检测方法，属于光学测量领域。
背景技术
移相干涉测量技术已经被越来越多的应用于光学元件面形测量及系统波像差的测量中。其原理是通过一定手段在干涉仪的两相干光程间引入有序的相移，在此过程中计算机对干涉图进行采样，将光强数字化后存储在存储器中，然后按照一定的数学模型根据光强的变化计算出相应的相位分布。移相过程中两相干光束形成的干涉场的光强分布可表示为：
其中(x，y)是干涉图的图像坐标系中任意一像素点的坐标，t是时间，I0是背景光强，V是对比度，
就是该点上的被测初始相位，Δδ(t)为移相量。当设定Δδ(t)分别等于0，π/2，π，3π/2时，可通过下式计算出原始相位值：
(2)式即为四步法的计算公式，之后发展起来的五步法、七步法、十三步法等方法均基于类似的原理。
以四步法为代表的这类时域的处理方法要求高精度的线性等步长移相才能满足高精度的测量要求，环境振动及气流扰动等因素对此类移相干涉测量系统有严重的影响，因此就不免要增加复杂的反馈控制系统或采用特殊的系统结构，这无疑会在一定程度上限制系统的应用范围和增加系统的复杂程度。
随着计算机技术以及数据存储技术的发展，数据处理过程中干涉图帧数的增加和信息量的增大已经不再成为问题，于是出现了基于干涉场时序光强的频谱分析的移相干涉技术。连续移相过程中，干涉场中某点的时序光强I(x，y，t)在时域中具有余弦形式的变化规律，可用下式表示：
其中，f0是连续移相的移相速率。
对光场中各点的时序光强信号进行采样并做傅里叶变换，然后通过提取频谱中基频峰值点的相位就可以计算出该点所对应的初始相位。这种单点提取相位的频域处理方法理论上确实具有一定的抗噪声能力，但是实际的测量过程中不仅会存在环境振动及随机噪声等扰动因素，而且移相速率的线性度也很难保证，此时某像素点的时序光强不再是标准的余弦形式，而是如下式所示：
其中，I0是光强均值，V是对比度，f0是移相速率的线性部分，t是时间，
是干涉图中坐标(x，y)处的被测初始相位，δ(x，y，t)是环境振动和移相非线性引起的相位畸变，δr(x，y，t)是随机扰动引起的相位畸变的等效值。
当对(4)式所示的时序光强采样并作傅里叶变换后，频谱中会因各种干扰因素的存在出现多个峰值频率，而且这些干扰频率的能量很容易超过基频的峰值甚至将基频峰值淹没，从而导致基频峰值点提取失败。可见，这种单点的频域处理方法要保证好的鲁棒性同样需要保证严格的线性连续移相才能满足高精度的测量要求，因此复杂的反馈控制系统同样必不可少。
发明内容
本发明为了解决定步长移相干涉术和单点频谱处理方法抗干扰能力差的问题，提出一种抗扰动的时频域波前检测方法。本方法通过对光场中时序光强序列的频谱的宽带滤波和时域相位统计序列的线性拟合，实现抗扰动的高精度测量。该方法简便易行，抗扰动效果好，可以应用在机械移相干涉仪和波长调谐移相干涉仪中。
本发明的一种抗扰动的时频域波前检测方法，具体实现步骤如下：
1)采集一系列的干涉图
以大于f0二倍的采样频率fs采集一系列分辨率相同的干涉图并顺序存储在计算机中；此时每张干涉图中任意相同位置的像素点(x，y)上的光强值组成一个时序的光强序列，通过下式表示：
其中，I0是光强均值，V是对比度，f0是移相速率的线性部分，tn是采样时间序列，
是干涉图中坐标(x，y)处被测的初始相位，δ(x，y，tn)是环境振动和移相非线性引起的相位畸变，δr(x，y，tn)是随机扰动引起的相位畸变的等效值，随机扰动包括CCD随机噪声和气流扰动等因素。
同时，建立一个与干涉图同维数的空白矩阵，用以保存第6)步计算所得到的各像素点上初始相位的计算值。
2)减去光强均值
选取第1)步得到的每张干涉图中任意一个相同位置的像素点(x0，y0)，将该像素点上的光强序列减去该序列的均值以去除傅里叶变换后频谱中的零频成分，得到的新序列记为S(x0，y0，tn)，可用下式表示：
3)光强序列的快速傅里叶变换
对于第2)步得到的减均值处理后的光强序列作快速傅里叶变换，得到减均值处理后的光强序列的频谱序列。对于傅里叶变换前的光强序列，无需考虑时域窗口的具体形式，采用最简单易行的矩形窗即可。经快速傅里叶变换后频谱的正负半轴部分就分别代表(6)式中正半部分和负半部分的傅里叶变换。
4)滤波和移频
采用宽度为0～fs/2的矩形窗对第3)步得到的频谱序列进行宽带滤波，将频谱序列的整个正半轴部分滤出。
滤波后，对滤出的频谱序列沿频率坐标轴移频-f0得到新的频谱序列。
5)快速傅里叶逆变换和提取相位统计序列
对第4)步得到的频谱序列作快速傅里叶逆变换，再对结果求对数，提取结果的虚部并解包裹后得到一个时序的相位统计序列Φ(x0，y0，tn)，通过下式表示：
由于采用了宽带滤波，此时得到的时序相位包含了所有的相位信息。
6)线性拟合并取出结果的第一个数存入矩阵
对第5)步得到的相位统计序列Φ(x0，y0，tn)作线性拟合，过程如下式表示：
其中，N为干涉图帧数。第4)步的实际操作过程中移频量与基频真值f0间通常会存在微小偏差Δf，从而(7)式和(8)式中会由此而增加一个斜率项Δftn，因此在线性拟合后我们将tn＝0时的
，即拟合所得序列的第一个值存入第1)步建立的矩阵中(x0，y0)处，作为当前坐标(x0，y0)处初始相位的计算值，此时斜率项Δftn被自动消除。该计算值可用下式表示：
其中，
代表非随机性扰动对计算结果的影响，代表随机噪声对计算结果的影响。在采样帧数足够多的情况下，经线性拟合后的随机噪声项
的值趋于零，而非随机性扰动的影响被归结为一个与像素点的坐标无关的值
，这个值在所有像素点上都是相等的，因此不影响被测波前的形状，只对整个波前增加一个活塞式的平移。
7)空域解包裹重建波前图
选取其他像素点(x′，y′)上的光强序列，重复第2)步到第6)步的操作，直到得到所有像素点的初始相位计算值并将其存入矩阵中与选取像素点相同的坐标位置，最终得到一个带有包裹的相位分布矩阵；对得到的矩阵经空域解包裹重建波前原始相位的分布。
有益效果
1.频谱的宽带滤波不需要精确定位基频位置，免受环境振动的影响，既不需要高度精确的移相步长，也不需要严格线性的移相过程，可使干涉仪脱离特殊的减振设备独立工作；
2.采集的干涉图帧数足够多的情况下，通过线性拟合将随机噪声均化为零，并将非随机的扰动归结为一个与像素坐标无关的值，因而具有很好的抗扰动效果；
3.本发明中的时频域波前检测方法对光源光强随时间的漂移也有一定免疫能力，因此可结合变频移相，使干涉仪体积更小，重量更轻，系统更为简单，为研制可用于现场测试的小型化、便携式干涉仪奠定了良好的基础，在机械移相和波长调谐移相的泰曼-格林干涉仪和斐索干涉仪中均可得到应用。
附图说明
图1为本发明抗扰的时频域波前检测方法计算流程图；
图2为采集到干涉图后像素点上得到的光强序列的时序分布图；
图3为光强序列减均值后得到的新序列的时序分布图；
图4为减均值后的光强序列经快速傅里叶变换后得到的频谱分布图及滤波窗口；
图5为经滤波和移频后得到的新的频谱图；
图6为快速傅里叶逆变换和对数操作后提取到的时序的相位统计序列以及对相位统计序列对线性拟合所得的结果图；
图7为经时频域波前检测方法计算所重建的相位分布图。
具体实施方式
下面结合附图和实施例对本发明做进一步说明：
实施例
1)采集一系列的干涉图
模拟的干涉图分辨率为128*128，干涉场的光强均值I0＝120，对比度V＝5/6，移相速率f0＝20Hz，环境振动为频率2Hz、振幅0.85λ的余弦，随机噪声引起的光强波动σ＝1。采样频率fs＝100Hz，采集干涉图帧数为100帧，即采样时间为1s。
以fs＝100Hz的采样频率采集一系列的干涉图并存储在计算机中。
此时每张干涉图中任意相同位置的像素点(x，y)上得到的时序的光强序列，通过下式表示：
同时建立一个128*128的空白矩阵A(128，128)，用于存放第6)步得到的每个像素点初始相位的计算值。
2)光强序列减去序列均值
选取第1)步采集的每张干涉图中坐标为(1，1)的像素点上的光强值，得到一个如(5)式所表示的时序的光强序列，如图2所示；将该像素点的光强序列减去该序列的均值后得到新的光强序列，如图3所示。
3)减均值后的光强序列的快速傅里叶变换
对于第2)步得到的减均值处理后的光强序列作快速傅里叶变换，得到的新序列，如图4所示。
4)滤波和移频
采用宽度为0～50的矩形窗对图4所表示的频谱序列进行宽带滤波，滤波窗口如图4中虚线框所示。滤波后，对滤出的部分频谱序列沿频率坐标轴向负方向平移20Hz得到新的频谱，新的频谱如图5所示。
5)快速傅里叶逆变换和提取相位统计序列
对图5中所示的新的频谱序列作快速傅里叶逆变换，然后对结果求对数并提取虚部得到一个带包裹的相位序列，对此相位序列解包裹后得到一个时序的相位统计序列，如图6中曲线所示。
6)线性拟合并取出结果的第一个数存入矩阵
对第5)步得到的时序相位统计序列进行线性拟合，其结果如图6虚线表示。将拟合所得序列中的第一个数取出作为坐标为(1，1)处的初始相位计算值存入第1)步建立的矩阵A中的(1，1)处。
7)空域解包裹重建波前图
选取其他像素点(x′，y′)上的光强序列，重复第2)步到第6)步的操作。当所有像素点计算完毕后，最终得到的A为一个带有包裹的相位分布矩阵；对得到的矩阵A经空域解包裹计算重建波前原始相位的分布，如图7所示。最终面形测量误差PMS＝0.0005λ，PV＝0.0038λ。
Claims (1)
1.一种抗扰动的时频域波前检测方法，具体实现步骤如下：
1)采集一系列的干涉图
以大于f0二倍的采样频率fs采集一系列分辨率相同的干涉图并顺序存储在计算机中；此时每张干涉图中任意相同位置的像素点(x，y)上的光强值组成一个时序的光强序列，通过下式表示：
其中，I0是光强均值，V是对比度，f0是移相速率的线性部分，tn是采样时间序列，是干涉图中坐标(x，y)处被测的初始相位，δ(x，y，tn)是环境振动和移相非线性引起的相位畸变，δr(x，y，tn)是随机扰动引起的相位畸变的等效值，随机扰动包括CCD随机噪声和气流扰动等因素；
同时，建立一个与干涉图同维数的空白矩阵，用以保存第6)步计算所得到的各像素点上初始相位的计算值；
2)减去光强均值
选取第1)步得到的每张干涉图中任意一个相同位置的像素点(x0，y0)，将该像素点上的光强序列减去该序列的均值以去除傅里叶变换后频谱中的零频成分，得到的新序列记为S(x0，y0，tn)，可用下式表示：
3)光强序列的快速傅里叶变换
对于第2)步得到的减均值处理后的光强序列作快速傅里叶变换，得到减均值处理后的光强序列的频谱序列；
4)滤波和移频
采用宽度为0～fs/2的矩形窗对第3)步得到的频谱序列进行宽带滤波，将频谱序列的整个正半轴部分滤出；
滤波后，对滤出的频谱序列沿频率坐标轴移频-f0得到新的频谱序列；
5)快速傅里叶逆变换和提取相位统计序列
对第4)步得到的频谱序列作快速傅里叶逆变换，再对结果求对数，提取结果的虚部并解包裹后得到一个时序的相位统计序列Φ(x0，y0，tn)，通过下式表示：
6)线性拟合并取出结果的第一个数存入矩阵
对第5)步得到的相位统计序列Φ(x0，y0，tn)作线性拟合，过程如下式表示：
其中，N为干涉图帧数，将tn＝0时的
即拟合所得序列的第一个值存入第1)步建立的矩阵中(x0，y0)处，作为当前坐标(x0，y0)处初始相位的计算值，此时斜率项Δftn被自动消除；该计算值可用下式表示：
7)空域解包裹重建波前图
选取其他像素点(x′，y′)上的光强序列，重复第2)步到第6)步的操作，直到得到所有像素点的初始相位计算值并将其存入矩阵中与选取像素点相同的坐标位置，最终得到一个带有包裹的相位分布矩阵
；对得到的矩阵经空域解包裹重建波前原始相位的分布。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN2010102702138A CN101968383B (zh) | 2010-09-02 | 2010-09-02 | 一种抗扰动的时频域波前检测方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN2010102702138A CN101968383B (zh) | 2010-09-02 | 2010-09-02 | 一种抗扰动的时频域波前检测方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN101968383A true CN101968383A (zh) | 2011-02-09 |
| CN101968383B CN101968383B (zh) | 2012-01-25 |
Family
ID=43547570
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN2010102702138A Expired - Fee Related CN101968383B (zh) | 2010-09-02 | 2010-09-02 | 一种抗扰动的时频域波前检测方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN101968383B (zh) |
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106383512A (zh) * | 2016-10-10 | 2017-02-08 | 中国科学技术大学 | 一种dsp控制系统的电磁防护加固方法 |
| CN106840418A (zh) * | 2017-01-22 | 2017-06-13 | 中国工程物理研究院机械制造工艺研究所 | 一种移相干涉仪的抗振动方法 |
| CN112066909A (zh) * | 2020-08-24 | 2020-12-11 | 南京理工大学 | 一种基于倾斜平面高精度提取的抗振动干涉测量方法 |
| CN114964181A (zh) * | 2022-05-27 | 2022-08-30 | 哈尔滨工业大学 | 基于波前零差干涉的高精度双轴激光水平仪及测量方法 |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1904569A (zh) * | 2006-08-07 | 2007-01-31 | 中国科学院光电技术研究所 | 一种基于线性相位反演的波前测量方法 |
| US7175278B2 (en) * | 2003-06-20 | 2007-02-13 | Visx, Inc. | Wavefront reconstruction using fourier transformation and direct integration |
| US20070222948A1 (en) * | 2006-03-23 | 2007-09-27 | Visx, Incorporated | Systems and methods for wavefront reconstruction for aperture with arbitrary shape |
| CN101387555A (zh) * | 2007-09-14 | 2009-03-18 | 富士能株式会社 | 光学拾波器用波前测量装置 |
| CN101726366A (zh) * | 2009-12-02 | 2010-06-09 | 山东师范大学 | 一种基于多针孔板的波前测量方法和装置 |
-
2010
- 2010-09-02 CN CN2010102702138A patent/CN101968383B/zh not_active Expired - Fee Related
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US7175278B2 (en) * | 2003-06-20 | 2007-02-13 | Visx, Inc. | Wavefront reconstruction using fourier transformation and direct integration |
| US20070222948A1 (en) * | 2006-03-23 | 2007-09-27 | Visx, Incorporated | Systems and methods for wavefront reconstruction for aperture with arbitrary shape |
| CN1904569A (zh) * | 2006-08-07 | 2007-01-31 | 中国科学院光电技术研究所 | 一种基于线性相位反演的波前测量方法 |
| CN101387555A (zh) * | 2007-09-14 | 2009-03-18 | 富士能株式会社 | 光学拾波器用波前测量装置 |
| CN101726366A (zh) * | 2009-12-02 | 2010-06-09 | 山东师范大学 | 一种基于多针孔板的波前测量方法和装置 |
Cited By (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106383512A (zh) * | 2016-10-10 | 2017-02-08 | 中国科学技术大学 | 一种dsp控制系统的电磁防护加固方法 |
| CN106383512B (zh) * | 2016-10-10 | 2019-02-01 | 中国科学技术大学 | 一种dsp控制系统的电磁防护加固方法 |
| CN106840418A (zh) * | 2017-01-22 | 2017-06-13 | 中国工程物理研究院机械制造工艺研究所 | 一种移相干涉仪的抗振动方法 |
| CN106840418B (zh) * | 2017-01-22 | 2019-04-16 | 中国工程物理研究院机械制造工艺研究所 | 一种移相干涉仪的抗振动方法 |
| CN112066909A (zh) * | 2020-08-24 | 2020-12-11 | 南京理工大学 | 一种基于倾斜平面高精度提取的抗振动干涉测量方法 |
| CN112066909B (zh) * | 2020-08-24 | 2022-04-08 | 南京理工大学 | 一种基于倾斜平面高精度提取的抗振动干涉测量方法 |
| CN114964181A (zh) * | 2022-05-27 | 2022-08-30 | 哈尔滨工业大学 | 基于波前零差干涉的高精度双轴激光水平仪及测量方法 |
| CN114964181B (zh) * | 2022-05-27 | 2023-04-25 | 哈尔滨工业大学 | 基于波前零差干涉的高精度双轴激光水平仪及测量方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN101968383B (zh) | 2012-01-25 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| KR101916499B1 (ko) | 간섭 감지 시스템에서 움직임을 보상하기 위한 방법 및 장치 | |
| CN101968383B (zh) | 一种抗扰动的时频域波前检测方法 | |
| TW200604494A (en) | Small displacement measuring method and instrument | |
| CN108872153B (zh) | 基于非均匀傅里叶变换的平行平板光学均匀性的测量方法 | |
| CN107917676A (zh) | 一种基于条纹图像频谱分析的干涉测量方法 | |
| CN112597947A (zh) | 一种基于傅里叶域光学相干层析成像技术的色散补偿方法 | |
| Wang et al. | An efficient and adaptive approach for noise filtering of SAR interferometric phase images | |
| CN109507704A (zh) | 一种基于互模糊函数的双星定位频差估计方法 | |
| Li et al. | Maximum likelihood estimation of optical path length in spectral interferometry | |
| CN104155011B (zh) | 一种二维干涉图的相位提取方法 | |
| CN111521112A (zh) | 一种傅里叶及窗口傅里叶变换的联合相位重构算法 | |
| CN106482633B (zh) | 一种基于π/4相移的多光束干涉相位提取方法 | |
| Kimbrough et al. | The spatial frequency response and resolution limitations of pixelated mask spatial carrier based phase shifting interferometry | |
| CN105865371B (zh) | 一种基于互相关计算的白光干涉显微轮廓复原方法 | |
| Cheng et al. | Additive-to-multiplicative moiré fringe transition in simultaneous dual-wavelength interferometry | |
| Kulkarni et al. | Multiple phase estimation in digital holographic interferometry using product cubic phase function | |
| CN104614083B (zh) | 一种恢复相移干涉图相位分布、及获取两幅图间相移量的方法 | |
| CN112883318A (zh) | 相减策略的多频衰减信号参数估计算法 | |
| Chong et al. | Spectral effects of dual wavelength low coherence light source in white light interferometry | |
| Guo et al. | Order-crossing removal in Gabor order tracking by independent component analysis | |
| Ganotra et al. | Profilometry for the measurement of three-dimensional object shape using radial basis function, and multi-layer perceptron neural networks | |
| CN106441082A (zh) | 一种相位恢复方法及装置 | |
| CN105181300B (zh) | 低相干频域干涉图的自适应干涉项提取方法 | |
| CN113674370A (zh) | 基于卷积神经网络的单帧干涉图解调方法 | |
| Lan et al. | A carrier removal method based on frequency domain self-filtering for interferogram analysis |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| C14 | Grant of patent or utility model | ||
| GR01 | Patent grant | ||
| C17 | Cessation of patent right | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20120125 Termination date: 20130902 |