CN104034285A - 整数线性规划搜索法的双频正弦光栅绝对相位解包裹方法 - Google Patents

Abstract

本发明提出一种整数线性规划搜索法的双频正弦光栅绝对相位解包裹方法，根据传统的双频解包裹原理，提出了关于相对相位、条纹级数的整数线性规划法，建立对应的目标函数和约束条件；以变形的正弦光栅图像为计算对象，采用该方法对其进行遍历搜索，可以快速准确地找到条纹级数的整数对，从而达到解包裹的目的。本发明不要求高频与低频之间成一定的比例关系，扩展了方法的使用范围；不分先后顺序，可以实现对高频与低频的同时解包裹，增强了方法的可操作性；不引入任何中间过度量和划分子区间，简化了计算过程，大大减少了计算量，提高了运算效率。

Description

整数线性规划搜索法的双频正弦光栅绝对相位解包裹方法

技术领域

本发明属于光学三维测量相关领域，涉及一种绝对相位展开的时域解包裹方法，尤其是涉及一种整数线性规划搜索法的绝对相位解包裹方法，该方法可以对任意互为质数的双频正弦结构光栅条纹实现快速准确地绝对相位展开。

背景技术

数字光栅投影技术具有非接触、低成本、高精度和高效率等优势，因此被广泛的应用在三维形貌测量中。该技术的关键之一是对变形的光栅条纹图像的快速准确解包裹(绝对相位展开)，即将相位值分布在[-π,π]区间的截断相位(相对相位)展开到[-fπ,fπ](f为光栅的频率)区间的连续相位(绝对相位)。

在实际工程测量中，为了获得可靠性强、精度高的绝对相位，时域解相位方法通常被广泛应用。一般利用光栅的模型和频率的组合实现解包裹。已有的光栅模型：正弦结构光栅、线性结构光栅、彩色结构光栅以及格雷码等；光栅的频率：单频、双频和多频。其中，对双频正弦结构光栅而言，绝对相位获取的数学模型为：

式中：

Φ—连续的相位分布，即绝对相位值；

—截断的相位分布，即相对相位值；

n—截断点处对应的条纹级数，取整数；

(x,y)—图像中水平方向和垂直方向上的像素序列；

下表low,high—低频f_low和高频f_high。

上式中相对相位值通过相移技术易求得，而关键在确定条纹级数n。其传统的解包裹方法分为两步：

第一步：利用截断相位容易展开的低频f_low正弦结构光栅投影，确定出n_low(x,y)。但是测量精度低，不能反映物体型面的细节。

第二步：采用满足测量精度的高频f_high正弦结构光栅投影。利用低频的测量结果，结合高低频光栅的关系来确定出n_high(x,y)，进而计算出高频光栅的绝对相位值Φ_high(x,y)。

该方法要求高频与低频之间成一定的比例关系，而且在计算的过程中，必须先计算出低频的绝对相位信息，然后才能确定出足测量精度的高频绝对相位信息。

相关学者对其做了改进，采用频率互为质数的正弦光栅投影与计算，填补了传统方法的缺陷。但是：

(1)引入了中间过度量Φ₀，即频率为1时的正弦光栅绝对相位值。

(2)增加了对相对相位值分布区间[-π,π]的划分，即利用高低频率将该区间划分成若干个子区间，然后在每一个子区间中可以同时找到高频与低频对应的条纹级数n_low(x,y)和n_high(x,y)，组成一个模板，应用到实际测量中，可同时求解Φ_low(x,y)和Φ_high(x,y)。

但在(2)中，并没有从真实变形的正弦光栅图像中来确定n_low(x,y)和n_high(x,y)，这就忽略了条纹级数在真实变形光栅图像中的实际值，与确定的模板存在一定的差异，这样会给解包裹带来偏差。

发明内容

