大面陣航攝影像重建高效像素插值算法研究

胡海彥，余旭初，楊韞瀾，方 勇，江振治

(1.信息工程大學,河南 鄭州 450001；2.西安測繪研究所，陜西 西安 710054)

復合型大面陣或超大面陣航攝相機提供給作業用戶的最終影像，其每一個像素都是從原始子影像重采樣而來，就計算耗時而言，根據經驗圖像重采樣所耗費的時間至少要占到整個大面陣航攝影像預處理的1/3以上，所以對此環節進行耗時優化，必然會極大提升大幅面影像生成效率[1-3]。圖像內插重建是基于一組等間距離散采樣點恢復連續空間圖像函數的過程，其是圖像縮放、子像素平移(非整像素)、圖像旋轉、圖像形變等基本圖像操作的技術基礎。當這些圖像操作所在位置缺少既有圖像灰度值，即需要圖像插值進行灰度填補。而且對于大幅面數字航攝相機而言，在影像復合過程中可通過高精度像素重采樣一并解決鏡頭畸變引起的圖像幾何失真、色彩配準等圖像重建問題。由采樣理論可知，sinc函數是最為理想的圖像插值內核函數(見圖1)。然而，實際應用中并不能直接使用該函數進行插值[4]，需要有限截斷并實現離散化計算，同時要在計算速度與數學精度之間取得平衡，即在盡可能接近sinc函數理論精度的前提下，設計出一個內核函數對sinc進行函數逼近。通常可以通過n階對稱分段多項式進行內核函數逼近，如常用的雙線性卷積核函數等即為其中一種特例，文獻[5]對此類內核函數進行了定量評估，相關研究也比較豐富[6～9]，但往往給出的插值算法精度接近像素級，也沒有挖掘此類算法的潛在效率，僅適用于圖像融合等某些常規應用。本文從一維信號重構原理著手，推導出n階對稱分段多項式內核函數的數學形式，并給出多項式定義域變換后的加速算法實現，最后對所提算法實作性能進行試驗分析。算法在盡可能不引入額外像素插值誤差，從而避免損失或破壞影像自身固有幾何量測精度前提下，大幅度降低數據處理計算復雜度，滿足大幅面數字航攝影像重建對效率的要求。

圖1 空間域中的sinc函數

1 像素插值算法

1.1 內核函數的多項式逼近

對于一維空間，近似sinc內核函數的n階對稱分段多項式在區間[-(n-1),(n-1)]上的定義為

(1)

其中：i=0,1,…,m-1，m表示整數截斷范圍，滿足n=2m-1。這樣總共會有(n+1)m個aij系數。由于內核函數具有對稱性及維度的可分離性，即二維內核函數為兩個一維內核函數的乘積，即h(x,y)=h(x)h(y)，所以為了描述的簡潔性，本文均在一維情況下予以討論。

為了保證節點處函數擬合的連續性與光滑性，除保證多項式與sinc函數在節點上的取值相同外，還需要對式(1)附加導數連續性約束條件。以下將3階對稱分段多項式內核函數作為樣例，給出推導過程，更高階多項式類同。在區間[-2,2]上sinc函數的近似3階對稱多項式內核函數為

hC(x)=

(2)

由約束條件(2)總共可列出7個方程，其中含8個未知數系數。進一步令a13=α作為調整參數，則未知數個數減1，使得方程組可解。有

a03=(α+2)，a02=-(α+3)，a01=0，a00=1

a13=α,a12=-5α，a11=8α，a10=-4α.

表1 平坦性約束下的3～9階多項式α參數

可以看出不同的α值，對應不同的插值計算，常用的3次卷積對應α=-1/2，也稱之為基斯(Keys)插值內核函數[10](見圖2)。

hkeys=

(3)

圖2 基斯3次卷積內核函數

1.2 多項式定義域轉換

內插權重計算即是利用分段多項式求得采樣點集上對應的多項式函數值，以此求得卷積核模板，與采樣點進行卷積運算便完成插值點函數插值。所以需要先行計算插值點到采樣點的相應距離ξ。在3次卷積核的情況下，這些距離值分別為ξ，1-ξ，1+ξ及2-ξ，其中ξ∈[0,1]。考慮計算效率與便捷性，可將所有多項式的定義域區間變換為區間[0,1]，利用多項式的分段對稱性(圖2(a)已給出圖示)，可對式(1)中的多項式p0(x)在ξ=x,ξ=1-x及p1(x)在ξ=x+1,ξ=2-x上進行分段處理。對于更高的n階多項式內核函數，實現的策略是相同的，這樣可以得到n+1個相同定義域的分段多項式，對于偶數i,即(0,2,4,6…)。令k=i/2，即(0,1,2,3…)。則有

(4)

其中，

(5)

對圖2進行3次多項式定義域變換后對應的函數圖形如圖3所示。

圖3 定義域變換為[0,1]的Keys 3次卷積核

為了節約計算時間，給出關于變換后多項式系數屬性的5條規則。圖4為圖形化說明(以5階插值內核函數系數表為例，α=3/64；表2為3階插值內核函數系數，α=-1/2)。

圖4 5階插值多項式內核函數系數圖示

規則1：僅第一個多項式f0(x)的最低階項系數為1，其余多項式fi(x)最低階項系數ai0均為0；

規則2：多項式fi(x)和fi+1(x)，最高階項系數ain成對相等、符號相反；

規則3：所有j≤-2階項系數成對出現，且奇數階項大小相等符號相反、偶數階項大小相等且符號相同；

規則4：多項式fn(x)的所有j≤n-2階數項的系數anj為0；

規則5：多項式f0(x)的所有j≤n-2(j為奇數)階數項的系數a0j為0。

