劉世杰，晏飛，3，王衛(wèi)安，李榮興

（1.同濟(jì)大學(xué) 空間信息科學(xué)及可持續(xù)發(fā)展應(yīng)用中心，上海200092；2.同濟(jì)大學(xué) 測(cè)繪與地理信息學(xué)院，上海200092；3.重慶市地理信息中心，重慶401121）

衛(wèi)星遙感影像立體匹配是航天攝影測(cè)量中的關(guān)鍵技術(shù)和核心問(wèn)題［1］，許多學(xué)者提出了各種匹配方法，包括最小二乘匹配［2］、歸一化互相關(guān)（normalized cross correlation，NCC）［3］等灰度區(qū)域匹配、全局／半全局匹配［4］以及相位相關(guān)等頻域匹配算法［5］.各種匹配方法差異反映了影像立體匹配的復(fù)雜性和高難度.在圖像獲取的過(guò)程中，三維空間投影到二維的平面圖像上，導(dǎo)致許多信息丟失，因此圖像匹配屬于逆向工程類問(wèn)題，這是一個(gè)病態(tài)問(wèn)題［6］.除了因遮擋造成匹配失敗，還有可能由于重復(fù)紋理導(dǎo)致誤匹配，以及由于圖像存在噪聲和幾何變形使得匹配算法不穩(wěn)定等問(wèn)題.針對(duì)這些問(wèn)題，有學(xué)者提出和研究了基于尺度不變特征轉(zhuǎn)換（SIFT）的立體匹配算法［7］，以克服立體影像尺度差異和幾何變形造成的匹配困難，但由于特征比較稀疏且分布容易不均，不利于生成較高分辨率的數(shù)字表面模型（digital surface model，DSM）；針對(duì)DSM生成，有學(xué)者研究了基于物方的匹配方法如鉛垂線軌跡法（VLL）［8］，該方法以地面為采樣計(jì)算單元，通過(guò)高程值的改變計(jì)算像方最大相關(guān)系數(shù)，從而在成功匹配的同時(shí)獲取地面DSM，該方法速度快，但仍面臨影像畸變和遮擋帶來(lái)的問(wèn)題.

本文綜合運(yùn)用SIFT算子和NCC算子，并考慮核線和視差約束，提出了一種改進(jìn)的基于三角網(wǎng)視差傳遞約束的立體影像分層匹配方法，并利用上海崇明區(qū)域的WorldView-1影像和浙江舟山區(qū)域的資源三號(hào)衛(wèi)星影像進(jìn)行實(shí)驗(yàn)，驗(yàn)證了本文方法的有效性，獲得了高精度的DSM.

1 影像成像模型

有理函數(shù)模型（rational function model，RFM）是目前高分辨率衛(wèi)星廣泛采用的成像模型，衛(wèi)星影像供應(yīng)商直接提供有理函數(shù)模型參數(shù).有理函數(shù)模型是將像點(diǎn)坐標(biāo)表示為以相應(yīng)地面點(diǎn)空間坐標(biāo)為自變量的多項(xiàng)式的比值，如式（1）所示［9］.

式中：（rn，cn）是像平面坐標(biāo)；（Xn，Yn，Zn）是地面坐標(biāo)；它們均是經(jīng)過(guò)標(biāo)準(zhǔn)化后的坐標(biāo)值，以減小計(jì)算過(guò)程中由于數(shù)量級(jí)差別過(guò)大引入的誤差；對(duì)于標(biāo)準(zhǔn)有理函數(shù)模型，pi（i=1，2，3，4）是一個(gè)20項(xiàng)三階多項(xiàng)式.

使用原始影像提供的有理函數(shù)模型參數(shù)（RPCs）直接進(jìn)行定位具有明顯的系統(tǒng)誤差，需利用地面控制點(diǎn)進(jìn)行誤差補(bǔ)償，提高定位精度.一般說(shuō)有兩種補(bǔ)償方案［10］，一種是在物方直接對(duì)定位結(jié)果進(jìn)行誤差糾正，另一種是先糾正像平面坐標(biāo)的誤差，再進(jìn)行立體定位.后一種方案理論上嚴(yán)密，對(duì)于高分辨率衛(wèi)星影像，一般仿射變換模型能滿足精度要求，且較穩(wěn)健，故本文采用像方仿射變換誤差補(bǔ)償模型.