针对上述技术问题，本发明公开了一种整数线性规划搜索法的双频正弦光栅绝对相位解包裹方法，可以对任意互为质数的双频正弦结构光栅条纹实现快速准确地绝对相位展开。

本发明的技术方案为：

所述一种整数线性规划搜索法的双频正弦光栅绝对相位解包裹方法，其特征在于：采用以下步骤：

步骤1：设计编码，经过N步相移，生成高频f_high、低频f_low互为质数的标准双频正弦光栅条纹图像2N幅，高频正弦光栅条纹图像低频正弦光栅条纹图像各N幅；(x,y)为图像中水平方向和垂直方向上的像素序列；

步骤2：利用投影仪将步骤1中生成的正弦光栅图像和依次投影到被测物体表面，并用摄像机采集被物体表面调制变形的正弦光栅图像和 $I_n^{C\_\mathrm{high}}(x,y);$

步骤3：对步骤2中变形正弦光栅图像和分别进行高斯平滑滤波，并通过建立调制灰度函数，剔除变形的正弦光栅图像中的背景无效像素信息，保留变形正弦光栅图像的有效像素信息；

步骤4：将步骤3处理得到的变形正弦光栅图像代入N步相移公式：

分别计算得到低频f_low和高频f_high对应的相对相位值和N步相移公式中根据计算的低频或高频对应取步骤3处理得到的变形正弦光栅图像；

步骤5：建立关于频率f_low、f_high，相对相位值和条纹级数n_low(x,y)、n_high(x,y)的目标函数以及约束条件的整数线性规划法则：

$\begin{array}{cc}\min & |n_{\mathrm{high}}(x,y)-\frac{f_{\mathrm{high}}}{f_{\mathrm{low}}}{n}_{\mathrm{low}}(x,y)-\frac{\mathrm{INt}}{f_{\mathrm{low}}}|\end{array}$

其中floor()为函数符号，floor(f_low/2)表示不大于f_low/2的最大整数，INt为图像矩阵的整数集合；INt采用以下步骤得到：

步骤5.1：建立变形的正弦光栅图像整数矩阵：

其中round表示四舍五入函数；

步骤5.2：对整数矩阵INT进行统计：INt＝tabulate(INT)，其中tabulate表示统计函数，得到整数矩阵INT中出现的元素值；

步骤6：采用步骤5中整数线性规划法则对步骤4中的相对相位值图像进行遍历搜索，确定出条纹级数n_low(x,y)、n_high(x,y)，并计算出绝对相位值Φ_low(x,y)和Φ_high(x,y)，完成绝对相位展开。

进一步的优选方案：

所述一种整数线性规划搜索法的双频正弦光栅绝对相位解包裹方法，其特征在于：步骤3中通过建立调制灰度函数，剔除变形的正弦光栅图像中的背景无效像素信息，保留变形正弦光栅图像的有效像素信息的过程为：

步骤3.1：建立高低频率的正弦光栅条纹图像的调制灰度函数：

$B^{C\_\mathrm{low}}(x,y)=\frac{2}{N}\sqrt{{\left[\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{I}_n^{C\_\mathrm{low}}(x,y)\sin (2\mathrm{\πn}/N)\right]}^2+{\left[\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{I}_n^{C\_\mathrm{low}}(x,y)\cos (2\mathrm{\πn}/N)\right]}^2}$

$B^{C\_\mathrm{high}}(x,y)=\frac{2}{N}\sqrt{{\left[\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{I}_n^{C\_\mathrm{high}}(x,y)\sin (2\mathrm{\πn}/N)\right]}^2+{\left[\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{I}_n^{C\_\mathrm{high}}(x,y)\cos (2\mathrm{\πn}/N)\right]}^2}$

步骤3.2：遍历整个变形的正弦光栅图像，如果B^C_low(x,y)＞θ且B^C_high(x,y)＞θ，则保留该像素信息，否则剔除该像素信息，θ为阈值。

有益效果

本发明的有益效果与上述技术相比在于：

