相位解包裹中欠采樣問題的實驗研究

郭 媛，陳小天，毛 琦

(齊齊哈爾大學(xué)計算機與控制工程學(xué)院，齊齊哈爾161006)

式中

引 言

數(shù)字散斑干涉測量技術(shù)是激光技術(shù)與數(shù)字圖像處理相結(jié)合的光測方法。它具有精度高、靈敏度高、非接觸無損測量等優(yōu)點［1-3］。在數(shù)字散斑測量中，待測光場的相位分布是通過求反正切運算得到的，得到的包裹相位值被限制在(-π，π］之間，這就需要通過相位解包裹技術(shù)得到真實相位。近年來，隨著光學(xué)測量技術(shù)的發(fā)展，很多相位解包裹方法被提出，如路徑跟蹤算法、最優(yōu)估計算法、最小范數(shù)算法等［4-6］。其中最常用的算法是基于最小范數(shù)的最小二乘解包裹算法，最小二乘方法最后得到泊松方程，求解泊松方程的方法有迭代法［7］和直接法［8-9］。同時，最小二乘相位解包裹又分為非加權(quán)最小二乘相位解包裹和加權(quán)最小二乘相位解包裹［10-11］。

欠采樣是由于干涉條紋較密集，而現(xiàn)有的CCD的有限空間帶寬限制，使得計算機的取樣間隔過大造成信息丟失造成的。CCD的空間寬帶積越小，欠采樣現(xiàn)象越嚴(yán)重。應(yīng)用傳統(tǒng)的相位解包裹算法將不能從欠采樣包裹相位圖中解出正確的連續(xù)相位分布，致使解包裹失敗，這個問題一直是相位解包裹中的一個難點。

剪切干涉測量是光學(xué)測量中的一種重要的測量方法，它通過同一光波與它的一個小錯位后的光波之間的干涉完成測量，根據(jù)剪切方向可分為橫向剪切、徑向剪切、旋轉(zhuǎn)剪切和方向剪切等。有學(xué)者將橫向剪切的理論引入欠采樣相位解包裹中，能夠有效地解決欠采樣問題［12-15］。橫向剪切干涉能夠有效地去除欠采樣現(xiàn)象的原因在于它利用波前相位與其自身小的錯位之間的偏差往往很小的特點，然后通過對偏差進(jìn)行處理，進(jìn)而實現(xiàn)對欠采樣包裹相位的恢復(fù)。但是，橫向剪切只是在一定的范圍內(nèi)能夠解決欠采樣問題，當(dāng)欠采樣較嚴(yán)重、相位變化過快時，現(xiàn)有的橫向剪切算法就會失效，不能正確地恢復(fù)正確的相位分布。本文中在橫向剪切最小二乘相位解包裹算法(least square algorithm based on lateral shearing interferometry，LSBLS)的基礎(chǔ)上提出新的相位解包裹算法，解決了欠采樣帶來的誤差，提高了原有的相位解包裹算法的抗采樣干擾能力，并通過理論分析和實驗驗證了其正確性和有效性。

1 橫向剪切最小二乘算法原理

數(shù)字散斑干涉測量中，從相干圖像提取的包裹相位為 ψi，j(- π≤ψi，j≤π)，與之對應(yīng)的真實相位為 φi，j，并有:式中，ki，j為整數(shù)，(i，j)是 M × N 圖片中像素點的坐標(biāo)，0≤i≤M-1，0≤j≤N-1。

定義x方向和y方向上的包裹相位差分別為Δxi，j和 Δyi，j:

式中，W為包裹算子，其目的是對包裹相位偏導(dǎo)數(shù)進(jìn)行加減 2π，從而確保 Δxi，j和 Δyi，j位于(- π，π］之間。

最小二乘相位解包裹算法的基本思想是尋求真實相位的偏導(dǎo)數(shù)和包裹相位差的差分的值最小，即:

對(3)式φi，j求導(dǎo)數(shù)并令其為0，整理得泊松方程:

式中

因為相位微分在包裹相位圖像邊緣是無效的，故泊松方程的Neumann邊界條件為:

顯然最小二乘相位解包裹算法轉(zhuǎn)化為對泊松方程的求解，常用的求泊松方程的方法有迭代法和離散余弦變換(discrete cosine transform，DCT)或者快速傅里葉變換(fast Fourier transform，F(xiàn)FT)。

傳統(tǒng)的方法是用上述方法求出包裹相位，當(dāng)相位變化快而使解包裹失敗。為解決這個問題，從而LSBLS算法被提出，其不同于傳統(tǒng)的解包裹方法是引入剪切干涉的原理，在(2)式的基礎(chǔ)上建立一個等效的2維復(fù)光場:

將兩光場相除，即可得到新的光場:

式中，Im()和Re()分別表示取復(fù)函數(shù)的虛部和實部運算。同理，在y方向上有:

將(10)式和(11)式代入(5)式中，通過DCT解包裹即可求出連續(xù)相位。

2 LSBLS算法的改進(jìn)

盡管LSBLS方法可以有效地克服欠采樣帶來的麻煩，與傳統(tǒng)最小二相位解包裹相比存在優(yōu)勢，但當(dāng)欠采樣非常嚴(yán)重時，就會導(dǎo)致LSBLS方法失效，真實相位而不能解出。

為解決現(xiàn)有的LSBLS方法的缺點，作者在現(xiàn)有的LSBLS方法的基礎(chǔ)上進(jìn)行改進(jìn)。同樣將剪切干涉的原理引入到相位解包裹中，在(5)式的基礎(chǔ)上構(gòu)建一個等效的 2 維復(fù)光場:

同樣沿x方法作1pixel的平移，創(chuàng)建新的光場:

將兩光場相除有:

則:

同理在y方向上有:

則(5)式變?yōu)?

故用常規(guī)的求解泊松方程的方法即可求出連續(xù)相位。與原有的LSBLS算法相比，改進(jìn)后的LSBLS算法并不會改變原有的LSBLS算法的計算量，在理論上，兩種算法的計算時間是相同的。經(jīng)過實驗研究，作者發(fā)現(xiàn)改進(jìn)后的LSBLS算法沒有影響原有的LSBLS算法的計算時間，但在解決欠采樣問題上，改進(jìn)后的LSBLS算法有明顯的優(yōu)勢。

盡管DCT能夠滿足常用的最小二乘相位解包裹，但是，當(dāng)殘差點較多時DCT非加權(quán)權(quán)最小二乘相位解包裹算法的計算結(jié)果比較平滑，從而有學(xué)者引入權(quán)重，抑制誤差的傳播，稱為加權(quán)離散余弦變換(weighted discrete cosine transform，WDCT)。為補償其平滑作用，采用加權(quán)來對最小二乘方法進(jìn)行改進(jìn)。權(quán)重定義

式中，pi，j和 qi，j分別為 x 方向和 y 方向相鄰包裹相位差分之差，unit()為歸一化處理，filt［］為均值濾波，求出權(quán)重 ui，j后，用 ui，j對 ρi，j進(jìn)行加權(quán)處理:

式中，ρui，j表示對 ρi，j加權(quán)處理的結(jié)果，然后利用(4)式泊松方程的求解即可求出連續(xù)相位。鑒于WDCT算法的優(yōu)勢，本文實驗中采用WDCT算法解泊松方程。

3 實驗與結(jié)果分析

為了驗證本文中改進(jìn)后的LSBLS方法的有效性，用MATLAB編程，以100倍的peaks函數(shù)構(gòu)建一個2維相位分布(512pixel×512pixel)進(jìn)行實驗驗證，該2維相位分布的相位最大值和相位最小值分布為-655.1120rad和810.6041rad。用LSBLS算法解該2維相位如圖1所示。

Fig.1 Results of unimproved LSBLSa—original phase b—wrapped phase c—unwrapped phase by LSBLS d—comparison between unwrapped phase by LSBLS and original phase

圖1a是原始相位，圖1b是原始包裹相位，圖1c是用LSBLS算法解包裹的結(jié)果，圖1c和圖1a相比可以直觀看出LSBLS解包裹誤差很大，圖1d是坐標(biāo)y為256pixel時的切面輪廓圖，從圖1d可以看出，LSBLS解包裹誤差極大。綜合圖1可以看出，LSBLS算法在原始相位變換太快時，欠采樣帶來的影響較為嚴(yán)重，相位解包裹失敗。

下面采用本文中改進(jìn)的LSBLS相位解包裹方法，為了對比，同樣采用100倍的peaks函數(shù)進(jìn)行實驗，改進(jìn)后的LSBLS算法解包裹實驗如圖2所示。

Fig.2 Results of improved LSBLSa—unwrapped phase by improved LSBLS b—comparison between unwrapped phase by improved LSBLS and original phase

從圖2可以看出，改進(jìn)后的LSBLS算法不僅能夠解出其包裹相位，而且改進(jìn)后的LSBLS算法誤差較小。另外，改進(jìn)后的LSBLS算法在運行時間上對原有的LSBLS算法運行時間沒有影響，改進(jìn)前后的LSBLS算法的運行時間都為3s左右。

在現(xiàn)實情況中由相移技術(shù)得到的包裹相位都存在一定的噪聲，本文中根據(jù)現(xiàn)實情況對有噪聲的包裹相位進(jìn)行了研究，對有噪聲的欠采樣相位解包裹問題的研究如圖3所示。

圖3是對有噪聲的包裹相位，通過窗口傅里葉濾波［16-18］后用改進(jìn)的LSBLS方法解包裹的結(jié)果。圖3c和圖3d、圖3b相比可以看出，原有的LSBLS算法解包裹失敗了，而改進(jìn)后的LSBLS算法能夠還原原始相位。圖3e是坐標(biāo)x=256pixel時的剖面輪廓圖，圖3e更是明確地說明了改進(jìn)后的LSBLS算法對欠采樣較嚴(yán)重的包裹相位能夠很好地解出真實相位。通過一系列的實驗可知，本文中改進(jìn)的LSBLS算法對欠采樣問題的解決有較大的實用性。

Fig.3 Results of phase unwrapping of the undersampling wrapped phase with noisea—wrapped phase with noise b—original phase c—unwrapped phase by LSBLS d—unwrapped phase by improved LSBLS e—comparison between two LSBLS algorithms and original phase

4 結(jié)論

相位變化太快而發(fā)生欠采樣，導(dǎo)致不能順利完成相位解包裹。欠采樣問題一直是相位解包裹中的一個難點，盡管LSBLS算法有一定的抗欠采樣能力，但當(dāng)欠采樣嚴(yán)重時，該算法同樣不能得到滿意的結(jié)果，甚至?xí)獍?。本文中依?jù)現(xiàn)有的LSBLS算法，對LSBLS算法加以改進(jìn)，并通過實驗驗證了改進(jìn)后的LSBLS算法的正確性和有效性。