2 影像立體匹配和DSM生成

本文提出的影像立體匹配生成DSM的流程如圖1所示，包括影像預(yù)處理、金字塔逐層約束匹配獲取密集匹配點(diǎn)、立體定位和點(diǎn)云內(nèi)插生成DSM等.

2.1 影像預(yù)處理

預(yù)處理包括影像去噪、增強(qiáng)和影像金字塔構(gòu)建.數(shù)字影像獲取不可避免會(huì)存在噪聲，因此需要進(jìn)行去噪處理.對(duì)于圖像平滑去噪算子而言，不僅需要過(guò)濾掉噪聲，更重要的是要保持邊緣信息.本文采用Saint-Marc提出的自適應(yīng)平滑算子［11］，屬于非線性濾波，通過(guò)一個(gè)迭代的過(guò)程，達(dá)到去除噪聲并保留邊緣的目的.影像增強(qiáng)采用 Wallis濾波［12］，用于增強(qiáng)紋理并提高信噪比，便于特征提取和匹配.在建立金字塔影像過(guò)程中，先利用高斯模糊算子對(duì)當(dāng)前層影像進(jìn)行低通濾波，然后進(jìn)行降采樣，得到上一層影像，最終建立包含4層的影像金字塔.

圖1 高分辨率立體影像生成DSM方法流程圖Fig.1 Workflow of DSM generation from high resolution satellite stereo images

2.2 影像匹配

2.2.1 SIFT初匹配

影像初始匹配點(diǎn)將作為下層加密匹配的約束條件，對(duì)準(zhǔn)確性的要求非常高.SIFT特征具有對(duì)尺度、旋轉(zhuǎn)、亮度變化保持不變的穩(wěn)定特性，故采用SIFT算子進(jìn)行初匹配，結(jié)合RANSAC（隨機(jī)抽樣一致）檢驗(yàn)算法［13］，檢驗(yàn)并剔除可能的誤匹配，保證初匹配的準(zhǔn)確性.

2.2.2 特征點(diǎn)和格網(wǎng)點(diǎn)匹配

本文采用F?rstner特征算子［14］，它能給出特征點(diǎn)的類型且精度較高.其基本思想是：對(duì)于角點(diǎn)，對(duì)最佳窗口內(nèi)通過(guò)每個(gè)像元的邊緣直線（垂直于梯度方向）進(jìn)行加權(quán)中心化，得到角點(diǎn)的坐標(biāo)；對(duì)于圓狀點(diǎn)，對(duì)最佳窗口內(nèi)通過(guò)每個(gè)像元的梯度直線進(jìn)行加權(quán)中心化，得到圓心的坐標(biāo).

利用特征算子提取的特征點(diǎn)主要分布在影像灰度存在變化具有明顯特征的地方，其分布不均勻且數(shù)量有限.為了生成DSM，需要更多密集的匹配點(diǎn).因此，除特征點(diǎn)外，在影像上按一定間隔劃分格網(wǎng)，并去除灰度方差低于一定閾值的格網(wǎng)點(diǎn)，然后對(duì)保留的格網(wǎng)點(diǎn)逐個(gè)在待匹配影像上進(jìn)行搜索匹配.

