龍騰飛，焦偉利，王 威

1.中國科學院 遙感與數(shù)字地球研究所，北京 100094；2.中國科學院大學 北京 100049

1 引 言

傳統(tǒng)的遙感圖像幾何校正模型是基于控制點建立的，當能夠準確獲取地面控制點時，采用地面控制點進行幾何校正能夠較好地改善圖像的定位精度；然而在很多困難地區(qū)的衛(wèi)星影像上，如沙漠、山區(qū)等，很難辨認出準確的物理點特征，因此難以獲得足夠數(shù)量和足夠準確的地面控制點，而線特征、面特征則更易獲得（如已知的道路、水體等）。另一方面，在以不同分辨率影像作為參考數(shù)據(jù)源時，點的位置坐標也難以準確確定。事實上，就特征提取而言，提取有意義的面特征比提取有意義的點特征更容易，而且，無論自然環(huán)境還是人工環(huán)境中都存在豐富的面狀地物，如建筑物、運動場、公園、湖泊等。參考影像、數(shù)字線劃地圖（DLG）、GIS矢量數(shù)據(jù)中存在大量的面特征，使用面特征進行幾何校正能夠充分地利用這些數(shù)據(jù)。

長期以來，人們主要關注基于點［1－3］、直線［4－6］和曲線［7］的圖像配準與幾何校正方法，近年來國內外的學者針對面狀特征的提取［8－9］、面狀特征的相 似性度 量［10－11］以 及 面 狀 特 征 的 匹 配 策 略［12］等方面進行了一些研究，基于面狀地物的圖像配準方法已經(jīng)應用于光學影像與SAR影像［13－14］、遙感影像與GIS數(shù)據(jù)的配準［15］，但是這些方法主要側重于面狀地物的提取和配準，很少涉及利用面特征來建立校正模型，進行圖像幾何校正。目前，已有的方法往往將面特征轉化為點特征使用（如提取多邊形的重心作為控制點［13，15］），最終基于點特征進行幾何校正，并未充分利用面特征的信息。由于傳感器的成像方式、側視角、地形等的影響，影像上的面狀地物會發(fā)生不同程度的變形，重心的位置也會發(fā)生改變，因此使用面狀地物的重心作為控制點的方法存在缺陷。本文從面狀特征的距離度量出發(fā)，直接利用面狀特征的輪廓信息，提出一種基于面特征的圖像幾何校正方法。

2 基于面特征的幾何校正模型

本文方法的關鍵技術是利用面特征建立和求解誤差方程，其中誤差方程的建立需要計算控制面多邊形與像平面多邊形之間的距離，多邊形之間的距離以多邊形上的點到多邊形的距離為基礎。因此，下面依次介紹點到多邊形的距離計算方法、多邊形到多邊形的距離計算方法以及基于面特征的幾何校正模型建立方法。

2.1 點到多邊形的距離計算方法

首先給出點到多邊形距離的定義：同一平面內，點p到多邊形A的距離為多邊形A所有邊界點到點p的最小距離，特別的，當點p在多邊形A內部時，距離為0。

記多邊形的邊界點集合為?A，任意點p到多邊形A的距離可表示為

欲求點p到多邊形A邊界的最小距離，首先分別計算點p到多邊形A各條邊的距離，然后得到其最小值。其中，點p到多邊形A某條邊l的距離分為兩種情況：過點p作線段l的垂線，若垂足在線段l上，則點p到垂足的距離為點p到l的距離；若垂足不在線段l上，則點p到線段l的兩個端點的距離的較小值為點p到l的距離。圖1為不同情況下點到線段的示意圖。圖1（a）中，垂足m在線段l上，則Vpm為點p到線段l的距離；圖1（b）中，垂足m在線段l外，且Vpa＜Vpb，則Vpa為點p到線段l的距離。

圖1 不同情況下點到線段的距離Fig.1 Distance from point to line segment in different cases

圖2 不同情況下點到多邊形的距離Fig.2 Distance from point to polygon in different cases

距離為標量，用于平差時缺乏方向性，故引入距離矢量的概念：若式（1）中，距離值在x＝a處取得，則定義點p到點a的矢量為點p到多邊形A的距離矢量，記為ρ（p，A）。特別的，如果ρ（p，A）＝0，則ρ（p，A）＝0。

點p到多邊形A的距離矢量算法如下。

步驟1：令ρ為一個很大的值max，令l為多邊形A的第1條邊。