利用這5個規則，可以明顯提升計算效率。首先對于所有多項式,距離x具有相同的值，因此可以對x冪進行預先計算。而且成對相同的系數及與x冪相乘的積僅需計算1次，同時可跳過具有零系數的計算項。

表2 多項式變換后的3次卷積核系數

所以，式(3)給出的Keys 3次卷積核，根據上述理論進行多項式轉換，可等價得到[0,1]空間上的4個多項式為

(6)

2 算法性能分析及試驗結果

2.1 性能分析

表3給出一維3階多項式卷積核插值權重的計算次數。對于一維3階卷積，這里提出的方法總共只需16個計算操作，若替代傳統方法的26個計算操作，可使計算加速1.625。表4為n階多項式卷積核運算次數，隨著多項式階數的增加，加速倍率逐漸降低(例如9階多項式的加速倍率為1.52)，對于無窮階多項式，加速倍率穩定為4/3。由于分段n階卷積是可分離的(卷積級聯)，所以高維N-D卷積的計算耗時是N×1D卷積時間。

當然，要對整個插值操作進行評估時，除需要考慮插值多項式權重計算外，還有采樣點與這些權重的卷積計算。該卷積計算的時間復雜度由(n+1)N個加法、乘法和計算機內存訪問耗時共同決定。這也意味著當進行高維插值(N>2)時，卷積計算量占主導地位。但通過將多項式定義域變換為[0，1]區間后，對任意維度下的插值都具加速效果，只是此方法對于一、二維低維度插值而言更為高效。

表3 3階多項式一維卷積核運算次數

表4 n階多項式卷積核最大運算次數估計

2.2 試驗結果

用C++程序語言對上述高階多項式內核卷積進行算法實現。在配置為Intel Core i3-370M(2.4 GHz)的計算機上完成試驗測試。為了進一步加速高階多項式的計算效率，編程中同時采用查找表(LUT)數據結構對內核函數進行預先計算，為保證插值精度，查找表中將單位長度設置為10 000個網格點(即插值像素點位可精確到0.000 1)。在3階多項式情況下，由于內存訪問相對多項式直接計算耗時，所以采用查找表比直接計算反而略顯低效——構建查找表耗時18 982 ms，整個插值過程耗時51 193 ms，但更高階內核插值試驗表明，采用查找表仍是必要的，所以在編碼實現時建議分3階以下/以上多項式兩種情況進行程序設計。為了排除計算的偶然性，獲取插值時長的平均值，對指定影像總共反復執行100次插值處理，二維插值的試驗結果如表5所示。

表5 試驗結果 ms

可以看出，新的多項式計算耗時相比常規計算提速比為1.73倍，接近理論值1.625。整個二維3次卷積計算與傳統方法相比提速比為1.15倍。這里所用影像數據為西安測繪研究所研制的DMZ大面陣航攝相機實攝航空影像，raw數據為12個子影像，內插合成后的影像像幅為24 K×24 K(見圖5)。

圖5 DMZ大面陣航攝相機實攝航空影像內插合成

3 結 論

本文提出分段多項式核函數高效插值算法，數學形式簡潔統一。相比較于常規插值方法，新的3階多項式算法耗時提速比為1.625倍,對于任意高階多項式，耗時提速比上限為1.3倍。這說明新插值算法對常規插值算法具有一定的替代性。

對于大幅面數字航空攝像機而言，因為用于圖像重建的插值計算量十分巨大，因此，本文方法對更為快速高效地進行原始攝影數據預處理意義非凡，已對西安測繪研究所研制的國產DMZ大面陣航測相機進行相關試驗，由于航測相機要求的量測精度極高(最低子像元級)，利用文中的高階多項式內核函數，可保證所引入的插值誤差被忽略不計(實踐表明采用5階多項式已足夠有效)，這使得大幅面影像生成兼顧效率與精度，進而保障相機的使用效能。還需說明的是，所提方法不限于大幅面數字航空相機的子影像復合拼接應用，需要高精度、高效率的圖像插值應用場合，該算法都可發揮一定的作用。當然，可將此算法進一步結合當代GPU并行編程技術予以實現，使得實際作業效率提升更為顯著。

參考文獻：

[1] 郭海青.信息化航空攝影測量生產與管理系統建設[J].測繪通報,2014(11):53-72.

[2] 劉翠.大面陣CCD成像系統關鍵技術[J].科學與財富,2015(2):127-128.

[3] 葛明,許永森,沈宏海,等.航空相機大面陣CCD多自由度拼接方法研究[J].紅外與激光工程,2015(3):923-928.

[5] MEIJERING E,NIESSEN W,VIERGEVER M.Quantitative evaluationofconvolution-basedmethodsformedicalimageinterpolation[J].Medical Image Analysis,2001,5:111-126.

[6] 鐘寶江,陸志芳,季家歡.圖像插值技術綜述[J].數據采集與處理,2016(6):1083-1096.

[7] 楊文波.航空圖像超分辨率重構技術研究[D].長春：中國科學院長春光學精密機械與物理研究所,2014.

[8] 喬少華,李潤鑫,劉輝,等.基于統計量的加權函數圖像重建方法[J].傳感器與微系統,2017(9):53-56.

[9] 趙秀影,王洪玉,從福仲.亞像元幾何圖像成像與超分辨恢復算法[J].光學技術,2010,36(5):795-798.

[10] KEYS R G.Cubic convolution interpolation for digital image processing[C].IEEE Transactions on Acoustics,Speech,and Signal Processing,1981,29(6):1153-1160.

[11] Weisstein E W.Horner’s rule.Eric Weisstein’s World of Mathematics.2003c,http://mathworld.wolfram.com/HornersRule.html