特征點(diǎn)和格網(wǎng)點(diǎn)匹配采用NCC算法.在傳統(tǒng)的匹配方法中，一般采用矩形窗口計(jì)算相關(guān)系數(shù)，這種方法適用于幾何變形較小的影像匹配，而對(duì)于幾何變形較大的情況，矩形窗口不再適用.因此，本文提出一種基于視差幾何糾正的相關(guān)系數(shù)計(jì)算方法.在該方法中，對(duì)于待匹配點(diǎn)，根據(jù)附近的視差三角形格網(wǎng)點(diǎn)上的視差計(jì)算幾何變換系數(shù)，糾正待匹配窗口并重采樣，使左右影像匹配窗口內(nèi)容一致.另外，由于采用固定模板進(jìn)行匹配會(huì)降低匹配成功率，本文采用一種自動(dòng)改變匹配窗口大小和形狀的匹配策略，它通過(guò)把匹配窗口向四個(gè)方向擴(kuò)展，直到相關(guān)系數(shù)不再增大為止，從而提高匹配的準(zhǔn)確性和成功率.2.2.3 特征線匹配

特征線是影像和地面上的重要特征，反映了影像灰度或地面高程的不連續(xù)性.通過(guò)對(duì)特征線的提取和匹配，可以實(shí)現(xiàn)對(duì)地面高程不連續(xù)的精確定位，對(duì)精確DSM的構(gòu)建具有重要意義.對(duì)于特征線段的提取，本文采用Akinlar提出的EDLines算法［15］.該方法無(wú)需對(duì)參數(shù)進(jìn)行調(diào)整而可得到精確的結(jié)果，且速度快.該檢測(cè)算子利用亥姆霍茲規(guī)則進(jìn)行線段驗(yàn)證，剔除錯(cuò)誤的線段.

對(duì)線段的匹配，除幾何位置、長(zhǎng)度和方向（與核線夾角）外，還考慮線段周圍的灰度信息，采用平行于線段的矩形滑動(dòng)窗口計(jì)算灰度相關(guān)系數(shù).在匹配搜索過(guò)程中，綜合利用已匹配點(diǎn)所形成的三角網(wǎng)約束和核線約束來(lái)確定候選匹配線段.

2.2.4 匹配搜索約束

（1）三角網(wǎng)約束

對(duì)于灘涂、島礁等特征不明顯、紋理重復(fù)或缺乏的區(qū)域，直接進(jìn)行影像匹配很容易造成誤匹配，因此本文采用三角網(wǎng)約束的策略進(jìn)行匹配［16］.點(diǎn)的匹配首先在最上層低分辨率的影像上進(jìn)行，匹配成功的點(diǎn)傳遞到較高分辨率的下一層重新進(jìn)行匹配以獲得更高的匹配精度，并把匹配上的點(diǎn)建立一個(gè)視差表面的狄洛尼三角網(wǎng)，用作加密匹配點(diǎn)的搜索約束.這個(gè)過(guò)程不斷重復(fù)，直到最下層原始分辨率影像.在每個(gè)影像層上，同名點(diǎn)搜索區(qū)域根據(jù)上一層的視差確定.從影像最高層到最低層，匹配上的同名點(diǎn)越來(lái)越多，獲得的減小搜索范圍的視差信息就越精確.

（2）核線約束

線陣CCD（電荷耦合器件）推掃式影像遵循一種“動(dòng)態(tài)”的行中心投影成像方式，因此無(wú)法像框幅式中心投影影像那樣基于成像幾何關(guān)系給出嚴(yán)格直觀的核線定義.目前，在衛(wèi)星影像近似核線幾何關(guān)系的各種描述中，基于投影軌跡法的核線定義建立在成像的幾何約束條件之上，在理論上最為嚴(yán)密.線陣CCD推掃式影像核線類似于雙曲線，但在影像范圍內(nèi)可以近視看作直線［17］.

根據(jù)三角網(wǎng)視差約束條件，得到左影像上某個(gè)待匹配點(diǎn)在右影像上的估計(jì)位置.然后以估計(jì)的位置為中心，建立矩形搜索窗口.再由左影像上的待匹配點(diǎn)計(jì)算在右影像上的核線.匹配搜索的時(shí)候，就只需在這個(gè)搜索窗口內(nèi)的核線段上搜索即可.由于計(jì)算的核線會(huì)存在一定的誤差，因此實(shí)際情況下會(huì)把核線左右兩側(cè)擴(kuò)寬一定像素，然后沿著這條具有一定寬度的核線段搜索.這樣不僅可以大大提高搜索效率，而且有效地降低了匹配的不確定性.

