DEM粗差剔除方法研究

1、第 頁 （共 18頁）DEM粗差剔除方法研究1引言1.1 選題意義數(shù)字高程模型(Digital Elevation Model,簡稱DEM)是以數(shù)字的形式按一定結(jié)構(gòu)組織在一起，表示實(shí)際地形特征空間分布的模型，也是地形形狀大小和起伏的數(shù)字描述，由一系列地面點(diǎn)X，Y位置及其相聯(lián)系的高程Z組成，是各種信息的載體，是地理信息系統(tǒng)(Geographic Information System，GIS)技術(shù)中最重要的內(nèi)容，是空間數(shù)據(jù)基礎(chǔ)設(shè)施的重要組成部分，在生產(chǎn)中具有很高的利用價值。首先，它能夠反映區(qū)域內(nèi)的地形條件，為各用圖部門提供地形基礎(chǔ)，其次，用它制作沙盤，具有快速、簡便、精確的優(yōu)點(diǎn)，可用于軍事指揮和模

2、型演示；還可以用于農(nóng)業(yè)部門的農(nóng)田水利規(guī)劃，水利部門的洪水淹沒損失估算及水利建設(shè)的土方量計算，交通、建筑等各建設(shè)部門的選址、規(guī)劃，通訊部門信號覆蓋范圍的規(guī)劃與計算，地址、勘探部門地形的分析，各旅游景點(diǎn)的規(guī)劃及土地資源調(diào)查；也能用于二維地理空間上連續(xù)分布并逐漸變化的各種非高程屬性數(shù)據(jù)的建模與分析上?？梢哉fDEM數(shù)據(jù)具有廣泛的應(yīng)用潛力。目前，由于地理信息系統(tǒng)(GIS)的普及以及空間數(shù)據(jù)基礎(chǔ)設(shè)施的發(fā)展和建設(shè)，DEM作為數(shù)字線化圖(DLG)、數(shù)字高程模型(DEM)、數(shù)字正射影像圖(DOM)和數(shù)字柵格圖(DRG)即“4D”產(chǎn)品之一，己經(jīng)成為國家空間數(shù)據(jù)基礎(chǔ)設(shè)施(National Spatial Data

3、Infrastructure，NSDI)的基本內(nèi)容之一，并被納入數(shù)字化空間數(shù)據(jù)框架(DGDF)進(jìn)行規(guī)?；a(chǎn)。今天，DEM已作為一個獨(dú)立的產(chǎn)品而存在，并越來越廣泛地被用來代替?zhèn)鹘y(tǒng)地圖中等高線對地形的描述，成為地理信息系統(tǒng)的核心數(shù)據(jù)庫以及地學(xué)分析的基礎(chǔ)數(shù)據(jù)。與其他工業(yè)產(chǎn)品一樣，DEM產(chǎn)品也必須有質(zhì)量管理和質(zhì)量控制。由于影響DEM質(zhì)量的因素是多種多樣的，無論采用哪種途徑都不能很好地解決所有的問題。從誤差理論的角度分析，可以劃分為系統(tǒng)誤差、隨機(jī)誤差和粗差。其中，系統(tǒng)誤差總是與特定的觀測系統(tǒng)聯(lián)系在一起，它們可能是由于立體像對的不正確定向造成的，也可能是攝影測量生產(chǎn)過程中忽略了棱鏡變形而導(dǎo)致。從統(tǒng)計意義

4、上說，系統(tǒng)誤差給觀測值帶來的偏差在空間上并不是一個常數(shù)，但當(dāng)它們己知時，可用一定的函數(shù)關(guān)系式來表示。隨機(jī)誤差在大小和符號上表現(xiàn)出偶然性，假如存在足夠的觀測值，隨機(jī)誤差服從下列特征:(1)正負(fù)誤差發(fā)生的概率大致相同;(2)小誤差發(fā)生的概率高于大誤差;(3)大誤差很少發(fā)生。隨機(jī)誤差既不可以通過反復(fù)觀測加以剔除，也不能用一定的函數(shù)關(guān)系式表示。在空間上，它們呈弱相關(guān)性，并在整個DEM中服從高斯分布。粗差是由于觀測者在操作儀器時的粗心，或讀數(shù)、記錄時的錯誤造成的，從統(tǒng)計的角度分析，這類觀測值與其他觀測值并不屬于同一分布。與系統(tǒng)誤差和隨機(jī)誤差相比較，粗差對DEM數(shù)據(jù)所反映的空間變化的扭曲更為嚴(yán)重，特別是現(xiàn)

5、代觀測數(shù)據(jù)量大的情形。因此，在現(xiàn)代測量數(shù)據(jù)處理中如何消除粗差的影響，特別是自動化觀測水平高，數(shù)據(jù)量大的情況下，就顯得越來越重要。粗差它的存在會導(dǎo)致DEM及產(chǎn)品的嚴(yán)重失真，甚至完全不能使用。因此，設(shè)計一些算法檢測并剔除DEM數(shù)據(jù)中的粗差是完全必要的。因此探測和剔除DEM數(shù)據(jù)中存在的粗差顯得尤為重要。本課題是在對現(xiàn)有的DEM粗差探測法進(jìn)行分析的基礎(chǔ)上，重點(diǎn)分析基于趨勢面的探測DEM粗差的方法及粗差的剔除，并證明其可行性。1.2 研究現(xiàn)狀 DEM是地理空間定位的數(shù)字?jǐn)?shù)據(jù)集合，它最初是美國麻省理工學(xué)院Miller教授為高速公路的自動設(shè)計于1956年提出來的，隨著各種相關(guān)技術(shù)的發(fā)展,特別是計算機(jī)技術(shù)在測

