廖宇鋮，伍世虔,4，鄧高旭，陳 彬,4

（1.武漢科技大學(xué) 信息科學(xué)與工程學(xué)院，湖北 武漢 430081；2.武漢科技大學(xué) 冶金裝備及其控制教育部重點(diǎn)實(shí)驗(yàn)室，湖北 武漢 430081；3.武漢科技大學(xué) 冶金自動(dòng)化與檢測技術(shù)教育部工程研究中心，湖北 武漢 430081；4.武漢科技大學(xué) 機(jī)器人與智能系統(tǒng)研究院，湖北 武漢 430081）

引言

相位展開是從(-π,+π]的包裹相位中獲得連續(xù)相位的過程，廣泛應(yīng)用于干涉合成孔徑雷達(dá)、干涉合成孔徑聲納、光學(xué)干涉和相位測量輪廓術(shù)[1-3]。在理想情況下，相鄰2 個(gè)像素點(diǎn)的包裹相位差不會(huì)超過2 π。這時(shí)，展開相位可以通過任意路徑的積分所獲得。事實(shí)上，展開相位的正確性始終會(huì)受到相位的不連續(xù)、噪聲和陰影等因素的影響。至今為止，相位展開仍然是一個(gè)具有挑戰(zhàn)性的問題[4-5]。

相位展開可分為時(shí)間相位展開法和空間相位展開法[6-9]。時(shí)間相位展開法可以有效地展開具有較大不連續(xù)性和孤立區(qū)域的包裹相位，但是需要在不同時(shí)刻投影多幅不同頻率的條紋圖案，因此會(huì)花費(fèi)很多時(shí)間[10-12]?？臻g相位展開法根據(jù)包裹相位圖的一些空間特征和2 π相位跳變等信息就能進(jìn)行相位展開，其速度快、應(yīng)用范圍廣[13-14]。空間相位展開算法大致可分為路徑跟蹤算法和最小范數(shù)法。路徑跟蹤算法高效且魯棒，它依賴不同的方法來保證展開路徑的可靠性。質(zhì)量引導(dǎo)算法是路徑跟蹤算法中具有代表性的一種算法，它根據(jù)相位的質(zhì)量圖來引導(dǎo)展開路徑，其中幾種常見的質(zhì)量圖有：偽相關(guān)系數(shù)圖、相位導(dǎo)數(shù)方差圖、最大相位梯度圖等。由于該方法是通過路徑展開相位，一旦1 個(gè)像素展開出錯(cuò)，誤差則會(huì)沿著路徑傳播，從而導(dǎo)致相位展開失敗。最小二乘法是最小范數(shù)法中經(jīng)典的算法，其原理是在全局范圍內(nèi)最小化代價(jià)函數(shù)，從而獲得連續(xù)的展開相位[15-17]。Ghiglia等[18]提出了一種最小二乘相位展開算法，并通過快速余弦變換展開相位。Lu 和Wang[19]提出了相位導(dǎo)數(shù)方差相關(guān)圖，該圖能夠很好地反映包裹相位的質(zhì)量，并將其作為加權(quán)系數(shù)帶入最小二乘法，得到了平滑的相位。Xia 等[20]提出了相位校準(zhǔn)的迭代最小二乘算法，可以抵抗較大的散斑噪聲。然而在實(shí)際應(yīng)用中，由于物體的高度或者其他因素，相位會(huì)突然變化從而產(chǎn)生相位的不連續(xù)。而這些算法在處理不連續(xù)相位時(shí)總會(huì)在不連續(xù)的區(qū)域上出現(xiàn)錯(cuò)誤。

針對(duì)不連續(xù)的問題，本文提出了快速、簡單的算法。該算法利用結(jié)構(gòu)張量計(jì)算出每個(gè)像素點(diǎn)的方向角，將其映射到三維空間。通過計(jì)算找到不連續(xù)的區(qū)域，并將其作為加權(quán)系數(shù)帶入到最小二乘法來展開相位。實(shí)驗(yàn)證明該方法能快速地識(shí)別不連續(xù)位置，能準(zhǔn)確地獲得正確的展開相位。

1 相關(guān)工作

1.1 加權(quán)最小二乘相位展開算法

包裹相位 ψi,j和展開相位 ?i,j之間的關(guān)系為

