基于單次傅里葉變換的分段衍射算法

胡 琪，王 喆，劉洪順，薛智文，鄧家春

(天津理工大學(xué) 理學(xué)院，天津 300384)

1 引 言

計算全息是在實驗全息術(shù)的基礎(chǔ)上發(fā)展起來的，利用惠更斯-菲涅爾衍射原理，通過對光場衍射行為的數(shù)值計算，可模擬光場在各種光學(xué)系統(tǒng)中的衍射行為。由于不受實驗環(huán)境因素的影響，這種分析方法具有較大的靈活性，經(jīng)常用來做復(fù)雜系統(tǒng)實驗的仿真模擬。隨著光場測量技術(shù)的不斷進步，計算全息在光場調(diào)控、光通訊、激光技術(shù)、衍射元件設(shè)計、光學(xué)系統(tǒng)開發(fā)等領(lǐng)域迸發(fā)出了巨大的活力[1-7]。在波前傳感器[8-11]取得巨大進步的今天，人們開始嘗試利用對物光波前的相位數(shù)據(jù)進行采集，并計算其通過光學(xué)成像系統(tǒng)后的成像結(jié)果，從而獲得層析圖像、深度圖像或者三維立體圖像，這些工作的開展都離不開較長距離衍射的精確計算[12-17]。因此，優(yōu)化較長距離的衍射算法，具有很大的研究價值。本文依據(jù)采樣定理[18]中對傳輸距離和衍射面畫幅尺寸的關(guān)系，提出了一種分段衍射算法，解決傳輸距離對衍射面畫幅尺寸和有效像素分辨率的影響，使算法可以靈活調(diào)整衍射面畫幅尺寸。數(shù)據(jù)顯示，計算結(jié)果與傳統(tǒng)算法結(jié)果一致，而分辨率可提高2～3個數(shù)量級。

2 基本理論

在利用菲涅爾衍射積分計算光場衍射行為的過程中，單次快速傅里葉變換算法(Single Fast Fourier Transform Algorithm，S-FFT)在衍射距離較長的情況下表現(xiàn)出失真率低、計算速度快、抗欠采樣性較強、內(nèi)存使用少等諸多優(yōu)勢[19]。但是，由于采樣定理的限制，此算法計算得到的衍射面畫幅尺寸受衍射距離影響較大，特別是在傳輸距離較長的情況下出現(xiàn)畫幅很大的情況，導(dǎo)致有效內(nèi)容的像素分辨率很低，不利于分析有效數(shù)據(jù)和應(yīng)用。本文基于采樣定理分析了上述問題的產(chǎn)生原因，并提出了通過分割衍射距離的方法解決上述問題。模擬計算表明，該方法行之有效，并保留了S-FFT算法的優(yōu)點。

根據(jù)惠更斯-菲涅爾原理，如圖1所示，衍射面上P0處的光場U(P0)為從孔徑S上各點發(fā)出的次級發(fā)散球面波exp(ikr01)/r01的相干疊加[9]。

(1)

圖1 在傍軸情況下衍射坐標示意圖 Fig.1 Schematic diagram of diffraction coordinates in paraxial

利用標量衍射理論[20]，通過求解麥克斯韋方程[21-22]，并運用傍軸近似將積分簡化，可以得到菲涅爾衍射積分的傅里葉變換形式[18]：

(2)

利用菲涅爾衍射積分，如圖1所示，我們可以計算得到初始光場U(x0,y0,0)傳播一段距離d之后等到的光場U(x,y,d)。

式(2)中設(shè)U(x0,y0)為物光平面波復(fù)振幅，U(x,y)為衍射面光波復(fù)振幅，傳播距離為d，則菲涅爾衍射積分可以表示為：

(3)

(4)

該算法稱為S-FFT算法(single fast Fourier transform algorithm)[19-20]。由于僅運用了一次傅里葉變換，其計算速度和失真率均較低。

另外計算光場衍射行為較為常見的還有D-FFT算法(double fast Fourier transform algorithm)和T-FFT(triple fast Fourier transform algorithm)算法。經(jīng)過實驗比較，三種算法各有所長。其中D-FFT算法和T-FFT算法適用于較短距離衍射范圍，其衍射面畫幅尺寸等于初始畫幅尺寸，二者之中D-FFT算法具有更強的抗欠采樣性，而T-FFT算法則具有更強的衍射圖案細節(jié)描繪能力。S-FFT算法更加適用于較長距離的衍射范圍，彌補了另外兩種算法在長距離上面的不足[19，23]。而本論文進一步優(yōu)化了S-FFT算法在長距離衍射方面的表現(xiàn)出的衍射面畫幅尺寸相比衍射圖案過大、衍射圖案分辨率低等缺點，為能夠完滿地獲得整個菲涅爾衍射區(qū)域的衍射場做出了一點貢獻。

3 采樣定理對S-FFT算法中衍射面畫幅尺寸的影響