步驟2：計算點p到線段l的距離矢量d及距離d。

步驟3：如果d＜ρ，則令ρ＝d，ρ＝d。

步驟4：如果l是多邊形A的最后一條邊，則轉到步驟5；否則，令l為多邊形A的下一條邊，并轉到步驟2。

步驟5：輸出點p到多邊形A的距離矢量ρ。

2.2 多邊形到多邊形的距離計算方法

在點到多邊形距離的基礎上給出多邊形到多邊形距離的定義：多邊形A到多邊形B的距離為多邊形A所有邊界點到多邊形B的距離中最大者。

特別說明，多邊形的邊是由連續(xù)的邊界點組成的，而多邊形的頂點是多邊形邊的端點。圖3說明了多邊形的邊界點和頂點的區(qū)別。

圖3 多邊形的邊界點和頂點Fig.3 Boundary points and vertices of a polygon

多邊形A到多邊形B的距離可以表示為

式中，p為多邊形A的邊界點。

利用該定義計算多邊形之間的距離矢量時間復雜度較大，下面給出一個命題，可以簡化多邊形距離的計算。

命題1：多邊形A到多邊形B的距離為多邊形A所有頂點到多邊形B的距離中的最大值。

借助命題1，如果記多邊形的頂點集合為~A，則多邊形A到多邊形B的距離可以表示為

式中，p為多邊形A的頂點；x為多邊形B的邊界點。

式（3）中，若多邊形A到多邊形B的距離值在p＝a，x＝b處取得，則定義點a到點b的矢量為多邊形A到多邊形B的距離矢量，記為ρ（A，B）。特別的，如果ρ（A，B）＝0，則ρ（A，B）＝0。

多邊形A到多邊形B的距離矢量算法如下。

步驟1：令ρ＝0，令p為多邊形A的第一個頂點。

步驟2：計算點p到多邊形B的距離矢量d及距離d。

步驟3：如果d＞ρ，則令ρ＝d，ρ＝d。

步驟4：如果p是多邊形A的最后一個頂點，則轉到步驟5，否則，令p為多邊形A的下一個頂點，并轉到步驟2。

步驟5：輸出多邊形A到多邊形B的距離矢量ρ。

2.3 基于面特征的幾何校正模型建立方法

經(jīng)典的成像幾何模型用來建立地面點三維空間坐標與相應像點二維平面坐標之間的關系，一般的，通用模型可表示為

式中，t＝［t1，t2，…，tn］T表示傳感器幾何校正模型的n個參數(shù)；（X，Y，Z）表示控制點的地面坐標；（x，y）表示控制點在影像上的量測坐標。目前，常用的幾何校正模型有嚴格成像模型、仿射變換模型、多項式模型、有理函數(shù)模型等，這些模型均可用式（4）來表示［6］，本文的方法可用于各種幾何校正模型。下面說明利用控制面特征建立誤差方程的過程。

式（5）表示多邊形A到多邊形A0的距離矢量為0。由上節(jié)的定義可知，多邊形A到多邊形A0的距離必在多邊形A的某一頂點處取得，不妨記為ak（xk，yk）∈~A，其中，k為整數(shù)，且1≤k≤m，則

式中，A0為常量，矢量ρ（ak，A0）具有水平方向和垂直方向的兩個分量，可記ρ（ak，A0）＝（ρx，ρy），且根據(jù)式（4）有

式中，（Xk，Yk，Zk）為地面多邊形A′的第k個頂點坐標，則式（6）等價于式（8）

對式（8）進行線性化可以得到控制面特征的誤差方程

式中，vx和vy為隨機誤差；lx＝－ρx（fx，fy），ly＝－ρy（fx，fy）；Δt＝［Δt1，Δt2，…，Δtn］T表示t的改正向量。

這里ρx和ρy沒有解析形式，可采用數(shù)值微分［16］的方法來逼近各函數(shù)偏導，如

一個控制面特征可以導出如式（9）的兩個誤差方程，因此利用不少于n／2個控制面特征即可求解出模型中的n個未知參數(shù)t，參數(shù)的求解可使用 LM 算法［17－18］。

3 試驗與分析