2.2.5 誤匹配剔除

（1）雙向匹配

雙向匹配算法的基本思想是在第一次匹配時(shí)，以左影像某個(gè)特征點(diǎn)或格網(wǎng)點(diǎn)作為待匹配點(diǎn)，在右影像上一定視差搜索范圍內(nèi)計(jì)算所有候選點(diǎn)的相關(guān)系數(shù)，取相關(guān)系數(shù)最大值點(diǎn)為有效匹配點(diǎn)，即進(jìn)行所謂的正匹配.然后將上一次匹配的結(jié)果（右影像上的點(diǎn)）作為待匹配點(diǎn)，在左影像上一定視差范圍內(nèi)進(jìn)行搜索，以同樣的方法進(jìn)行二次匹配，即所謂的逆匹配.如果兩次結(jié)果的差異小于設(shè)定的閾值（如1個(gè)像素），則接受該匹配結(jié)果，否則舍棄，達(dá)到剔除誤匹配點(diǎn)的目的.

（2）視差平滑約束

盡管在匹配過(guò)程中采用了三角網(wǎng)視差約束和雙向匹配的算法，仍然不能保證所有的匹配點(diǎn)都是正確的.在本文的匹配算法中，引入了視差平滑約束條件.對(duì)于每一對(duì)同名點(diǎn)，使用它們的鄰域內(nèi)的同名點(diǎn)擬合一個(gè)視差平面（不跨越邊線）.如果該點(diǎn)視差與視差平面的距離超過(guò)閾值（3倍中誤差），認(rèn)為該點(diǎn)是誤匹配并剔除.實(shí)驗(yàn)表明，該方法可有效剔除絕大多數(shù)錯(cuò)誤匹配.

2.3 DSM構(gòu)建

成功匹配的同名點(diǎn)通過(guò)成像模型立體定位得到對(duì)應(yīng)地面點(diǎn)的三維坐標(biāo)，但影像匹配定位得到的三維點(diǎn)的排列是不規(guī)則的，因此需要通過(guò)內(nèi)插獲得規(guī)則格網(wǎng)的DSM.綜合考慮時(shí)間復(fù)雜度和插值精度，本文采用反距離加權(quán)插值方法構(gòu)建DSM.

3 實(shí)驗(yàn)與分析

3.1 崇明WorldView-1衛(wèi)星影像生成DSM

（1）實(shí)驗(yàn)數(shù)據(jù)

該實(shí)驗(yàn)使用上海市崇明島東灘區(qū)域的WorldView-1影像，分辨率0.5m，為不同時(shí)期（2012／01／25 和 2012／04／06）異 軌 立 體，交 會(huì) 角28.7°.實(shí)驗(yàn)區(qū)域大小約為18km×18km，見(jiàn)圖2a和圖2b，區(qū)域內(nèi)地勢(shì)平坦，主要為農(nóng)田和灘涂；該區(qū)域位于相鄰兩景影像內(nèi)（供應(yīng)商按一定分幅規(guī)格進(jìn)行分幅，且相鄰影像具有一定重疊度），兩期共四幅影像.為了糾正影像幾何定位誤差，使用了17個(gè)地面控制點(diǎn)，其地面坐標(biāo)由GPS（全球定位系統(tǒng)）測(cè)得；同時(shí)，在重疊區(qū)域選取了9個(gè)連接點(diǎn)參與平差解算，以保證上下兩幅影像生成的DSM高程連續(xù).

圖2 崇明WorldView-1影像及控制點(diǎn)與連接點(diǎn)分布（三角點(diǎn)：控制點(diǎn)，圓點(diǎn)：連接點(diǎn)）和生成的DSMFig.2 Chongming WorldView-1image and distribution of control points and tie points （triangle：control point，circle：tie point）and generated DSM