6、繪方面的應(yīng)用使得測繪學(xué)科逐漸向數(shù)字化、實(shí)時處理與多用途的方向發(fā)展1。質(zhì)量控制是數(shù)字高程模型(Digital Elevation Model- DEM) 生產(chǎn)的重要環(huán)節(jié)之一。然而粗差對DEM 數(shù)據(jù)所造成的空間扭曲往往最為嚴(yán)重,有時能導(dǎo)致DEM及其產(chǎn)品嚴(yán)重失真,甚至完全不能使用,因此很有必要研究一些方法對DEM的粗差進(jìn)行探測和修正。要有效地探測DEM數(shù)據(jù)中存在的粗差,必須在數(shù)據(jù)采集后形成的原始DEM 數(shù)據(jù)中進(jìn)行。原始的DEM數(shù)據(jù)結(jié)構(gòu)有規(guī)則和不規(guī)則兩類。從實(shí)用的角度來看,在規(guī)則格網(wǎng)的DEM數(shù)據(jù)中探測粗差相對簡單一些,因而研究成果也相對豐富,如Hannah2(1981)的基于坡度信息算法,FolicS

7、imon (1994)的統(tǒng)計參數(shù)法, 以及Lper3(1997) 的主成分分析法等等。如果原始DEM 為不規(guī)則數(shù)據(jù),要轉(zhuǎn)成規(guī)則格網(wǎng)點(diǎn),則需經(jīng)過數(shù)學(xué)內(nèi)插的方法處理。在此過程中,原始數(shù)據(jù)中粗差點(diǎn)會影響到轉(zhuǎn)換后的多個格網(wǎng)點(diǎn),從而增加格網(wǎng)點(diǎn)粗差檢測的難度。本文是在探討常用的各種算法的基礎(chǔ)上進(jìn)一步描述基于趨勢面的DEM粗差探測與剔除方法。2 DEM粗差探測與剔除方法研究2.1 基于三維可視化的DEM粗差探測與剔除通過地形表面的三維可視化建模來審查DEM中可疑數(shù)據(jù)點(diǎn)，從而剔除嚴(yán)重影響數(shù)據(jù)質(zhì)量的粗差或者說錯誤。DEM有著非常適宜于建立3維可視化的特點(diǎn)，采用DEM3維可視化技術(shù)，該方法可以交互式的來檢查DEM

8、中出現(xiàn)的可疑數(shù)據(jù)，剔除嚴(yán)重影響DEM數(shù)據(jù)質(zhì)量的粗差4。一般對于一個特定的研究區(qū)域，在三維透視圖上可疑點(diǎn)是否表現(xiàn)為粗差非常直觀，很容易據(jù)此作出正確地判斷。實(shí)際上，由于DEM有著非常適宜于建立三維可視化的特點(diǎn)，所以可以首先通過目視效果對粗差進(jìn)行檢測。通常粗差的地形很不自然，因此在實(shí)際應(yīng)用中，可以首先通過目視進(jìn)行粗差的檢測這種方法需要高效可靠的構(gòu)網(wǎng)技術(shù)、快速的交互相應(yīng)效率以及對異常值敏感的視化圖形，如線框透視圖、暈渲圖等常用的可視化圖形，在技術(shù)層面上，操作的經(jīng)驗(yàn)和工作態(tài)度對結(jié)果也會有相當(dāng)大的影響。三維可視化的前提是要建立數(shù)字地面模型，為了保證所分析都基于原始數(shù)據(jù),可選的方法是直接利用原始數(shù)據(jù)建立不規(guī)

9、則三角網(wǎng)絡(luò)模型(Tin)。該方法不利因素在于，一方面需要高效可靠的建模技術(shù)以及可視化處理的策略;另一面，它僅適用于較大粗差的判釋，對于中小粗差并不敏感。2.2 基于高程信息的不規(guī)則分布數(shù)據(jù)粗差檢測與剔除呈散亂分布的數(shù)據(jù)點(diǎn)粗差探測技術(shù)在原理上與規(guī)則格網(wǎng)比較類似，但由于散亂分布的數(shù)據(jù)點(diǎn)的分布特征，坡度信息獲取比較困難，具體實(shí)現(xiàn)上有兩點(diǎn)不同:第一，窗口確定，在規(guī)則格網(wǎng)上采用3*3局部窗口是適宜的，但不規(guī)則分布點(diǎn)的鄰域范圍要進(jìn)行指定，一般可采用窗口尺寸或窗口區(qū)域的采樣點(diǎn)數(shù)量兩種方式確定。第二，一致性標(biāo)準(zhǔn)確定，規(guī)則格網(wǎng)上比較容易獲取坡度信息，而不規(guī)則分布上獲取坡度信息比較困難，因此，由于高程和坡度同是刻

