童小華，葉 真，劉世杰，謝 歡，金雁敏，陳 鵬

高分辨率遙感衛(wèi)星在軍用和民用方面都具有極其重要的價值，其發(fā)展和應(yīng)用得到了各國的高度重視。衛(wèi)星平臺的姿態(tài)穩(wěn)定性，是決定高分辨率衛(wèi)星幾何定位精度和輻射質(zhì)量的關(guān)鍵因素之一，然而衛(wèi)星在軌運行期間，外部空間環(huán)境的變化和內(nèi)部機械的振動等會導致姿態(tài)存在不穩(wěn)定和周期性微小振動現(xiàn)象，即衛(wèi)星平臺顫振現(xiàn)象，如果不進行有效處理會大幅影響衛(wèi)星的測圖精度和成像質(zhì)量[1]。隨著成像傳感器分辨率的提高，平臺顫振對高分衛(wèi)星影像幾何精度的影響愈發(fā)顯著。因此，衛(wèi)星平臺顫振是高分辨率對地觀測系統(tǒng)和智能化攝影測量處理亟須解決的難題之一。

衛(wèi)星平臺震顫是高分辨率衛(wèi)星在軌運行普遍存在的復雜現(xiàn)象，國內(nèi)外學者在衛(wèi)星平臺顫振分析方面已開展了相關(guān)研究。目前，采用衛(wèi)星影像或產(chǎn)品數(shù)據(jù)來反演顫振信息已成為了熱門方法，與傳統(tǒng)的基于地面密集控制和姿態(tài)數(shù)據(jù)的探測方法相比，該方法具有便捷高效和探測能力強等優(yōu)點[2]?；谟跋耖g密集視差變化的視差法顫振探測已經(jīng)成功應(yīng)用于ASTER、PLEIADES-HR、UK-DMC、天繪一號和資源三號等國內(nèi)外衛(wèi)星和傳感器的顫振處理分析[3-7]。實現(xiàn)高精度穩(wěn)健的密集匹配是視差法顫振探測的基礎(chǔ)，選擇有效的匹配方法來獲取影像間顫振導致的相對視差變化是決定衛(wèi)星平臺顫振探測質(zhì)量的關(guān)鍵之一。相位相關(guān)匹配算法因其高精度、高計算效率和對比度不敏感等特點，已廣泛應(yīng)用于視差法顫振探測框架[5-7]。

相位相關(guān)是一種基于傅里葉變換的區(qū)域影像匹配方法，利用頻率域的影像信息和操作來精確估計影像間的偏移量，與歸一化互相關(guān)等其他常用影像相關(guān)方法相比，能夠保證更高的精度和穩(wěn)定度[8-9]，已經(jīng)在遙感影像配準[10]、窄基線立體重建[11]和地表動態(tài)監(jiān)測[12]等智能處理方面得到了廣泛關(guān)注和應(yīng)用。相位相關(guān)的理論基礎(chǔ)是傅里葉變換的平移特性，亞像素的相位相關(guān)算法研究主要分為兩類，第1類算法通過確定互功率譜逆傅里葉變換后精確的峰值位置來獲取偏移量，第2類算法通過直接估計互功率譜的線性相位差來獲取偏移量[12-13]。在第1類算法的研究方面，文獻[13]假設(shè)影像間的亞像素偏移通過存在整像素偏移的高分辨率圖像降采樣獲得，利用sinc函數(shù)近似Dirichlet核，理論上推導了亞像素計算的解析形式。文獻[14]根據(jù)一維sinc函數(shù)擬合推導出peak evaluation formula(PEF)方法，采用主相關(guān)峰及其兩側(cè)不同距離的位置組成的多組三元組，根據(jù)推導公式，通過非迭代的最小二乘法來估計亞像素偏移。在第2類算法的研究方面，針對非纏繞的相位角理論上是一個二維平面這一特性，文獻[15]提出利用最小二乘估計來擬合平面參數(shù)，平面的斜率對應(yīng)亞像素的偏移值。文獻[16]提出利用奇異值分解對互功率譜進行子空間分解，通過擬合奇異值向量對應(yīng)的相位角斜率來獲取亞像素偏移值。

