WO2013107076A1

WO2013107076A1 - Adaptive window fourier phase extraction method in optical three-dimensional measurement

Info

Publication number: WO2013107076A1
Application number: PCT/CN2012/071731
Authority: WO
Inventors: 达飞鹏; 王辰星
Original assignee: 东南大学
Priority date: 2012-01-19
Filing date: 2012-02-28
Publication date: 2013-07-25
Also published as: CN102628676B; KR20130119910A; CN102628676A; KR101475382B1

Abstract

An adaptive window Fourier phase extraction method in optical three-dimensional measurement comprises the following steps: (1) projecting a gray-scale sinusoidal fringe pattern on an object to be measured, and a CCD collecting a formed deformed fringe pattern; (2) performing empirical mode decomposition on the deformed fringe pattern by row to obtain a series of intrinsic mode functions and a residual component as a decomposition result for each row; (3) performing Hilbert transform on each intrinsic mode function of signals by row to obtain a corresponding instantaneous frequency and a marginal spectrum, and determining, through analysis of the marginal spectrum, the instantaneous frequency capable of accurately describing change of the fringe pattern; (4) determining a background component of a fringe signal by row; (5) adaptively locating a local stable region of the fringe signal by row using the determined instantaneous frequency; (6) removing the background component from the original fringe pattern by row to obtain a new signal; (7) in each local stable region, performing adaptive window Fourier transform on the new signal by row, to obtain distribution of a wrapped phase of the fringe signal; and (8) unwrapping the wrapped phase to obtain an absolute phase and converting the absolute phase into height information, thereby implementing three-dimensional measurement of the object. The phase extraction method has strong adaptability and high robustness, thereby greatly improving the measurement precision of an object having a complex surface or abrupt edges.

Description

一种光学三维测量中的自适应窗口傅里叶相位提取法技术领域 Adaptive window Fourier phase extraction method in optical three-dimensional measurement

本发明属于三维信息重构的领域。基于灰度正弦光栅投影，采用希尔伯特-黄变换法分析变形条纹图像并求得条纹的瞬时频率及背景分量，用所求的瞬时频率自适应计算窗口傅里叶变换的窗口尺度因子，通过对去除背景的条纹信号直接做傅里叶变换或自适应窗口傅里叶变换从而得到精确的包裹相位，这种方法有利于提高对形貌复杂物体或携带陡峭边缘物体的测量精度。 The invention belongs to the field of three-dimensional information reconstruction. Based on the gray scale sinusoidal grating projection, the deformation fringe image is analyzed by Hilbert-Huang transform method and the instantaneous frequency and background component of the fringe are obtained. The window scale factor of the window Fourier transform is adaptively calculated by the instantaneous frequency. By directly performing Fourier transform or adaptive window Fourier transform on the stripe signal of the background to obtain accurate wrap phase, this method is beneficial to improve the measurement accuracy of complex objects with steep shapes or objects with steep edges.

背景技术 Background technique

光学三维测量技术能够准确获得物体的三维面形数据，可用于三维模型重建、物体表面轮廓测量、工业环境中的尺寸和形位参数的检测等，因此它在虚拟现实、影视特技、医学整形和美容、工业产品的外观设计、艺术雕塑和文物保护等领域都有广阔的应用前景。 Optical 3D measurement technology can accurately obtain 3D shape data of objects, which can be used for 3D model reconstruction, object surface contour measurement, size and shape parameter detection in industrial environment, etc., so it is in virtual reality, film and television special effects, medical plastic surgery and Cosmetics, industrial products, design, art sculpture and cultural relics protection have broad application prospects.

光栅投影法是一种重要的三维测量技术，通过向物体表面投射正弦光栅，将物体的高度信息以相位的形式隐含在光栅中，利用 CCD获得物体表面的光栅条纹图像，并使用特定算法对条纹图像进行处理，提取其中的相位，从而建立物体的三维信息。 Grating projection method is an important three-dimensional measurement technique. By projecting a sinusoidal grating onto the surface of the object, the height information of the object is implicitly contained in the grating in the form of phase. The raster image of the surface of the object is obtained by CCD, and a specific algorithm is used. The fringe image is processed to extract the phase therein, thereby establishing three-dimensional information of the object.

常用的求解条纹图像相位的方法有基于时间域的方法和基于变换域的方法。基于变换域的方法由于只需采集一幅变形条纹图即可完成相位的测量，有利于实现物体的动态测量，因此而被广泛研究和应用，其中包括傅里叶变换法、小波变换法、 S变换法、窗口傅里叶变换法等。傅里叶变换是一种全局的信号分析工具，不能提取局部信号的特征，只适用于平稳信号。近年来，各种时频分析技术被广泛研究以期能够准确获取条纹图像的细节相位信息。连续小波变换由于具有多分辨率分析的特性被引入光学三维测量领域，通过检测小波变换最大脊处的相位，得到变形条纹图的包裹相位。但是这种方法有个前提条件，即被检测的相位必须是近似线性变化且变化缓慢，否则该方法的理论推导是不成立的。同时，小波母函数及相关参数的选取没有成熟的理论依据，也为小波变换轮廓术的广泛应用带来了限制。 S变换与小波变换求解条纹图包裹相位的原理高度相似，故它也不可避免的受到与小波变换同样的条件限制。窗口傅里叶变换也具有较好的时频分析性能，但窗口尺寸即窗口尺度因子的自适应选取一直是研究的热点。现有的选取窗口尺度因子的方法大都是先通过连续小波变换或 S变换检测最大脊处的尺度因子，将该尺度因子直接作为窗口傅里叶的尺度因子，或将该尺度因子取倒数作为条纹信号的瞬时频率，再通过相关算法计算窗口的尺度因子。这些方法需要事先设定基函数及经验值，因此自适应性不好，同时也同样受到被检测相位前提条件的限制，使求解的尺度因子或瞬时频率过平滑化，而不能很好的描述条纹信号的局部特征，从而导致这些方法只能测量表面变化相对平滑的物体，而对表面变化复杂或含有陡峭边缘物体的测量会受到较大限制。发明内容 Common methods for solving the phase of a fringe image are a time domain based method and a transform domain based method. The transform domain-based method can complete the phase measurement by collecting only one deformed fringe pattern, which is beneficial to realize the dynamic measurement of the object. Therefore, it is widely studied and applied, including Fourier transform method, wavelet transform method, S Transformation method, window Fourier transform method, etc. Fourier transform is a global signal analysis tool that cannot extract the characteristics of local signals and is only suitable for stationary signals. In recent years, various time-frequency analysis techniques have been extensively studied in order to accurately acquire detailed phase information of a fringe image. The continuous wavelet transform is introduced into the field of optical three-dimensional measurement due to the characteristics of multi-resolution analysis. By detecting the phase of the maximum ridge of the wavelet transform, the wrap phase of the deformed fringe pattern is obtained. However, this method has a precondition that the phase to be detected must be approximately linear and change slowly, otherwise the theoretical derivation of the method is not valid. At the same time, the selection of wavelet master functions and related parameters has no mature theoretical basis, and it also imposes limitations on the wide application of wavelet transform profilometry. The principle of S-transformation and wavelet transform for solving the phase of the fringe pattern is highly similar, so it is inevitably limited by the same conditions as the wavelet transform. The window Fourier transform also has better time-frequency analysis performance, but the adaptive selection of window size, ie window scale factor, has been a hot topic. The existing methods for selecting the window scale factor are mostly to first measure the scale factor at the maximum ridge by continuous wavelet transform or S transform, and directly use the scale factor as the scale factor of the window Fourier, or take the scale factor as the fringe as the stripe. The instantaneous frequency of the signal, and then the scale factor of the window is calculated by the correlation algorithm. These methods require a base function and an empirical value to be set in advance, so The adaptability is not good, and it is also limited by the detected phase preconditions, so that the scale factor or instantaneous frequency of the solution is over-smoothed, and the local features of the fringe signal cannot be well described, so that these methods can only measure surface changes. Relatively smooth objects, and measurements on surfaces with complex changes or objects with steep edges are subject to greater limitations. Summary of the invention

技术问题：针对光学三维测量中窗口傅里叶变换求取条纹图包裹相位的精确性和自适应性等问题，本发明在自适应进行窗口傅里叶变换及提高包裹相位的求解精度等方面给出了解决办法。本方法借助希尔伯特-黄变换，精确且自适应地求解能准确描述条纹图变化情况的瞬时频率，从而自适应计算窗口傅里叶变换的窗口尺度，同时，在不需要额外计算的情况下可以有效去除条纹图的背景分量从而大大减少窗口傅里叶变换处理过程中零频频谱的干扰。本方法能大大提高对复杂面形或含有陡峭边缘物体的测量精度。 Technical Problem: For the accuracy and adaptability of the fringe pattern wrapping phase in the window Fourier transform in optical three-dimensional measurement, the present invention gives an adaptive window luminance Fourier transform and improved package phase accuracy. A solution has been made. The method uses the Hilbert-Huang transform to accurately and adaptively solve the instantaneous frequency which can accurately describe the variation of the fringe pattern, thereby adaptively calculating the window scale of the window Fourier transform, and at the same time, without additional calculations. The background component of the fringe pattern can be effectively removed to greatly reduce the interference of the zero-frequency spectrum during the window Fourier transform process. This method can greatly improve the measurement accuracy of complex shapes or objects with steep edges.

技术方案：一种光学三维测量中的自适应窗口傅里叶相位提取法，步骤如下： Technical Solution: An adaptive window Fourier phase extraction method in optical three-dimensional measurement, the steps are as follows:

步骤 1: 将灰度正弦条纹图投影到被测物体表面，用 CCD对被测物体表面进行拍摄，得到一幅宽度为^ 高度为 /的变形条纹图像 g(y,

+ φ(γ, x)]，其中，是背景分量，是物体表面反射率， /。是正弦条纹基频， ^(y,x)是待求的包裹相位分布，（y, c)表示变形条纹图像中各像素点的二维坐标，取值范围分别为： l< y<l, l<x<c , β()^)( )8[2;Γ/。;《; + (}Μ)]为基频分量，此处将每个像素点看作一个信号； Step 1: Project the gray sine fringe pattern onto the surface of the measured object, and use the CCD to shoot the surface of the measured object to obtain a deformed fringe image g (y, with a width of ^).

+ φ(γ, x)], where is the background component, the surface reflectivity of the object, /. Is the sinusoidal fringe fundamental frequency, ^(y, x) is the parcel phase distribution to be sought, and (y, c) is the two-dimensional coordinate of each pixel in the deformed fringe image, and the value ranges are: l<y<l,l<x<c , β()^)( )8[2;Γ/. ; "; + (}Μ)] is the fundamental frequency component, where each pixel is treated as a signal;

步骤 2: 令 y = l， n = l, 用经验模式分解法即 EMD对进行分解，方法如下：步骤 2.1: 将记为一行信号 g )，其中 X仍然满足 l≤ c≤c，该行信号的相位为 φ(χ、 , 寻找 g )的极大值点和极小值点，分别采用公知的三次样条插值法对这些极大值点和极小值点进行插值，然后将这些值连接得到极大值包络线 g_max (； c)和极小值包络线 Step 2: Let y = l, n = l, decompose using the empirical mode decomposition method, ie EMD pair, as follows: Step 2.1: will be recorded as a line of signal g), where X still satisfies l ≤ c ≤ c, the line signal The phase is the maximum and minimum points of φ(χ, , looking for g), and these maxima and minima points are interpolated by well-known cubic spline interpolation, respectively, and then these values are connected. Obtain the maximal envelope g _max (; c) and the minimum envelope

8^ ^χ)；步骤 2.2: 用原始信号减去极大值包络线 g_max ( 与极小值包络线 ^皿 ( c)的平均值，得到 h(x) = g( )-g_max ) + g_mn ) . 步骤 2.3: 分别计算的平均幅度 mi¾z/i_izm/7( )和包络幅度 em _izm/7( ) 8^ ^χ ); Step 2.2: Subtract the maximum value envelope g _max (with the minimum value of the envelope ^ c (c) from the original signal to obtain h(x) = g( )-g _m ax ) + g _m n ) . Step 2.3: Calculated average amplitude mi3⁄4z/i_izm/7( ) and envelope amplitude em _izm/7( )

mean _ amp (x)

Mean _ amp (x)

步骤 2.4:如果 A )同时满足以下三个条件，则得到一个本征模函数 IMF,记为 ( ， K_Sn(x) = h(x)，同时将 s_n(x)从原信号 g(x)中分离，得到新的信号 = g(x)-s_n(x)，否则，直接将从原信号 g )中分离，得到新的信号 g'( ) = g( )-/z( )，所述的三个条件是： Step 2.4: If A) satisfies the following three conditions at the same time, an eigenmode function IMF is obtained, denoted as ( , K _Sn (x) = h(x), and s _n (x) from the original signal g(x) In the separation, a new signal = g(x)-s _n (x) is obtained, otherwise, it will be directly separated from the original signal g) to obtain a new signal g'( ) = g( )-/z( ), The three conditions stated are:

条件 1： mean _amp(x) / env _amp(x) < 0.5，条件 2: 满足不等式《¾6¾^_" ( /6^_" (； <0.5的像素个数占同一整行像素总数的比例小于 0.05， Condition 1: mean _amp(x) / env _amp(x) < 0.5, Condition 2: Satisfy the inequality "3⁄463⁄4^_" ( /6^_" (; <0.5 The number of pixels in the total number of pixels in the entire row is less than 0.05,

条件 3: 极大值和极小值的个数之和与过零点的个数相等或至多相差 1个，步骤 2.5：若步骤 2.4得到的 g \x)也同时满足以上三个条件，令 g ( ) = g

, η = η + 1 , 返回步骤 2.1; 否则分解停止，并将记为得到最后的分解结果为： g(x) = _js_n(x) + res(x)，其中 /1为一个 IMF即的序数， <η<Ν , Λ/为 IMF的总个数； Condition 3: The sum of the maximum value and the minimum value is equal to or at most one difference from the zero crossing point. Step 2.5: If g \x) obtained in step 2.4 also satisfies the above three conditions, let g ( ) = g

, η = η + 1 , return to step 2.1; otherwise the decomposition stops, and will be recorded as the final decomposition result: g(x) = _j s _n (x) + res(x), where /1 is an IMF The ordinal number, <η<Ν , Λ/ is the total number of IMFs;

步骤 3: 通过希尔伯特-黄变换，确定能够准确描述一行条纹信号变化规律的瞬时频率，具体过程如下： Step 3: Determine the instantaneous frequency that accurately describes the variation of a stripe signal by Hilbert-Huang transform. The specific process is as follows:

步骤 3.1: 对第个 11\1「即 ( 做希尔伯特变换，得到： )=^)*丄=丄厂^^，其中 " *"为卷积运算符， Γ为积分变量，为希尔伯特变换的结果，分别构造 W的解析信号 ζ„( )： z_n(x) = s ) + = λ_η(χ)β^{ιφ χ)} , 其中 i 为虚部单位， λ (χ) = _η ²(χ) + 为该解析信号 ζ„( ) 的模值， φ_η (X) = arctan(^s» ^{χ/⁾ , )为解析信号 ζ_η (χ)的相位；步骤 3.2: 对第 w个 IMF求其瞬时频率/ Step 3.1: For the first 11\1 "that is (do Hilbert transform, get: ) = ^) * 丄 = 丄 factory ^ ^, where "*" is the convolution operator, Γ is the integral variable, for As a result of the Albert transform, the analytical signal of W is constructed separately ζ„( ): z _n (x) = s ) + = λ _η (χ)β ^{ιφ χ)} , Where i is the imaginary unit, λ (χ) = _η ² (χ) + is the modulus of the analytical signal ζ„( ), φ _η (X) = arctan( ^s » ^{χ / ⁾ , ) is the analytical signalζ Phase of _η (χ); Step 3.2: Find the instantaneous frequency of the wth IMF /

1 άφ_η (χ) 1 άφ _η (χ)

2π dx 2π dx

对第 w个 IMF求其边际谱 M„ (/)

Find the marginal spectrum M„ (/) for the wth IMF

步骤 3.3: 确定含基频分量最多的 IMF的序数 K: Step 3.3: Determine the ordinal number of the IMF with the most fundamental component K:

K = E{\f_n ^M--f₀\}-, 其中 /„^Mm"为第《个 IMF的边际谱最大值 M_max所对应的频率值， E{ }为取 I f_n ^M -/。 I的最小值所对应的 IMF的序数， K = E{\f _n ^M --f ₀ \}-, where /„ ^Mm " is the frequency value corresponding to the maximum marginal M _{max of} the IMF, and E{ } is I f _n ^M −/ . The ordinal number of the IMF corresponding to the minimum value of I,

如果 IMF所对应的 l/„^M™ -f₀1最小值个数为 1个以上，则从具有 l/„^M™ -f₀1最小值的 If the minimum number of l/„ ^M TM −f ₀ 1 corresponding to the IMF is one or more, then the minimum value of l/„ ^M TM −f ₀ 1 is obtained.

IMF中选择出具有最大的边际谱最大值 M_max的 IMF, 并以具有最大的边际谱最大值 M_max 的 IMF所对应的序数为所求的值， The IMF with the largest marginal spectral maximum value M _{max is selected} in the IMF, and the ordinal corresponding to the IMF having the largest marginal spectral maximum value M _max is the obtained value.

选取 f_K 作为能够准确描述该行条纹信号变化规律的瞬时频率； Select f _K as the instantaneous frequency that can accurately describe the variation law of the stripe signal of the line;

步骤 4:确定一行条纹信号的背景分量，具体过程如下：根据步骤 3.2已得的每个 IMF 的瞬时频率，求出每个 IMF的瞬时频率均值，并从中找到最小的瞬时频率均值，然后确定最小的瞬时频率均值所对应的 IMF序数，记为 K_b, 则以步骤 2.5所得的第个 IMF 到第 N个 IMF以及残余分量之和即 t ( ) + r^( )为该行条纹信号的背景分量组合；步骤 5: 根据一个局部平稳区域内的瞬时频率最大值小于 2倍最小值这一属性，用步骤 3.3确定的瞬时频率 f_K(x)自适应定位一行条纹信号的局部平稳区域，具体过程如下：步骤 5.1: 瞬时频率向量 ^( c)为 [/^(1)， f_K(2 …… f_K(c)], 其中元素全部为不小 .(2) Step 4: Determine the background component of the stripe signal. The specific process is as follows: According to the instantaneous frequency of each IMF obtained in step 3.2, find the instantaneous frequency mean of each IMF, and find the minimum instantaneous frequency mean, and then determine the minimum. The IMF ordinal corresponding to the instantaneous frequency mean, denoted as K _b , then the sum of the first IMF to the Nth IMF and the residual component obtained in step 2.5, ie t ( ) + r^( ), is the background of the stripe signal of the line. Component combination; Step 5: According to the property that the maximum value of the instantaneous frequency in a local stationary region is less than 2 times the minimum value, adaptively locate the local stationary region of the stripe signal with the instantaneous frequency f _K (x) determined in step 3.3. The process is as follows: Step 5.1: The instantaneous frequency vector ^(c) is [/^(1), f _K (2 ...... f _K (c)], where all the elements are not small. (2)

于 0的瞬时频率，对^ ( c)求转置向量，得到尸 W为 fAO 步骤 5.2: 令/ ^TW中每个元素分别减去整个 2x W向量，得到一个 _Cxc的方阵 At the instantaneous frequency of 0, find the transpose vector for ^(c), and get the corpse W as fAO. Step 5.2: Decrease the entire 2x W vector for each element in / ^T W to get a square matrix of _C xc

令方阵内所有负值为 0，并找到方阵中全部元素为 0的零对角线，得到： Let all the negative values in the square matrix be 0, and find the zero diagonal of all the elements in the square matrix, and get:

F = ；

F = ;

步骤 5.3:对方阵 F建立坐标系，坐标原点为 F左上角第一个元素，其坐标值记为 (1,1)，横坐标的坐标方向为方阵的行方向，横坐标取值范围为 l ~ c，纵坐标的坐标方向为方阵的列方向，纵坐标的取值范围仍为 l~c， Step 5.3: The coordinate system is established by the partner array F. The coordinate origin is the first element in the upper left corner of F. The coordinate value is recorded as (1, 1), the coordinate direction of the abscissa is the row direction of the square matrix, and the horizontal coordinate value range. For l ~ c, the coordinate direction of the ordinate is the column direction of the square matrix, and the value of the ordinate is still l~c.

在方阵 F中寻找出以方阵 F的零对角线的一部分或全部为对角线的最大全零子方阵，从方阵 f 的坐标原点开始，沿零对角线方向，依次将这些最大全 0子方阵的零对角线上最后一个元素相对于方阵 f的坐标记录下来，依次用（_Α,_Ρι)， (ρ₂,ρ₂) , ……， (p_r,p_r Find a maximum all-zero square matrix with a part or all of the zero diagonal of the square matrix F as the diagonal in the square matrix F, starting from the coordinate origin of the square matrix f, along the zero diagonal direction, in turn The last element of the zero-diagonal line of these largest all-zero square matrices is recorded relative to the coordinates of the square matrix f, followed by ( _Α , _Ρι ), (ρ ₂ , ρ ₂ ), ..., (p _r , p _r

(c, c)表示，其中（_Ρι,_Ρι)为第一个最大全 0子方阵的零对角线上最后一个元素相对于方阵 f 的坐标，依次排列，（Α,/ 为倒数第二个最大全 0子方阵的零对角线上最后一个元素相对于方阵 F的坐标，（_C,c)为最后一个最大全 0子方阵的零对角线上最后一个元素相对于方阵 _f的坐标，最后，该行条纹信号的平稳区域划分情况为：（1,……， A)，（A+1,……， Ρ₂)， ……， (p_r + ， c); 步骤 6:将原条纹信号 g )减去步骤 4中确定的背景分量，得到步骤 1中去除后的基频分量 g。 (x)，即为 g。 (x)=B(x, y) cos[2^/^ + (χ, y)]；步骤 7: 判断步骤 5确定的局部平稳区域个数是否为 1，若 "是"，则直接对 g。( c)做快速傅里叶变换得到傅里叶频谱，记作 F_f，若 "否"，则对 g。( c)进行自适应高斯窗口傅里叶变换，具体过程为：步骤 7.1: 对 g。( c)做自适应高斯窗口傅里叶变换： (c, c) means that ( _Ρι , _Ρι ) is the coordinate of the last element on the zero diagonal of the first largest all-zero square matrix relative to the square matrix f, in order, (Α, / is the last number The coordinates of the last element on the zero-diagonal line of the two largest all-zero square matrices relative to the square matrix F, ( _C , c) being the last element on the zero-diagonal line of the last largest all-zero square matrix relative to The coordinates of the square matrix _f , finally, the smooth region division of the stripe signal of the line is: (1, ..., A), (A+1, ..., Ρ ₂ ), ..., (p _r + , c) Step 6: Subtract the original fringe signal g) from the background component determined in step 4 to obtain the removed fundamental frequency component g in step 1. (x), which is g. (x)=B(x, y) cos[2^/^ + (χ, y)]; Step 7: Determine whether the number of local stationary regions determined in step 5 is 1, if "yes", then directly to g . (c) Do a fast Fourier transform to obtain the Fourier spectrum, denoted as F _f , and if "No", then for g. (c) Perform an adaptive Gaussian window Fourier transform, the specific process is: Step 7.1: For g. (c) Do an adaptive Gaussian window Fourier transform:

WF_a ^b[g₀ (x)] =厂 [ ( ^»邛(_ ] ， WF _a ^b [g ₀ (x)] = factory [ ( ^»邛(_ ] ,

J J

其中 b为平移因子， b依次取值为 1、 2、 3、 …、 c α为每一个局部平稳区域内每个像素点对应的高斯窗口的尺度因子， a = U 4 , 其中 _为相应的局部平稳区域的长度值，单位为像素， )²]为高斯窗函数，Where b is the translation factor, and b is taken as 1, 2, 3, ..., c α as the scale factor of the Gaussian window corresponding to each pixel in each local stationary region, a = U 4 , where _ is the corresponding The length value of the local stationary region, in pixels, ² ] is a Gaussian window function.

对一行条纹信号所有的信号依次进行自适应高斯窗口傅里叶变换后，得到由 WF_fl ^fc[g。(x)]构成的一个二维复数方阵，大小为 c x c，二维复数方阵每一行的元素为每一个窗口内信号的频谱，二维复数方阵一共有 c行，代表 b从 1取值到 c即一共有 c个窗口的信号的频谱； After all the signals of one line of the stripe signal are sequentially subjected to adaptive Gaussian window Fourier transform, the result is obtained by WF _fl ^fc [g. (x)] A two-dimensional complex square matrix, the size of which is cxc, the element of each row of the two-dimensional complex square matrix is the spectrum of the signal in each window, and the two-dimensional complex matrix has a total of c rows, representing b from 1 The value to c is the spectrum of the signal of a total of c windows;

步骤 7.2: 将步骤 7.1所得的二维复数方阵按列叠加，得到总频谱，记作 f_/; Step 7.2: The two-dimensional complex square matrix obtained in step 7.1 is superimposed in columns to obtain a total spectrum, which is denoted as f _/;

步骤 8: 根据 /^确定基频频谱的范围，并提取出来，记为 F_{f 0 ;} 对求傅里叶逆变换，并根据对 f/o求傅里叶逆变换所得的结果，求得相位角即包裹在 [-π, π]之间的相位，令) /=y+l，返回步骤 2，若>^ > /，则进入步骤 9; 步骤 9: 对包裹相位进行展开得到绝对相位，根据经典光栅投影的相位到高度的转换公式，最终求得测量物体的三维信息。有益效果： Step 8: Determine the range of the fundamental frequency spectrum according to /^, and extract it, denoted as F _{f 0 ;} to find the inverse Fourier transform, and obtain the phase based on the inverse Fourier transform of f/o The angle is wrapped in the phase between [-π, π], let) /=y+l, return to step 2, if >^ > /, go to step 9; Step 9: Unfold the parcel phase to get the absolute phase, According to the phase-to-height conversion formula of the classical grating projection, the three-dimensional information of the measured object is finally obtained. Beneficial effects:

相比相移法，本方法不会受到设备关于正弦性、相移的精度和速度等条件限制，且仅需要一幅变形条纹图即可完成包裹相位的提取，有利于实现动态测量。与其他基于变换域的现有技术相比，本发明具有如下优点： Compared with the phase shift method, the method is not limited by the sinusoidality, phase shift accuracy and speed of the device, and only requires a deformed fringe pattern to complete the extraction of the package phase, which is beneficial for dynamic measurement. Compared with other prior art based on transform domain, the present invention has the following advantages:

首先，本发明采用希尔伯特-黄变换法求解条纹信号的瞬时频率，相比于现有常用的小波变换最大脊法或者 S变换最大脊法，避免了被测相位必须线性逼近且变化缓慢的条件限制，得到的瞬时频率更加真实的描述了变形条纹信号的变化，通过本方法所得的瞬时频率确定的窗口尺度因子能真实的根据信号的局部特征进行自适应变化，从而达到窗口傅里叶精确提取局部细节相位的目的； Firstly, the present invention uses the Hilbert-Huang transform method to solve the instantaneous frequency of the fringe signal. Compared with the conventional wavelet transform maximum ridge method or S transform maximum ridge method, the measured phase must be linearly approximated and slowly changed. The conditional limit, the obtained instantaneous frequency more realistically describes the change of the deformed fringe signal, and the window scale factor determined by the instantaneous frequency obtained by the method can be adaptively changed according to the local characteristics of the signal, thereby achieving the window Fourier The purpose of accurately extracting the local detail phase;

其次，本发明提出的自适应定位局部平稳区域的方法，不需要事先确定任何经验值，也不需要迭代计算，效率高、速度快，大大提高了整个自适应窗口傅里叶相位提取的处理效率；最后，本发明在求解条纹图的瞬时频率时，可以很方便的 "顺便 "去除条纹图的背景分量，大大减少了提取基频分量时背景分量所带来的零频频谱干扰，当确定的窗口尺度因子较小时，这一干扰将会非常严重，因此本发明这一处理会有很大的优势。 Secondly, the method for adaptively positioning the local stationary region proposed by the invention does not need to determine any empirical value in advance, and does not need iterative calculation, and has high efficiency and high speed, and greatly improves the entire adaptive window Fourier phase extraction. Finally, the present invention can conveniently remove the background component of the fringe pattern "by the way" when solving the instantaneous frequency of the fringe pattern, and greatly reduce the zero-frequency spectrum interference caused by the background component when extracting the fundamental frequency component. When the determined window scale factor is small, the interference will be very serious, so the process of the present invention has great advantages.