10、畫地形曲面連續(xù)性的指標(biāo)，在散亂數(shù)據(jù)分布的區(qū)域上，高程信息取代坡度成為一致性標(biāo)準(zhǔn)。在每一個窗口中，用高程信息計算統(tǒng)計指標(biāo)以及確定閡值，方法與規(guī)則格網(wǎng)類似。2.3 檢測規(guī)則分布數(shù)據(jù)中粗差的算法2.3.1 基于坡度信息的粗差檢測坡度是地表的固有屬性，在局部連續(xù)空間的漸變模型上，坡度變化也是連的，因此可采用采樣點(diǎn)與周圍點(diǎn)的坡度變化是否一致來檢測是否含有粗差，通以局部3*3窗口對每一采樣點(diǎn)進(jìn)行判斷。其基本思想是對每個表面上的點(diǎn)，在坡度上，高程或突變量引起的形狀不連續(xù)，可能被懷疑有誤差，通過坡度上每個點(diǎn)，應(yīng)用坡度逼近或改變量來計算，考慮坡度變化的相對值，并以這些相對值計算一個統(tǒng)計值為判斷該點(diǎn)合法性的閩值，

11、使計算結(jié)果更為可靠。如表1所示，P點(diǎn)在高程矩陣中的行列號為(I，J)，它的相臨8個點(diǎn)1，2，3，4，6，7，8，9的行列號分別為:1點(diǎn)(I+1，J-l)，2點(diǎn)(I+1，J)，3點(diǎn)(I+，J+1)，4點(diǎn)(I，J-l)，6點(diǎn)(I，J+1)，7點(diǎn)(I-l，J-l)，8點(diǎn)(I-1，J)，9點(diǎn)(I-l，J+1)。 表1 P點(diǎn)坡度計算1 （I+1，J-1）2（I+1，J）3（I+1，J+1）4（I，J-1）P（I，J）6（I，J+1）7（I-1，J-1）8（I-1,J）9（I-1，J+1）該算法分三步進(jìn)行，首先以檢測點(diǎn)P的8個鄰域點(diǎn)分別計算I，J方向的坡度值，然后計算各個方向的坡度變化值DSC，根據(jù)每個數(shù)

12、據(jù)點(diǎn)在同一方向上的兩個DSC值相加，其值用于計算均方根差(RMSE)，如果坡度變化一致的話，則同一點(diǎn)在同一方向上的兩個DSC值和的絕對值將是很小的值(接近0)，反之如果坡度變化不一致的話，這個值將比較大。若某一數(shù)據(jù)點(diǎn)在行列方向上的DSC值都大于閉值(閉值為RMSE的K倍)，則可確信它含有粗差。對于K值，不同情況下，可以使用不同的值，如果DSC值分布比較均勻(此時RMSE值比較小)，K可以取大一些的值;反之(RMSE值較大)，則K值以取小一些的值。對判斷含粗差的點(diǎn)進(jìn)行改正，以保證數(shù)據(jù)質(zhì)量的提高5。2.3.2 基于參數(shù)統(tǒng)計的粗差檢測假設(shè)超限誤差只是局部相關(guān)的，F(xiàn)elicisimo算法研究的對象是某

13、點(diǎn)高程值，和其鄰域點(diǎn)內(nèi)插出的高程值之，之間的差值民，如果使用雙線性內(nèi)插而言，只要使用鄰域的四個點(diǎn)就夠了，在這種情況下，對矩陣第i行，第j列交叉處點(diǎn)高程值可按下式計算6:高程估值與DEM中高程值之差為： 如果此過程應(yīng)用于DEM中所有點(diǎn)，可以得到高程差值的代數(shù)均值和標(biāo)準(zhǔn)偏差，假設(shè)服從均值為，標(biāo)準(zhǔn)偏差為（均從這個采樣中獲得）的高斯分布 ，可以用雙尾檢驗(yàn)來驗(yàn)證是否屬于該分布的集合。通過引進(jìn)t統(tǒng)計量，可以看做一個標(biāo)準(zhǔn)化殘差，因從模型中獲得的數(shù)據(jù)量很大，故可假設(shè)服從分布，對于置信水平=0.001，統(tǒng)計量的臨界值為3.219，由此可進(jìn)行兩個假設(shè)的檢驗(yàn)，其中零假設(shè)為，而被選假設(shè)為，任何殘差，使的數(shù)據(jù)點(diǎn)都被懷疑

14、為含有粗差，但事實(shí)上，大的值并不能指出粗差，而僅僅是一個警告符號。表2 3*3窗口計算單元粗差的修正應(yīng)該緊緊伴隨著粗差的探測過程，在DEM的柵格矩陣中，由于不允許空格點(diǎn)的存在，一旦某個粗差被檢測就應(yīng)該加以改正，這里最簡單的方法就是用臨近點(diǎn)的高程均值來代替可疑點(diǎn)的高程值，需要指出的是任何可疑數(shù)據(jù)點(diǎn)都不應(yīng)參與高程估值的計算。上述過程可以迭代進(jìn)行，其計算方法如表2 以3*3窗口為計算單元，每一次迭代中統(tǒng)計量和都會變化，直到?jīng)]有超過臨界值的殘差出現(xiàn)。2.3.3 基于主成分分析的粗差檢測主成分分析是把多個指標(biāo)化為少數(shù)幾個綜合指標(biāo)的一種統(tǒng)計方法7。在實(shí)際的研究中，為了全面分析問題，往往使用眾多有關(guān)的變量。

15、但是，變量太多不但會增加計算的復(fù)雜性，而且也給合理地分析問題和解釋問題帶來困難。一般來說，雖然每個變量都提供了一定的信息，但其重要性有所不同。實(shí)際上，在很多情況下，眾多變量間有一定的相關(guān)關(guān)系，人們希望利用這種相關(guān)性對這些變量加以“改造”，用為數(shù)較少的新變量來反映原變量所提供的大部分信息，通過對新變量的分析達(dá)到解決問題的目的。主成分分析便是在這種降維的思想下產(chǎn)生的處理高維數(shù)據(jù)的統(tǒng)計方法。主成分分析的基本方法是通過構(gòu)造原變量的適當(dāng)?shù)木€性組合，以產(chǎn)生一系列互不相關(guān)的新變量，從中選出少數(shù)幾個新變量并使它們含有盡可能多的原變量帶有的信息，從而使得用這幾個新變量代替原變量分析問題和解決問題成為可能。當(dāng)研究