1)不要求高频与低频之间成一定的比例关系，扩展了方法的使用范围；

2)不分先后顺序，可以实现对高频与低频的同时解包裹，增强了方法的可操作性；

3)不引入任何中间过度量和划分子区间，简化了计算过程，大大减少了计算量，提高了运算效率；

4)基于真实的变形正弦光栅图像来确定条纹级数n_low(x,y)、n_high(x,y)，增强了该方法的科学性，提高了相位精度。

附图说明

图1整数线性规划搜索法的技术路线流程图

图2双频正弦光栅图像

图3变形的双频正弦光栅图像

图4相对相位值曲线图像

其中，(a)相对相位值曲线(第167行)图像；(b)相对相位灰度图像

图5整数线性规划搜索法的模型示意图

图6绝对相位展开图像

其中，(a)绝对相位值曲线(第167行)图像；(b)绝对相位灰度图像；(c)绝对相位灰度的高度图像。

具体实施方式

下面结合附图、表格对本发明作进一步的详细说明。

Step1：生成标准的双频正弦光栅图像

设计编码，经过N步相移，生成高、低频率(低频f_low、高频f_high)互为质数的标准双频正弦光栅条纹图像2N幅，高频、低频各N幅；

上述Step1的具体实现方法是：

1.1选择任意互为质数的一对高低频率f_high、f_low；

1.2基于相移技术和正弦光栅的数学模型，生成标准的双频正弦光栅条纹图像，其高、低频率图像可用光强分布表示为

$I_n^{\mathrm{low}}(x,y)=a_0+b_0\cos (2\mathrm{\π}{f}_{\mathrm{low}}\frac{x}{\mathrm{width}}+\frac{2\mathrm{\π}(n-1)}{N})---\left(2\right)$

$I_n^{\mathrm{high}}(x,y)=a_0+b_0\cos (2\mathrm{\π}{f}_{\mathrm{high}}\frac{x}{\mathrm{width}}+\frac{2\mathrm{\π}(n-1)}{N})---\left(3\right)$

式中：

a₀—光栅图像的时域DC值；b₀—光栅图像的时域AC值；

f_low—光栅的低频率；f_high—光栅的低频率；

x—光栅图像的水平像素序列，x＝1,2,...,width；

width—光栅图像的水平尺寸；

n—第n步相移，n＝1,2,...,N；

N—相移总步数。

本实施例中，设置生成光栅的参数，见表1，并将表1中的参数，分别代入式(2)和式(3)，在Matlab平台上设计编码，生成低频为21，高频为22的20幅标准正弦光栅图像，高低频各10幅，见附图2所示。

表1 标准正弦光栅图像的参数一览表

Step2：投影并采集变形的双频正弦光栅图像

选取DLP LightCrafter^TM 3000投影仪和DMK 21BU04.H CCD摄像机的参数，设置参数见表2。

表2 投影仪与摄像机的参数设置一览表

利用该投影仪将Step1中生成的标准正弦光栅图像和依次投影到喜洋洋石膏像的表面，并用上述的摄像机采集经石膏像表面调制变形的正弦光栅图像 $I_n^{C\_\mathrm{low}}(x,y)$ 和 $I_n^{C\_\mathrm{high}}(x,y),$ 见附图3。

Step3：变形的双频正弦光栅图像的预处理

对Step2中变形正弦光栅图像和分别进行高斯平滑滤波，并建立调制灰度函数，设置一定的阈值范围，剔除变形的正弦光栅图像中背景等无效像素信息，保留变形正弦光栅图像的有效像素信息；

上述Step3的具体实现方法是：

3.1在Matlab平台中，根据高斯平滑算法编写代码，对Step2中变形正弦光栅图像和一一进行高斯平滑滤波；

3.2建立高、低频率的正弦光栅条纹图像的调制灰度函数