視差法顫振探測需要選擇兩幅或多幅不同時間成像的影像來進行相互匹配，選擇的影像間由于成像時間差異或譜段差異不可避免地存在復雜的輻射差異。因此，保證匹配算法的輻射穩(wěn)健性對衛(wèi)星平臺顫振探測的質(zhì)量相當重要。相位相關(guān)匹配在頻域進行了歸一化處理并只利用了影像的相位信息，因此對影像強度變化特別是亮度和對比度的全局線性變化并不非常敏感。然而，關(guān)于局部和非線性等復雜輻射差異對亞像素相位相關(guān)匹配精度影響和相應(yīng)有效改進的研究非常少。因此，為了進一步提高亞像素相位相關(guān)算法對影像間輻射差異的穩(wěn)健性，本文從影像表達方式入手，提出一種基于局部頻率分析的改進亞像素相位相關(guān)算法，結(jié)合相位一致性和局部平均相位生成輻射不變的影像表達，從而保證亞像素相位相關(guān)算法在復雜輻射差異情況下的精度優(yōu)勢。將提出的亞像素相位相關(guān)算法融入視差法顫振探測框架，可以在存在噪聲和復雜輻射差異等干擾情況下精確獲取影像間的顫振視差，從而提高顫振探測的精度，擴展用于顫振探測影像的選擇，增加從影像中反演顫振信息的可信度。

1 基于亞像素相位相關(guān)的視差法顫振探測

1.1 基于局部頻率分析的亞像素相位相關(guān)算法

1.1.1 相位相關(guān)理論基礎(chǔ)

相位相關(guān)的理論基礎(chǔ)是傅里葉變換的平移特性，即空域下的平移在頻域下顯示為線性相位差。假設(shè)兩幅影像f(x,y)和g(x,y)存在水平方向偏移x0和垂直方向偏移y0，其對應(yīng)傅里葉變換表示為G(u,v)和F(u,v)，傅里葉變換的平移特性表示為

G(u,v)=F(u,v)exp-i(ux0+vy0)

(1)

式中，i表示虛數(shù)。計算兩幅影像間的歸一化互功率譜矩陣來提取相位信息

expi(ux0+vy0)

(2)

式中，*表示復共軛。在實際離散影像情況下，采用快速傅里葉變換，歸一化互功率譜進行逆傅里葉變換后理論上可以得到一個二維脈沖函數(shù)，通過獲取峰值的位置可以求取整像素級的匹配結(jié)果[13]。

1.1.2 輻射不變的影像表達

為了減少輻射差異對影像匹配和配準結(jié)果的影響，局部結(jié)構(gòu)特征常被采納加入影像匹配配準框架[17]。通過局部頻率分析計算局部相位信息被證實符合人類對影像特征的視覺感知[18]，相當于提取影像內(nèi)在的局部結(jié)構(gòu)特性，在不同光照條件和不同模態(tài)下保持一致，與梯度和局部熵等其他結(jié)構(gòu)特征相比，受影像強度變化和噪聲的影響相對更小[19]。因此，本文借鑒文獻[19—20]的方法通過局部頻率分析提取相位一致性和局部平均相位信息來構(gòu)造輻射穩(wěn)健的結(jié)構(gòu)影像表達。

(3)

在特定尺度n下，局部幅值和相位為

(4)

局部平均相位可以計算為

(5)

考慮噪聲補償、頻率寬度和相位差，相位一致性計算公式定義為[21]

(6)

為了分析影像的局部頻率信息，需要從一維情況擴展到二維?？梢岳脝窝菪盘柪碚搧韺崿F(xiàn)二維相位一致性的計算[19-20]。

結(jié)合局部平均相位和相位一致性可以得到輻射不變的影像表達為

(7)