為驗證上述方法的有效性，本文選用Landsat TM、ALOS PRISM和QuickBird 3種影像進行了試驗。Landsat影像為安徽地區(qū)影像，攝于2008年4月，影像大小為6856像素×5733像素，分辨率為30m，影像范圍內最大高程差為400m左右；ALOS影像為吉林地區(qū)影像，攝于2010年9月，影像大小為29 493像素×16 000像素，分辨率為2.5m，影像范圍內最大高程差為400m左右；QuickBird影像為西寧地區(qū)的影像，攝于2010年10月，影像大小為24 523像素×33 549像素，分辨率為0.6m，影像范圍內海拔均在2000m以上，最大高程差為1000m左右。

對每種影像分別進行4組試驗。第1組試驗中，均勻地選取若干個控制面特征，用本文方法進行幾何校正；第2組試驗中，在第1組試驗所選的控制面特征上各選取一對對應的頂點作為控制點，用這些控制點進行幾何校正；在1、2組試驗的控制面和控制點中人為地加入兩個粗差（5～15個像素），對每種影像的控制面和對應的控制點所加入的粗差相同；第3組試驗用含有粗差的控制面進行幾何校正；第4組試驗用含有粗差的控制點進行幾何校正。利用相同的檢查點進行四組試驗的精度評價。Landsat影像的試驗采用其衛(wèi)星軌道模型［19］，ALOS影像和QuickBird影像的試驗采用基于像方補償?shù)挠欣砗瘮?shù)模型［20］。面狀特征的提取采用基于曲線演化的水平集方法［21］。

圖4～圖6分別顯示了Landsat影像、ALOS影像和QuickBird影像試驗中控制點、控制面和檢查點的分布及示意圖。各圖中左邊為控制點、控制面和檢查點在整景影像中的分布圖，其中青色三角形表示檢查點，紅色圓點表示控制點或控制面（圖像縮小后控制面用點代替顯示，與控制點重合）；右圖為局部放大后的控制面和控制點示意圖，控制面用紅色輪廓表示，控制點用紅色圓點表示。

圖4 Landsat影像控制點、控制面及檢查點示意圖Fig.4 GCPs，GCAs and check points of Landsat image

圖5 ALOS影像控制點、控制面及檢查點示意圖Fig.5 GCPs，GCAs and check points of ALOS image

圖6 QuickBird影像控制點、控制面及檢查點示意圖Fig.6 GCPs，GCAs and check points of QuickBird image

表1列出了3種影像的4組試驗結果，包括控制面的擬合精度、控制點的擬合精度以及檢查點的精度。

表1 試驗結果對比Tab.1 Result comparison of tests

通過對表1結果的分析，可以得出以下結論。

（1）基于控制面的幾何校正方法可以用于不同的衛(wèi)星影像和不同的成像模型。

（2）當控制資料不含粗差時，使用控制面和控制點建立校正模型，圖像校正精度基本一致。如表1所示，在Landsat、ALOS及QuickBird影像無粗差的試驗中，使用控制面和使用控制點建立校正模型，進行幾何校正，控制點和檢查點的精度均能達到一個像素以內，兩種方式的校正精度基本一致。

（3）試驗中控制面的殘差為根據(jù)幾何模型計算出的多邊形到參考多邊形的距離，由多邊形到多邊形距離的定義可知，殘差為計算多邊形所有頂點到參考多邊形距離的最大值，因此控制面的殘差一般大于控制點的殘差。

（4）當控制資料中含有粗差時，基于控制面的幾何校正模型比基于控制點的幾何校正模型具有更強的容錯能力。如表1所示，在Landsat和QuickBird影像含有粗差的試驗中，用控制點進行幾何校正后，檢查點的最大殘差超過了3像素，而用控制面進行幾何校正后，檢查點的最大殘差在2像素以內；在ALOS影像含有粗差的試驗中，用控制點進行幾何校正后，檢查點的最大殘差超過了3像素，而用控制面進行幾何校正后，檢查點的最大殘差在1像素以內。其主要原因為：在Landsat和QuickBird的面特征校正試驗中，加入粗差后，粗差點參與模型的解算會導致模型參數(shù)偏離真值，因此，該粗差對基于面特征的校正模型精度影響力受到限制；而在ALOS的面特征校正試驗中，加入粗差后，粗差點落到了參考面多邊形內部，其到參考面多邊形的距離為0，對多邊形距離的計算沒有貢獻，因此，該粗差對基于面特征的校正模型精度幾乎沒有影響。