$B^{C\_\mathrm{low}}(x,y)=\frac{2}{N}\sqrt{{\left[\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{I}_n^{C\_\mathrm{low}}(x,y)\sin (2\mathrm{\πn}/N)\right]}^2+{\left[\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{I}_n^{C\_\mathrm{low}}(x,y)\cos (2\mathrm{\πn}/N)\right]}^2}---\left(4\right)$

$B^{C\_\mathrm{high}}(x,y)=\frac{2}{N}\sqrt{{\left[\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{I}_n^{C\_\mathrm{high}}(x,y)\sin (2\mathrm{\πn}/N)\right]}^2+{\left[\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{I}_n^{C\_\mathrm{high}}(x,y)\cos (2\mathrm{\πn}/N)\right]}^2}---\left(5\right)$

由式(4)和式(5)分别计算出变形正弦光栅图像的调制灰度值B^C_low(x,y)和B^C_high(x,y)，高频与低频的局部调制灰度值见表3。

表3 高频与低频的局部调制灰度值

3.3设置阈值θ＝3.0，遍历每幅图像，水平像素序列x:1→608，垂直像素序列y:1→684，如果B^C_low(x,y)＞3.0，B^C_high(x,y)＞3.0，则保留该有效像素信息，否则剔除该无效像素信息。

Step4：计算相对相位值

将Step3处理得到的变形正弦光栅图像代入N步相移公式，计算得到低频f_low和高频f_high对应的相对相位值和

N步相移公式为：N步相移公式中根据计算的低频或高频对应取Step3处理得到的变形正弦光栅图像。本实施例中由相移公式(6)，计算出高低频率变形正弦光栅图像对应的相对相位值和见附图4。

Step5：遍历搜索进行绝对相位展开

建立关于频率f_low、f_high，相对相位值和条纹级数n_low(x,y)、n_high(x,y)的目标函数以及约束条件的整数线性规划法则：

$\begin{array}{cc}\min & |n_{\mathrm{high}}(x,y)-\frac{f_{\mathrm{high}}}{f_{\mathrm{low}}}{n}_{\mathrm{low}}(x,y)-\frac{\mathrm{INt}}{f_{\mathrm{low}}}|\end{array}---\left(7\right)$

其中floor()为函数符号，floor(f_low/2)表示不大于f_low/2的最大整数，INt为图像矩阵的整数集合；

上述Step5中的目标函数具体实现方法是：

5.1相位测量中，两种频率分量测得的物体表面同一点的高度h(x,y)是相同的，用公式表示为：

h(x,y)＝K_lowΦ_low(x,y)＝K_highΦ_high(x,y) (9)式中：K_high＝L/2πf_highD；K_low＝L/2πf_lowD；L—摄像机中心到参考平面的距离；D—摄像机中心与投影仪中心的基线长度。

5.2联立式(1)和(9)，可得

上式等号右边是整数运算，结果必为整数；而等号左边理论上也应该是整数，为了使等号左边满足该条件，对其进行四舍五入运算：

式中：round—Matlab中的四舍五入函数。

5.3对于一幅变形的正弦光栅图像而言，遍历操作，INT就为图像尺寸大小(height×width)的一个整数矩阵，统计该数组中出现的各元素，即进行唯一化运算：

INt＝tabulate(INT) (12)

tabulate—Matlab中的统计函数，得到整数矩阵INT中出现的元素值。

5.4联立式(10)和式(12)，整数线性规划的目标函数式(7)。

上述Step5中的约束条件具体实现方法是：

5.5相对相位值的范围为：

5.6联立式(1)和式(13)，采用放大缩小原理得到整数线性规划的约束条件式(8)。

本实施例中应用式(7)、(8)、(11)、(12)对Step4中的相对相位图像进行遍历搜索，建立条纹级数n_low(x,y)、n_high(x,y)和整数集合INt的查找表，见表4。

表4 条纹级数n_low(x,y)、n_high(x,y)和整数集合INt的查找表