综上，本发明仅用一幅变形条纹图，便能够准确地获取被测物体的包裹相位，且克服了传统方法对具有复杂面形或陡峭边缘的物体测量不高的缺点。 In summary, the present invention can accurately acquire the wrapping phase of the object to be measured by using only one deformed fringe pattern, and overcomes the disadvantage that the conventional method does not measure the object having a complicated surface or a steep edge.

附图说明 DRAWINGS

图 1是本发明整个过程的流程图。 BRIEF DESCRIPTION OF THE DRAWINGS Figure 1 is a flow chart of the entire process of the present invention.

图 2是 CCD采集到的以塑料泡沫为被测物体的变形条纹图像。 Figure 2 is a deformed fringe image of a measured object with a plastic foam as the CCD.

图 3是步骤 2中对一行条纹信号进行经验模式分解的具体过程流程图。 FIG. 3 is a specific process flow chart for performing empirical mode decomposition on a line of stripe signals in step 2.

图 4是图 2中直线所代表的一行条纹信号及其 EMD分解结果。 Figure 4 is a line of stripe signals represented by the straight line in Figure 2 and its EMD decomposition results.

图 5是图 4中相应 IMF的瞬时频率。 Figure 5 is the instantaneous frequency of the corresponding IMF in Figure 4.

图 6是图 4中相应 IMF的边际谱。 Figure 6 is the marginal spectrum of the corresponding IMF in Figure 4.

图 7是步骤 5的示意图。 Figure 7 is a schematic diagram of step 5.

图 8是对图 2中直线所在行条纹信号所定位的局部平稳区域。 Figure 8 is a partial plateau region positioned for the stripe signal of the line in Figure 2.

图 9是去掉背景分量的变形条纹图。 Figure 9 is a deformed fringe diagram with background components removed.

图 10是所求的尺度因子分布图。 Figure 10 is a scale factor distribution diagram obtained.

图 11是所求的包裹相位图。 Figure 11 is a plot of the package phase sought.

图 12是最终恢复的相位图。 Figure 12 is a phase diagram of the final recovery.

具体实施方式 detailed description

下面结合附图对本发明的具体实施方式做进一步说明。在 Windows操作***下选用 VC++6.0作为编程工具，对 CCD采集到的变形条纹图像进行处理。该实例采用塑料泡沫作为被测物体，最终得到比较精确的含有三维信息的全场包裹相位分布和。 The specific embodiments of the present invention are further described below in conjunction with the accompanying drawings. In the Windows operating system, VC++6.0 is used as a programming tool to process the deformed fringe image collected by the CCD. This example uses plastic foam as the object to be measured, and finally obtains a more accurate full-field package phase distribution containing three-dimensional information.

图 1是本发明的整个过程的流程图。 BRIEF DESCRIPTION OF THE DRAWINGS Figure 1 is a flow chart of the entire process of the present invention.

光学三维测量中条纹图进行包裹相位提取时，针对现有的时频分析方法对表面复杂或含有陡峭边缘的物体的测量精度较低等问题，本发明提出基于希尔伯特-黄变换的多尺度窗口傅里叶相位提取法。首先采用经验模态分解法逐行分解采集到的变形条纹图像并进行希尔伯特-黄变换处理，逐行分析确定能够描述变形条纹变化情况的瞬时频率，同时确定条纹图的背景分量。将原条纹图去除背景分量。根据确定的瞬时频率逐行定位每一行信号的局部平稳区域，若任意一行信号的局部平稳区域个数为 1，则对该行信号直接进行傅里叶变换然后提取基频分量的频谱。若任意一行信号的局部平稳区域个数不为 1，根据所确定的局部平稳区域的长度，计算每一个局部平稳区域内每个信号相应的高斯窗函数的尺度因子。由不同的尺度因子引导不同尺寸的高斯窗，对相应的信号逐点进行由高斯窗函数引导的窗口傅里叶变换，每个像素点的变换结果都是一个一维复数数组即傅里叶频谱，将基频频谱提取出来。对一行信号逐个信号进行相同的处理，所有的处理结果依次按行排列，得到一个二维复数矩阵。根据高斯窗函数的时频统计特性，将所有的提取出的基频频谱按列相加，得到最终的一行信号的基频频谱。对相加后的一行信后的基频频谱进行傅里叶反变换，得到一行信号的基频分量的相位。对条纹图逐行进行处理后，最终可得整幅条纹图的全场包裹相位分布。由于所采用的方法不会使所求的瞬时频率被局部平滑化，故根据其所确定的窗口尺寸能够较好的反映条纹信号的局部信息，从而能够精确地测量条纹图细节处的相位。同时，由于背景分量被去除，故在基频提取的过程中极大地减少了零频频谱的干扰，有利于基频频谱的精确提取。 When the fringe pattern is extracted in the optical three-dimensional measurement, the present invention proposes a multi-Hilbert-yellow transform based on the problem that the existing time-frequency analysis method has low measurement accuracy for objects with complex surfaces or steep edges. Scale window Fourier phase extraction method. Firstly, the acquired deformation fringe image is decomposed by the empirical mode decomposition method and Hilbert-Yellow transform processing is performed. The line-by-line analysis determines the instantaneous frequency which can describe the variation of the deformation fringe and determines the background component of the fringe pattern. Remove the background component from the original fringe pattern. Position each line line by line according to the determined instantaneous frequency In the local stationary region of the signal, if the number of local stationary regions of any one line signal is 1, the Fourier transform is directly performed on the row signal and then the spectrum of the fundamental frequency component is extracted. If the number of local stationary regions of any one of the signals is not 1, the scale factor of the corresponding Gaussian window function of each signal in each local stationary region is calculated according to the determined length of the local stationary region. The Gaussian window of different sizes is guided by different scale factors, and the window Fourier transform guided by the Gaussian window function is performed point by point on the corresponding signal. The transformation result of each pixel is a one-dimensional complex array, ie, Fourier spectrum. , extract the fundamental frequency spectrum. The same processing is performed on one line of signals one by one, and all the processing results are sequentially arranged in rows to obtain a two-dimensional complex matrix. According to the time-frequency statistical characteristics of the Gaussian window function, all the extracted fundamental frequency spectra are added in columns to obtain the fundamental frequency spectrum of the final line of signals. The inverse Fourier transform is performed on the fundamental frequency spectrum of the added one-line signal to obtain the phase of the fundamental frequency component of the one-line signal. After the stripe pattern is processed line by line, the whole field wrapping phase distribution of the entire stripe pattern is finally obtained. Since the method used does not cause the instantaneous frequency to be locally smoothed, the local information of the fringe signal can be better reflected according to the determined window size, so that the phase at the detail of the fringe pattern can be accurately measured. At the same time, since the background component is removed, the interference of the zero-frequency spectrum is greatly reduced in the process of the fundamental frequency extraction, which is advantageous for the accurate extraction of the fundamental frequency spectrum.

本发明的具体实施步骤如下： The specific implementation steps of the present invention are as follows:

+ (γ, χ)]，其中，是背景分量，是物体表面反射率， /。是正弦条纹基频， (y,x)是待求的包裹相位分布，（y, c)表示变形条纹图像中各像素点的二维坐标，取值范围分别为： l< y<l, l<x<c , β()^)( 8[2;Γ/。;ϊ + ^(}Μ)]为基频分量，此处将每个像素点看作一个信号。 Step 1: Project the gray sine fringe pattern onto the surface of the measured object, and use the CCD to shoot the surface of the measured object to obtain a deformed fringe image g (y, with a width of ^).

+ (γ, χ)], where is the background component, the surface reflectivity of the object, /. Is the sinusoidal fringe fundamental frequency, (y, x) is the parcel phase distribution to be sought, and (y, c) is the two-dimensional coordinate of each pixel in the deformed fringe image, and the value ranges are: l<y<l, l <x<c , β()^)( 8[2;Γ/.;ϊ + ^(}Μ)] is the fundamental frequency component, where each pixel is treated as a signal.

图 2为所得的变形条纹图像，大小为 1020X 1020，图中直线为 y=170时的一行示例条纹信号 g(x)。 Fig. 2 is an example of a stripe signal g(x) of a stripe image obtained by a deformed fringe image having a size of 1020X 1020 and a straight line of y=170.

步骤 2: 令 y = l， n = l, 用经验模式分解法即 EMD对进行分解，方法如下：步骤 2.1: 将记为一行信号 g )，其中 X仍然满足 l≤ c≤c，该行信号的相位为 φ{χ) , 寻找 g )的极大值点和极小值点，分别采用公知的三次样条插值法对这些极大值点和极小值点进行插值，然后将这些值连接得到极大值包络线 g_max )和极小值包络线 g^ )。此处也可选择其他插值法，具体的选取办法可参考文献 "N. E. Huang, Z. Shen, and Step 2: Let y = l, n = l, decompose using the empirical mode decomposition method, ie EMD pair, as follows: Step 2.1: will be recorded as a line of signal g), where X still satisfies l ≤ c ≤ c, the line signal The phase is φ{χ), looking for the maxima and minima points of g), interpolating these maxima and minima points using well-known cubic spline interpolation, respectively, and then connecting these values Obtain the maximum envelope g _max ) and the minimum envelope g^). Other interpolation methods can also be selected here. For specific selection methods, please refer to the literature "NE Huang, Z. Shen, and