［1］ XU X，WANG K F，GU G Q，et al.Measurement of internal material flaws based on out-of-plane displacement digital speckle pattern interferometry［J］.Laser Technology，2012，36(4):548-552(in Chinese).

［2］ XIE M M，TAO B Q，JIANG J H.Nondestructive inspection technology based on grating big shearing digital speckle interference［J］.Laser Technology，2001，25(4):279-282(in Chinese).

［3］ ZHANG Y P，WANG K F.Application of LabVIEW and MATLAB in ESPI image processing［J］.Laser Technology，2009，33(6):582-589(in Chinese).

［4］ QIAN X F，RAO F，LI X H，et al.Accurate least-squares phase unwrapping algorithm［J］.Chinese Journal of Lasers，2012，39(2):0209001(in Chinese).

［5］ WEI Z Q，XU F，JIN Y Q.Phase unwrapping for SAR interferometry based on an ant colony optimization algorithm［J］.International Journal of Remote Sensing，2008，29(3):711-725.

［6］ QIAN X F，ZHANG Y A，LI X Y，et al.Phase unwrapping algorithm based on mask and least-squares iteration［J］.Acta Optica Sinica，2010，30(2):440-444(in Chinese).

［7］ GUO Y，CHEN X T.Study of improved phase unwrapping algorithm based on least squares［J］.Chinese Journal of Lasers，2014，41(5):0508005(in Chinese).