（2）DSM 和精度分析

首先利用地面控制點(diǎn)和影像連接點(diǎn)糾正影像成像模型誤差（采用像方仿射模型），平差后成像模型誤差在0.5像素以內(nèi)，最后生成的拼接后的DSM如圖2c所示.從生成的DSM上可清晰看出地面紋理，特別是線性特征，如道路和不同地表覆蓋的邊線等，反映了線特征提取和匹配的作用.利用該區(qū)域內(nèi)均勻分布的100個(gè)高程注記點(diǎn)作為檢核來(lái)評(píng)價(jià)DSM的精度，計(jì)算得DSM的高程中誤差約為0.7m.

3.2 舟山資源三號(hào)衛(wèi)星影像生成DSM

（1）實(shí)驗(yàn)數(shù)據(jù)

實(shí)驗(yàn)影像數(shù)據(jù)為中國(guó)資源三號(hào)衛(wèi)星2013年1月24日拍攝的全色前、后視立體影像，前、后視夾角為44°，地面分辨率為3.6m，區(qū)域大小約為58km×58km，如圖3a所示.該區(qū)域包含眾多島嶼，高程變化在海拔550m范圍內(nèi).

（2）DSM和精度分析

由于沒(méi)有該區(qū)域的地面控制點(diǎn)，故直接采用影像附帶有理函數(shù)模型進(jìn)行地面定位，生成的DSM如圖3b所示.為了檢驗(yàn)DSM的精度，利用該區(qū)的ASTER GDEM（先進(jìn)星載熱輻射反射探測(cè)儀全球數(shù)字高程模型）作為參考，先將生成的DSM與ASTER GDEM進(jìn)行配準(zhǔn)，然后比較并統(tǒng)計(jì)高程差異.考慮到舟山群島區(qū)域地形復(fù)雜，按山地（高程大于20m的山地區(qū)域）、城區(qū)（城鎮(zhèn)區(qū)域，主要分布在沿海）和平地（高程小于20m的沿海空曠區(qū)域）三類地形檢驗(yàn)DSM精度.結(jié)果顯示，山地區(qū)域DSM與ASTER GDEM的偏差均方根為15.5m，城區(qū)為8.9m，平地為4.6m，兩者吻合較好.該差異一部分是因?yàn)锳STER GDEM 本身存在誤差，另一部分是因?yàn)锳STER和資源三號(hào)影像拍攝相隔10余年（ASTER GDEM是利用1999發(fā)射的ASTER傳感器所獲取的15m分辨率影像生產(chǎn)得到，資源三號(hào)影像為2013年拍攝），由于人類活動(dòng)而造成了地面變化（如山體開(kāi)挖、城區(qū)建設(shè)等）.

4 結(jié)語(yǔ)

本文綜合運(yùn)用SIFT算子的可靠性和灰度相關(guān)算子的密集匹配特性，并考慮核線和視差約束，提出了一種改進(jìn)的基于三角網(wǎng)視差傳遞約束的衛(wèi)星影像金字塔分層匹配方法，實(shí)現(xiàn)特征點(diǎn)、格網(wǎng)點(diǎn)和特征線的匹配，用于從高分辨率衛(wèi)星立體影像生成DSM.該方法在金字塔分層匹配中引入三角網(wǎng)視差傳遞約束和核線約束以提高匹配效率和準(zhǔn)確性，同時(shí)采用自適應(yīng)窗口改進(jìn)來(lái)提高匹配成功率，對(duì)于誤匹配點(diǎn)通過(guò)雙向匹配和視差平滑約束進(jìn)行剔除，最后將匹配成功的點(diǎn)經(jīng)立體幾何模型定位生成點(diǎn)云，并內(nèi)插獲得DSM.由于采用了自適應(yīng)窗口和多種約束及誤匹配檢測(cè)算法，該DSM生成方法具有較好的魯棒性.實(shí)驗(yàn)證明，該方法能提高匹配的準(zhǔn)確性和成功率，獲得高精度的DSM.

