CN110108643A - 一种用于光声检测的干涉条纹相位提取方法 - Google Patents
一种用于光声检测的干涉条纹相位提取方法 Download PDFInfo
- Publication number
- CN110108643A CN110108643A CN201910305054.1A CN201910305054A CN110108643A CN 110108643 A CN110108643 A CN 110108643A CN 201910305054 A CN201910305054 A CN 201910305054A CN 110108643 A CN110108643 A CN 110108643A
- Authority
- CN
- China
- Prior art keywords
- phase
- point
- irradiation level
- interference fringe
- indicate
- 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
- 238000000605 extraction Methods 0.000 title claims abstract description 31
- 238000001514 detection method Methods 0.000 title claims abstract description 26
- 238000000034 method Methods 0.000 claims abstract description 53
- 230000001186 cumulative effect Effects 0.000 claims abstract description 27
- 238000001914 filtration Methods 0.000 claims abstract description 10
- 238000005305 interferometry Methods 0.000 claims abstract description 8
- 238000006243 chemical reaction Methods 0.000 claims abstract description 7
- 230000010354 integration Effects 0.000 claims abstract description 5
- 238000004364 calculation method Methods 0.000 claims description 3
- 238000005303 weighing Methods 0.000 claims 1
- 238000004422 calculation algorithm Methods 0.000 abstract description 11
- 239000000284 extract Substances 0.000 abstract description 2
- 230000003287 optical effect Effects 0.000 description 5
- 238000004458 analytical method Methods 0.000 description 2
- 238000006073 displacement reaction Methods 0.000 description 2
- 238000005259 measurement Methods 0.000 description 2
- 230000010363 phase shift Effects 0.000 description 2
- 238000004867 photoacoustic spectroscopy Methods 0.000 description 2
- 238000004088 simulation Methods 0.000 description 2
- 238000012360 testing method Methods 0.000 description 2
- 238000010586 diagram Methods 0.000 description 1
- 239000006185 dispersion Substances 0.000 description 1
- 238000002986 genetic algorithm method Methods 0.000 description 1
- 230000002068 genetic effect Effects 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 238000001228 spectrum Methods 0.000 description 1
- 230000003068 static effect Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01J—MEASUREMENT OF INTENSITY, VELOCITY, SPECTRAL CONTENT, POLARISATION, PHASE OR PULSE CHARACTERISTICS OF INFRARED, VISIBLE OR ULTRAVIOLET LIGHT; COLORIMETRY; RADIATION PYROMETRY
- G01J9/00—Measuring optical phase difference; Determining degree of coherence; Measuring optical wavelength
- G01J9/02—Measuring optical phase difference; Determining degree of coherence; Measuring optical wavelength by interferometric methods
- G01J9/0246—Measuring optical wavelength
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N21/00—Investigating or analysing materials by the use of optical means, i.e. using sub-millimetre waves, infrared, visible or ultraviolet light
- G01N21/17—Systems in which incident light is modified in accordance with the properties of the material investigated
- G01N21/1702—Systems in which incident light is modified in accordance with the properties of the material investigated with opto-acoustic detection, e.g. for gases or analysing solids
Landscapes
- Physics & Mathematics (AREA)
- Spectroscopy & Molecular Physics (AREA)
- General Physics & Mathematics (AREA)
- General Health & Medical Sciences (AREA)
- Analytical Chemistry (AREA)
- Biochemistry (AREA)
- Chemical & Material Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Immunology (AREA)
- Pathology (AREA)
- Health & Medical Sciences (AREA)
- Image Processing (AREA)
- Length Measuring Devices By Optical Means (AREA)
Abstract
本发明公开了一种用于光声检测的干涉条纹相位提取方法，首先利用极值—中心取点法确定单幅干涉条纹中需要提取辐照度的五个位置点；然后设置线阵图像传感器的驱动频率和积分时间；利用高速模数转换器获取一幅干涉条纹图并进行中值滤波；提取五个位置点的辐照度并计算干涉条纹相位值；计算前后相邻干涉条纹相位差；对相位差做累加并对累加和结果做低通滤波；连续获取多幅干涉条纹图，依次重复上述过程即可得到反应光声变化情况的连续相位曲线。该方法可自动实现干涉条纹的相位提取，减轻噪声和随机误差的影响，并且算法实时性好，资源开销低。
Description
技术领域
本发明涉及一种用于光声检测的干涉条纹相位提取方法，适用于光声光谱痕量气体检测，也可用于其他应用干涉仪且需要提取干涉条纹相位的场合，属于光学检测领域。
背景技术
光学干涉检测以光学干涉原理为基础，是光学检测中精度和灵敏度较高的检测方法。其原理是利用相互干涉的两束光在感光器件上产生干涉条纹，当测量光束相对于参考光束发生变化时，感光器件上就会产生左右移动的干涉条纹。通过高速高精度采集装置获取变化的干涉条纹图，经过理论计算，提取干涉条纹相位变化，即可得到测量光束的位移量。在光声检测中，按照一定调制规律驱动激光器，并控制其在一定波长范围内扫描，这样就能得到该波长范围内每个扫描波长处的干涉条纹相对位移量，对该位移量做处理，与标准谱库对比，就能得到待测气体中含有的气体种类及浓度信息。由此可见，干涉条纹相位提取方法是光学检测中最基础和重要的部分。
目前已有的干涉条纹相位提取算法可根据提取相位所需干涉条纹图幅数，大致分为单幅和多幅干涉条纹相位提取算法。单幅干涉条纹相位提取算法是指仅需要一幅干涉条纹图即可提取相位信息，适用于动态和对实时性要求较高的测量上，包括同步法、锁相法、傅里叶变换法、正则化相位跟踪法和基于遗传算法的相位分析方法等，同步法通常需要有载波生成器、移相模块、两路乘法器和两路低通滤波器等模块；锁相法需要乘法器、低通滤波器和压控振荡器等模块；傅里叶变换法、正则化相位跟踪法和基于遗传算法的相位分析方法等都更适于在计算机上实现；总之该类算法在可编程逻辑器件上实现时硬件和逻辑资源开销大，性价比低；多幅干涉条纹相位提取算法是指提取相位信息需要多幅干涉条纹图，适用于静态和对实时性要求不高的测量上，主要包括相移法等。
由于光声光谱痕量气体检测要求检测连续，实时性高；考虑成本，对器件资源敏感。已有的干涉条纹相位提取算法均存在不足，无法高效的应用到光声检测中。
发明内容
本发明解决的技术问题是：克服现有技术不足，提供一种用于光声检测的干涉条纹相位提取方法，适用于实时性要求高、对资源敏感的场合，能够自动实现相位检测，避免人为相移误差。
本发明的技术解决方案为：
一种用于光声检测的干涉条纹相位提取方法，包括如下步骤：
(1)利用极值—中心取点法确定单幅干涉条纹中需要提取辐照度的五个位置点；
(2)设置线阵图像传感器的驱动频率、积分时间；
(3)利用高速模数转换器采集线阵图像传感器的输出模拟信号，所述线阵图像传感器的输出模拟信号为单幅干涉条纹的所有像素点模拟量，高速模数转换器对采集的模拟信号进行模数转换后，得到包含一幅干涉条纹各个像素点处辐照度的数组，即一幅数字化的干涉条纹；
(4)从干涉条纹中提取步骤(1)确定的五个位置点处的辐照度，据此计算干涉条纹的相位值；
(5)继续获取下一幅干涉条纹；
(6)重复步骤(4)，得到新获取的干涉条纹的相位值；
(7)判断步骤(6)得到的本次干涉条纹相位值与上次干涉条纹相位值的大小，若本次比上次大π，则将本次结果减去2π后作为本次干涉条纹相位值的最新结果；若本次比上次小π，则将本次结果加上2π后作为本次干涉条纹相位值最新结果，将本次干涉条纹与上次干涉条纹的相位差保存；
(8)对步骤(7)得到的干涉条纹相位差做累加，并对累加和结果做低通滤波；
(9)连续获取多幅干涉条纹，每次获取一幅干涉条纹，即重复步骤(6)-(8)；
(10)对累加和滤波后的曲线进行锁相处理，得到曲线的幅值，该幅值即为对应的光声变化量。
所述步骤(1)中，利用极值—中心取点法确定单幅干涉条纹中五个位置点的具体实现方式如下：
(2.1)以单幅干涉条纹上辐照度最高的像素点为相位π点；
(2.2)相位π点左边最近，且为辐照度极小值的点为相位0点；
(2.3)相位π点右边最近，且为辐照度极小值的点为相位2π点；
(2.4)根据相位0点和相位π点的位置，采用取像素坐标中值方法或取平衡点的方法获得相位π/2点的位置，根据相位2π点和相位π点的位置，采用取像素坐标中值方法或取平衡点的方法获得相位3π/2点的位置。
所述步骤(2.4)中，取像素坐标中值方法适用于相位0点和相位π点间及相位2π点和相位π点间均为奇数个点的情况，具体实现方式如下：
设相位π点的像素坐标为(x3,y3)，相位0点的像素坐标为(x1,y1)，相位2π点为(x5,y5)，则相位π/2点的像素坐标(x2,y2)满足：
相位3π/2点的像素坐标(x4,y4)满足：
所述步骤(2.4)中，取平衡点方法适用于相位0点和相位π点间或相位2π点和相位π点间有偶数个点的情况，具体实现方式如下：
若相位0点和相位π点间有奇数个点，而相位π点和相位2π点之间有偶数个点，则采用取像素坐标中值方法获得相位π/2点的像素坐标，若相位π/2点的辐照度大于相位0点和相位π点间的辐照度中值，则取相位π点和相位2π点中间两个点中右边的点作为相位3π/2点，若相位π/2点的辐照度小于相位0点和相位π点间的辐照度中值，则取相位π点和相位2π点中间两个点中左边的点作为相位3π/2点；
若相位π点和相位2π点间有奇数个点，而相位0点和相位π点之间有偶数个点，则采用取像素坐标中值方法获得相位3π/2点的像素坐标，若相位3π/2点的辐照度大于相位π点和相位2π点间辐照度中值，则取相位0点和相位π点中间两个点中左边的点作为相位π/2点，若相位3π/2点的辐照度小于相位π点和相位2π点间辐照度中值，则取相位0点和相位π点中间两个点中右边的点作为相位π/2点；
若相位0点和相位π点间及相位π点和相位2π点间均为偶数个点，则取相位0点和相位π点间距离辐照度中值最近的点作为相位π/2点，若相位π/2点的辐照度大于相位0点和相位π点间辐照度中值，则取相位π点和相位2π点中间两个点中右边的点作为相位3π/2点，若相位π/2点的辐照度小于相位0点和相位π点间辐照度中值，则取相位π点和相位2π点中间两个点中左边的点作为相位3π/2点。
所述辐照度中值计算方式如下：
设相位0点处的辐照度为I1(x1,y1)，相位π点处的辐照度为I3(x3,y3)，则相位0点和相位π点间辐照度中值Im(xm,ym)满足：
所述步骤(4)中，利用如下公式计算干涉条纹相位值φ(x,y)：
I1(x1,y1)为相位0点的辐照度，I2(x2,y2)为相位π/2点的辐照度，I3(x3,y3)为相位π点的辐照度，I4(x4,y4)为相位3π/2点的辐照度，I5(x5,y5)为相位2π点的辐照度。
五个位置点的辐照度满足如下公式：
I1(x1,y1)＝Idc(x1,y1)+Iac(x1,y1)cos[φ(x,y)]
I2(x2,y2)＝Idc(x2,y2)-Iac(x2,y2)sin[φ(x,y)]
I3(x3,y3)＝Idc(x3,y3)-Iac(x3,y3)cos[φ(x,y)]
I4(x4,y4)＝Idc(x4,y4)+Iac(x4,y4)sin[φ(x,y)]
I5(x5,y5)＝Idc(x5,y5)+Iac(x5,y5)cos[φ(x,y)]
I1(x1,y1)表示相位0点的辐照度，Idc(x1,y1)表示相位0点处辐照度的直流分量，Iac(x1,y1)表示相位0点处辐照度的交流分量幅值，I2(x2,y2)表示相位π/2点的辐照度，Idc(x2,y2)表示相位π/2点处辐照度的直流分量，Iac(x2,y2)表示相位π/2点处辐照度的交流分量幅值，I3(x3,y3)表示相位π点的辐照度，Idc(x3,y3)表示相位π点处辐照度的直流分量，Iac(x3,y3)表示相位π点处辐照度的交流分量幅值，I4(x4,y4)表示为相位3π/2点的辐照度，Idc(x4,y4)表示相位3π/2点处辐照度的直流分量，Iac(x4,y4)表示相位3π/2点处辐照度的交流分量幅值，I5(x5,y5)表示相位2π点的辐照度，Idc(x5,y5)表示相位2π点处辐照度的直流分量，Iac(x5,y5)表示相位2π点处辐照度的交流分量幅值。
高速模数转换器对干涉条纹每个像素点的模拟信号均采集三次，然后取中间值，并存储到FPGA中；相邻干涉条纹数据在FPGA中做乒乓缓存。
所述步骤(8)中，采用二阶IIR滤波器做低通滤波，用x(n)表示当前累加和结果，x(n-1)表示上一次累加和结果，x(n-2)表示再上一次累加和结果，y(n)表示当前滤波结果，y(n-1)表示上一次滤波结果，y(n-2)表示再上一次滤波结果，b1、b2、b3、a1、a2和a3分别表示各个时期累加和及滤波结果前的系数，二阶IIR滤波器公式如下，初始滤波没有的值，按0计算：
a1·y(n)＝b1·x(n)+b2·x(n-1)+b3·x(n-2)-a2·y(n-1)-a3·y(n-2)
二阶IIR滤波器公式中的系数按如下方法选取：
a1′＝1,a2′＝-1.9704,a3′＝0.9736,b1′＝0.0132,b2′＝0,b3′＝-0.0132
a1＝a1′左移14位
a2＝a2′左移14位
a3＝a3′左移14位。
b1＝b1′左移14位
b2＝b2′左移14位
b3＝b3′左移14位
本发明与现有技术相比的优点在于：
(1)本发明针对单幅干涉条纹提出了极值—中心取点法，利用该方法确定五个位置点，得到了单幅干涉条纹相位获取算法，通过连续获取多幅干涉条纹图，对前后相邻干涉条纹相位差做累加并对累加和结果做低通滤波，得到反应光声变化情况的连续相位曲线。本发明可自动实现干涉条纹的相位提取，适于在可编程逻辑器件上运行，实时性好，资源开销低。
(2)本发明对原始干涉条纹进行了取中值处理，对累加和结果进行了滤波处理，从而能够有效降低噪声和随机误差。
附图说明
图1为本发明实现流程图；
图2为本发明利用高速模数转换器采集CMOS线阵图像传感器并滤波后的512个像素干涉条纹图像；
图3为本发明采用极值—中心取点法获取的五点在干涉条纹上的位置；
图4为本发明应用于干涉仪中悬臂梁振动时的相位累加和结果图(采用极值—中心取点法获取五点，计算相位)；
图5为本发明应用于干涉仪中悬臂梁振动时的相位累加和的锁相结果图；
图6为本发明在XC2V3000型FPGA上运行一个完整的相位提取周期耗时情况的算法仿真图；
图7为本发明在XC2V3000型FPGA上实现时的资源开销情况。
具体实施方式
下面结合附图对本发明作进一步说明。
如图1所示，本发明的实现步骤如下：
(1)利用极值—中心取点法确定单幅干涉条纹中需要提取辐照度的五个位置点，具体方法如下：
(1.1)以单幅干涉条纹上辐照度最高的像素点为相位π点；
(1.2)相位π点左边最近，且为辐照度极小值的点为相位0点；
(1.3)相位π点右边最近，且为辐照度极小值的点为相位2π点；
(1.4)根据相位0点和相位π点的位置，采用取像素坐标中值方法或取平衡点的方法获得相位π/2点的位置，根据相位2π点和相位π点的位置，采用取像素坐标中值方法或取平衡点的方法获得相位3π/2点的位置。
取像素坐标中值方法适用于相位0点和相位π点间及相位2π点和相位π点间均为奇数个点的情况，具体实现方式如下(下式中的像素坐标y分量均为0)：
设相位π点的像素坐标为(x3,y3)，相位0点的像素坐标为(x1,y1)，相位2π点为(x5,y5)，则相位π/2点的像素坐标(x2,y2)满足：
相位3π/2点的像素坐标(x4,y4)满足：
取平衡点方法适用于相位0点和相位π点间或相位2π点和相位π点间有偶数个点的情况，具体实现方式如下：
若相位0点和相位π点间有奇数个点，而相位π点和相位2π点之间有偶数个点，则取相位0点和相位π点间像素坐标中值作为相位π/2点的像素坐标，若相位π/2点的辐照度大于相位0点和相位π点间辐照度中值，则取相位π点和相位2π点中间两个点中右边的点作为相位3π/2点，若相位π/2点的辐照度小于相位0点和相位π点间辐照度中值，则取相位π点和相位2π点中间两个点中左边的点作为相位3π/2点。
若相位π点和相位2π点间有奇数个点，而相位0点和相位π点之间有偶数个点，则取相位π点和相位2π点间像素坐标中值作为相位3π/2点的像素坐标，若相位3π/2点的辐照度大于相位π点和相位2π点间辐照度中值，则取相位0点和相位π点中间两个点中左边的点作为相位π/2点，若相位3π/2点的辐照度小于相位π点和相位2π点间辐照度中值，则取相位0点和相位π点中间两个点中右边的点作为相位π/2点。
若相位0点和相位π点间及相位π点和相位2π点间均为偶数个点，则取相位0点和相位π点间距离辐照度中值最近的点作为相位π/2点的像素坐标，若相位π/2点的辐照度大于相位0点和相位π点间辐照度中值，则取相位π点和相位2π点中间两个点中右边的点作为相位3π/2点，若相位π/2点的辐照度小于相位0点和相位π点间辐照度中值，则取相位π点和相位2π点中间两个点中左边的点作为相位3π/2点。
图3展示了在一幅干涉条纹中采用极值—中心取点法获取的五点位置。图中，512点像素位置处对应的辐照度用空心圆圈表示，在五个相位点(0，π/2，π，3π/2和2π)上分别用“*”填充空心圆圈，并在旁边标明。
利用如下公式计算干涉条纹相位值φ(x,y)：
五个位置点的辐照度满足如下公式：
I1(x1,y1)＝Idc(x1,y1)+Iac(x1,y1)cos[φ(x,y)]
I2(x2,y2)＝Idc(x2,y2)-Iac(x2,y2)sin[φ(x,y)]
I3(x3,y3)＝Idc(x3,y3)-Iac(x3,y3)cos[φ(x,y)]
I4(x4,y4)＝Idc(x4,y4)+Iac(x4,y4)sin[φ(x,y)]
I5(x5,y5)＝Idc(x5,y5)+Iac(x5,y5)cos[φ(x,y)]
I1(x1,y1)表示相位0点的辐照度，Idc(x1,y1)表示相位0点处辐照度的直流分量，Iac(x1,y1)表示相位0点处辐照度的交流分量幅值，I2(x2,y2)表示相位π/2点的辐照度，Idc(x2,y2)表示相位π/2点处辐照度的直流分量，Iac(x2,y2)表示相位π/2点处辐照度的交流分量幅值，I3(x3,y3)表示相位π点的辐照度，Idc(x3,y3)表示相位π点处辐照度的直流分量，Iac(x3,y3)表示相位π点处辐照度的交流分量幅值，I4(x4,y4)表示为相位3π/2点的辐照度，Idc(x4,y4)表示相位3π/2点处辐照度的直流分量，Iac(x4,y4)表示相位3π/2点处辐照度的交流分量幅值，I5(x5,y5)表示相位2π点的辐照度，Idc(x5,y5)表示相位2π点处辐照度的直流分量，Iac(x5,y5)表示相位2π点处辐照度的交流分量幅值。
(2)设置线阵图像传感器的驱动频率、积分时间；
以CMOS线阵图像传感器S9227为例，S9227具有512个像素点，最大图像数据速率为5MHz。将其驱动频率设置为2MHz，积分时间设置为101.5us。
(3)利用高速模数转换器采集线阵图像传感器的输出模拟信号，所述线阵图像传感器的输出模拟信号为单幅干涉条纹的所有像素点模拟量，高速模数转换器对采集的模拟信号进行模数转换后，得到包含一幅干涉条纹各个像素点处辐照度的数组，即一幅数字化的干涉条纹；
如高速模数转换器选用AD9240，是正5伏单电源供电，转换位数14位，转换速率达到100纳秒的模数转换器。将AD9240的采样速率设置为6MHz，保证每个像素点可被采集3个数据。对每个像素点得到的采样结果取中值，并存储到FPGA内部开辟的随机存储器(RAM)空间中，该RAM空间在FPGA中开辟两份，每份位宽为14位，深度512个。相邻干涉条纹数据做乒乓缓存。高速模数转换后得到的一幅512点干涉条纹的平滑连线图如图2所示。
采用乒乓缓存方式进行存储是为了提高处理速度，即可编程逻辑器件内部开辟两块存储空间A和B，也可将外部存储器开辟两块存储空间A和B，第一次干涉条纹数据存储在A中，第二次干涉条纹数据存储在B中，第三次干涉条纹数据存储于A中，第四次干涉条纹数据存储于B中，依此规律进行存储操作。
可编程逻辑器件采用现场可编程逻辑门阵列FPGA(也可采用复杂可编程逻辑器件CPLD)，外部存储器采用非易失性的磁性随机存储器MRAM。
(4)从干涉条纹中提取步骤(1)确定的五个位置点处的辐照度，据此计算干涉条纹的相位值；
(5)继续获取下一幅干涉条纹；
(6)重复步骤(4)，得到新获取的干涉条纹的相位值；
(7)判断步骤(6)得到的本次干涉条纹相位值与上次干涉条纹相位值的大小，若本次比上次大π，则将本次结果减去2π后作为本次干涉条纹相位值的最新结果；若本次比上次小π，则将本次结果加上2π后作为本次干涉条纹相位值最新结果，将本次干涉条纹与上次干涉条纹的相位差保存；
(8)对步骤(7)得到的干涉条纹相位差做累加，并对累加和结果做低通滤波；
具体地，采用二阶IIR滤波器做低通滤波，用x(n)表示当前累加和结果，x(n-1)表示上一次累加和结果，x(n-2)表示再上一次累加和结果，y(n)表示当前滤波结果，y(n-1)表示上一次滤波结果，y(n-2)表示再上一次滤波结果，b1、b2、b3分别为当前累加和结果、上一次累加和结果以及再上一次累加和结果前的系数、a1、a2和a3分别表示当前滤波结果、上一次滤波结果以及再上一次滤波结果前的系数，二阶IIR滤波器公式如下，初始滤波没有的值，按0计算：
a1·y(n)＝b1·x(n)+b2·x(n-1)+b3·x(n-2)-a2·y(n-1)-a3·y(n-2)
二阶IIR滤波器公式中的系数按如下方法选取：
a1′＝1,a2′＝-1.9704,a3′＝0.9736,b1′＝0.0132,b2′＝0,b3′＝-0.0132
a1＝a1′左移14位
a2＝a2′左移14位
a3＝a3′左移14位。
b1＝b1′左移14位
b2＝b2′左移14位
b3＝b3′左移14位
(9)连续获取多幅干涉条纹，每次获取一幅干涉条纹，即重复步骤(6)-(8)；
(10)经过步骤(8)后得到的累加和滤波后曲线即为反应光声变化情况的连续相位曲线，其幅值对应光声变化量，需要获取该幅值。对曲线进行锁相处理，即可得到该幅值。
设定光声气体检测仪持续固定扫描波数1035cm-1时，按照上述方法，连续获取多幅干涉条纹图，采用极值—中心取点法获取五点，按照上述步骤计算得到相位累加和结果如图4所示。
对图4所示的相位累加和结果进行锁相处理，得到悬臂梁的振动幅值情况，如图5所示。该相位幅值变化的标准差为11.167，离散度为0.0107，即相位变化约为1％，表明相位数据很稳定。
在XC2V3000型FPGA上对本发明方法进行仿真测试，运行一个完整的相位提取周期耗时情况如图6所示。图6中左边的时间轴位于8106700ns时刻，此时信号“square_result_ok”所代表的相位提取方法刚刚运行完毕，得到最新的计算结果，即左边时间轴位于“square_result_ok”信号的下降沿；右边的时间轴位于15256700ns时刻，此时信号“square_result_ok”所代表的相位提取方法再次运行完毕，又一次得到最新的计算结果，即右边时间轴同样位于“square_result_ok”信号的下降沿。左右时间轴间的时间差7150000ns，即相位提取算法实际运行时间为7.15ms。
在XC2V3000型FPGA(300万门)上实现本发明的一种用于光声检测的干涉条纹相位提取方法时的资源开销情况如图7所示。触发器(Flip-Flop)占用39％，四输入查找表(LUT)总的占用为36％，内部块存储器(Block RAM)占用4％，乘法器(MULT18X18s)占用8％。说明本发明方法资源开销较低。
本发明未详细说明部分属于本领域公知技术。
Claims (10)
1.一种用于光声检测的干涉条纹相位提取方法，其特征在于包括如下步骤：
(1)利用极值—中心取点法确定单幅干涉条纹中需要提取辐照度的五个位置点；
(2)设置线阵图像传感器的驱动频率、积分时间；
(3)利用高速模数转换器采集线阵图像传感器的输出模拟信号，所述线阵图像传感器的输出模拟信号为单幅干涉条纹的所有像素点模拟量，高速模数转换器对采集的模拟信号进行模数转换后，得到包含一幅干涉条纹各个像素点处辐照度的数组，即一幅数字化的干涉条纹；
(4)从干涉条纹中提取步骤(1)确定的五个位置点处的辐照度，据此计算干涉条纹的相位值；
(5)继续获取下一幅干涉条纹；
(6)重复步骤(4)，得到新获取的干涉条纹的相位值；
(7)判断步骤(6)得到的本次干涉条纹相位值与上次干涉条纹相位值的大小，若本次比上次大π，则将本次结果减去2π后作为本次干涉条纹相位值的最新结果；若本次比上次小π，则将本次结果加上2π后作为本次干涉条纹相位值最新结果，将本次干涉条纹与上次干涉条纹的相位差保存；
(8)对步骤(7)得到的干涉条纹相位差做累加，并对累加和结果做低通滤波；
(9)连续获取多幅干涉条纹，每次获取一幅干涉条纹，即重复步骤(6)-(8)；
(10)对累加和滤波后的曲线进行锁相处理，得到曲线的幅值，该幅值即为对应的光声变化量。
2.根据权利要求1所述的一种用于光声检测的干涉条纹相位提取方法，其特征在于：所述步骤(1)中，利用极值—中心取点法确定单幅干涉条纹中五个位置点的具体实现方式如下：
(2.1)以单幅干涉条纹上辐照度最高的像素点为相位π点；
(2.2)相位π点左边最近，且为辐照度极小值的点为相位0点；
(2.3)相位π点右边最近，且为辐照度极小值的点为相位2π点；
(2.4)根据相位0点和相位π点的位置，采用取像素坐标中值方法或取平衡点的方法获得相位π/2点的位置，根据相位2π点和相位π点的位置，采用取像素坐标中值方法或取平衡点的方法获得相位3π/2点的位置。
3.根据权利要求2所述的一种用于光声检测的干涉条纹相位提取方法，其特征在于：所述步骤(2.4)中，取像素坐标中值方法适用于相位0点和相位π点间及相位2π点和相位π点间均为奇数个点的情况，具体实现方式如下：
设相位π点的像素坐标为(x3,y3)，相位0点的像素坐标为(x1,y1)，相位2π点为(x5,y5)，则相位π/2点的像素坐标(x2,y2)满足：
相位3π/2点的像素坐标(x4,y4)满足：
4.根据权利要求2所述的一种用于光声检测的干涉条纹相位提取方法，其特征在于：所述步骤(2.4)中，取平衡点方法适用于相位0点和相位π点间或相位2π点和相位π点间有偶数个点的情况，具体实现方式如下：
若相位0点和相位π点间有奇数个点，而相位π点和相位2π点之间有偶数个点，则采用取像素坐标中值方法获得相位π/2点的像素坐标，若相位π/2点的辐照度大于相位0点和相位π点间的辐照度中值，则取相位π点和相位2π点中间两个点中右边的点作为相位3π/2点，若相位π/2点的辐照度小于相位0点和相位π点间的辐照度中值，则取相位π点和相位2π点中间两个点中左边的点作为相位3π/2点；
若相位π点和相位2π点间有奇数个点，而相位0点和相位π点之间有偶数个点，则采用取像素坐标中值方法获得相位3π/2点的像素坐标，若相位3π/2点的辐照度大于相位π点和相位2π点间辐照度中值，则取相位0点和相位π点中间两个点中左边的点作为相位π/2点，若相位3π/2点的辐照度小于相位π点和相位2π点间辐照度中值，则取相位0点和相位π点中间两个点中右边的点作为相位π/2点；
若相位0点和相位π点间及相位π点和相位2π点间均为偶数个点，则取相位0点和相位π点间距离辐照度中值最近的点作为相位π/2点，若相位π/2点的辐照度大于相位0点和相位π点间辐照度中值，则取相位π点和相位2π点中间两个点中右边的点作为相位3π/2点，若相位π/2点的辐照度小于相位0点和相位π点间辐照度中值，则取相位π点和相位2π点中间两个点中左边的点作为相位3π/2点。
5.根据权利要求4所述的一种用于光声检测的干涉条纹相位提取方法，其特征在于：所述辐照度中值计算方式如下：
设相位0点处的辐照度为I1(x1,y1)，相位π点处的辐照度为I3(x3,y3)，则相位0点和相位π点间辐照度中值Im(xm,ym)满足：
。
6.根据权利要求3或4所述的一种用于光声检测的干涉条纹相位提取方法，其特征在于：所述步骤(4)中，利用如下公式计算干涉条纹相位值φ(x,y)：
I1(x1,y1)为相位0点的辐照度，I2(x2,y2)为相位π/2点的辐照度，I3(x3,y3)为相位π点的辐照度，I4(x4,y4)为相位3π/2点的辐照度，I5(x5,y5)为相位2π点的辐照度。
7.根据权利要求6所述的一种用于光声检测的干涉条纹相位提取方法，其特征在于：五个位置点的辐照度满足如下公式：
I1(x1,y1)＝Idc(x1,y1)+Iac(x1,y1)cos[φ(x,y)]
I2(x2,y2)＝Idc(x2,y2)-Iac(x2,y2)sin[φ(x,y)]
I3(x3,y3)＝Idc(x3,y3)-Iac(x3,y3)cos[φ(x,y)]
I4(x4,y4)＝Idc(x4,y4)+Iac(x4,y4)sin[φ(x,y)]
I5(x5,y5)＝Idc(x5,y5)+Iac(x5,y5)cos[φ(x,y)]
I1(x1,y1)表示相位0点的辐照度，Idc(x1,y1)表示相位0点处辐照度的直流分量，Iac(x1,y1)表示相位0点处辐照度的交流分量幅值，I2(x2,y2)表示相位π/2点的辐照度，Idc(x2,y2)表示相位π/2点处辐照度的直流分量，Iac(x2,y2)表示相位π/2点处辐照度的交流分量幅值，I3(x3,y3)表示相位π点的辐照度，Idc(x3,y3)表示相位π点处辐照度的直流分量，Iac(x3,y3)表示相位π点处辐照度的交流分量幅值，I4(x4,y4)表示为相位3π/2点的辐照度，Idc(x4,y4)表示相位3π/2点处辐照度的直流分量，Iac(x4,y4)表示相位3π/2点处辐照度的交流分量幅值，I5(x5,y5)表示相位2π点的辐照度，Idc(x5,y5)表示相位2π点处辐照度的直流分量，Iac(x5,y5)表示相位2π点处辐照度的交流分量幅值。
8.根据权利要求1所述的一种用于光声检测的干涉条纹相位提取方法，其特征在于：高速模数转换器对干涉条纹每个像素点的模拟信号均采集三次，然后取中间值，并存储到FPGA中；相邻干涉条纹数据在FPGA中做乒乓缓存。
9.根据权利要求1所述的一种用于光声检测的干涉条纹相位提取方法，其特征在于：所述步骤(8)中，采用二阶IIR滤波器做低通滤波，用x(n)表示当前累加和结果，x(n-1)表示上一次累加和结果，x(n-2)表示再上一次累加和结果，y(n)表示当前滤波结果，y(n-1)表示上一次滤波结果，y(n-2)表示再上一次滤波结果，b1、b2、b3、a1、a2和a3分别表示各个时期累加和及滤波结果前的系数，二阶IIR滤波器公式如下，初始滤波没有的值，按0计算：
a1·y(n)＝b1·x(n)+b2·x(n-1)+b3·x(n-2)-a2·y(n-1)-a3·y(n-2) 。
10.根据权利要求9所述的一种用于光声检测的干涉条纹相位提取方法，其特征在于：二阶IIR滤波器公式中的系数按如下方法选取：
a1′＝1,a2′＝-1.9704,a3′＝0.9736,b1′＝0.0132,b2′＝0,b3′＝-0.0132
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910305054.1A CN110108643B (zh) | 2019-04-16 | 2019-04-16 | 一种用于光声检测的干涉条纹相位提取方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910305054.1A CN110108643B (zh) | 2019-04-16 | 2019-04-16 | 一种用于光声检测的干涉条纹相位提取方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110108643A true CN110108643A (zh) | 2019-08-09 |
| CN110108643B CN110108643B (zh) | 2021-12-07 |
Family
ID=67485502
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910305054.1A Active CN110108643B (zh) | 2019-04-16 | 2019-04-16 | 一种用于光声检测的干涉条纹相位提取方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110108643B (zh) |
Citations (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US5479259A (en) * | 1991-05-20 | 1995-12-26 | Hitachi, Ltd. | Method and apparatus for detecting photoacoustic signal |
| CN101140186A (zh) * | 2007-10-19 | 2008-03-12 | 中国科学院上海光学精密机械研究所 | 偏振移相双剪切干涉波面测量仪及其检测方法 |
| CN101639383A (zh) * | 2008-12-25 | 2010-02-03 | 北京理工大学 | 一种主被动结合的移相干涉波面测量方法 |
| CN102306331A (zh) * | 2011-08-17 | 2012-01-04 | 广东外语外贸大学 | 一种干涉条纹智能计数方法 |
| CN102425988A (zh) * | 2011-11-20 | 2012-04-25 | 中国科学院光电技术研究所 | 一种用于移相干涉条纹图的相位提取方法 |
| CN103328921A (zh) * | 2011-01-25 | 2013-09-25 | 麻省理工学院 | 单镜头全视场反射相显微镜 |
| CN104111120A (zh) * | 2014-07-25 | 2014-10-22 | 中国科学院上海光学精密机械研究所 | 基于朗奇剪切干涉仪的相位提取方法 |
| CN107462581A (zh) * | 2016-06-02 | 2017-12-12 | 株式会社岛津制作所 | 缺陷检测方法以及缺陷检测装置 |
| CN107917676A (zh) * | 2017-10-24 | 2018-04-17 | 南京理工大学 | 一种基于条纹图像频谱分析的干涉测量方法 |
| EP3376522A1 (en) * | 2017-03-14 | 2018-09-19 | Technische Universität Berlin | Method and apparatus for carrying out a time-resolved interferometric measurement |
-
2019
- 2019-04-16 CN CN201910305054.1A patent/CN110108643B/zh active Active
Patent Citations (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US5479259A (en) * | 1991-05-20 | 1995-12-26 | Hitachi, Ltd. | Method and apparatus for detecting photoacoustic signal |
| CN101140186A (zh) * | 2007-10-19 | 2008-03-12 | 中国科学院上海光学精密机械研究所 | 偏振移相双剪切干涉波面测量仪及其检测方法 |
| CN101639383A (zh) * | 2008-12-25 | 2010-02-03 | 北京理工大学 | 一种主被动结合的移相干涉波面测量方法 |
| CN103328921A (zh) * | 2011-01-25 | 2013-09-25 | 麻省理工学院 | 单镜头全视场反射相显微镜 |
| CN102306331A (zh) * | 2011-08-17 | 2012-01-04 | 广东外语外贸大学 | 一种干涉条纹智能计数方法 |
| CN102425988A (zh) * | 2011-11-20 | 2012-04-25 | 中国科学院光电技术研究所 | 一种用于移相干涉条纹图的相位提取方法 |
| CN104111120A (zh) * | 2014-07-25 | 2014-10-22 | 中国科学院上海光学精密机械研究所 | 基于朗奇剪切干涉仪的相位提取方法 |
| CN107462581A (zh) * | 2016-06-02 | 2017-12-12 | 株式会社岛津制作所 | 缺陷检测方法以及缺陷检测装置 |
| EP3376522A1 (en) * | 2017-03-14 | 2018-09-19 | Technische Universität Berlin | Method and apparatus for carrying out a time-resolved interferometric measurement |
| CN107917676A (zh) * | 2017-10-24 | 2018-04-17 | 南京理工大学 | 一种基于条纹图像频谱分析的干涉测量方法 |
Non-Patent Citations (5)
| Title |
|---|
| DAVID BADDELEY ET AL.: "4Pi microscopy deconvolution with a variable point-spread function", 《APPLIED OPTICS》 * |
| JING-TAO DONG ET AL.: "A five-point stencil based algorithm used for phase shifting low-coherence interference microscopy", 《OPTICS AND LASERS IN ENGINEERING》 * |
| R ZHU ET AL.: "Advanced carrier squeezing interferometry: a technique to suppress gamma distortion in fringe projection profilometry", 《APPLIED OPTICS》 * |
| 王大鹏等: "基于激光干涉的电容位移传感器非线性误差标定方法", 《激光与光电子学进展》 * |
| 郭媛等: "基于剪切干涉的单幅干涉条纹相位恢复算法", 《中国激光》 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN110108643B (zh) | 2021-12-07 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110596785B (zh) | 一种适用于原子干涉重力仪的振动噪声修正补偿方法及便携式装置 | |
| EP2736403B1 (en) | Apparatus for quantitative phase tomography through linear scanning with coherent and non-coherent detection | |
| Nashed et al. | Distributed automatic differentiation for ptychography | |
| Koniar et al. | Segmentation of motion regions for biomechanical systems | |
| DE19645953A1 (de) | Korrelationsfunktionsmeßverfahren und -gerät, und Verfahren und Gerät zur Sichtbarmachung eines Wellenquellenbildes basierend auf einer Korrelationsfunktionsmessung | |
| Fung et al. | Multigrid optimization for large-scale ptychographic phase retrieval | |
| Huang et al. | Fast digital lossy compression for X-ray ptychographic data | |
| CN110108643A (zh) | 一种用于光声检测的干涉条纹相位提取方法 | |
| CN102620674B (zh) | 一种简易电子散斑干涉的实时相移方法 | |
| Lopez et al. | Monitoring of oil tank filling with spaceborne SAR using coherent scatterers | |
| CN110211091A (zh) | 一种高分辨图像重构方法、装置及隐裂无损检测系统 | |
| US8520216B2 (en) | Shape measuring apparatus | |
| Merc et al. | Asteroseismology of the heartbeat star KIC 5006817 | |
| CN111964586A (zh) | 一种基于随机噪声修正的白光干涉信号处理方法 | |
| Viskum et al. | A search for delta Scuti stars in northern open clusters. I. CCD photometry of NGC 7245, NGC 7062, NGC 7226 and NGC 7654 | |
| CN106482633A (zh) | 一种基于π/4相移的多光束干涉相位提取方法 | |
| Shough | Beyond fringe analysis | |
| CN111829954B (zh) | 一种提高全场扫频光学相干层析测量量程的系统及方法 | |
| Alarousu et al. | Optical coherence tomography evaluation of internal random structure of wood fiber tissue | |
| Styk et al. | Smart pixel camera based signal processing in an interferometric test station for massive parallel inspection of MEMS and MOEMS | |
| Montgomery et al. | High-speed, on-line 4D microscopy using continuously scanning white light interferometry with a high-speed camera and real-time FPGA image processing | |
| Strub et al. | Global analysis of lisa data with galactic binaries and massive black hole binaries | |
| Volynskiĭ et al. | Recursion algorithms for processing video information in optical-coherence-tomography systems | |
| CN109557035A (zh) | 一种正弦相位差分解调的多普勒优化成像方法 | |
| EP4354178A1 (en) | Method and system for surface wear inspection using millimeter wave radar |
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 |