張金強(qiáng)，索志勇，李真芳，保 錚

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

0 引 言

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

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

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

1 CSAR子孔徑圖像序列分析

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

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

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

(1)

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

(2)

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

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

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

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

2 聯(lián)合相關(guān)DEM提取方法

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

圖3 聯(lián)合相關(guān)DEM提取方法流程圖Fig.3 Flowchart of joint cross-correlation DEM extraction method

聯(lián)合相關(guān)DEM提取方法以CSAR子孔徑圖像序列作為輸入。

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

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

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

聯(lián)合相關(guān)系數(shù)用公式表示為

(3)

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

(4)

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

ρ1,m=

(5)

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

ρ1,m

(6)

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

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

3 實測數(shù)據(jù)處理及分析

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

圖5 AFRL CSAR系統(tǒng)的觀測場景Fig.5 Imaged area of AFRL CSAR system

3.1 測度函數(shù)靈敏度分析

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

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

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

3.2 DEM提取結(jié)果及分析

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

圖7 和相關(guān)法提取DEMFig.7 DEM extracted by accumulation cross-correlation

圖8 改進(jìn)和相關(guān)法提取DEMFig.8 DEM extracted by improved accumulation cross-correlation

圖9 聯(lián)合相關(guān)法提取DEMFig.9 DEM extracted by joint cross-correlation

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

(5)

圖10 停車場中汽車的圖片F(xiàn)ig.10 Photos of vehicles in parking lot

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

m

表2 汽車的高程估計結(jié)果

4 結(jié) 論

本文通過對CSAR子孔徑圖像序列相關(guān)性的分析發(fā)現(xiàn),子孔徑圖像間的相關(guān)性隨著方位夾角的增大而降低?；谏鲜龇治?本文提出利用子孔徑圖像序列間的聯(lián)合相關(guān)系數(shù)作為測度函數(shù)提取觀測場景DEM。該測度函數(shù)利用了保持相關(guān)性的一段圓弧內(nèi)所有子孔徑圖像之間的相關(guān)性信息,提高了測度函數(shù)靈敏度和DEM提取精度。實測數(shù)據(jù)處理結(jié)果和分析驗證了本文算法的有效性和精確性。