（5）在QuickBird影像無粗差試驗中，所選控制面特征的影像多邊形和地面多邊形的邊界有一定的偏差，因此，控制多邊形的擬合精度稍差，超過了2像素，但使用控制面進行幾何校正的精度與使用控制點進行幾何校正的精度相差不大，可以達到1像素左右。這也表明控制面具有較好的容錯性。

（6）QuickBird影像范圍海拔均在2000m以上，最大高程差為1000m左右，表1中QuickBird影像試驗結果表明，基于控制面特征可用于高山地區(qū)影像的幾何校正。

4 結 論

在實際應用中，足夠數(shù)量和足夠準確的地面控制點有時難以獲得，而參考影像、數(shù)字線劃地圖、GIS矢量數(shù)據(jù)中存在大量的面特征，借助本文基于面特征的幾何校正方法能夠充分地利用這些數(shù)據(jù)。

（1）本文方法的本質仍然是點與點之間的計算，但其不要求地面多邊形和影像多邊形的頂點間具有一一對應的關系，因此對控制面的約束要求要比傳統(tǒng)的控制點弱的多。尤其是對于圓弧狀的地物，通常情況下無法找到對應的特征點，但將它作為一個整體的控制面，就可以用來進行圖像的幾何校正。

（2）本文方法不受具體成像模型的限制，可以根據(jù)影像自身的特點使用衛(wèi)星軌道模型、有理函數(shù)模型、仿射變換模型等，具有通用性。

（3）由于面特征是由許多點組成的，個別點的誤差對于控制面本身來說影響并不大，因此，基于面特征的幾何校正方法比基于點特征的幾何校正方法具有更強的容錯能力。試驗結果表明本文的方法可用于不同的衛(wèi)星影像和成像模型，當控制資料不含粗差時，控制面和控制點的校正精度基本一致，可達到1像素以內；當控制資料中含有粗差時，用控制面校正的精度明顯優(yōu)于用控制點校正的精度。

本文方法是基于已有的控制面特征進行的，而控制面特征提取和匹配的自動化程度較低、耗時較多，提高控制資料準備過程的自動化程度是后續(xù)工作的一個方向。此外，由于成像時間、成像角度、數(shù)據(jù)源等條件的不同，面特征自動提取的方法對不同影像所提取的同名面特征往往并不具有完全一致的輪廓，這種不一致對幾何校正結果的影響還需要進一步的研究和分析。

［1］ JI Shunping，YUAN Xiuxiao.Automatic Matching of High Resolution Satellite Images Based on RFM ［J］.Acta Geodatetica et Cartographica Sinica，2010，39 （6）：592－598.（季順平，袁修孝.基于RFM的高分辨率衛(wèi)星遙感影像自動匹配研究［J］.測繪學報，2010，39（6）：592－598.）

［2］ YANG Huachao，ZHANG Lei，YAO Guobiao，et al.An Automated Image Registration Method with High Accuracy Based on Local Homography Constraint［J］.Acta Geodatetica et Cartographica Sinica，2012，41（3）：401－408.（楊化超，張磊，姚國標，等.局部單應約束的高精度圖像自動配準方法［J］.測繪學報，2012，41（3）：401－408.）

［3］ ZHANG Guo，CHEN Tan，PAN Hongbo，et al.Patch－based Least Squares Image Matching Based on Rational Polynomial Coefficients Model［J］.Acta Geodatetica et Cartographica Sinica，2011，40（5）：592－597.（張過，陳鉭，潘紅播，等.基于有理多項式模型的物方面元最小二乘匹配［J］.測繪學報，2011，40 （5）：592－597.）

［4］ LI Cailin，GUO Baoyun，LI Chang.The High－accurate Extraction of Line Features of Object Contour［J］.Acta Geodatetica et Cartographica Sinica，2011，40（1）：66－70.（李彩林，郭寶云，李暢.目標輪廓直線特征的高精度提?。跩］.測繪學報，2011，40（1）：66－70.）

［5］ LONG Tengfei，JIAO Weili，WANG Wei.Geometric Rectification Using Feature Points Supplied by Straight－lines［J］.Procedia Environmental Sciences，2011，11：200－207.

［6］ LONG Tengfei，JIAO Weili，WANG Wei.Block Adjustment for Generic Geometric Model Using Points and Straight－Lines［J］.Advanced Materials Research，2011，268－270：584－589.

［7］ TANG Jinjun，CAO Kai.An Adaptive Trajectory Curves Map－matching Algorithm ［J］.Acta Geodaetica et Cartographica Sinica，2008，37（3）：308－315.（唐進君，曹凱.一種自適應軌跡曲線地圖匹配算法［J］.測繪學報，2008，37（3）：308－315.）