16、的問題確定之后，變量中所含的“信息”的大小通常用該變量的方差來度量。主成分分析法進(jìn)行DEM粗差的探測包括以下幾個步驟:(l)對于給定的大小為n*m的DEM，將其分別劃分為列主方向和行主方向度為w的帶狀區(qū)域，實(shí)施下列步驟直到滿足一定的準(zhǔn)則。(2)處理列主方向的帶狀區(qū)域。a.確定可能含有粗差的列;b.在每一列中確定差的位置，從而獲取備選粗差數(shù)組1。(3)處理行主方向的帶狀區(qū)域。a.確定可能含有粗差的行;b.在每一行中確定差的位置，從而獲取備選粗差數(shù)組2。(4)將兩套備選粗差進(jìn)行比較。(5)提供計算準(zhǔn)則，改正所有的粗差。(6)結(jié)束。由于粗差數(shù)據(jù)影響統(tǒng)計量的計算，因而上述過程需要迭代進(jìn)行，在每一次代中

17、，都可以得到一組“粗差候選序列”。如果確定為真正的粗差，就對相應(yīng)數(shù)據(jù)點(diǎn)加以修正，并從“候選序列”中去除，接著進(jìn)行下一步迭代過程，并設(shè)一定的條件來控制程序的進(jìn)程。2.4 檢測不規(guī)則分布數(shù)據(jù)中粗差的算法規(guī)則格網(wǎng)DEM具有很多優(yōu)點(diǎn)，它的數(shù)據(jù)結(jié)構(gòu)簡單，便于存儲和處理;但有不足之處，它對地表的描述沒有非規(guī)則DEM的精確，也無法顧及地形變化特征點(diǎn)和線。除此以外，在實(shí)際的生產(chǎn)中，規(guī)則格網(wǎng)DEM通常都是由非規(guī)則網(wǎng)DEM通過內(nèi)插而得到的。如果原始DEM為非規(guī)則數(shù)據(jù)，要轉(zhuǎn)成規(guī)則格網(wǎng)點(diǎn)在此過程中，原始數(shù)據(jù)中粗差點(diǎn)會影響到轉(zhuǎn)換后的多個格網(wǎng)點(diǎn)，從而增加格網(wǎng)粗差檢測的難度。因此，很有必要討論和驗(yàn)證基于非規(guī)則格DEM的粗差

18、探算法。目前，最典型的基于非規(guī)則DEM的粗差探測算法是李志林提出的點(diǎn)方算法8。2.4.1 基于點(diǎn)方式的算法首先，確定待定點(diǎn)P周圍的鄰域點(diǎn)范圍，然后計算窗口范圍內(nèi)所有點(diǎn)的平高程(或加權(quán)平均值)作為P點(diǎn)的估值，最后計算P點(diǎn)高程值與估值的高程差，如果高程差值大于閩值，則認(rèn)為P點(diǎn)含有誤差。(l)鄰域點(diǎn)的范圍確定待定點(diǎn)P周圍的鄰域點(diǎn)范圍，可根據(jù)以P為中心的窗口指定，窗口的確定有三種方法，一種是定義窗口的尺寸，另一種是定義窗口覆蓋區(qū)域內(nèi)高程點(diǎn)的數(shù)量，還有一種方法是同時使用上述兩種方法來確定窗口的大小，通過計算區(qū)域內(nèi)點(diǎn)的數(shù)量和坐標(biāo)范圍確定一平均窗口。(2)代表值的計算在點(diǎn)方式算法中，把待測點(diǎn)鄰域點(diǎn)的平均高程

19、作為該點(diǎn)的代表值。有兩種方法可計算鄰域點(diǎn)的平均高程，一種是簡單的計算高程值的算術(shù)平均值，另一種是對每一個鄰域點(diǎn)賦以不同的權(quán)值。如果P的鄰域點(diǎn)包含粗差，用簡單算術(shù)平均值法更加可信，計算速度也比較快。(3)計算閉值假設(shè)城是以第i個點(diǎn)為中心的鄰域點(diǎn)的算術(shù)平均值，耳為從與第i個點(diǎn)的高程值的差值，即：對DEM中所有的點(diǎn)，可以得到一系列的值，計算均值和標(biāo)準(zhǔn)偏差9：其中，為由中心點(diǎn)的窗口范圍內(nèi)數(shù)據(jù)點(diǎn)計算的高程估值，為點(diǎn)的高程測量值。假設(shè)服從均值為，標(biāo)準(zhǔn)偏差為的正態(tài)分布，則檢測粗差的閾值為的倍（為常數(shù)）。閾值確定后，對DEM任一數(shù)據(jù)點(diǎn)，如果,則認(rèn)為含有粗差。一旦數(shù)據(jù)點(diǎn)被檢測超限，用其估值來代替可疑的高程值，迭