S. . Long et al, "The empirical mode decomposition and the Hilbert spectrum for nonliner and non-stationary time series analysis," J. Proc. Lond. A. Great Britain, 454, 903-995(1998)",及" G. Rilling, P. Flandrin and P. Goncalves, "On empirical mode decomposition and its algorithm," in lEEEEURASIP Workshop on Nonlinear Signal and Image Processing, NSIP-03, GRADO(I) (IEEE, 2003)"; S. . Long et al, "The empirical mode decomposition and the Hilbert spectrum for nonliner and non-stationary time series analysis," J. Proc. Lond. A. Great Britain, 454, 903-995 (1998) ", and" G. Rilling, P. Flandrin and P. Goncalves, "On empirical mode decomposition and its algorithm," in lEEEEURASIP Workshop on Nonlinear Signal and Image Processing, NSIP-03, GRADO (I) (IEEE, 2003)";

步骤 2.2: 用原始信号减去极大值包络线 g_max ( 与极小值包络线 ( c)的平均值，得到

Step 2.2: Subtract the maximum value envelope g _max (with the minimum value envelope (c) from the original signal to obtain

步骤 2.3: 分别计算的平均幅度 mi¾z/i_izm/7(x)和包络幅度 e/w_izm/7(x)： mean _ amp (x) > max v ^ nun Step 2.3: Calculated average amplitude mi3⁄4z/i_izm/7(x) and envelope amplitude e/w_izm/7(x): mean _ amp (x) > max v ^ nun

" 2

步骤 2.4:如果 A )同时满足以下三个条件，则得到一个本征模函数 IMF,记为 s„ )， s_n(x) = h(x)，同时将 s_n(x)从原信号 g(x)中分离，得到新的信号 = g(x)-s_n(x)，否则，直接将从原信号 g )中分离，得到新的信号 g'( ) = g( )-/z( )，所述的三个条件是： Step 2.4: If A) satisfies the following three conditions at the same time, an eigenmode function IMF is obtained, denoted as s„), s _n (x) = h(x), and s _n (x) is from the original signal g Separate in (x) to obtain a new signal = g(x)-s _n (x), otherwise, it will be directly separated from the original signal g) to obtain a new signal g'( ) = g( )-/z( ), the three conditions mentioned are:

条件 1： mean _amp(x) / env _amp(x) < 0.5；条件 2: 满足不等式《¾6¾^_" ( /6^_" (； <0.5的像素个数占同一整行像素总数的比例小于 0.05; Condition 1: mean _amp(x) / env _amp(x) < 0.5; Condition 2: Satisfy the inequality "3⁄463⁄4^_" ( /6^_" (; <0.5 The number of pixels in the total number of pixels in the entire row is less than 0.05;

条件 3: 极大值和极小值的个数之和与过零点的个数相等或至多相差 1个；步骤 2.5:若步骤 2.4得到的 g '( )也同时满足以上三个条件，令 = g '( )，《 = « + 1，返回步骤 2.1; 否则分解停止，并将记为得到最后的分解结果为： Condition 3: The sum of the maximum value and the minimum value is equal to or at most one difference from the zero crossing point; Step 2.5: If g '( ) obtained in step 2.4 also satisfies the above three conditions, let = g '( ), " = « + 1, return to step 2.1; otherwise the decomposition stops and will be recorded as the final decomposition result:

其中 n为一个 IMF即 5(x)的序数， 1<η<Ν , Λ/为 IMF的总个数。图 3 为一行信号分解过程的流程图。这里以变形条纹图的一行变形信号为例，即对如图 2中直线所代表信号 g )进行处理。在处理过程中，三个条件的阈值 0.5、 0.05和 Where n is the ordinal number of an IMF, ie 5(x), 1<η<Ν , Λ/ is the total number of IMFs. Figure 3 is a flow chart of a signal decomposition process. Here, a row of the deformation signal of the deformed fringe pattern is taken as an example, that is, the signal g) represented by the straight line in FIG. 2 is processed. During processing, the thresholds for the three conditions are 0.5, 0.05 and

0.05 为建议阈值，实际应用过程中可根据需要适当微调整，调整的前提是使分解彻底，即剩余分量不能再分离出 IMF, 同时又不能过分解。图 4为原信号及其被分解后的结果，从中可看见， IMFs随着尺度的变化由高到底依次排列，第一个 IMF为高频噪声分量， res(x) 的尺度最大，趋近于 0，它能够描述一行信号的整体趋势。 0.05 is the recommended threshold. In the actual application process, it can be adjusted as needed. The premise of adjustment is to make the decomposition completely, that is, the residual component can no longer separate the IMF, and at the same time can not be decomposed. Figure 4 shows the original signal and its decomposed results. It can be seen that the IMFs are arranged in order from high to low. The first IMF is the high-frequency noise component, and the res(x) scale is the largest, approaching 0, it can describe the overall trend of a line of signals.

步骤 3.1: 对第个 IMF即 ( 做希尔伯特变换，得到： s_n (x) = s_n(x)*— =― -άτ， Step 3.1: For the first IMF, (do the Hilbert transform, get: s _n (x) = s _n (x) * - = ― -άτ,

πχ π "^J χ - 其中 " *"为卷积运算符， r为积分变量，为希尔伯特变换的结果，分别构造 s„ W的解析信号 _Z„(x)： z_n(x) = s ) +

= λ_η(χ)β^{ιφ χ)} , 其中 /' 为虚部单位， Λ_η(χ) = _η ²(χ) + [^ωΐ为该解析信号 z„( ) 的模值， φ„ (x) = arctan(^s« ^{x)/ )为解析信号 z_n (x)的相位；步骤 3.2: 对第 w个 IMF求其瞬时频率/ χ π π " ^J χ - where " * " is the convolution operator and r is the integral variable, which is the result of the Hilbert transform, respectively constructing the analytical signal _Z „( W): z _n (x) = s ) +

= λ _η (χ)β ^{ιφ χ)} , where /' is the imaginary unit, Λ _η (χ) = _η ² (χ) + [^ωΐ is the modulus of the analytical signal z„( ), φ„ (x ) = arctan( ^s « ^{x) / ) is the phase of the resolved signal z _n (x); Step 3.2: Find the instantaneous frequency of the wth IMF /

1 άφ_η (χ) 1 άφ _η (χ)

2π dx 2π dx

对第 w个 IMF求其边际谱 M /)：

Find the marginal spectrum M /) for the wth IMF:

步骤 3.3: 确定含基频分量最多的 IMF的序数 K:

其中 /„^M皿为第《个 IMF的边际谱最大值 M_max所对应的频率值， }为取 I I的最小值所对应的 IMF的序数， Step 3.3: Determine the ordinal number K of the IMF with the most fundamental component:

Wherein / ^M is the frequency value corresponding to the marginal maximum M _{max of} the IMF, and } is the ordinal number of the IMF corresponding to the minimum value of II.

如果 IMF所对应的 I f_n ^M -/。 I最小值个数为 1个以上，则从具有 I -/。 I最小值的 IMF中选择出具有最大的边际谱最大值 M_max的 IMF, 并以具有最大的边际谱最大值 M_max 的 IMF所对应的序数为所求的值。 If the IMF corresponds to I f _n ^M -/. When the number of I minimum values is one or more, the slave has I - /. I minimum The IMF having the largest marginal spectral maximum value M _{max is selected} in the IMF, and the number corresponding to the IMF having the largest marginal spectral maximum value M _max is the obtained value.

选取 f_K 作为能够准确描述该行条纹信号变化规律的瞬时频率。 f _{K is} chosen as the instantaneous frequency that can accurately describe the variation of the stripe signal of the line.

如图 5为图 4中每个 IMF所得的瞬时频率，图 6为相应的边际谱。图 6中，由于条纹的基频是 0.05，可见 IMF2的边际谱最大值所对应的频率值与基频 0.05的差值最小，故 IMF2为含有基频分量最多的分量， 11\^-/₂则被选定为能够准确描述该行条纹信号变化规律的瞬时频率。 Figure 5 shows the instantaneous frequency obtained by each IMF in Figure 4, and Figure 6 shows the corresponding marginal spectrum. In Fig. 6, since the fundamental frequency of the fringe is 0.05, it can be seen that the difference between the frequency value corresponding to the maximum marginal spectrum of IMF2 and the fundamental frequency of 0.05 is the smallest, so IMF2 is the component containing the most fundamental frequency component, 11\^-/ ₂ It is then selected as the instantaneous frequency that accurately describes the variation of the line fringe signal.

步骤 4:确定一行条纹信号的背景分量，具体过程如下：根据步骤 3.2已得的每个 IMF 的瞬时频率，求出每个 IMF的瞬时频率均值，并从中找到最小的瞬时频率均值，然后确定最小的瞬时频率均值所对应的 IMF序数，记为 K_b, 则以步骤 2.5所得的第个 IMF 到第 N个 IMF以及残余分量之和即

) + r^(;c)为该行条纹信号的背景分量组合。步骤 5: 根据一个局部平稳区域内的瞬时频率最大值小于 2倍最小值这一属性，用步骤 3.3确定的瞬时频率 f_K(x)自适应定位一行条纹信号的局部平稳区域，具体过程如下：步骤 5.1: 瞬时频率向量 ^( c)为 [/ 1)， f_K(2) , …… f_K(c)], 其中元素全部为不小 .(2) Step 4: Determine the background component of the stripe signal. The specific process is as follows: According to the instantaneous frequency of each IMF obtained in step 3.2, find the instantaneous frequency mean of each IMF, and find the minimum instantaneous frequency mean, and then determine the minimum. The IMF ordinal corresponding to the instantaneous frequency mean, denoted as K _b , then the sum of the first IMF to the Nth IMF and the residual component obtained in step 2.5

) + r^(;c) is the background component combination of the stripe signal of the line. Step 5: According to the property that the maximum value of the instantaneous frequency in a local stationary region is less than 2 times the minimum value, the local stationary region of the stripe signal is adaptively located by using the instantaneous frequency f _K (x) determined in step 3.3. The specific process is as follows: Step 5.1: The instantaneous frequency vector ^(c) is [/ 1), f _K (2) , ... f _K (c)], where the elements are all not small. (2)