［8］ ZHANG Ying，ZHANG Jingxiong.Edge Detection and Areal Feature Extraction Based on Improved LOG Operator［J］.Journal of Geomatics，2011，36（5）：8－10.（張盈，張景雄.基于LOG改進算子的邊緣檢測與面狀特征提?。跩］.測繪信息與工程，2011，36（5）：8－10.）

［9］ XIN Liang，ZHANG Jingxiong.Fast Extraction of Conjugated Area Features and Accurate Registration of Remote Sensing Image［J］.Geomatics and Information Science of Wuhan University，2011，36（6）：678－682.（辛亮，張景雄.共軛面狀特征的快速提取與遙感影像精確配準［J］.武漢大學學報：信息科學版，2011，36（6）：678－682.）

［10］ WANG Xin.Research on Identical Entity Geometric Matching in Multi－source Spatial Data［D］.Zhengzhou：Information Engineering University，2008.（王馨.多源空間數(shù)據(jù)同名實體幾何匹配方法研究［D］.鄭州：信息工程大學，2008.）

［11］ HAO Yanling，TANG Wenjing，ZHAO Yuxin，et al.Areal Feature Matching Algorithm Based on Spatial Similarity［J］.Acta Geodaetica et Cartographica Sinica，2008，37（4）：501－506.（郝燕玲，唐文靜，趙玉新，等.基于空間相似性的面實體匹配算法研究［J］.測繪學報，2008，37（4）：501－506.）

［12］ DONG Xiaohua，DENG Susu，SHI Wenzhong.A Probabilistic Theory－based Matching Method ［J］.Acta Geodaetica et Cartographica Sinica，2007，36（2）：210－217.（董小華，鄧愫愫，史文中.基于概率的地圖實體匹配方法［J］.測繪學報，2007，36（2）：210－217.）

［13］ DARE P，DOWMANR I.An Improved Model for Automatic Feature－based Registration of SAR and SPOT Images［J］.ISPRS Journal of Photogrammetry ＆ Remote Sensing，2001，56（1）：13－28.

［14］ ZHANG Dengrong，YU Le，CAI Zhigang.A Region Feature Based Automatic Matching for Optical and SAR Images［J］.Journal of China University of Mining ＆ Technology.2007，36（6）：843－847.（張登榮，俞樂，蔡志剛.基于面特征的光學與SAR影像自動匹配方法［J］.中國礦業(yè)大學學報，2007，36（6）：843－847.）

［15］ ZHANG Xiaodong，LI Deren，GONG Jianya，et al.A Matching Method of Remote Sensing Image and GIS Data Based on Area Feature［J］.Journal of Remote Sensing，2006，10（3）：373－380.（張曉東，李德仁，龔健雅，等.一種基于面特征的遙感影像與GIS數(shù)據(jù)配準方法［J］.遙感學報，2006，10（3）：373－380.）

［16］ BURDEN R L，F(xiàn)AIRES J D.Numerical Analysis［M］.Ninth Edition.Boston：Brooks／Cole，2011：174－179.

［17］ MORE J J.The Levenberg－Marquardt Algorithm：Imple－Mentation and Theory［J］.Numerical Analysis：Lecture Notes in Mathematics，1978，630：105－116.

［18］ LOURAKIS M I A.A Brief Description of the Levenberg－Marquardt Algorithm Implemented by Levmar［EB／OL］.2005－02－11［2012－01－02］.http：∥ www.ics.forth.gr／lourakis／levmar／.

［19］ CHEN Pengshan，JIAO Weili，JIA Xiupeng.Robust LM Algorithm and Its Application on Rigorous Physical Model［J］.Science Technology and Engineering，2009，9（16）：4614－4618.（陳朋山，焦偉利，賈秀鵬.抗差LM算法求解遙感影像嚴格物理模型［J］.科學技術與工程，2009，9（16）：4614－4618.）

［20］ GRODECKI J，DIAL G.Block Adjustment of Highresolution Satellite Images Described by Rational Polynomials［J］.Photogrammetric Engineering ＆ Remote Sensing，2003，69（1）：59－68.

［21］ SHI Y G，KARL W C.A Real－time Algorithm for the Approximation of Level－set Based Curve Evolution ［J］.IEEE Transaction on Image Processing，2008，17 （5）：645－656.