圓跡SAR子孔徑圖像序列聯合相關DEM提取方法

張金強，索志勇，李真芳，保 錚

(西安電子科技大學雷達信號處理國家重點實驗室，陜西西安710071)

0 引 言

合成孔徑雷達(synthetic aperture radar, SAR)獲取觀測場景地形高程模型(digital elevation model, DEM)主要通過兩種技術途徑:雷達立體攝影技術[1-2]和雷達干涉測高技術[3-5]。基于不同雷達下視角下獲取的兩幅SAR圖像,雷達立體攝影技術利用視差提取目標高程,而雷達干涉測高技術利用相位差提取目標高程。由于單一照射方向獲取的SAR圖像存在疊掩或陰影區域,導致利用上述技術途徑無法獲取觀測場景全方位DEM[6-7]。

圓跡SAR(circular SAR, CSAR)隨雷達平臺以觀測場景為中心作360°圓周運動,波束始終照射同一地面場景,實現對其全方位觀測[8-11]。利用CSAR可以提取觀測場景全方位DEM[12-14]。文獻[12-14]提出了基于CSAR子孔徑圖像序列的DEM提取方法。首先,將360°圓環劃分為多段圓弧,每段圓弧劃分為多個子孔徑,對各子孔徑回波進行聚焦成像處理,獲取子孔徑圖像序列。然后,針對每段圓弧,利用圓弧內部分子孔徑圖像序列間的相關性估計觀測場景DEM。文獻[12]利用圓弧內各子孔徑圖像與中心子孔徑圖像間相關系數的和值作為測度函數,每段圓弧提取一幅DEM。文獻[13-14]利用圓弧內各子孔徑圖像與中心子孔徑圖像間的相關系數分別作為測度函數,獲取多幅DEM后進行平均,最終獲取一幅DEM。相比于文獻[13],文獻[14]在對圓弧內各子孔徑圖像與中心子孔徑圖像進行配準時,采用了金字塔分層搜索的方法,提高了配準的穩健性。由于疊掩或陰影的影響,不同段圓弧提取的DEM中存在不同的無效高程值區域。最后,融合不同段圓弧DEM獲取全方位DEM。

針對每段圓弧,文獻[12-14]只利用了各子孔徑圖像與中心子孔徑圖像之間的相關性信息,非中心子孔徑圖像之間的相關性信息沒有得到利用。為了解決上述問題,本文提出了聯合相關法。該方法將圓弧內所有子孔徑圖像之間的聯合相關系數作為測度函數。該測度函數利用了圓弧內全部子孔徑圖像之間的相關性信息,能夠提高測度函數靈敏度和DEM提取精度。同時,將子孔徑圖像序列向三維空間投影以校正其幾何形變,提高子孔徑圖像序列之間的相關性。本文首先對CSAR子孔徑圖像序列進行分析,然后給出聯合相關法處理流程,最后利用CSAR實測數據驗證所提算法的有效性和精確性。

1 CSAR子孔徑圖像序列分析

CSAR成像幾何模型如圖1所示。根據分辨率要求確定子孔徑方位角寬度,然后將360°圓環劃分為多個子孔徑,利用后向投影(back projection, BP)算法[11]對各子孔徑回波進行聚焦成像處理,獲取CSAR子孔徑圖像序列。不失一般性,假設各子孔徑圖像成像于地面坐標系。

圖1 CSAR成像幾何模型Fig.1 CSAR imaging geometry

成像平面高程與目標實際高程不一致,將導致目標在子孔徑圖像上的平面位置與實際不一致,引起子孔徑圖像發生幾何形變[12]。如圖1所示,假設子孔徑圖像#i對應的雷達軌跡中心為Ci,雷達下視角為θi,方位角為φi,觀測場景中一目標點P實際高程與聚焦平面高程之差為Δh,則其在子孔徑圖像#i上的平面位置相對于實際位置偏移量為

(1)

根據式(1)可知,目標實際高程與成像平面高程差值越大,子孔徑圖像幾何形變越大。不同子孔徑圖像幾何形變不同,與其雷達下視角和方位角有關。如圖1所示,假設子孔徑圖像#j對應的雷達軌跡中心為Cj,雷達下視角為θj,方位角為φj,目標點P在子孔徑圖像#j上的平面位置相對于其在子孔徑圖像#i上的位置偏移量為

(2)