式中，μ為加權(quán)參數(shù)，通過計算影像的熵來自適應(yīng)確定[19]，較低的圖像熵意味較少的紋理信息，因此需要較大的加權(quán)參數(shù)。前一部分通過計算相位一致性生成邊緣保持表達EPM，提供結(jié)構(gòu)信息，后一部分歸一化局部平均相位得到紋理保持表達TPM，提供紋理信息。如圖1所示，不同曝光設(shè)置拍攝的影像間存在顯著的輻射差異，而通過局部頻率分析生成的影像表達明顯減少了輻射差異。

圖1 不同曝光影像生成的輻射不變影像表達Fig.1 Illumination-invariant image representations from images with different exposures

1.1.3 改進的亞像素相位相關(guān)算法

將基于局部頻率分析生成的結(jié)構(gòu)影像表達融入亞像素相位相關(guān)匹配框架可以得到輻射穩(wěn)健的亞像素相位相關(guān)算法。生成新影像表達的改進方式是通用的，從本質(zhì)上減少輻射差異對亞像素偏移值估計的影響。因其高效率與高可信度，本文結(jié)合PEF算法[14]提出改進的亞像素相位相關(guān)算法。

PEF算法根據(jù)多組三元組擬合一維sinc函數(shù)來確定亞像素偏移值。根據(jù)文獻[13]的推導，亞像素的相位相關(guān)生成了一個二維Dirichlet核，分解為一維情況，可以用一維sinc函數(shù)來近似

(8)

式中，β表示峰值大??；r(x)為逆傅里葉變換后的相關(guān)函數(shù)。對于一個三元組r(p-d)、r(p)和r(p+d)滿足

(p-d+dx)·r(p-d)+(p+d+x0)·r(p+d)=

2(p+x0)·cos(πd)·r(p)

(9)

調(diào)整后可得

(10)

式中，p和d取值確定三元組的位置，根據(jù)相關(guān)函數(shù)峰值整像素位置附近的多組三元組利用最小二乘法擬合可估計出亞像素偏移值。

通過結(jié)合基于局部頻率分析的影像表達和相位相關(guān)匹配，本文方法充分利用了影像的局部和全局相位信息來確保在噪聲和輻射差異等干擾因素影響下亞像素偏移值估計的精度和穩(wěn)定性。提出的亞像素相位相關(guān)算法流程為：

(1) 對于輸入影像，根據(jù)式(7)計算輻射不變影像表達，附加Tukey窗減弱邊緣效應(yīng)的影響。

(2) 對生成的影像表達進行傅里葉變換，根據(jù)式(2)計算歸一化互功率譜矩陣。

(3) 利用高斯函數(shù)對互功率譜進行頻域濾波，減弱噪聲和混疊的影響，并進行逆傅里葉變換，計算二維相關(guān)函數(shù)。

(4) 采用奇異值分解計算秩一近似將相關(guān)函數(shù)分離成一維情況。

(5) 根據(jù)式(10)估計亞像素偏移值。

1.2 視差法顫振探測方法

理論上，如果空間分辨率和成像質(zhì)量足夠高，衛(wèi)星平臺任何微小的振動都會影響CCD傳感器的成像并記錄在生成的影像中。從影像中反演衛(wèi)星平臺顫振信息，可以采用視差法顫振探測方法，根據(jù)短時間相隔成像的線陣CCD間的視差變化來估計顫振引起的隨時間變化的像移[3-7]。視差法顫振探測方法大致可以分為密集匹配，視差顫振轉(zhuǎn)換和頻譜分析3個步驟。其中，高精度的密集匹配是視差法顫振探測的基礎(chǔ)，匹配精度決定顫振探測的能力。為了提高存在輻射差異和噪聲等影響情況下的匹配精度，保證衛(wèi)星平臺顫振探測的可信度和可用性，將提出的局部頻率分析亞像素相位相關(guān)算法融入常規(guī)的視差法顫振探測框架，來估計影像間的顫振視差。

改進的視差法顫振探測方法具體流程如下：

(1) 對 影像進行預(yù)處理，包括影像增強和初始配準。初始配準目的是為了消除了兩幅影像間較大的變形和差異，減少影像密集匹配時誤匹配點的可能性。