20、代進(jìn)行，直到?jīng)]有超過閾值的高程較差出現(xiàn)。2.4.2 基于粗差簇群的算法基于點(diǎn)方式的粗差檢測算法是針對檢測數(shù)據(jù)中僅存在單個粗差的情況，而事實(shí)上DEM中的粗差經(jīng)常以一種排列緊湊、數(shù)據(jù)巨大的簇群方式存在這在自動相關(guān)技術(shù)獲取的數(shù)據(jù)中經(jīng)常存在。與點(diǎn)方式算法相同，首先需要定義一以為中心的窗口，將窗口中的第一點(diǎn)從窗口中移去，從窗口中剩余的點(diǎn)計算新的“代表值”即平均值，然后計算并記錄這個平均值與移去的數(shù)據(jù)點(diǎn)的值之差，此過程迭代進(jìn)行，直到窗口中所有的點(diǎn)都通過檢驗(yàn)10。假設(shè)在窗口中有M個點(diǎn)，那么通過下式可計算M個差值： 式中，是窗口所有剩余數(shù)據(jù)點(diǎn)的平均值，P是窗口中所有數(shù)據(jù)點(diǎn)的平均值。余下的過程與點(diǎn)方式檢測粗差的

21、算法相同，也是M個值用來計算一個統(tǒng)計值，并使用該統(tǒng)計值生成閩值，如果某一差值超過了這個閩值，則認(rèn)為這個數(shù)據(jù)點(diǎn)含有粗差而需要將其進(jìn)行剔除。2.5 基于等高線拓?fù)潢P(guān)系的粗差檢測與剔除眾所周知，相鄰等高線的高程值之間的關(guān)系有且僅有三種:遞增、遞減或相等，根據(jù)這些關(guān)系，可對等高線的高程值是否有錯作一判斷。在等高線地形圖上由于存在等高線密集、注記的壓蓋、斷崖地形等情況，常常造成等高線的不連續(xù)有時丟失的情形，因此不能僅僅依靠等高距來確定可疑處是否錯誤。在對所有的可疑處自動檢測后，應(yīng)當(dāng)對每個可疑處根據(jù)等高線的關(guān)系由人工進(jìn)行校驗(yàn)并進(jìn)行修改，剔除粗差。以上所提出的幾種探測方法各有千秋,而我們應(yīng)在實(shí)踐中根據(jù)自己的

22、DEM的數(shù)據(jù)源、數(shù)據(jù)特點(diǎn)、現(xiàn)有的條件選擇合適的粗差檢測與剔除方法。3 基于趨勢面的粗差剔除方法研究 3.1 獲取原始DEM數(shù)據(jù)DEM的數(shù)據(jù)來源方式有地面測量（利用自動記錄的測距經(jīng)緯儀在野外實(shí)測）、現(xiàn)有地圖數(shù)字化（利用數(shù)字化儀對已有地圖上的信息如等高線、地形線等進(jìn)行數(shù)字化，目前常用的數(shù)字化儀有手扶跟蹤數(shù)字化儀與掃描數(shù)字化儀）、空間傳感器（利用GPS等進(jìn)行數(shù)據(jù)采集）、數(shù)字?jǐn)z影測量方法（這是DEM數(shù)據(jù)采集最常用的一種方法）。DEM數(shù)據(jù)采集的方式有以上幾種方法，在實(shí)際生產(chǎn)中，究竟要采用哪種方法和生產(chǎn)工藝，主要是取決于DEM的分辨率、精度要求、應(yīng)用范圍的大小以及成本和速度快慢等綜合考慮?，F(xiàn)把常用的DEM

23、數(shù)據(jù)采集方法以及各自特性的比較列于表3。表3 DEM采集方法及各自特性的比較獲取方式DEM的精度速度成本更新程度應(yīng)用范圍地面測量非常高耗時很高很困難小區(qū)域范圍立體遙感低比較快比較高周期性全球范圍內(nèi)的數(shù)據(jù)收集GPS比較高很快比較高容易小范圍地形圖手扶跟蹤數(shù)字化比較低比較耗時低周期性國家范圍內(nèi)以及軍事上的數(shù)據(jù)采集，中小比例尺地形圖的數(shù)據(jù)獲取地形圖屏幕數(shù)字化比較低非?？毂容^低周期性同上激光雷達(dá)、干涉雷達(dá)非常高很快非常高容易高分辨率、各種范圍攝影測量比較高比較快比較高周期性國家范圍內(nèi)的數(shù)據(jù)收集3.2 DEM 數(shù)據(jù)預(yù)處理DEM數(shù)據(jù)預(yù)處理是DEM 內(nèi)插之前的準(zhǔn)備工作，它是整個數(shù)據(jù)處理的一部分，一般包括數(shù)據(jù)

24、格式的轉(zhuǎn)換、坐標(biāo)系統(tǒng)的變換、數(shù)據(jù)的編輯、柵格數(shù)據(jù)的矢量化及數(shù)據(jù)分塊等內(nèi)容。3.2.1 格式轉(zhuǎn)換由于數(shù)據(jù)采集的硬、軟件系統(tǒng)各不相同，因而數(shù)據(jù)的格式可能也不相同。常用的代碼有ASCII碼、BCD碼及二進(jìn)制碼，每一次記錄的各項內(nèi)容及每項內(nèi)容的類型位數(shù)也可能各不相同，要根據(jù)相應(yīng)的DEM內(nèi)插軟件的要求，將各種數(shù)據(jù)轉(zhuǎn)換成該軟件所要求的數(shù)據(jù)格式11。3.2.2 坐標(biāo)變換若采集的數(shù)據(jù)一般要轉(zhuǎn)換到地面坐標(biāo)系。地面坐標(biāo)系一般采用國家坐標(biāo)系，也可采用局部坐標(biāo)系。3.2.3 數(shù)據(jù)編輯將采集的數(shù)據(jù)用圖形方式顯示在計算機(jī)屏幕上（或展繪在數(shù)控測圖儀上），作業(yè)人員采用圖形交互方式進(jìn)行數(shù)據(jù)編輯。包括剔除錯誤的、過密的與重復(fù)的點(diǎn)