根據式(2)可知,子孔徑圖像間方位夾角增大,導致子孔徑圖像間幾何形變變化增大,引起子孔徑圖像間相關性降低。同時,目標散射特性隨方位角的變化,也將導致子孔徑圖像間相關性降低[13]。因此,子孔徑圖像間相關性隨著方位夾角的增大而降低。下面利用美國空軍研究實驗室(air force research laboratory, AFRL)公開發布的CSAR實測數據(gotcha volumetric SAR data set, Version 1.0)驗證該結論。

子孔徑圖像序列生成過程與本文第3節相同。在子孔徑圖像序列中根據方位夾角隨機選擇兩幅圖像,逐像素估計其相關系數,相關系數估計窗口大小為15像素×15像素。每一方位夾角均隨機選擇5對子孔徑圖像,計算所有像素相關系數的平均值,其隨方位夾角的變化如圖2所示,方位夾角變化范圍為3°～42°。根據圖2可知,子孔徑圖像對間相關系數隨著方位夾角的增大而減小。

圖2 CSAR 子孔徑圖像對間相關系數隨方位夾角的變化Fig.2 Cross-correlation coefficients between subaperture imagepairs of CSAR versus azimuth angle intervals

CSAR子孔徑圖像序列的相關性將影響DEM提取精度,因此對CSAR子孔徑圖像序列相關性的分析結果對于下文聯合相關DEM提取方法中圓環的劃分具有重要參考價值。

2 聯合相關DEM提取方法

CSAR具有360°全方位觀測能力,但是對于不同的方位觀測角,子孔徑圖像中存在不同的疊掩或陰影區域。為了獲取觀測場景全方位DEM,將整個圓環劃分為多段圓弧,利用每段圓弧內的子孔徑圖像序列提取觀測場景DEM。不同段圓弧提取的DEM中由于疊掩或陰影導致的無效高程值區域不同,融合不同段圓弧提取的DEM即可獲取全方位DEM。聯合相關DEM提取方法流程圖如圖3所示,下面對主要步驟進行詳細介紹。

圖3 聯合相關DEM提取方法流程圖Fig.3 Flowchart of joint cross-correlation DEM extraction method

聯合相關DEM提取方法以CSAR子孔徑圖像序列作為輸入。

為了消除CSAR子孔徑圖像間幾何形變變化對相關性的影響,將子孔徑圖像向三維空間投影,如圖4所示。首先,根據觀測場景的先驗平面位置和目標高程范圍,在地面坐標系下建立三維格網;然后,針對每一格網點根據式(1)進行SAR反向定位,獲取格網點在子孔徑圖像中的像素坐標;最后,通過插值提取該像素的幅值賦予格網點。當格網點對應的高程與目標高程一致時,子孔徑圖像幾何形變得到校正,子孔徑圖像間相關性最高。

根據第1節分析可知,子孔徑圖像間方位夾角越大相關性越低,相關性降低將導致DEM提取精度降低。根據式(2)可知,子孔徑圖像間方位夾角較小時,目標平面位置偏移量對目標高程不敏感,DEM提取精度較低。因此,在進行子孔徑圖像序列劃分時需要根據其相關性選擇圓弧方位角寬度。

圖4 CSAR子孔徑圖像向三維空間投影示意圖Fig.4 Project subaperture images of CSAR on 3D space

聯合相關系數用公式表示為

(3)

式中,M為圓弧內子孔徑圖像數量;sm(l,k)為所選窗口內圖像#m的像素幅度值;μm為相應的幅度平均值;(2L+1)×(2K+1)為所選窗口大小。所選窗口為以待估計像素為中心的矩形窗,并且窗口內其他像素與中心待估計像素的高程值相同。對于每一平面坐標(x,y),沿高程方向(即圖4中Z軸方向)計算JC隨h的變化。然后選擇使JC取最大值時的h作為估計高程,即

(4)

文獻[13-14]利用圓弧內各子孔徑圖像與中心子孔徑圖像間的相關系數作為測度函數,獲取多幅DEM后進行平均,最終獲取一幅DEM。兩幅子孔徑圖像間的相關系數用公式表示為

ρ1,m=

(5)

式中,下標1表示中心子孔徑圖像編號;下標m∈[2,3,…,M]表示其他各子孔徑圖像編號。文獻[12]利用圓弧內各子孔徑圖像與中心子孔徑圖像間相關系數的和值作為測度函數,估計得到一幅DEM。

ρ1,m

(6)

對比式(3)、式(5)和式(6)可知,文獻[12-14]僅利用了圓弧內部分子孔徑圖像間的相關性信息,而本文方法利用了圓弧內全部子孔徑圖像間的相關性信息。

將所有圓弧提取到的DEM進行融合[12],獲取全方位DEM并輸出。

3 實測數據處理及分析

