遙感影像GCP選取方法研究

王圣堯，劉善磊，石善球，王光輝

(1. 江蘇省基礎地理信息中心，江蘇 南京 210013； 2. 國家測繪地理信息局衛星測繪應用中心，北京 100830)

一、引 言

地面控制點(ground control point，GCP)的選取是遙感影像幾何校正的重要步驟之一。從目前發展來看，傳統的利用參考圖目視獲取GCP的方式還將繼續存在，因此有必要研究一些新的GCP選取方式，以輔助人工選點，減少人為因素的影響；此外，GCP的分布、數量等指標直接影響遙感影像幾何精度，因此選取分布均勻、數量足夠的地面控制點非常重要。

GCP的選取是近年來遙感圖像處理領域研究的熱點之一。林乘逸等[1]通過數學形態學方法選取道路交叉點作為GCP；王永明等[2]綜合SIFT(scale invariant feature transform)局部特征和Harris角點檢測的優點實現了同尺度圖像間的GCP自動選取；袁修孝等[3]采用邊緣特征作為配準要素，以2期圖像上的邊緣匹配點應滿足邊緣梯度角之差與邊緣點間距離的加權和最小的條件作為匹配策略；張萍[4]采用歸一化相關系數進行同名點的精確定位，給出了精確的匹配結果，實現了SAR圖像GCP的自動匹配；楊小岡等[5]設計了一種基于影像方差的檢測模板，提出了相應的檢測準則及搜索方案，并在畸變圖像中利用歸一化積相關影像匹配算法實現了與基準控制點相對應的畸變控制點的精確定位；馬靈霞等[6]利用FAST(features from accelerated segment test)算子在基準底圖上快速選取均勻分布的候選特征點，通過圖像自身攜帶的地理定位信息確定初始同名點對，經處理后獲取亞像素的GCP；周家香等[7]將Voronoi圖應用于GCP的選取中，減少了幾何校正中GCP選取的主觀性。此外，有的學者提出建立GCP影像數據庫，通過GCP影像在待校正影像上匹配獲取控制點[8]。

本文在眾多學者的研究基礎上，以0.5 m分辨率的航空影像作為基準影像、2.5 m分辨率的ALOS影像作為待校正影像，研究不同尺度、異源影像間的GCP選取。針對遙感影像中的道路交叉口及具有典型幾何特征的永久性地物，在基準影像上采用GrabCut分割的方法輔助人工選點，并將Voronoi圖應用于GCP的選取中，減少GCP選取的主觀性。

二、研究方法

遙感影像GCP選取方法的基本思想是：在人工參與下根據基準影像和待校正影像自帶的地理定位信息，獲取道路交叉點和典型地物重心為目標創建GCP選取視窗。如果已有Voronoi圖，則根據Voronoi圖自動創建GCP選取視窗。在視窗內根據GrabCut分割算法獲取目標區域，經過二值化、數學形態學預處理后根據不同的目標分別采用兩種方法選取GCP：① 細化(骨架化)和加權鄰域檢測方法選取道路交叉點作為GCP；② 輪廓檢測和橢圓擬合方法選取區域重心作為GCP。當誤差在精度要求范圍內且GCP數量及分布達到最優，則形成GCP，同時更新Voronoi圖，否則重新選擇GCP選取視窗獲取GCP。其流程如圖1所示。

三、基于GrabCut分割的GCP選取

1. GrabCut分割

將圖割理論應用于計算機視覺領域是近年來發展起來的一種新方法。GrabCut 算法是對Graph cuts算法的改進：一是用高斯混合模型(Gaussian mixture model，GMM)取代灰度直方圖，實現了彩色圖像的分割；二是采用多次迭代方式取代一次估計GMM參數，完成能量最小化，提高分割精度；三是通過非完全標號降低用戶的交互工作量[9]。

圖1 GCP選取流程圖

GrabCut分割的基本過程是：首先在目標區域周圍畫一個矩形框，初始化三元圖T，以非完全標號的方式將矩形框內標為未知區域Tu、矩形框外標為背景區域TB，將Tu內的像素值標為1、TB內的像素值標為0，同時設目標區域TF為空，根據用戶標定的像素標號值為0和1的集合，分別進行目標和背景的GMM初始化；然后將Tu劃分為目標和背景2類，對新劃分的目標和背景像素點進行切割，更新GMM參數，在新的GMM參數下繼續對Tu迭代，直到滿足迭代終止條件，確定GMM參數；最后根據確定的GMM對Tu進行一次切割，得到目標區域[10]。

2. 加權鄰域檢測