下面討論一下衍射面觀察屏的最大坐標，以x軸方向為例。圖2中繪制了單條光線在x-z平面中傳播的過程，其傳播方向為波矢k的方向，與x軸夾角為α。圖中，x0軸為次級波源所在x-y平面的橫軸坐標，x軸為衍射面觀察屏所在平面的橫軸坐標；虛線為光場的等陣面；光線從次級波源O點(x0=0,z=0)發(fā)出，傳輸?shù)窖苌涿嬗^察屏上的P點(x=x,z=d)。根據(jù)等陣面可知，光場在x方向的空間周期為X。

由圖2可以得到光傳播方向與空間頻率的關(guān)系

(5)

(6)

圖3 衍射計算中衍射面觀察屏再現(xiàn)范圍示意圖 Fig.3 Sketch map of viewing screen reproduction range in diffraction calculation

根據(jù)式(6)可以知，基于S-FFT算法計算衍射光路時，衍射觀察屏的大小與采樣數(shù)、光的波長、衍射距離以及初始衍射面大小(即衍射圖案大小)4個因素相關(guān)。在實際中，由于波長、采樣數(shù)、初始衍射面大小均由實驗設(shè)備決定，衍射面觀察屏的大小只與衍射距離有關(guān)，無法改變。當平行光照射物體，衍射距離過長時，衍射面畫幅尺寸過大，會得到一個很小區(qū)域的衍射圖像；當球面波照射物體，并使用凸透鏡成像時，衍射距離過長會導(dǎo)致衍射面畫幅尺寸過小，只能得到圖像的某一個區(qū)域。這會使該算法在實際運用中，受到極大限制并嚴重影響圖像中有效畫面的像素質(zhì)量。圖5模擬計算了波長為532 nm的均勻平面光波垂直照射一個邊長為5 mm的初始衍射面，中央開有邊長為2.5 mm的正方形透光孔(如圖4所示)后，經(jīng)過不同距離d后在衍射面觀察屏上的衍射光場及其光強分布。

圖4 程序模擬的透光孔 Fig.4 Transparent hole simulated by program

圖5 在不同衍射距離下的衍射圖樣 Fig.5 Diffraction patterns at different diffraction distances

由圖5可以看出，在默認初始衍射圖案大小變化不大的情況下，隨著衍射距離d的變化，衍射面畫幅尺寸的區(qū)域變化很大。當衍射距離d過長時，衍射再現(xiàn)結(jié)果相比較衍射面畫幅尺寸來說過小，很難觀測。

為了解決上述問題，本文提出一種分段衍射算法。該方法是將衍射距離d劃分為兩段，分別記為d1、d2。這樣，在每次衍射計算中，衍射面畫幅尺寸在x和y方向的尺寸(Lx1、Ly1、Lx2、Ly2)都會因波長λ、采樣數(shù)(M、N)、衍射距離(d1、d2)和初始衍射面大小(Lx0、Ly0、Lx1、Ly1)而改變。其中，衍射d1距離衍射面畫幅尺寸為:

(7)

繼續(xù)衍射d2距離的衍射面畫幅尺寸為:

(8)

不難看出，衍射面畫幅尺寸的決定因素就由d1，變成了d2/d1。

因此，該方法可有效、靈活的通過調(diào)整d1、d2的比值來控制最終衍射面畫幅尺寸的大小，獲得更高分辨率的有效像素數(shù)據(jù)，并保留S-FFT算法的高精度及低失真率。

4 模擬計算結(jié)果與分析

圖6和圖7分別計算了波長為532 nm的均勻平面光波垂直照射一個邊長為5 mm中央開有邊長為2.5 mm的正方形透光孔(如圖4所示)的初始衍射面，經(jīng)過d=1 m(圖6)和d=10 m(圖7)后，衍射面觀察屏上的衍射光場及其光強分布。圖中分別按照直接計算(圖6(a)、圖7(a))和分段衍射算法(圖6(b)～(d)、圖7(b)～7(d))計算得到的觀察屏光強分布。

圖6 衍射1 m后的衍射圖像 Fig.6 Diffraction patterns after 1 meter of diffraction

圖7 衍射10 m后的衍射圖像 Fig.7 Diffraction patterns after 10 meters of diffraction

可以看出，分段衍射算法可以靈活有效的調(diào)整衍射面畫幅尺寸大小，有效數(shù)據(jù)的分辨率得到了調(diào)控。在傳輸距離為1 m時，由于式(7)衍射面畫幅大小為2 970.25 mm2，衍射圖案占據(jù)的范圍不到畫幅范圍的0.26%(以下簡稱“占比”)，有效內(nèi)容像素只有676個。而利用分段衍射算法(以d1=d2=0.5 m為例)，根據(jù)式(8)可將畫幅調(diào)整為25 mm2，使得衍射圖案占據(jù)總畫幅的43.58%，有效內(nèi)容像素為114 244個，數(shù)據(jù)量增加了169倍。同樣，在衍射距離為10 m時，傳統(tǒng)S-FFT算法得到的衍射面畫幅為296 807.04 mm2，其中有效像素僅為9個(占比0.34%)，而通過分段衍射算法(以d1=d2=5 m為例)實現(xiàn)的S-FFT算法得到的畫幅為25 mm2，有效像素仍保持在154 449個(占比58.77%，數(shù)據(jù)量增加了17 285倍)，極大的提高了有效像素的數(shù)量。同時，通過表1和表2可知，其計算結(jié)果與傳統(tǒng)算法保持一致。