于 0的瞬时频率，对^ ( c)求转置向量，得到尸 W为 fAO At the instantaneous frequency of 0, find the transpose vector for ^ ( c) and get the corpse W as fAO

步骤 5.2: 令/ ^TW中每个元素分别减去整个 2x W向量，得到一个 _Cxc的方阵 f_K(l)-2xf_K(l) f_K(l)-2xf_K(2) f_K(l)-2xf_K(c) Step 5.2: Decrease the entire 2x W vector for each element in / ^T W to obtain a square matrix of _C xc f _K (l)-2xf _K (l) f _K (l)-2xf _K (2) f _K (l)-2xf _K (c)

f_K(2)-2xf_K(l) f_K(2)-2xf_K(2) f_K(2)-2xf_K(c) f _K (2)-2xf _K (l) f _K (2)-2xf _K (2) f _K (2)-2xf _K (c)

F f_K{c)-2 f_K{\) f_K(c)-2xf_K(2) …… f_K(c)-2xf_K(c) F f _K {c)-2 f _K {\) f _K (c)-2xf _K (2) ...... f _K (c)-2xf _K (c)

令方阵内所有负值为 0，并找到方阵中全部元素为 0的零对角线，得到: Let all the negative values in the square matrix be 0, and find the zero diagonal of all the elements in the square matrix, and get:

0 f_K(D- xf_K(2) f_K(D-2xf_K(c) 0 f _K (D- xf _K (2) f _K (D-2xf _K (c)

f_K(2)-2xf_K(l) 0 f_K(2)-2xf_K(c)f _K (2)-2xf _K (l) 0 f _K (2)-2xf _K (c)

F f_K{c)-2 f_K{\) f_K(c)-2xf_K(D …… 0 F f _K {c)-2 f _K {\) f _K (c)-2xf _K (D ...... 0

步骤 5.3:对方阵 F建立坐标系，坐标原点为 F左上角第一个元素，其坐标值记为 (1,1)，横坐标的坐标方向为方阵的行方向，横坐标取值范围为 l ~ c，纵坐标的坐标方向为方阵的列方向，纵坐标的取值范围仍为 l~c， Step 5.3: The coordinate system is established by the partner array F, and the coordinate origin is the first element in the upper left corner of F, and the coordinate value is recorded as (1, 1). The coordinate direction of the abscissa is the row direction of the square matrix, the abscissa is in the range of l ~ c, the coordinate direction of the ordinate is the column direction of the square matrix, and the ordinate is still in the range l~c.

(c, c)表示，其中（_Ρι,_Ρι)为第一个最大全 0子方阵的零对角线上最后一个元素相对于方阵 f 的坐标，依次排列，（Α,/ 为倒数第二个最大全 0子方阵的零对角线上最后一个元素相对于方阵 F的坐标，（_C,c)为最后一个最大全 0子方阵的零对角线上最后一个元素相对于方阵 _f的坐标，最后，该行条纹信号的平稳区域划分情况为：（1,……， A)，（A+1,……， Ρ₂)，……， (p_r + ， c)。 (c, c) means that ( _Ρι , _Ρι ) is the coordinate of the last element on the zero diagonal of the first largest all-zero square matrix relative to the square matrix f, in order, (Α, / is the last number The coordinates of the last element on the zero-diagonal line of the two largest all-zero square matrices relative to the square matrix F, ( _C , c) being the last element on the zero-diagonal line of the last largest all-zero square matrix relative to The coordinates of the square matrix _f , finally, the smooth region division of the stripe signal of the line is: (1, ..., A), (A+1, ..., Ρ ₂ ), ..., (p _r + , c) .

图 7为自适应定位一行条纹信号局部平稳区域的直观描述，图 8为图 3中对直线代表信号局部平稳区域的划分。 Fig. 7 is a visual depiction of a local stationary region of a line of stripe signals adaptively positioned, and Fig. 8 is a division of a local stationary region of a line representative signal in Fig. 3.

步骤 6:将原条纹信号 g )减去步骤 4中确定的背景分量，得到步骤 1中去除后的基频分量 g。 (x)，即为 g。 (x)=B(x, y) cos[2^/^ + (χ, y)]。 Step 6: Subtract the original fringe signal g) from the background component determined in step 4 to obtain the removed fundamental frequency component g in step 1. (x), which is g. (x)=B(x, y) cos[2^/^ + (χ, y)].

图 9为去除背景分量后的条纹图。 Figure 9 is a stripe diagram after removing the background component.

步骤 7: 判断步骤 5确定的局部平稳区域个数是否为 1，若 "是"，则直接对 g。( c)做快速傅里叶变换得到傅里叶频谱，记作 F_f，若 "否"，则对 g。( c)进行自适应高斯窗口傅里叶变换，具体过程为： Step 7: Determine whether the number of local stationary regions determined in step 5 is 1, and if yes, directly to g. (c) Do a fast Fourier transform to obtain the Fourier spectrum, denoted as F _f , and if "No", then for g. (c) Perform an adaptive Gaussian window Fourier transform, the specific process is:

步骤 7.1: 对 g。( c)做自适应高斯窗口傅里叶变换： Step 7.1: For g. (c) Do an adaptive Gaussian window Fourier transform:

WF_a ^b[g₀(x)] = J厂 [ ( ( 6邛(_ ] ， WF _a ^b [g ₀ (x)] = J factory [ ( ( 6邛(_ ] ,

其中 b为平移因子， b依次取值为 1、 2、 3、…、 c α为每一个局部平稳区域内每个像素点对应的高斯窗口的尺度因子， a = LI4 , 其中 _为相应的局部平稳区域的长度值，单位为像素， ^( =丄1^(^)=~^6邛[-丄(^)²]为高斯窗函数， Where b is the translation factor, and b is taken as 1, 2, 3, ..., c α is the scale factor of the Gaussian window corresponding to each pixel in each local stationary region, a = LI4, where _ is the corresponding local The length of the stationary region, in pixels, ^( =丄1^(^)=~^6邛[-丄(^) ² ] is a Gaussian window function.

a a j2na 2 a 对一行条纹信号所有的信号依次进行自适应高斯窗口傅里叶变换后，得到由 WF_fl ^fc[g。(x)]构成的一个二维复数方阵，大小为 c x c，二维复数方阵每一行的元素为每一个窗口内信号的频谱，二维复数方阵一共有 c行，代表 b从 1取值到 c即一共有 c个窗口的信号的频谱； Aa j2na 2 a After all the signals of one line of the stripe signal are sequentially subjected to adaptive Gaussian window Fourier transform, the result is obtained by WF _fl ^fc [g. (x)] A two-dimensional complex square matrix, the size of which is cxc, the element of each row of the two-dimensional complex square matrix is the spectrum of the signal in each window, and the two-dimensional complex matrix has a total of c rows, representing b from 1 The value to c is the spectrum of the signal of a total of c windows;

步骤 7.2: 将步骤 7.1所得的二维复数方阵按列叠加，得到总频谱，记作 ^。 Step 7.2: Superimpose the two-dimensional complex square matrix obtained in step 7.1 by column to obtain the total spectrum, which is denoted as ^.

步骤 8: 根据 /^确定基频频谱的范围，并提取出来，记为 F_{f 0 ;} 对求傅里叶逆变换，并根据对 f/o求傅里叶逆变换所得的结果，求得相位角即包裹在 [_π, π]之间的相位，令) /=y+l，返回步骤 2，若>^ > /，则进入步骤 9。步骤 9: 对包裹相位进行展开得到绝对相位，此处相位展开的方法采用质量图相位展开法。最后，根据经典光栅投影的相位到高度的转换公式求得测量物体的三维信息。高度转换公式如下： Step 8: Determine the range of the fundamental frequency spectrum according to /^, and extract it, denoted as F _{f 0 ;} to find the inverse Fourier transform, and obtain the phase based on the inverse Fourier transform of f/o The angle is wrapped in the phase between [_π, π], let) /=y+l, return to step 2, if >^ > /, proceed to step 9. Step 9: Unfold the parcel phase to obtain the absolute phase. Here, the phase unwrapping method uses the mass map phase expansion method. Finally, the three-dimensional information of the measured object is obtained according to the phase-to-height conversion formula of the classical grating projection. The height conversion formula is as follows:

ΙΑφ ΙΑφ

h x, y) = ~， h x, y) = ~,

Αφ - ω₀ά 其中， /、 d是测量***的几何参数， /是投影仪到测量平面的距离， d是 CCD摄像头到投影仪的距离， Δ = (υ)-% (·^, } 表示相位变化量， (·¾, } 为展开相位结果，。（为初始相位结果，由测量参考面决定，《。为投影光栅的角频率，由***标定可得。 Αφ - ω ₀ ά where /, d is the geometric parameter of the measurement system, / is the distance from the projector to the measurement plane, d is the distance from the CCD camera to the projector, Δ = (υ)-% (·^, } The amount of phase change, (·3⁄4, } is the unwrapped phase result, (for the initial phase result, determined by the measurement reference plane, “. is the angular frequency of the projected grating, which is available from the system calibration.

如图 10为尺度因子的全场分布图，可见用本方法能够较好地描述物体的变化规律。图 11为对一整幅变形条纹图进行本发明方法处理后所得的全场包裹相位分布图。图 12 为所测物体的相位变化分布图，即先将全场包裹相位用质量图法进行相位展开，然后用同样方法得到没有被物体调制的参考平面条纹图的绝对相位图，将包含物体三维信息的变形条纹图的展开相位与参考平面条纹图的展开相位相减即可得到相位变化分布图。用本方法所测的包裹相位，在洪水法进行相位展开时，由于物体边缘相位测量的精度较高，基本没有误差传递，而其他方法会有大片的误差传递现象。故本方法对复杂面形或含有陡峭边缘物体的测量时精度会有具有较大的提高。 Figure 10 shows the full-field distribution of the scale factor. It can be seen that this method can better describe the variation law of the object. Figure 11 is a diagram showing the phase distribution of the full field package obtained by performing the method of the present invention on a whole deformed fringe pattern. Figure 12 is a phase change distribution diagram of the measured object, that is, the phase expansion of the full-field wrapping phase is performed by the mass map method, and then the absolute phase map of the reference plane fringe pattern not modulated by the object is obtained by the same method, and the object is three-dimensionally included. The phase change profile is obtained by subtracting the unwrapped phase of the deformed fringe pattern of the information from the unwrapped phase of the reference plane fringe pattern. When the phase of the package measured by the method is used for phase unwrapping in the flood method, since the accuracy of the phase measurement of the edge of the object is high, there is basically no error transmission, and other methods have a large error transmission phenomenon. Therefore, the method has a great improvement in the measurement accuracy of a complex surface or an object having a steep edge.

Claims

权利要求书 Claim

1. 一种光学三维测量中的自适应窗口傅里叶相位提取法，具体步骤如下： 1. An adaptive window Fourier phase extraction method in optical three-dimensional measurement, the specific steps are as follows:

步骤 1: 将灰度正弦条纹图投影到被测物体表面，用 CCD对被测物体表面进行拍摄，得到一幅宽度为^、高度为 /的变形条纹图像 Step 1: Project the gray sine fringe pattern onto the surface of the measured object, and use the CCD to shoot the surface of the measured object to obtain a deformed fringe image with a width of ^ and a height of /.

g(y,

+ (γ, χ)]， g(y,

+ (γ, χ)],

其中，是背景分量，是物体表面反射率， /。是正弦条纹基频， (y,x)是待求的包裹相位分布，（y, c)表示变形条纹图像中各像素点的二维坐标，取值范围分别为： l< y<l, l<x<c , £(}^)(^[2；2"/。_1 + (}^)]为基频分量，此处将每个像素点看作一个信号； Among them, is the background component, which is the surface reflectivity of the object, /. Is the sinusoidal fringe fundamental frequency, (y, x) is the parcel phase distribution to be sought, and (y, c) is the two-dimensional coordinate of each pixel in the deformed fringe image, and the value ranges are: l< y<l, l <x<c , £(}^)(^[2;2"/._1 + (}^)] is the fundamental frequency component, where each pixel is treated as a signal;

步骤 2: 令 y = l， n = l , 用经验模式分解法即 EMD对进行分解，方法如下：步骤 2.1: 将记为一行信号 g )，其中 c仍然满足 l≤ c≤c，该行信号的相位 φ χ、, 寻找 g )的极大值点和极小值点，分别采用公知的三次样条插值法对这些极大值点和极小值点进行插值，然后将这些值连接得到极大值包络线 g_max )和极小值包络线 Step 2: Let y = l, n = l, decompose using the empirical mode decomposition method, ie EMD pair, as follows: Step 2.1: Will be recorded as a line of signals g), where c still satisfies l ≤ c ≤ c, the line signal The phase φ χ , , find the maxima and minima points of g), interpolate these maxima and minima points by well-known cubic spline interpolation, respectively, and then connect these values to the pole Large value envelope g _max ) and minimum envelope