本文選擇AFRL公開發布的CSAR實測數據中航過1、HH極化方式的數據驗證所提算法的有效性和精確性。為了對比分析所提算法的性能,同時利用文獻[12]中的算法(稱為和相關法)和文獻[14]中的算法(稱為改進和相關法)進行DEM提取。上述實測數據由X波段、640 MHz的雷達系統錄取。美國地理學會(US geography society, USGS)獲取的觀測場景的光學正射影像如圖5(a)所示。將360°圓環劃分為120個子孔徑,每個子孔徑對應的方位角寬度為3°。利用BP算法對每個子孔徑的回波數據在地面坐標系下進行聚焦成像,成像平面高程為0 m,格網間距為0.2 m,圖像大小為501像素×501像素。對所有子孔徑圖像進行非相干疊加[15]獲取的SAR圖像如圖5(b)所示。

圖5 AFRL CSAR系統的觀測場景Fig.5 Imaged area of AFRL CSAR system

3.1 測度函數靈敏度分析

利用測度函數主瓣3 dB寬度(稱之為主瓣寬度)和第一旁瓣峰值與主瓣峰值的比值(稱之為峰值旁瓣比)評價其靈敏度[12]。主瓣寬度和峰值旁瓣比越小,表明測度函數靈敏度越高,DEM提取精度越高。

根據方位角寬度在子孔徑圖像序列中隨機選擇若干幅方位角連續的子孔徑圖像組成一段圓弧。針對每一平面坐標(x,y)計算和相關系數和聯合相關系數隨高程h的變化,然后分別估計其主瓣寬度和峰值旁瓣比。每一方位角寬度均隨機選擇5組子孔徑圖像,窗口大小選為5像素×5像素和11像素×11像素。圖6給出了和相關函數和聯合相關函數的主瓣寬度和峰值旁瓣比均值隨圓弧方位角寬度的變化,圓弧方位角寬度的變化范圍為21°～87°。由圖6可知,與和相關函數法相比,本文方法測度函數的主瓣寬度和峰值旁瓣比更小,靈敏度更高。

圖6 測度函數靈敏度隨圓弧方位角寬度的變化Fig.6 Sensitivities of measure function versus azimuthangle intervals of arcs

3.2 DEM提取結果及分析

根據圖1和圖2的分析結果,選擇圓弧方位角寬度為60°,將360°圓環劃分為相互重疊的24個圓弧[12]。選擇大小為5像素×5像素和11像素×11像素的兩種窗口計算相關系數隨目標高程的變化,將得到的兩組相關系數在對應高程處相乘,然后基于此估計目標高程。上述處理流程能夠提高高程估計性能[12]。觀測場景內目標高程范圍為-2～2 m,格網間隔選擇0.2 m。和相關法提取DEM如圖7所示,改進和相關法提取DEM如圖8所示,聯合相關法提取DEM如圖9所示。對比圖7～圖9的結果可知,本文方法提取的DEM中汽車輪廓更清晰。由于觀測場景內大部分地物后向散射系數較小、子孔徑圖像序列間相關性較低,導致DEM提取精度較低。

圖7 和相關法提取DEMFig.7 DEM extracted by accumulation cross-correlation

圖8 改進和相關法提取DEMFig.8 DEM extracted by improved accumulation cross-correlation

圖9 聯合相關法提取DEMFig.9 DEM extracted by joint cross-correlation

為了定量評估算法性能,選擇圖5(b)中所標注的7輛汽車評估所提取的DEM精度,圖中箭頭指向汽車的車頭方向。7輛汽車的圖片如圖10所示,7輛汽車的實際長、寬和高如表1所示。汽車#F的DEM提取結果如圖11所示,選擇圖11(a)中內框所示的區域評估DEM提取精度。統計所提取目標高程的均值和均方根誤差[16],其中均方根誤差的計算公式為

(5)

圖10 停車場中汽車的圖片Fig.10 Photos of vehicles in parking lot

圖11 汽車#F的DEM提取結果Fig.11 Extracted DEM of car #F

m

表2 汽車的高程估計結果

4 結 論

本文通過對CSAR子孔徑圖像序列相關性的分析發現,子孔徑圖像間的相關性隨著方位夾角的增大而降低。基于上述分析,本文提出利用子孔徑圖像序列間的聯合相關系數作為測度函數提取觀測場景DEM。該測度函數利用了保持相關性的一段圓弧內所有子孔徑圖像之間的相關性信息,提高了測度函數靈敏度和DEM提取精度。實測數據處理結果和分析驗證了本文算法的有效性和精確性。