表1和表2分別分析并統(tǒng)計了直接計算和分段衍射算法兩種方法獲得的衍射面畫幅尺寸相應(yīng)比例圖像的失真度。其方法是(以d=1 m，d1=d2=0.5 m為例)首先根據(jù)衍射面畫幅尺寸(直接計算的衍射屏邊長為La=54.5 mm，分段衍射算法計算的衍射屏邊長為Lb=5 mm)比例截取其中主要能量部分并作歸一化(所有像素點數(shù)值除以像素中最大數(shù)值)，所得到圖像作為對比圖，把直接計算結(jié)果作為參考圖，將對比圖中的對應(yīng)像素點做加權(quán)平均后形成一張與參考圖像素個數(shù)相同(x像素數(shù)為M，y像素數(shù)為N)的比照圖，然后用比照圖fe(x,y)和參考圖f(x,y)中的對應(yīng)像素點數(shù)值對比計算均方根誤差erms和峰值信噪比Rsn[24]。其計算方法為：

(9)

(10)

表1 衍射距離d=1 m時，直接計算與分段衍射算法計算圖像的erms和Rsn

由表1可見，兩幅圖的均方根誤差均小于0.1，而峰值信噪比均大于20，兩幅圖的保真度較高。其誤差產(chǎn)生的主要原因，除了對應(yīng)像素點個數(shù)不是整數(shù)之外，還應(yīng)考慮到分辨率較高的對比圖本身細節(jié)增多，其平均強度值與參考圖中對應(yīng)像素的數(shù)值理應(yīng)存在一定偏差。這種偏差其實也反映了圖像分辨率增加后導(dǎo)致圖案細節(jié)增多。因此在一定程度上，本問題中分析失真度時，誤差并非越小越好。

同樣的，運用上述方法在表2中總結(jié)了衍射距離為10 m時，兩種計算方法獲得結(jié)果的均方根誤差和峰值信噪比。從中可以看出，與d=1 m時結(jié)果類似，對照圖和參考圖相似度也很高。而由于在較長衍射距離的情況下，直接計算得到的衍射面畫幅尺寸過大，導(dǎo)致有效內(nèi)容的像素點僅有9個，致使誤差比前種情況大。根據(jù)前面的分析，這其中有較大部分原因是由于對照圖中分辨率的提升導(dǎo)致圖案細節(jié)更加豐富，而在對比時使用了對比圖中對應(yīng)像素點的平均值，制作與參考圖像素個數(shù)一致的對照圖，因此背景噪聲有所增加，在超過一定分辨率放大率的情況下，均方根誤差應(yīng)該有所升高而峰值信噪比應(yīng)該下降，這是提高分辨率增加有效像素數(shù)量的一種體現(xiàn)。

表2 衍射距離d=10 m時，直接計算與分段衍射算法計算圖像的erms和Rsn

5 結(jié) 論

S-FFT算法中，由于采樣定理的約束，衍射距離將嚴重限制衍射面畫幅尺寸，在較長衍射距離上有效數(shù)據(jù)分辨率很低。由于實際應(yīng)用中采樣率、波長、原始衍射面尺寸等因素均需要根據(jù)實驗確定，因此直接計算無法解決上述問題。本文模擬計算結(jié)果表明，分段衍射算法可利用采樣定理實現(xiàn)對衍射面畫幅尺寸的靈活調(diào)整，從而使上述問題得以解決。如傳輸距離為1 m時，畫幅大小為2 970.25 mm2，占比約為0.26%，有效像素只有676個。而利用分段衍射算法(以d1=d2=0.5 m為例)，畫幅調(diào)整為25 mm2，占比提高到43.58%，有效像素為114 244個，數(shù)據(jù)量增加了169倍。造成誤差的一個主要原因來自有效數(shù)據(jù)分辨率提高后，細節(jié)分布與低分辨率像素值之間的差別。在圖像細節(jié)較豐富時，其差別較大。因此這種差別應(yīng)視為優(yōu)于直接計算的一種結(jié)果。波前采集進而對立體物體進行層析成像計算是課題組接下來研究的重點方向，而提高圖像細節(jié)的有效數(shù)據(jù)分辨率，靈活控制衍射面畫幅尺寸，對該工作的研究具有重要的基礎(chǔ)作用。結(jié)果表明，分段衍射算法對于提高優(yōu)秀算法的適用衍射距離，具有較高的推廣價值。