步骤 2.2: 用原始信号 g )减去极大值包络线 g_maxW与极小值包络线的平均值，得到

Step 2.2: subtracting the average value of the maximum value envelope g _max W from the minimum envelope by using the original signal g)

步骤 2.3: 分别计算 ( 的平均幅度" ^im_amp( )和包络幅度 ew_amp( )： mean _ amp (x) Step 2.3: Calculate the average magnitude of ( ^im_amp( ) and the envelope magnitude ew_amp( ): mean _ amp (x)

2

步骤 2.4: 如果 A )同时满足以下三个条件，则得到一个本征模函数 IMF, 记为 s_n(x) , 且 W = A ) ，同时将从原信号 g ) 中分离，得到新的信号 g\x) = g(x)-s_n(x) , 否则，直接将 h(x) 从原信号中分离，得到新的信号 g\x) = g(x)-h(x) , 所述的三个条件是： Step 2.4: If A) satisfies the following three conditions at the same time, an eigenmode function IMF is obtained, denoted as s _n (x), and W = A), and will be separated from the original signal g) to obtain a new signal. g\x) = g(x)-s _n (x) , otherwise, directly separate h(x) from the original signal to obtain a new signal g\x) = g(x)-h(x) The three conditions stated are:

条件 1： mean _ amp(x) I env _ amp(x) < 0.5， Condition 1: mean _ amp(x) I env _ amp(x) < 0.5,

条件 2: 满足不等式《^ _" ( /6^_" (； <0.5的像素个数占同一整行像素总数的比例小于 0.05，条件 3: 极大值和极小值的个数之和与 A )过零点的个数相等或至多相差 1个，步骤 2.5: 若步骤 2.4 得到的也同时满足以上三个条件，令 g( ) = g'( )， n = n + 返回步骤 2.1; 否则分解停止，并将记为 r_es(;c)，得到最后的分解结果为：

Condition 2: Satisfy the inequality "^ _" ( /6^_"(;<0.5 of the number of pixels in the total number of pixels in the entire row is less than 0.05, Condition 3: The sum of the maximum and minimum values is equal to or equal to the number of A) zero crossings. Step 2.5: If the conditions obtained in step 2.4 also satisfy the above three conditions, let g( ) = g'( ), n = n + returns to step 2.1; otherwise the decomposition stops, and will be denoted as r _e s(;c), and the final decomposition result is:

其中《为一个 IMF即 ;)的序数， \<η<Ν， N为 IMF的总个数； The ordinal number of "IMF is ;", \<η<Ν, N is the total number of IMFs;

步骤 3.1: 对第个 IMF即 ( 做希尔伯特变换，得到： s_n x) = ( ) *丄 =丄 Step 3.1: For the first IMF ie (do the Hilbert transform, get: s _n x) = ( ) *丄=丄

πχ π厂 ^~άτ Χχ π factory ^~άτ

∞ χ— τ ∞ χ — τ

其中 " *"为卷积运算符， r为积分变量，为希尔伯特变换的结果， Where " *" is the convolution operator and r is the integral variable, which is the result of the Hilbert transform.

分别构造 ( C)的解析信号 _Ζ„( )： z_n(x) = s ) + (χ)β^{φ χ)} , 其中 i 为虚部单位， (χ) = _η ²(χ) +

为该解析信号 Z„( ) 的模值， φ„(χ) = arciSin(^S"^(x/⁾ )为解析信号 z„( )的相位；步骤 3.2: 对第 w个 IMF求其瞬时频率/ Each configuration (C) of the analytic signal _{_{Ζ "(): z n (}} x) = s) + (χ) β φ χ), where i is the imaginary _{^{unit, (χ) = η 2 (}} χ) +