25、，發(fā)現(xiàn)某些需要補(bǔ)測，對斷面掃面數(shù)據(jù)，還要進(jìn)行掃描的系統(tǒng)誤差的糾正。3.2.4 柵格數(shù)據(jù)轉(zhuǎn)換為矢量數(shù)據(jù)由地圖掃描數(shù)字化獲取的掃描影像是一灰度陣列，首先經(jīng)過二值化處理，再經(jīng)過濾波或形態(tài)處理（利用數(shù)學(xué)形態(tài)進(jìn)行各種運(yùn)算），并進(jìn)行邊緣跟蹤，獲得等高線上按順序排列的點(diǎn)坐標(biāo)，即矢量數(shù)據(jù)。3.2.5 數(shù)據(jù)分塊由于數(shù)據(jù)采集方式不同，數(shù)據(jù)的排列順序也不同，例如等高線是按各條等高線采集的先后順序排列的。但在內(nèi)插DEM時，待定點(diǎn)常常只與其周圍的數(shù)據(jù)點(diǎn)有關(guān)，為了能在大量的數(shù)據(jù)點(diǎn)中迅速地找的所需要的數(shù)據(jù)點(diǎn)，必須將其進(jìn)行分塊。在某些軟件中，需要將數(shù)據(jù)點(diǎn)劃分成計算單元，每個計算單元之間有一定的重疊度，以保證單元的連續(xù)性。分塊

26、的方法是先將整個區(qū)域分成等間隔得格網(wǎng)（通常比DEM格網(wǎng)大），然后將數(shù)據(jù)點(diǎn)按格網(wǎng)分成不同的類，通常有交換法和鏈指針法。3.2.6 子邊界的提取根據(jù)離散的數(shù)據(jù)點(diǎn)內(nèi)插規(guī)則格網(wǎng)DEM，通常是將地面看做是一個光滑的連續(xù)曲面，但是地面上存在著各種各樣的斷裂帶，如陡崖、絕壁以及各種人工地物，如路堤等，使地面并不光滑，這就需要將地面分成若干區(qū)域，即子區(qū)，是每一個子區(qū)的表面呈一連續(xù)光由特征線滑曲面。這些子區(qū)的邊界由特征線（如斷裂線）與區(qū)域的邊界線組成。確定每一子區(qū)的邊界可以采用專門的數(shù)據(jù)結(jié)構(gòu)或利用圖論等多種方法來解決。3.3 基于趨勢面的粗差檢測與剔除方法基于趨勢面的粗差探測基礎(chǔ)是地形所具有的自相關(guān)性，即地形變

27、化符合一定的自然趨勢，表現(xiàn)為連續(xù)空間的漸變模型，并且這種連續(xù)變化可以用一種平滑的數(shù)學(xué)曲面趨勢面來加以描述。對粗差的檢測，可以通過模型誤差即實(shí)際觀測值與趨勢面計算值(模型值)之差來判斷其是否屬于異常數(shù)據(jù)，因此趨勢面分析的一個典型應(yīng)用就是揭示研究區(qū)域中不同于總趨勢的最大偏離部分。因此，可以采用趨勢面分析找出偏離趨勢超過一定閩值的異常數(shù)據(jù)可疑點(diǎn)。應(yīng)用趨勢面分析進(jìn)行粗差探測需要注意兩點(diǎn)，一是趨勢面函數(shù)的確定，另一個是閾值問題，即當(dāng)觀測值和計算值相差多大時，該點(diǎn)可被懷疑為粗差點(diǎn)，理論上任何復(fù)雜的曲面都可用高階多項式去逼近，但用高階多項式的解本身不穩(wěn)定，同時多項式系數(shù)的物理意義也不清楚，可能會導(dǎo)致不符合實(shí)

28、際地形起伏的趨勢。關(guān)于閾值問題，通常按統(tǒng)計方法，假設(shè)誤差頻率呈正態(tài)分布，取三倍的中誤差為極限值，但畢竟設(shè)定閩值的人為因素很大，且現(xiàn)實(shí)世界中的變化也可導(dǎo)致不正確的判斷，造成粗差點(diǎn)的遺漏或錯選。趨勢面可有各種不同的形式，在現(xiàn)實(shí)生活中，我們通常選擇由下式構(gòu)成的最小二乘趨勢面:據(jù)處理區(qū)域的形味大小，可以靈活地選擇不同階次的多項式，對大而復(fù)雜的區(qū)域采用較高階次。根據(jù)統(tǒng)計規(guī)律，常用2倍或3倍中誤差作為極限誤差，即模型誤差大于極限誤差的觀測數(shù)據(jù)被認(rèn)為是粗差。趨勢面分析技術(shù)的特點(diǎn)是能將問題簡單化、局部化，能找出大部分的可疑數(shù)。3.4 實(shí)驗(yàn)分析 根據(jù)提供的 10 個數(shù)據(jù)點(diǎn)的坐標(biāo)和待求點(diǎn)的平面坐標(biāo)，利用移動二次曲