假設p為待檢測的交叉點，n(p)表示p的8鄰域內非0像素的個數，在通過GrabCut分割算法獲取目標區域的基礎上(如圖2(b)所示)，根據文獻[11]完成細化后的影像。要想檢測出交叉點，算法只需滿足n(p)≥3。雖然影像細化前經過了數學形態學的腐蝕運算和開運算[12]及去除碎小區域等預處理，但仍會出現多點檢測的問題(如圖2(c)所示)。

為了確保道路交叉點檢測結果的唯一性和正確性，考慮從3個方面約束交叉點，提出一種加權鄰域檢測的方法，取權重最大的點作為要獲取的交叉點。

(1) 交叉口類型約束

十字路口權重為0.6，丁字路口權重為0.4，可在一定程度上確保十字路口作為交叉點的優先性。

(2) 交叉點相關線段長度約束

(3) 交叉點距細化后影像重心的距離約束

假設重心位置為(xc,yc)，通過

圖2 道路交叉點選取

3. 橢圓擬合

一般情況下，影像中的物體不是一個點，而是采用其重心作為物體的位置。對于目標區域(如圖3(b)所示)，雖然經過前期的數學形態學腐蝕運算和開運算等預處理得到了比較平滑的邊界，但是由于圖像中物體的復雜性，實際區域并非規則區域(人眼處理時是基于規則區域的)，這樣獲取的區域重心和實際物體重心會有偏差，因此當根據文獻[10]完成輪廓檢測后，必須對它進行橢圓擬合(如圖3(c)所示)。一個區域的重要形狀結構信息是它的位置、尺寸和方向，而橢圓能夠有效地反映這3種信息，而且根據圖像矩的相關理論可知，橢圓中心即為區域重心，橢圓的2個軸分別通過重心，相對于2個軸的二階中心矩在長軸和短軸方向上分別達到最大和最小。

圖3 典型地物重心選取

四、基于Voronoi圖的GCP分布與數量分析

Voronoi圖又稱泰森多邊形，在二維空間中，對含有若干離散點的平面進行有序分割，用于表達離散點之間的鄰近關系與影響范圍等重要信息，是解決空間劃分、最優資源配置、最佳選址位置及目標點合理分布等問題的有效技術手段[13]。Voronoi圖將含有n個離散點的平面分割為n個多邊形，每個多邊形含有一個離散點，該多邊形內部的任一點到特征點的距離都小于到其他特征點的距離。

假設P={p1,p2,...,pn}為平面上的一個離散點集，即

式中，i=1,2，…，n。由V(pi)對該平面進行劃分，稱為以pi為生成元的Voronoi圖，簡稱Voronoi圖。其中，pi被稱作一個生成元，包含pi的多邊形被稱作pi生成元的Voronoi區，pi生成元與相鄰的pi+1生成元具有公共的Voronoi邊，d(p,pi)為點p與pi之間的歐氏距離。

Voronoi圖可以定義為由距該點比所有其他點更近的目標所組成的區域，其中每一個空間生長目標對應唯一一個Voronoi圖多邊形，凡落在其Voronoi圖多邊形范圍內的空間點距其最近。因此Voronoi圖多邊形在一定程度上反映了其影響范圍，或稱為勢力范圍。而GCP的均勻分布與否可由每個GCP的控制范圍來度量，因此GCP的影響范圍可由Voronoi圖來度量。由GCP生成Voronoi圖，計算各Voronoi圖多邊形的面積大小，如果各區域面積差異較大，則在面積較大的區域內增加符合單點控制精度要求的GCP，直到各分區面積相當；同時檢驗GCP總單點殘差(RMS)和幾何校正精度以確定GCP的數量，當兩者同時符合精度要求且趨于穩定時則完成GCP的選擇。

五、試驗結果與分析

試驗數據為位于東經118°25′—119°12′、北緯31°45′—32°32′的一幅2011年未校正的2.5 m分辨率的ALOS影像(GCS_WGS_1984)，以及一幅對應范圍內0.5 m分辨率的航空影像(GCS_Xian_1980)。試驗區以南京市為主，另包括鎮江市和揚州市部分區域，面積約為3 278.47 km2，擁有驗證幾何校正精度的14個GCP(GCS_Xian_1980)(如圖4所示)。

首先，在保證試驗區邊角優先選取GCP的條件下，以RMS不大于2個像素、總誤差在1個像素以內的標準，通過GrabCut分割算法選取GCP，其中表1和表2是分別利用GrabCut分割和人工方法對道路交叉點(如圖2(a)所示)和典型地物(如圖3(a)所示)分5次選取的單點精度；其次，在滿足去除邊角不完整的Voronoi圖后，各分區最大面積與最小面積比例小于2的條件下，基于Voronoi圖的GCP分布與數量進行試驗，數據見表3、表4，圖5為GCP數目為44的Voronoi圖。

圖4 試驗區及GCP分布 圖5 Voronoi圖(GCP數目:44)