(2) 通過等間隔匹配的方式，利用提出的局部頻率，分析亞像素相位相關(guān)算法對參考影像和待匹配影像進行密集匹配，獲取待匹配影像各像素相對于參考影像的顫振視差圖。

(3) 根據(jù)計算的相關(guān)系數(shù)約束，剔除誤匹配點，并以行為單位，對該行影像的所有有效視差取均值，獲取顫振視差曲線。

(4) 通過單頻率分量的空域轉(zhuǎn)換方式或多頻率分量的頻域轉(zhuǎn)換方式實現(xiàn)視差到顫振的轉(zhuǎn)換[7]。

(5) 利用傅里葉頻譜分析的方法分析獲取的平臺顫振的頻率成分，計算顫振頻率和幅值等特性參數(shù)。

2 試驗與分析

本文的試驗部分包括以下兩部分內(nèi)容：

(1) 通過存在實際輻射差異的已知真實偏移的數(shù)據(jù)集，開展亞像素相位相關(guān)算法對比，驗證本文局部頻率分析亞像素相位相關(guān)算法的可信度。

(2) 采用改進的視差法顫振探測方法估計資源三號衛(wèi)星顫振特性，驗證其有效性。

2.1 亞像素相位相關(guān)算法對比

通過兩組存在實際輻射差異的影像數(shù)據(jù)集來驗證提出的亞像素相位相關(guān)方法在輻射變化情況下的效果提升，模擬生成影像間的亞像素偏移，以便直接計算不同方法的匹配精度，開展光源變化和曝光變化情況下的亞像素相位相關(guān)算法對比。

模擬試驗數(shù)據(jù)以Middlebury Stereo Datasets[22]中的兩幅實際輻射差異影像為基礎(chǔ)影像生成，對于每幅影像，分別通過3種不同曝光設(shè)置和3種不同光源配置來進行拍攝。為了分析光源變化的影響，固定拍攝的曝光設(shè)置，變化光源配置，共有3種組合：I1/I2，I1/I3和I2/I3。為了分析曝光變化的影響，固定拍攝的光源配置，變化曝光設(shè)置，同樣共有3種組合：E1/E2，E1/E3和E2/E3。對于每個組合，采用文獻[9]相似的低通濾波和降采樣方式來生成128對已知亞像素偏移值的模擬影像對，并添加微量的零均值高斯噪聲以求更加符合實際影像情況和分析不同干擾因素與輻射變化的綜合影響。圖2為模擬試驗的示例影像，根據(jù)I1/I3組合影像和E1/E3組合影像以x和y方向0.5個像素的偏移值生成。

圖2 模擬試驗示例影像Fig.2 Sample images in the simulated test

根據(jù)上述試驗數(shù)據(jù)，將提出的基于局部頻率分析影像表達的算法與其他主流的亞像素相位相關(guān)算法對比，以驗證提出算法的可信度。對比算法包括Foroosh算法[13]、PEF算法[14]、Stone算法[15]和Hoge算法[16]。對每一光照和曝光組合，通過總共128對影像根據(jù)真實偏移值來計算各方法x方向和y方向綜合匹配誤差的平均值和均方根作為評價測度。

圖3繪制了不同光源或曝光設(shè)置組合下，各亞像素傅里葉影像相關(guān)方法計算的絕對誤差的均值和均方根結(jié)果。從圖3中可以發(fā)現(xiàn)，無論是光源變化引起的局部輻射差異還是曝光變化引起的全局輻射差異，對亞像素傅里葉影像相關(guān)匹配的精度都有一定程度的影響，而且隨著輻射差異的增大，各方法計算的匹配誤差整體趨勢上增大。相對比其他4種主流的相位相關(guān)算法，本文提出的基于局部頻率分析的亞像素算法在不同影像組合下取得了一致的更優(yōu)和更穩(wěn)定的效果，匹配誤差的均值和均方根最小，驗證了基于局部頻率分析提取的輻射不變影像表達的可信度，以及用于減小輻射差異對亞像素相位相關(guān)算法影響的可行性。即使存在明顯輻射差異和噪聲等的影響，本文亞像素相位相關(guān)算法仍能達到0.05個像素的匹配精度。