式中：k為整數(shù)；W{·}為包裹算子。

相位展開是為了將 ψi,j展開為 ?i,j，加權(quán)最小二乘的原理是在全局范圍內(nèi)使展開相位差值和包裹相位差值的差異最小，即：

式中：wi,j為二值權(quán)重系數(shù)，表示該點(diǎn)的可靠度。

對(duì)公式（2）求偏導(dǎo)，可以得到如公式（5）所示的離散泊松方程：

對(duì)公式（5）利用預(yù)處理共軛梯度算法（preconditioned conjugate gradient，PCG）進(jìn)行迭代相位展開。對(duì)于加權(quán)最小二乘法，權(quán)重系數(shù)圖的質(zhì)量直接決定了相位展開的精度和速度，本文提出了計(jì)算權(quán)重系數(shù)圖的新算法。

1.2 基于方向圖變換的權(quán)重系數(shù)圖求解

1.2.1 方向圖求取

為了減少傳感器和其他影響導(dǎo)致的圖像缺陷，將圖像進(jìn)行歸一化。G(x,y)表示點(diǎn)(x,y)歸一化后的灰度值，圖像歸一化如公式（7）所示：

式中：μ和 σ2分別表示圖像f的平均灰度值和方差；μd和分別表示期望的平均值和方差(一般來說，。

圖像的梯度 ?g對(duì)邊緣檢測有很大的作用，但不適合尋找平行結(jié)構(gòu)或者方向。為了拿到包裹相位圖的方向，使用公式（8）所示的結(jié)構(gòu)張量進(jìn)行求解：

式中：?gσ是經(jīng)過高斯平滑后的圖像梯度。

公式（8）得到2×2 的矩陣，可簡化為

該矩陣為圖像的結(jié)構(gòu)張量，通常用來分析和提取圖像的局部信息，即局部范圍內(nèi)所表示的方向以及方向的相關(guān)程度。為了計(jì)算出準(zhǔn)確的方向，對(duì)結(jié)構(gòu)張量進(jìn)行增強(qiáng)，如公式（10）所示：

式中：E為圖像的能量；G(wk,wl)為傅里葉變換后的圖像，其關(guān)系如公式（11）所示:

對(duì)于增強(qiáng)后的結(jié)構(gòu)張量，計(jì)算出特征值 μ1和μ2(μ1≥μ2)。μ1和 μ2反映了圖像的局部結(jié)構(gòu)特征，其中特征值 μ2對(duì)應(yīng)的特征向量 ω2表示邊緣方向的估計(jì)，即相位條紋方向，如公式（12）和公式（13）所示：

用特征值 μ2對(duì)應(yīng)的特征向量 ω2計(jì)算圖像中每個(gè)點(diǎn)的方向角θ(x,y)，如公式（14）所示：

計(jì)算出每個(gè)點(diǎn)的方向角θ(x,y)，最終獲得方向圖θ。

1.2.2 不連續(xù)位置識(shí)別

求出方向角之后，為了找出不連續(xù)位置，將方向圖θ按行生成一維向量 λ（其中 θi,j變換為 λk），對(duì)其中每個(gè)元素 λk進(jìn)行公式（15）所示的變換。然后求出相鄰2 點(diǎn)間的距離，如公式（16）所示：

式中：m和n分別表示方向圖的行和列。Δfx(λk)對(duì)方向的變化很敏感，其值隨著方向的變化而變化。求出Δfx(λk)之后，將其還原為矩陣得到Δxθi,j。然后將方向圖按列排列成向量，用同樣的方法求得Δyθi,j。在對(duì)方向圖進(jìn)行變化差分后，將設(shè)定一個(gè)閾值求得權(quán)重圖，如公式（17）所示：

式中：wi,j是圖像的權(quán)重因子；Δθi,j為與的和。值為1表示可靠度較高的像素，值為0表示可靠度較低的像素。在本文中，T的取值為之間。該權(quán)重因子即為包裹相位的權(quán)重系數(shù)，通過加權(quán)最小二乘迭代求得包裹圖像的展開解。

2 數(shù)值仿真結(jié)果

本節(jié)利用計(jì)算機(jī)仿真的條紋圖驗(yàn)證所提出算法的有效性，并和其他3 種算法進(jìn)行對(duì)比。分別是相位導(dǎo)數(shù)方差相關(guān)圖（derivative variance correlation map,DVCM）[18]、可靠性圖（reliability map,RM）[14]以及相位拉普拉斯導(dǎo)數(shù)方差（phase Laplace derivative variance,PLDV）[17]。如圖1所示，用計(jì)算機(jī)生成連續(xù)的傾斜相位圖像，其函數(shù)關(guān)系式為

圖1 計(jì)算機(jī)生成的連續(xù)相位圖Fig.1 Continuous phase diagrams generated by computer

式中：ampPhase為幅度；x和y分別為圖像像素的橫、縱坐標(biāo)；N表示圖像的大小。在本實(shí)驗(yàn)中，N為512，即圖像尺寸為512×512。生成連續(xù)的相位圖后，用MATLAB 的wrapToPi 函數(shù)將其包裹。

在本實(shí)驗(yàn)中，仿真實(shí)驗(yàn)將在無噪聲和低噪聲環(huán)境下進(jìn)行，然后使用4 個(gè)方法分別進(jìn)行不連續(xù)區(qū)域的檢測。

圖2（a1）和圖2（b1）分別是無噪聲和低噪聲條件下的包裹相位，圖2（a2）～圖2（a5）和圖2（b2）～圖2（b5）分別對(duì)應(yīng)無噪聲和低噪聲條件下RM、DVCM、PLDV 和本文所提出的方法檢測出的不連續(xù)點(diǎn)。

圖2 不連續(xù)區(qū)域的檢測Fig.2 Detection of discontinuous areas

從圖2 中可以看出，當(dāng)圖像不存在噪聲時(shí)，4 種算法均能檢測出不連續(xù)的區(qū)域。DVCM 和PLDV 算法檢測出的不是一條直線，而是離散的點(diǎn)或者線段。隨著噪聲的加強(qiáng)，不連續(xù)的邊緣也會(huì)被噪聲影響。RM 算法雖然檢測到了一條直線，但是其周圍出現(xiàn)了一些噪聲點(diǎn)（為了方便對(duì)比，將噪聲點(diǎn)以*標(biāo)記出來），會(huì)影響展開效果，從而導(dǎo)致解包的精度下降。DVCM 和PLDV 算法能檢測到的點(diǎn)逐漸稀疏，但大致上還能看出是一條直線。本文提出的方法在無噪聲和低噪聲條件下，可以準(zhǔn)確檢測出完整的不連續(xù)位置。

為了比較這4 個(gè)方法展開圖像的質(zhì)量，本文計(jì)算了不連續(xù)位置的檢測準(zhǔn)確度和展開圖像的均方誤差（RMSE）。對(duì)于不連續(xù)位置的檢測準(zhǔn)確度，提出了如公式（19）的平均誤差計(jì)算方法，用于計(jì)算檢測出的不連續(xù)位置和真實(shí)不連續(xù)位置的誤差。

式中：gi是檢測到的不連續(xù)點(diǎn)；pi是真實(shí)的不連續(xù)點(diǎn)；k為檢測出不連續(xù)點(diǎn)的總數(shù)。

RMSE表示計(jì)算所得展開圖像和真實(shí)展開圖像之間的差別，如公式（20）所示：

式中：?i為計(jì)算所得展開相位值；為理論展開相位值。RMSE的值越小，代表展開相位的質(zhì)量越好。在表1 中可以看到，所提出的方法在展開的相位圖上表現(xiàn)最好，在無噪聲和低噪聲的條件下RMSE都在0.4 以下。

表1 不同方法的測量誤差Table 1 Measurement error of different methods

圖3（a1）和圖3（b1）分別是無噪聲和低噪聲條件下的展開相位；圖3（a2）～圖3（a5）和圖3（b2）～圖3（b5）分別是RM、DVCM、PLDV 和本文所提出的方法在無噪聲和低噪聲條件下所檢測出的不連續(xù)點(diǎn)；圖3（a6）～圖3（a9）和圖3（b6）～圖3（b9）分別是對(duì)應(yīng)的展開相位。

圖3 不同算法的不連續(xù)點(diǎn)檢測和展開相位圖Fig.3 Detection of discontinuous points and unwrapping phase diagrams of different algorithms

無噪聲下檢測的不連續(xù)圖和展開效果圖如圖3（a1）～圖3（a9）所示，低噪聲下的不連續(xù)圖和展開效果圖如圖4（b1）～圖4（b9）所示。圖3（a3）、圖3（a4）分別是WLS-DVCM、WLS-PLDV 計(jì)算出的不連續(xù)點(diǎn)，如表1所示，其不連續(xù)區(qū)域的誤差比提出的方法要小很多。但是在實(shí)際展開的圖上，WLS-PLDV 和WLS-DVCM 的矩形區(qū)域相位跳變很小，灰度變化不明顯。由于最小二乘法旨在最小化目標(biāo)函數(shù)，因此矩形區(qū)域線條斷開的地方圖像的相位會(huì)受到影響，導(dǎo)致相位跳變過小，進(jìn)一步影響周圍的點(diǎn)。WLS-RM 和本文提出的方法都檢測出了完整的矩形。

圖4 是展開相位圖在第160 行所取得的相位值。其中黑色的曲線代表原始相位圖上的變化趨勢，而藍(lán)色的曲線表示本文的方法所展開的相位值，其值和變化趨勢都非常接近原始相位。而WLS-PLDV、WLS-DVCM 和WLS-RM 都無法正確地展開包裹相位。

圖4 第160 行的展開相位值Fig.4 Unwrapped phase value on line 160

3 實(shí)驗(yàn)過程

為了驗(yàn)證提出算法的有效性，使用DLP 投影儀（Acer K137i DLP projector）和NiKON D810 相機(jī)組成的結(jié)構(gòu)光對(duì)籃球和膠帶進(jìn)行拍照，其裝置圖如圖5所示。

圖5 結(jié)構(gòu)光裝置Fig.5 Structured light device

用計(jì)算機(jī)生成分辨率為800×1 280 pixels 的正弦光柵條紋。圖6（a）和圖6（c）為拍攝到的物體圖，圖6（b）和圖6（d）是計(jì)算得到的包裹相位圖。

圖6 物體的光柵圖片和包裹相位圖Fig.6 Raster images of objects and wrapped phase diagram

由于籃球和膠帶有一定的高度，條紋在投射到平面和物體上時(shí)，在籃球的曲面和膠帶側(cè)面上將會(huì)產(chǎn)生一定的彎曲，這樣就會(huì)造成相位的不連續(xù)。用本文提出的算法對(duì)其進(jìn)行相位展開，其展開效果圖如圖7（a）和圖7（c）所示。圖7 中可以看出，籃球所在區(qū)域的相位都會(huì)發(fā)生變化，整體有一個(gè)跳變。圖7（b）和圖7（d）是將其他算法和本算法計(jì)算出的展開相位的對(duì)比圖。

圖7 展開相位圖的比較Fig.7 Comparison of unwrapping phase diagrams

在圖7（b）和圖7（d）中，矩形中的區(qū)域表示了相位的不連續(xù)，在此區(qū)域中，相位的變化應(yīng)和包裹相位有差異。而圖7（b）和圖7（d）中WLS-RM、WLS-DVCM 和WLS-PLDV 的變化與包裹相位的變化一致，表示沒有檢測出該地方的不連續(xù)。圖7（d）中箭頭處說明WLS-PLDV 進(jìn)行了錯(cuò)誤的跳變，從圖7（d）中可以看出，只有本文算法檢測到由膠帶表面高度調(diào)制對(duì)相位所帶來的變化。表2 給出了檢測出不連續(xù)位置所用的時(shí)間，可以看出所提算法消耗的時(shí)間最短。

表2 計(jì)算不連續(xù)區(qū)域所用的時(shí)間Table 2 Time used to calculate discontinuous areas

4 結(jié)論

包裹相位圖中的相位不連續(xù)會(huì)導(dǎo)致相位在不連續(xù)位置的展開存在較大誤差，為了解決這一問題，本文提出了基于方向圖變換的快速、簡單的加權(quán)最小二乘相位展開算法。利用包裹相位圖的結(jié)構(gòu)張量信息，精確地檢測出不連續(xù)的位置，從而提高了相位展開的精確度。計(jì)算機(jī)仿真實(shí)驗(yàn)和實(shí)際實(shí)驗(yàn)結(jié)果表明，本文提出的算法在速度與精確度上都有較好的結(jié)果。