圖3 舟山資源三號(hào)影像和生成的DSMFig.3 Zhoushan ZY-3image and generated DSM

［1］ Gruen A.Development and status of image matching in photogrammetry［J］.The Photogrammetric Record，2012，27（137）：36.

［2］ Zhang L.Automatic digital surface model（DSM）generation from linear array images［D］，Zurich：Swiss Federal Institute of Technology Zurich，2005.

［3］ Tsai D M，Lin C T，Chen J F.The evaluation of normalized cross correlations for defect detection［J］.Pattern Recognition Letters，2003，24：2525.

［4］ Alobeid A，Jacobsen K，Heipke C.Comparison of matching algorithms for DSM generation in urban areas from ikonos imagery［J］.Photogrammetric Engineering &Remote Sensing，2010，76（9）：1041.

［5］ Keller Y，Averbuch A.A projection-based extension to phase correlation image alignment［J］.Signal Processing，2007，87（1）：124.

［6］ Terzopoulos D.Regularization of inverse visual problems involving discontinuities［J］.Pattern Analysis and Machine Intelligence，IEEE Transactions on，1986，4：413.

［7］ Lowe D G.Distinctive image features from scale-invariant keypoints［J］.International Journal of Computer Vision，2004，60（2）：91.

［8］ 鄭順義，張祖勛，張劍清.基于物方影像匹配和概率松弛的斷面自動(dòng)提?。跩］.測(cè)繪信息工程，2004，29（2）：26.ZHENG Shunyi，ZHANG Zuxun，ZHANG Jianqing.Automatic profile extraction based on object space image matching and probability relaxation algorithm［J］.Journal of Geomatics，2004，29（2）：26.

［9］ Madani M.Real-time sensor-independent positioning by rational functions［C］／／Proceedings of ISPRS Workshop on Direct Versus the HRSC data.Indirect Methods of Sensor Orientation.Barcelona：［s.n.］，1999：64-75.

［10］ 劉世杰.高分辨率衛(wèi)星遙感影像成像模型與定位技術(shù)研究［D］.上海：同濟(jì)大學(xué)，2008 LIU Shijie.Imaging models and geo-positioning technology for high resolution satellite imagery ［D］.Shanghai： Tongji University，2008.

［11］ Saint-Marc P，Chen J S，Medioni G.Adaptive smoothing：a general tool for early vision［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，1991，13（6）：514.

［12］ 張力，張祖勛，張劍清.Wallis濾波在影像匹配中的應(yīng)用［J］.武漢測(cè)繪科技大學(xué)學(xué)報(bào)，1999（1）：24.ZHANG Li，ZHANG Zuxun，ZHANG Jianqing.The image matching based on wallis filtering［J］.Journal of Wuhan Technical University of Surveying and Mapping，1999（1）：24.

［13］ Fischler M A，Bolles R C.Random sample consensus：a paradigm for model fitting with applications to image analysis and automated cartography［J］.Communications of the ACM，1981，24（6）：381.

［14］ Forstner W，Gulch E.A fast operator for detection and precise location of distinct points，corners and centres of circular features ［C ］／／Proceedings of ISPRS intercommission conference on fast processing of photogrammetric data，Interlaken：［s.n.］，1987：281-305.

［15］ Akinlar C，Topal C.EDLines：a real-time line segment detector with a false detection control［J］.Pattern Recognition Letters，2011，32（13）：1633.

［16］ Hwangbo J W.Integration of orbital and ground imagery for automation of rover localization ［D］.Ohio：Ohio State University，2010.

［17］ 胡芬，王密，李德仁，等.利用投影軌跡的衛(wèi)星影像近似核線重排快速算法［J］.武漢大學(xué)學(xué)報(bào)：信息科學(xué)版，2009，34（12）：1431.HU Fen，WANG Mi，LI Deren，etal.Approximate epipolar rearrangement algorithm of satellite stereo-imagery by projection trajectory［J］.Geomatics and Information Science of Wuhan University，2009，34（12）：1431.