For the modulus of the analytical signal Z„( ), φ„(χ) = arciSin( ^S " ^(x / ⁾ ) is the phase of the analytical signal z„( ); Step 3.2: Find the instantaneous frequency of the wth IMF/

2π dx 2π dx

对第 n个 IMF求其边际谱 M_n (f)

Find the marginal spectrum M _n (f) for the nth IMF

步骤 3.3: 确定含基频分量最多的 IMF的序数 Step 3.3: Determine the ordinal number of the IMF with the most fundamental component

K=E{\f_n ^M- -f₀ \}-, 其中 /„^Mm"为第 n个 IMF的边际谱最大值 M_max所对应的频率值， { }为取 I f_n ^M - -/ I的最小值所对应的 IMF的序数， K=E{\f _n ^M - -f ₀ \}-, where /„ ^Mm ′ is the frequency value corresponding to the marginal spectral maximum value M _max of the nth IMF, { } is taken as I f _n ^M − −/ The ordinal number of the IMF corresponding to the minimum value of I,

如果 IMF所对应的 -/。 I最小值个数为 1个以上，则从具有 l/„^M™ -f₀1最小值的 IMF中选择出具有最大的边际谱最大值 M 的 IMF, 并以具有最大的边际谱最大值 M_max 的 IMF所对应的序数为所求的值，选取^ ( 作为能够准确描述该行条纹信号变化规律的瞬时频率； If the IMF corresponds to -/. If the number of I minimum values is one or more, the IMF having the largest marginal spectral maximum value M is selected from the IMF having the minimum value of l/„ ^M TM −f ₀ 1 , and has the largest marginal spectral maximum value M The ordinal corresponding to the IMF of _max is the value sought, Select ^ (as the instantaneous frequency that can accurately describe the variation of the stripe signal of the line;

步骤 4: 确定一行条纹信号的背景分量，具体过程如下：根据步骤 3.2 已得的每个 IMF 的瞬时频率，求出每个 IMF 的瞬时频率均值，并从中找到最小的瞬时频率均值，然后确定最小的瞬时频率均值所对应的 IMF序数，记为 ¾，则以步骤 2.5所得的第 +1个 Step 4: Determine the background component of the stripe signal. The specific process is as follows: According to the instantaneous frequency of each IMF obtained in step 3.2, find the instantaneous frequency mean of each IMF, and find the minimum instantaneous frequency mean, then determine the minimum. The IMF ordinal corresponding to the instantaneous frequency mean is recorded as 3⁄4, then the +1th obtained in step 2.5

IMF 到第 N个 IMF 以及残余分量之和即 t ( + 为该行条纹信号的背景分量组合. The sum of the IMF to the Nth IMF and the residual component, t ( + is the background component combination of the stripe signal of the line.

步骤 5: 根据一个局部平稳区域内的瞬时频率最大值小于 2倍最小值这一属性，用步骤 3.3确定的瞬时频率 f_K 自适应定位一行条纹信号的局部平稳区域，具体过程如下：步骤 5.1: 瞬时频率向量 ^( c)为 [/ 1)， f_K(2 …… f_K(c)], 其中元素全部为不小 (2) Step 5: According to the property that the maximum value of the instantaneous frequency in a local stationary region is less than 2 times the minimum value, the local stationary region of the stripe signal is adaptively located by using the instantaneous frequency f _K determined in step 3.3. The specific process is as follows: Step 5.1: The instantaneous frequency vector ^(c) is [/ 1), f _K (2 ...... f _K (c)], where all elements are not small (2)

于 0的瞬时频率，对^ ( c)求转置向量，得到 W为

At the instantaneous frequency of 0, find the transpose vector for ^ ( c) and get W as

步骤 5.2: 令 /^T )中每个元素分别减去整个 2x (； c)向量，得到一个 _Cxc的方阵 F: f_K(l)-2xf_K(V) f_K(l)-2xf_K(2) f_K(l)-2xf_K(c) Step 5.2: Decrease the entire 2x (; c) vector for each element in / ^T ) to obtain a square matrix F of _C xc: f _K (l) - 2xf _K (V) f _K (l) - 2xf _K (2) f _K (l)-2xf _K (c)

F f_K(c)-2xf_K(l) f_K(c)-2xf_K(2) …… f_K(c)-2xf_K(c) 令方阵内所有负值为 0，并找到方阵中全部元素为 0的零对角线，得到: F f _K (c) -2xf _K (l) f _K (c) -2xf _K (2) ...... f _K (c) -2xf _K (c) Let all negative values in the square matrix be 0, and find the square matrix The zero diagonal of all elements in the 0 is:

0 f_K(l)-2xf_K(2) f_K(l)-2xf_K(c) 0 f _K (l)-2xf _K (2) f _K (l)-2xf _K (c)

f_K(2)-2xf_K(l) 0 f_K(2)-2xf_K(c) f _K (2)-2xf _K (l) 0 f _K (2)-2xf _K (c)

F f_K(c)-2xf_K(l) f_K(c)-2xf_K(l) …… 0 步骤 5.3: 对方阵 F建立坐标系，坐标原点为 F左上角第一个元素，其坐标值记为 (1,1)，横坐标的坐标方向为方阵的行方向，横坐标取值范围为 1 ~ c，纵坐标的坐标方向为方阵的列方向，纵坐标的取值范围仍为 l ~ c， F f _K (c) -2xf _K (l) f _K (c) -2xf _K (l) ...... 0 Step 5.3: The coordinate system is established by the partner array F, and the coordinate origin is the first element in the upper left corner of F, and its coordinates The value is recorded as (1,1), the coordinate direction of the abscissa is the row direction of the square matrix, the abscissa is in the range of 1 ~ c, the coordinate direction of the ordinate is the column direction of the square matrix, and the range of the ordinate is still For l ~ c,

在方阵 F中寻找出以方阵 F的零对角线的一部分或全部为对角线的最大全零子方阵，从方阵 F的坐标原点开始，沿零对角线方向，依次将这些最大全 0子方阵的零对角线上最后一个元素相对于方阵 F 的坐标记录下来，依次用（_Α,_Ρι)， (ρ₂,ρ₂) , ……，Find a maximum all-zero square matrix with a part or all of the zero diagonal of the square matrix F as the diagonal in the square matrix F, starting from the coordinate origin of the square matrix F, along the zero diagonal direction, in turn The last element of the zero-diagonal line of these largest all-zero square matrices is recorded relative to the coordinates of the square matrix F, using ( _Α , _Ρι ), (ρ ₂ , ρ ₂ ), ..., in turn.

(_Pr,p_r), (c， c)表示，其中第一个最大全 0子方阵的零对角线上最后一个元素相对于方阵 F的坐标，依次排列，（Α,/ 为倒数第二个最大全 0子方阵的零对角线上最后一个元素相对于方阵 F的坐标，（_C,c)为最后一个最大全 0子方阵的零对角线上最后一个元素相对于方阵的坐标， ( _Pr , p _r ), (c, c) represents the last elemental phase on the zero diagonal of the first largest all-zero square matrix For the coordinates of the square matrix F, arrange in order, (Α, / is the coordinate of the last element on the zero diagonal of the penultimate maximum all zero sub-matrix relative to the square matrix F, ( _C , c) is the last one The coordinates of the last element on the zero diagonal of the largest all-zero square matrix relative to the square matrix,

最后，该行条纹信号的平稳区域划分情况为：（1,……， _Pl) ，Finally, the smooth region of the stripe signal is divided into: (1, ..., _Pl ),

(Pi+ ， P₂ ，（P_r+1, ， c); (Pi+ , P ₂ , (P _r +1, , c);

步骤 6: 将原条纹信号 g )减去步骤 4中确定的背景分量，得到步骤 1中去除后的基频分量 _go (x)，即为 g。 (x)=B(x, y) cos[2^/₀^ + φ(χ, y)]； Step 6: Subtract the original fringe signal g) from the background component determined in step 4 to obtain the base frequency component _go (x) removed in step 1, which is g. (x)=B(x, y) cos[2^/ ₀ ^ + φ(χ, y)];

步骤 7: 判断步骤 5确定的局部平稳区域个数是否为 1，若 "是"，则直接对 g。( c)做快速傅里叶变换得到傅里叶频谱，记作，若 "否"，则对 g。( c)进行自适应高斯窗口傅里叶变换，具体过程为： Step 7: Determine whether the number of local stationary areas determined in step 5 is 1, if "yes", directly to g. (c) Do a fast Fourier transform to obtain the Fourier spectrum, denoted as "if", then for g. (c) Perform an adaptive Gaussian window Fourier transform, the specific process is:

WF_a ^b[g₀(x)] =厂 [g。(x)W_a ⁶ Wexp(— )] ， WF _a ^b [g ₀ (x)] = factory [g. (x)W _a ⁶ Wexp(- )] ,

J J

其中 b为平移因子， b依次取值为 1、 Where b is the translation factor and b is taken as 1,

2、 2,

3、 …、 c， α为每一个局部平稳区域内每个像素点对应的高斯窗口的尺度因子， a = U4 , 其中 L为相应的局部平稳区域的长度值，单位为像素， ( =丄 W(^)=~^exp [-丄 (^)²]为高斯窗函数， 3, ..., c, α is the scale factor of the Gaussian window corresponding to each pixel in each local stationary region, a = U4, where L is the length of the corresponding local stationary region, in pixels, (=丄W (^)=~^exp [-丄(^) ² ] is a Gaussian window function.

a a 2πα 2 a a a 2πα 2 a

对一行条纹信号所有的信号依次进行自适应高斯窗口傅里叶变换后，得到由 After all the signals of one line of the stripe signal are sequentially subjected to adaptive Gaussian window Fourier transform,

WF_a ⁶[g。( )]构成的一个二维复数方阵，大小为 cxc，二维复数方阵每一行的元素为每一个窗口内信号的频谱，二维复数方阵一共有 c行，代表 b从 1取值到 c即一共有 c个窗口的信号的频谱； WF _a ⁶ [g. ( )] constitutes a two-dimensional complex square matrix, the size is cxc, the element of each row of the two-dimensional complex square matrix is the spectrum of the signal in each window, and the two-dimensional complex matrix has a total of c rows, representing the value of b from 1 To c is the spectrum of the signal of a total of c windows;

步骤 7.2: 将步骤 7.1所得的二维复数方阵按列叠加，得到总频谱，记作 F_/; Step 7.2: Superimpose the two-dimensional complex square matrix obtained in step 7.1 by column to obtain the total spectrum, which is denoted as F _/;

步骤 8: 根据确定基频频谱的范围，并提取出来，记为 F_{f Q;} 对求傅里叶逆变换，并根据对求傅里叶逆变换所得的结果，求得相位角即包裹在 [-π,π]之间的相位 (χ) , 令 }^1，返回步骤 2， ^y>l, 则进入步骤 9; Step 8: According to determining the range of the fundamental frequency spectrum, and extracting it, denoted as F _{f Q;} for inverse Fourier transform, and based on the result obtained by inverse Fourier transform, the phase angle is wrapped in [ Phase (χ) between -π,π], let }^1, return to step 2, ^y>l, then proceed to step 9;

步骤 9: 对包裹相位进行展开得到绝对相位，根据经典光栅投影的相位到高度的转换公式，最终求得测量物体的三维信息。 Step 9: Unfold the parcel phase to obtain the absolute phase. According to the phase-to-height conversion formula of the classical grating projection, the three-dimensional information of the measured object is finally obtained.