29、面擬合法，由格網(wǎng)點(diǎn)周圍的 10 個已知點(diǎn)內(nèi)插出待求格網(wǎng)點(diǎn)的高程，編制相應(yīng)的程序并進(jìn)行調(diào)試，最后解算出格網(wǎng)點(diǎn)P的高程并提交源程序代碼12。操作過程：1、首先要導(dǎo)入三個程序矩陣模塊stdafx.hSMatrix.h SingleImageResection.h2、輸入已知點(diǎn)數(shù)據(jù)的坐標(biāo)3、選擇二次曲面作為擬合曲面：列出誤差方程： 4、組方程，求出6個系數(shù)：5、計算出待求點(diǎn)高程。表4 已知數(shù)據(jù)點(diǎn)坐標(biāo)點(diǎn)號XYZ1102110152109113183105115194103103175108105216105108157115104208118108159116113171011311822編程計算點(diǎn)（1

30、10,110）上的高程。運(yùn)用C+進(jìn)行DEM內(nèi)插法的編程程序的代碼:#include stdafx.h#include SMatrix.h#include SingleImageResection.hint main(int argc, char* argv)SMatrix X(10,1); 系數(shù)矩陣SMatrix Y(10,1); 系數(shù)矩陣SMatrix Z(10,1); 系數(shù)矩陣X00 = 102; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /X10 = 109; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /X20 = 105; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /X30 = 103; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /X40 = 108; /

31、 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /X50 = 105; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /X60 = 115; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /X70 = 118; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /X80 = 116; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /X90 = 113; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Y00 = 110; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Y10 = 113; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Y20 = 115; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Y30 = 103; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Y40 = 105; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Y50 = 108; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Y60 = 104; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /

32、Y70 = 108; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Y80 = 113; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Y90 = 118; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Z00 = 15; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Z10 = 18; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Z20 = 19; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Z30 = 17; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Z40 = 21; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Z50 = 15; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Z60 = 20; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Z70 = 15; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Z80 = 17; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo) /Z90 = 22; / 輸入已知數(shù)據(jù)點(diǎn)坐標(biāo)

33、 /for(int i = 0 ;i 10;i+)Xi0 = Xi0 -110; Yi0 = Yi0 -110;SMatrix M(10,6);for(i = 0 ;i 10;i+)Mi0 =Xi0 *Xi0 ;Mi1 =Xi0 *Yi0 ;Mi2 =Yi0 *Yi0 ;Mi3 =Xi0;Mi4 =Yi0 ;Mi5 = 1 ;SMatrix P(10,10);for(i = 0 ;i 10;i+)Pii = 1/(Xi0 *Xi0+Yi0 *Yi0);SMatrix x(6,1);x = (M.T() * P * M).Invert() * M.T() * P * Z;printf(待定點(diǎn)的高

34、程是：%fn,x50);return 0;運(yùn)行結(jié)果截圖按上述方法計算模型誤差即實(shí)際觀測值與趨勢面計算值(模型值)之差來判斷其是否屬于異常數(shù)據(jù)。采用趨勢面分析找出偏離趨勢超過一定閾值的異常數(shù)據(jù)可疑點(diǎn)。粗差的修正應(yīng)該緊緊伴隨著粗差的探測過程，在DEM的柵格矩陣中，由于不允許空格點(diǎn)的存在，一旦某個粗差被檢測就應(yīng)該加以改正，這里最簡單的方法就是用臨近點(diǎn)的高程均值來代替可疑點(diǎn)的高程值，需要指出的是任何可疑數(shù)據(jù)點(diǎn)都不應(yīng)參與高程估值的計算。對于不同的數(shù)據(jù)窗口，分別利用最小二乘法進(jìn)行移動曲面擬合，得到各個數(shù)據(jù)窗口內(nèi)高程殘差值列向量，以此為基礎(chǔ)計算驗(yàn)后殘差值方差式，并以k值作為檢測粗差的閉值，由于同一個DEM高

35、程點(diǎn)在不同的擬合區(qū)域中計算出不同的殘差值，而這些殘差值對同一DEM高程點(diǎn)是否為粗差的判斷并不一致，即在一個擬合面內(nèi)認(rèn)為該點(diǎn)是粗差，而在另一個擬合面內(nèi)認(rèn)為是正常點(diǎn)，因此統(tǒng)計每個點(diǎn)在相關(guān)的擬合面內(nèi)被懷疑為粗差的百分比（percent），當(dāng)percent大于某一限值，則認(rèn)為點(diǎn)i含有粗差。中位數(shù)算法是以殘差的中值作為檢驗(yàn)統(tǒng)計量來檢測粗差，均值算法是以殘差的平均值作為檢驗(yàn)統(tǒng)計量來檢測粗差，雖然它們都是對同一問題的不同解決方式，但中位數(shù)算法比均值算法更為簡單。由于粗差未知，平差仍從慣常的最小二乘法開始，但在每次平差后，根據(jù)其殘差和有關(guān)其他參數(shù)，按所選擇的權(quán)函數(shù)，計算每個觀測值在下步迭代平差中的權(quán)。如果權(quán)函