通過分析表1和表2中的平均值和均方差可知，采用本文的方法得到的GCP，即使由不同的作業員采集，彼此間的誤差一般也不會超過1個像素。而這1個像素以內的誤差主要是由于人工獲取待選取GCP的圖像塊時產生的。完全人工選取GCP與本文采用的先人工粗定位后算法選取的方法相比，受作業員的經驗、工作態度、判斷力和反應能力等因素的影響較大，容易出現較大誤差，需要重新定位，甚至重新選取，這在一定程度上影響了GCP的采集效率。通過分析表1和表2中的點位均方差可知，當待選取GCP地物面積較大時，更能凸顯本文方法的選取精度。

表1 道路交叉點數據 像素

表2 典型地物重心數據 像素

表3 不同數目GCP的精度 像素

表4 不同數目GCP的Voronoi圖面積 km2

通過分析表3和圖6可知，在單點RMS不大于2個像素的條件下，當GCP數目超過33個時，總RMS和幾何校正精度基本趨于穩定。通過分析表4和圖7可知，對試驗區GCP建立Voronoi圖后，區內邊角處總有不完整的Voronoi圖(如圖5所示)，若統計全部Voronoi圖面積，各分區最大面積與最小面積能控制在10∶1以內，當GCP數目超過44個時，面積均方差趨于穩定；若只取完整Voronoi圖且Voronoi圖各分區最大面積與最小面積控制在2∶1以內，當GCP數目超過33個時，面積均方差基本趨于穩定。

綜上可得，利用0.5 m分辨率的航空影像作為基準影像，可以為 2.5 m分辨率ALOS影像的幾何校正提供單點定位精度可靠的GCP，本區域需要33～44個GCP，單點控制范圍為87～130 km2。

圖6 GCP數目與精度關系

圖7 GCP數目與Voronoi圖面積均方差關系

六、結束語

針對不同尺度、異源影像間的GCP選取，本文利用0.5 m分辨率的航空影像作為基準影像，為 2.5 m分辨率ALOS影像的幾何校正提供GCP，將道路交叉口及具有典型幾何特征的永久性地物作為研究對象，采用GrabCut分割的方法輔助人工選點；為得到分布均勻的GCP，將Voronoi圖應用于GCP的選取中，以確定GCP的數目，同時檢驗GCP的總RMS和幾何校正精度。試驗結果表明，采用本文方法可以減少人為因素對GCP選取的影響，對2.5 m分辨率的衛星遙感影像的幾何校正有一定的指導意義，GCP單點控制范圍為87～130 km2。自動選取GCP及加權Voronoi圖在GCP選取中的應用是今后研究的重點。

參考文獻：

[1] 林乘逸，陳繼潘.高解析力衛星影像控制點之自動化萃取[J].航測及遙測學刊，2003，8(2)：1-11.

[2] 王永明，王貴錦.圖像局部不變性特征與描述[M].北京：國防工業出版社，2010：182-186.

[3] 袁修孝，宋妍.基于邊緣特征匹配的遙感影像變化檢測預處理方法[J].武漢大學學報：信息科學版，2007，32(5)：381-384.

[4] 張萍.SAR影像GCP自動匹配方法的實現[J].測繪通報，2008(1)：30-32.

[5] 楊小岡，繆棟，付光遠，等.影像幾何校正控制點的自動確定方法[J].紅外與激光工程，2002，31(3)：194-198.

[6] MA L X，ZOU T Y，XU J. Remote Sensing Image Geometric Accurate Rectification Based on Automatic Matching and Triangulation[J].Journal of Remote Sensing，2011，15(5)：927-932.

[7] 周家香，朱曉東，周曉光.基于Voronoi圖的遙感影像校正地面控制點選取[J].地理信息世界，2007，5(2)：67-70.

[8] 周剛，徐德偉，楊燕明，等.基于IDL語言與控制點影像數據庫的圖像自動配準技術研究[J].遙感技術與應用，2010，25(5)：627-631.

[9] 徐秋平，郭敏，王亞榮.基于分水嶺變換和圖割的彩色圖像快速分割[J].計算機工程，2009，35(19)：210-212.

[10] 王建青，郭敏，徐秋平.結合小波變換與圖像分割的快速目標選取[J].計算機工程與應用，2010，46(33)：215-217.

[11] 孫燮華.數字圖像處理原理與方法[M].北京：機械工業出版社，2012：169-175.

[12] SOILLE P，王小鵬.形態學圖像分析原理與應用[M].北京：清華大學出版社，2008：77-103.

[13] 王新生，李全，郭慶勝，等.Voronoi圖的擴展、生成及其應用于界定城市空間影像范圍[J].華中師范大學學報：自然科學版，2002，36(1)：107-111.