将条纹级数代入式(1)中，计算出Φ_low(x,y)和Φ_high(x,y)，即绝对相位展开，见附图6。

Claims (2)

1.一种整数线性规划搜索法的双频正弦光栅绝对相位解包裹方法，其特征在于：采用以下步骤：

步骤1：设计编码，经过N步相移，生成高频f_high、低频f_low互为质数的标准双频正弦光栅条纹图像2N幅，高频正弦光栅条纹图像低频正弦光栅条纹图像各N幅；(x,y)为图像中水平方向和垂直方向上的像素序列；

步骤2：利用投影仪将步骤1中生成的正弦光栅图像和依次投影到被测物体表面，并用摄像机采集被物体表面调制变形的正弦光栅图像和 $I_n^{C\_\mathrm{high}}(x,y);$

步骤3：对步骤2中变形正弦光栅图像和分别进行高斯平滑滤波，并通过建立调制灰度函数，剔除变形的正弦光栅图像中的背景无效像素信息，保留变形正弦光栅图像的有效像素信息；

步骤4：将步骤3处理得到的变形正弦光栅图像代入N步相移公式：

分别计算得到低频f_low和高频f_high对应的相对相位值和N步相移公式中根据计算的低频或高频对应取步骤3处理得到的变形正弦光栅图像；

步骤5：建立关于频率f_low、f_high，相对相位值和条纹级数n_low(x,y)、n_high(x,y)的目标函数以及约束条件的整数线性规划法则：

$\begin{array}{cc}\min & |n_{\mathrm{high}}(x,y)-\frac{f_{\mathrm{high}}}{f_{\mathrm{low}}}{n}_{\mathrm{low}}(x,y)-\frac{\mathrm{INt}}{f_{\mathrm{low}}}|\end{array}$

其中floor()为函数符号，floor(f_low/2)表示不大于f_low/2的最大整数，INt为图像矩阵的整数集合；INt采用以下步骤得到：

步骤5.1：建立变形的正弦光栅图像整数矩阵：

其中round表示四舍五入函数；

步骤5.2：对整数矩阵INT进行统计：INt＝tabulate(INT)，其中tabulate表示统计函数，得到整数矩阵INT中出现的元素值；

步骤6：采用步骤5中整数线性规划法则对步骤4中的相对相位值图像进行遍历搜索，确定出条纹级数n_low(x,y)、n_high(x,y)，并计算出绝对相位值Φ_low(x,y)和Φ_high(x,y)，完成绝对相位展开。

2.根据权利要求1所述一种整数线性规划搜索法的双频正弦光栅绝对相位解包裹方法，其特征在于：步骤3中通过建立调制灰度函数，剔除变形的正弦光栅图像中的背景无效像素信息，保留变形正弦光栅图像的有效像素信息的过程为：

步骤3.1：建立高低频率的正弦光栅条纹图像的调制灰度函数：

$B^{C\_\mathrm{low}}(x,y)=\frac{2}{N}\sqrt{{\left[\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{I}_n^{C\_\mathrm{low}}(x,y)\sin (2\mathrm{\πn}/N)\right]}^2+{\left[\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{I}_n^{C\_\mathrm{low}}(x,y)\cos (2\mathrm{\πn}/N)\right]}^2}$

$B^{C\_\mathrm{low}}(x,y)=\frac{2}{N}\sqrt{{\left[\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{I}_n^{C\_\mathrm{low}}(x,y)\sin (2\mathrm{\πn}/N)\right]}^2+{\left[\underset{n=1}{\overset{N}{\mathrm{\Σ}}}{I}_n^{C\_\mathrm{low}}(x,y)\cos (2\mathrm{\πn}/N)\right]}^2}$

步骤3.2：遍历整个变形的正弦光栅图像，如果B^C_low(x,y)＞θ且B^C_high(x,y)＞θ，则保留该像素信息，否则剔除该像素信息，θ为阈值。