36、數(shù)選擇得當(dāng)，且粗差可定位，則含粗差觀測值的權(quán)將愈來愈小，直到趨近于零，迭代中止時，相應(yīng)的殘差將直接指出粗差的值，而平差的結(jié)果將不受粗差的影響，從而實(shí)現(xiàn)粗差的定位剔除。4 總結(jié)與展望認(rèn)真總結(jié)所做的工作以及在論文撰寫中的遇到的困難，認(rèn)為在DEM的粗差探測中仍然有一些問題有待進(jìn)一步研究與解決:本文在數(shù)據(jù)處理方面較為不足，如果能夠?qū)B可視化編程語言和MATLAB工具相結(jié)合，充分利用MATLAB的運(yùn)算功能和VB開發(fā)界面方便的特點(diǎn)進(jìn)行混合編程，即用VB設(shè)計界面作為主程序，調(diào)用MATLAB子程序，數(shù)據(jù)處理將會簡單、方便。今天，DEM已作為一個獨(dú)立的產(chǎn)品而存在，并越來越廣泛地被用來代替?zhèn)鹘y(tǒng)地圖中等高線對地形

37、的描述，成為地理信息系統(tǒng)的核心數(shù)據(jù)庫以及地學(xué)分析的基礎(chǔ)數(shù)據(jù)。目前DEM粗差探測領(lǐng)域中，針對于規(guī)則格網(wǎng)DEM算法研究成果較為豐富，這與規(guī)則格網(wǎng)形式的DEM具有簡單的數(shù)據(jù)結(jié)構(gòu)和便利的存儲方式有很大的關(guān)系。相反，由于非規(guī)則形式的DEM數(shù)據(jù)結(jié)構(gòu)相對復(fù)雜，這方面的算法設(shè)計并不很多，較為實(shí)用的是基于點(diǎn)方式的算法。點(diǎn)方式算法是通過中心點(diǎn)的內(nèi)插值計算檢驗(yàn)量，然后進(jìn)行粗差的判釋。無論是規(guī)則格網(wǎng)DEM算法還是非規(guī)則格網(wǎng)DEM算法來剔除粗差，這些方法從理論上講都是將粗差歸入函數(shù)模型來實(shí)現(xiàn)粗差的探測。在現(xiàn)代測量數(shù)據(jù)處理中如何消除粗差的影響，特別是自動化觀測水平高，數(shù)據(jù)量大的情況下，就顯得越來越重要。參 考 文 獻(xiàn)1

38、柯正誼，何建邦，池天河.數(shù)字地面模型M.中國科學(xué)技術(shù)出版社，1993：103-215.2 李志林,朱慶.數(shù)字高程模型M.武漢:武漢大學(xué)出版社, 2001：77-83.3 陳再輝，路曉峰.基于自適應(yīng)抗差最小二乘的DEM數(shù)據(jù)粗差剔除J.海洋測繪，2006，26(6):15-17.4 ?？?，易思蓉，韓春華.基于坡度RMSE與3維可視化的格網(wǎng)DEM粗差檢測與剔除J.測繪通報，2007，6(9):23-25.5 韓玲.格網(wǎng)DEM中基于坡度信息的粗差檢測與剔除方法實(shí)驗(yàn)J.長安大學(xué)學(xué)報，2003，25(l):74-78.6 黃宏波，梁鑫，楊曉云，等.基于參數(shù)統(tǒng)計的DEM粗差探測算法J.測繪工程2008，8（

39、2）：12-15.7 王貴滿,王東華.DEM的粗差探測與生產(chǎn)質(zhì)量控制方法J.測繪與空間地理信息.2011，2（1）：1-7.8 楊曉云,顧利亞,岑敏儀.基于不同大小窗口的移動曲面擬合法探測不規(guī)則DEM粗差的一種方法J.測繪學(xué)報,2005，34 (2):14-15.9 楊元喜.自適應(yīng)抗差最小二乘估計J.測繪學(xué)報,1996,25(3): 206-211.10 黃幼才.數(shù)據(jù)探測與抗差估計M.北京:測繪出版社,1990：180-185.11 徐士良.數(shù)值方法與計算機(jī)實(shí)現(xiàn)M.北京:清華大學(xué)出版社,2006：33-56.12 傅祖蕓,趙梅娜,丁巖,等,譯. C語言數(shù)值算法程序大全M.北京:電子工業(yè)出版社,

40、1995：98-115.DEM Gross Error Picking Method ResearchAbstract: Digital elevation model (DEM) as one of 4 D products (DEM, DOM - numbers are projective like, DLG-digital line drawing, DRG - digital raster graph), is a kind of digital description and simulation of the surface of the earth. The earth is

41、an important part of spatial data infrastructure , it can be used to establish a terrain database and all kinds of GIS (geographic information system), the basis of quantitative analysis, etc. Digital elevation model is finding wider and wider application, it can be used to study analysis, the two-d

42、imensional geographical spatial continuous distribution and gradual change of all kinds of elevation modeling and analysis of attribute data. DEM data access usually are the ground measurement (using the automatic recording of ranging theodolite in field survey), existing digital map (using the digi

43、tizer on existing information such as contours on maps and terrain digital line, at present commonly used digitizer have manual tracking digital instrument and scan digitizer), space sensors (using GPS data acquisition etc.), digital photogrammetric method (which is one of the DEM data acquisition i

44、s the most commonly used method).In fact, regardless of the measuring method, measuring address data always contain all kinds of error; the error can be divided into three kinds, namely, systematic error, random error and gross error. Compared with systematic error, random error, gross error is, in

45、fact, a kind of error, the possibility that they appear in the measurement is compared commonly small. But compared with the former two kinds of error, gross error of digital elevation model which reflects the spatial variation of distortion is more serious, so gross error is one of the important factors affecting DEM quality. According to statistical law, commonly used in two or three times errors as the limiting error, the model error is greater than the limit e