［8］ GHIGLIA D C，ROMERO L A.Robust two-dimensional weighted and unweighted phase unwrapping that uses fast transforms and iterative methods［J］.Journal of the Optical Society of America，1994，A11(1):107-117.

［9］ PRITT M D，SHIPMAN F S.Least-squares two-dimensional phase unwrapping using FFT’s［J］.IEEE Transactions on Geoscience and Remote Sensing，1994，32(3):706-708.

［10］ LIU J F，LI Y Q，LIU K.Weighted least-square phase unwrapping algorithm for noisy phase-map［J］.Acta Optica Sinica，2008，34(5):643-646(in Chinese).

［11］ DONG Ch M，LIU G L，YU Sh W，et al.An improved algorithm of weighted least squares for phase unwrapping［J］.Journal of Geodesy and Geodynamics，2011，31(s1):90-93(in Chinese).

［12］ FAN Q，YANG H R，LI G P，et al.Method for phase recovery from a single undersampled wrapped phase map［J］.Acta Optica Sinica，2011，31(3):0310002(in Chinese).

［13］ ZHANG X，QIAN X F.An improved on least-squares phase unwrapping algorithm for undersampled interferofram［J］.Acta Photonica Sinica，2011，40(1):121-125(in Chinese).

［14］ WANG H Y，ZHANG Z H，ZHAO B Q，et al.Unwrapping algorithm of under-sampled wrapped-phase map［J］.High Power Laser and Particle Beams，2012，24(10):2311-2317(in Chinese).

［15］ LI B，QIAN X F，LI X H，et al.Phase-unwrapping algorithm based on radial shearing principle［J］.Laser Technology，2013，37(1):44-47(in Chinese).

［16］ LI K，QIAN K M.Dynamic phase retrieval in temporal speckle pattern interferometry using least squares method and windowed Fourier filtering［J］.Optics Express，2011，19(19):18058-18066.

［17］ QIAN K M.Two-dimensional windowed Fourier transform for fringe pattern analysis:Principles，applications and implementations［J］.Optics and Lasers in Engineering ，2007，45(2):304-317.

［18］ QIAN K M.Windowed Fourier transform for fringe pattern analysis［J］.Applied Optics，2004，43(13):2695-2702.