圖3 模擬試驗亞像素相位相關(guān)算法對比Fig.3 Comparison of different subpixel phase correlation methods in the simulated test

2.2 資源三號衛(wèi)星平臺顫振探測

本文采用資源三號衛(wèi)星多光譜影像來進行平臺顫振探測試驗，試驗影像為2012年拍攝的河南區(qū)域影像。顫振探測方法采用本文提出的局部頻率分析亞像素相位相關(guān)算法，適用于采用存在明顯輻射差異影像的顫振探測，例如采用多光譜可見光波段和近紅外波段影像。擴展顫振探測影像的譜段，一方面可以完善顫振信息的提取，提供相互驗證手段，另一方面有助于避免視差法顫振探測的盲點頻率[4]。資源三號衛(wèi)星的多光譜傳感器包含紅(B1)、綠(B2)、藍(B3)和近紅外(B4)4個波段，本文選取B1與B4以及B1與B2的波段組合。

密集匹配后逐行平均得到的垂軌方向顫振視差曲線結(jié)果如圖4所示，從中可以看出，兩組波段組合計算得到的顫振視差曲線都出現(xiàn)了由顫振引起的隨時間周期性變化的波動，并且具有相似的變化趨勢。

圖4 垂軌方向衛(wèi)星顫振視差曲線結(jié)果Fig.4 Jitter disparity curve in the cross-track direction

為了驗證獲取顫振探測結(jié)果的準確性，通過視差顫振轉(zhuǎn)換和頻譜分析獲取顫振的頻率與幅值特性[7]，并與從姿態(tài)數(shù)據(jù)反演的顫振特性進行對比分析。表1列出了利用兩組波段組合計算的顫振頻率和幅值，以及從姿態(tài)數(shù)據(jù)反演的顫振頻率與幅值，可以看出資源三號衛(wèi)星存在頻率約0.66 Hz的周期性平臺顫振，從不同波段組合影像計算出的顫振頻率和幅值與從姿態(tài)數(shù)據(jù)得到的結(jié)果相似，表明了基于局部頻率分析亞像素相位相關(guān)算法的視差法顫振探測方法的有效性。

表1 顫振頻率與幅值結(jié)果對比

3 結(jié) 論

衛(wèi)星平臺顫振的精密探測和反演對保證高分辨率衛(wèi)星的測圖精度和成像質(zhì)量至關(guān)重要，是攝影測量領(lǐng)域亟須解決的難題之一。本文結(jié)合高精度機器視覺，提出了基于局部頻率分析亞像素相位相關(guān)的高分辨率衛(wèi)星平臺顫振精密探測方法。針對從影像反演顫振對匹配算法精度和不受輻射差異影響的要求，提取相位一致性和局部平均相位構(gòu)造輻射不變的結(jié)構(gòu)影像表達來改進亞像素相位相關(guān)算法，結(jié)構(gòu)影像表達考慮了影像的內(nèi)在信息，在不同光照條件和不同模態(tài)下保持一致，并將改進的相位相關(guān)算法融入視差法顫振探測框架，實現(xiàn)穩(wěn)健高效的衛(wèi)星平臺顫振精密探測。利用實際輻射變化的數(shù)據(jù)集開展了相位相關(guān)算法對比試驗，本文亞像素相位相關(guān)算法的效果優(yōu)于其他4種主流亞像素算法，在輻射差異和噪聲等影響情況下，匹配精度達到0.05個像素。采用本文顫振探測方法對資源三號衛(wèi)星平臺顫振的特性進行了有效地反演，與姿態(tài)數(shù)據(jù)估計結(jié)果對比驗證了其可信度。筆者的后續(xù)工作將嘗試將提出的顫振精密探測方法應(yīng)用于更多高分辨率衛(wèi)星的顫振分析。