毕业论文——数字高程模型的粗差探测与剔除

1、数字高程模型的粗差探测与剔除摘要：数字高程模型作为地形表面的重要数字表达，其质量问题受到人们的普 遍关注。影响dem精度的因素多种多样，其中有失误引起的粗差，会造成dem 空间上的严重扭曲，甚至导致dem及其产品严重失真，因此探测并剔除粗差显 得尤为重要。传统的粗差处理都是基于平差原理的，但是，如果不存在平差问题，也就不 能在平差问题中对粗差进行定位。所以要检查dem数据中存在的错误，不能简 单用一般的平差方法，并且只分析单个独立的数据也是不够的，只有从整体或局 部对数据进行分析处理才能使问题得到解决。通过影像匹配建立起来的数字高程模型，粗差通常发生在自动影像相关时影 像的错误匹配。在dem的

2、生成上，和人工操作比较，通过匹配自动生成dem的 方法有很多优点，但是通过这种方法产生的dem可能包含各种来源的匹配错误 （特征丢失，纹理重复等），从而产生粗差。当立体图像的几何条件和辐射条件 不同时，这种错误会增加，这样就降低了 dem的精度，同时增加了后续人工编 辑的操作。立体影像的错误匹配常常产生粗差，甚至粗差簇，粗差的存在会造成数字高 程模型空间上的严重扭曲，因此需要预先剔除影像匹配获得的原始数据中存在的 粗差。该方法先利用dem对原始的立体影像对进行正射纠正，生成两幅正射影像，接着对其进行匹配，然后根据匹 配结果的视差对dem进行高程改正，重复进行该操作，直至匹配结果的视差小 于某阈

3、值。关键字:数字高程模型；粗差探测；影1匹配；视差；高程改正1绪论数字高程模型作为地理信息系统的重要数据來源，由于dem描 述的是地面高程信息，它在测绘、水文、气象、地貌、地质、土壤、 工程建设、通讯、军事等国民经济和国防建设以及人文和自然科学领 域有着广泛的应用。如在工程建设中，可用于如土方量的计算、通视 分析等；在防洪减灾方面，dem是进行水文分析如汇水区的分析、 水系网络分析、降雨分析、蓄洪计算、淹没分析等的基础；在无线通 讯上，可用于蜂窝电话的基站分析等。因此其质量问题成为人们日益 关注的焦点。一般来说，影响dem精度的主要原因是原始数据精度和高程精 度，而无论是数据的采集还是dem内

4、插的过程都有可能产生误差， 其中，粗差的影响最严重，会造成空间数据的严重扭曲。所以要设计 相关的算法，对粗差进行探测和剔除。1.2数字高程模型的粗差探测与剔除的研究现状dem数据有多种形式：规则格网，不规则格网，等高线，剖面图等， 它的原始数据可能以规则格网形式存在也可能以不规则形式存在。规 则格网形式数据有独特的特性，这些特性有助于粗差检测算法的实现。 所以说，适合于格网数据的粗差探测算法可能对不规则数据的粗差检 测毫无作用，故对不同类型的数据有必要设计不同的粗差探测方法。1.2.1基于趋势面按照口然地形地貌的成因，绝大多数口然地形表面都符合一定的 自然趋势，表现为连续的空间渐变模型，并且这

5、种连续变化可以用趋 势面来描述。对粗差的检测，可以通过模型误差即实际观测值与趋势 面计算值之差来判定其是否属于异常数据，由此可见，可以采用趋势 面分析找出偏离总趋势超过一定阈值的可疑数据。通过趋势面分析可 以找去大部分可疑数据，从而把问题局部化，简单化，但是趋势面分 析的一个缺点就是可以找出可疑数据，但不能确定数据是否为真止的 粗差。1.2.2基于坡度信息由于坡度是地面上点的一个基本属性，因此可以利用坡度的连续 性和一致性来检测格网数据中的粗差。如果在一个点的周围一定局部 区域内约束的允许坡度和允许坡度变化量大于给定的约束条件，就认 第2页（共16页）为该点可能存在粗差。1.2.3可视化采用d

6、em三维表面可视化技术，该方法可以交互式的来检查这 些可疑数据，并剔除严重影响数据质量的粗差或者错误。一般对于一 个特定的研究区域，在三维透视图上可疑点是否表现为粗差非常直观, 很容易作出判断。实际上，由于dem有着非常适宜于建立三维可视 化的特点，所以可以首先通过目视效果对粗差进行检测。通常有粗差 的地形是很不自然的。因此在实际应用中，可以首先通过目视对粗差 进行检测，然后在使用各种方法进行自动或半自动的粗差探测与剔除。2 dem的粗差探测方法的研究2.1 dem粗差探测方法研究的背景上述粗差剔除的算大多都是源于对数据本身的理解和研究，粗差在数 据屮可能孤立的分布，也可能成簇的存在。对于单个

7、粗差的存在，上 述算法的效果可能比较好,但对于以簇群形式存在的误差，效果可能 就不会很好。通过影像匹配口动生成的dem,作为一种口动化程度 较高的生产方式，在实际应用当中有着广泛的应用潜力。通过影像匹 配建立起來的数字高程模型，粗差经常发生在自动影像相关时影像错 误匹配的位置，在这种情况下，如果一个点含有很大的粗差，那么受 他的影响，在它周围有这个点参与其高程值内插的点都含有粗差，那 么这些点所建表面将会有很大的偏差，这就带來了以成簇的形式存在 的粗差。dem的编辑工作一般是由作业员在计算机或者测图仪上手工完 成，是一项劳动强度很大的工作。在摄影测量的领域，对自动匹配的 dem,急需一种误差口

8、动改止算法来减少人工的工作量。如何有效的 消除这种粗差形式，针对这种问题，在摄影测量领域展开了广泛的研 究。为了有效的完成dem的编辑工作，2989年，schenk21提出了利 用迭代的影像纠止来生成demo 1996年norvelll3】用这种方法对高程 数据进行改正。1996年lobonc从理论上证明了利用它进行高程改 正的可行性，并且针对航空影像的实验结果进行了说明。由于线阵推 扫式卫星影像的不同特性，尤其是对于异轨立体的卫星影像，由于儿 何条件和辐射条件都差异较大，因此匹配的错误的可能性较高。2. 2影像纠正法的主要思想本文针对如何有效的消除在匹配阶段造成的误差，引进了影像纠 正方法于

9、线阵推扫式粗差，从而减少后续编辑工作。利用中心投影的航空影像生成正射影像的过程，实际上是消除投 影差的过程，作为垂直投影的正射影像是没有投影误差的。但是对于 含有粗差的dem纠正生成的正射影像，不是完全意义上的正摄投影 或者垂直投影，因为投影差并没有完全被消除，x方向与y方向均有 残余投影差（中心投影的几何变形引起的两个方向的投影差），理论 上说这两个方向上的投影差与dem粗差相关。因此可以根据正射影 像的残余投影差来消除dem生成屮的粗差进而完成对dem的高程改 正。若通过匹配得到含有粗差的dem,那么根据dem和两幅影像的定 向参数可以生成两幅正射影像。然后对正射影像进行匹配，它的匹配 结

10、果所反映出来的视差就是残余投影差的一种表现，因此我们可以根 据视差这个线索来完成dem中粗差的剔除。2. 3影像纠正法的过程（1）立体影像、外部定向参数、dem数据的获取；（2）对两幅影像，利用同样的dem数据分别生成止射影像；（3）通过止射影像的立体匹配寻找同名像点；（4）根据同名像点计算视差，然后根据视差对dem数据进行改止；（5）重复2-4的操作，直至视差小于某给定的阈值2.4高程改正在影像纠正的过程中，最重要的就是利用投影差来改正高程。对 于正射影像的匹配结果，其视差分布与dem的粗差的关系可以由公 式（1）来表示，可以将视差转换为高程的改正值。对于有粗差的dem 纠正生成的正射影像，

11、由于其投影差并没有完全被消除，因此不是完 全意义上的正射投影或垂直投影，在x方向和y方向上均有投影差, 但是y方向的视差远小于x方向的视差，故只根据x方向的残余视差 来消除粗差进而完成对dem的高程的改正dh=dxh!h）（1）dx = x-x2这里的dh为高程改正值，dx为左右正摄影的残余视差，h/h为基高 比的倒数。3数字高程模型的粗差探测的过程根据上述所描述的算法，采用spot影像数据（覆盖范围为辽宁地区） 进行实验，表1为spot数据的具体参数。表1spot数据的具体参数左像右像彩像中心经度125.7268125.7396影像中心纬度43.854143.8541影像入射角14.6640

12、15.9988影响获取时间20010910200110173.1正射影像的生成在系统的主菜单中，选择“产品”一“生成正射影像”项，自动制作 当前模型的正射影像，屏幕显示计算提示界面。计算完后，自动生成 当前模型的正射影像。正射影像生成后，应显示其影像，检查正射影 像是否正确或完整。在系统的主菜单中，选择“显示”一“正射影像” 项，屏幕显示当前模型的正射影像。将光标移至影像屮。按鼠标 右键弹出菜单，供选择不同的比例，可对影像进行缩放。3. 2影像匹配在系统主菜单中，选择菜单“处理”-“影像匹配”项，出现影像 匹配的进程显示窗口，自动进行影像匹配。按照上述所说的理论，dem是否含有粗差可以由立体相

13、对生成的止 射影像来进行比较和检测，如果两幅止射影像差异较大就说明含有粗 差，通过对含有粗差的dem进行迭代改止，直至止射影像的差异小 于某一阈值，正射影像的差异利用匹配后的视差来描述。实验结果如下：图1 dem和左影像生成的正射影响图2 dem和又影像生成的正射影像图1是dem和左影像生成的正射影像，图2是dem和右影像生成 的止射影像。此时的dem是在helave数字摄影测量工作站上生成之 后经过人工编辑的产品。从止射影像的匹配结果看基本上消除了残余 视差。对上述dem数据人为引入的模拟粗差，实际上是将dem的整个表 面降低了 100米。图3含有粗差的dem和左影像生成的正射影像图4含有粗

14、差的dem和右蹈像生成的止射影像图3是含有粗非的dem和左影像生成的止射影像，图4是含有 粗差的dem和右影像生成的止射影像。两幅影像的匹配结果为：平 均视差为10.4个像素，平均视差的方差为0.86,两幅影像的基高比 为0.58,按照公式（1）计算的高程的改止量为107.6米，从中可以 看岀改止误差7.6米小于一个像素的视差所引起的高程改正量，在模 拟数据的实验过程中，由于各点的粗差是相同的，因此视差的分布是均匀的。图5冃接生成叮视化结果的图像图6经过迭代方法产生的可视化结果由于经过人工编辑的helave数字摄影测量工作站dtm产品已经不含 粗差，在不考虑dem精度情况下，基本上代表了地表的

15、起伏。先将 不同迭代次数的dtm产品与其进行比较，比较结果如表2所示：表2不同迭代次数的dtm产品比较结果迭代次数高程误差的最人值（m）高程的中误差（m）054.3737.15136.4523.23221.8317.67320.5617.24上述两幅图像是采用实际数据的可视化结果，图5是直接生成可 视化结果的图像，图6是经过采用迭代方法产生的dem可视化结果， 可以看出图5屮的明显粗差已经在图6屮得到了较好的消除。4数字高程模型的粗差剔除过程利用下面的程序剔除粗差：option explicitpublic sub matcherrorpro()'调用此函数剔除粗差dim rs as

16、new adodb. recordsetdim i as integer, x as single, constr as stringconstr 二 "provider二microsoft jet oledb 4. 0;data source二& datapath &demmdbdim t as singledim hh as doublet 二 timerscreen.mousepointer 二 11hh = mymatchmethod. interval - mymatchmethod. destarea frmwait showfrmwaitprogress

17、barl.visible 二 falsedim lx() as double, ly() as double, rx() as double, ry()as double, n as long, ptid() as long, delid() as longn 二 0rs.open "select * from tbldata where modelno二 &dblpicname &'， constr, adopenstatic, adlockoptimisticdo while not rs. eof()n 二 n + 1redim preserve lx(

18、1 to n), ly(1 to n), rx(1 to n),ry(1 to n), ptid(l to n)lx(n)二 rs. fields(lx)ly (n)二 rs. fields(ly)rx(n)二 rs. fields ("rx)ry (n)二 rs. fields (z，ryz，)ptid(n)二 rs. fields(id)rs.movenextlooprs. closeif n < 10 then exit subdim isused () as boolean, v () as long, vmin as double,index as integer,

19、m as long, k as integer, d as double, j aslongm 二 10dim cc as longfor k 二 1 to nredim isused(l to n), v(1 to m)for i 二 1 to mvmin = 9999999for j = 1 to nif j <> k and abs(ly(j) 一 ly (k) >= i iiithend 二(lx(j) - lx (k) * (lx(j)-lx (k)+ (ly(j) - ly (k) * (ly(j) - ly (k)if isused(j) = false and

20、 d < vminthenvmin 二 dindex 二 jend ifend ifnext jv(i)二 indexisused (index)二 truenext iif checkerr (lx (k) - rx (k), ly(k) - ry (k), lx(),ly(), rx(), ry(), v(), m) = false thenif checkerr3 (lx(k) - rx(k), ly(k) - ry(k), lx(),ly(), rx(), ry(), v() thencc 二 cc + 1redim preserve delid(1 to cc)del id (

21、cc) = ptid(k)end ifnext rs.open "select * from tbldata where modelno二 &dblpicname &' “，constr, adopenstatic, adlockoptimisticdo while not rs. eof()for i = 1 to ccif rs.fields(id)二 delid(i) thenrs. deleteexit forend ifnextrs. movenextlooprs. closeunload frmwaitscreen.mousepointer 二 0

22、msgbox，z共剔除粗差点数:cc, vblnformation,历时:& timer - t &秒end subprivate function checkerr(dx as double, dy as double, lx() as double, ly() as double, rx() as double, ry() as double, index () as long, m as long) as booleandim i as integer, id as integerdim avgdx as double, avgdy as doublefor i 二 1

23、to mid = index(i)avgdx = avgdx + lx (id) - rx(id)avgdy 二 avgdy + ly(id) ry (id)nextavgdx 二 avgdx / mavgdy 二 avgdy / mdim ddx as double, ddy as doublefor i = 1 to mid = index (i)ddx 二 ddx + (lx(id) - rx(id) - avgdx) * (lx (id) - rx(id) -avgdx)ddy 二 ddy + (ly(id) - ry (id) - avgdy) * (ly(id) - ry (id)

24、 -avgdy)nextddx 二 sqr (ddx / m)ddy 二 sqr (ddy / m)if ddx 二 0 then ddx 二 1if ddy 二 0 then ddy 二 1if abs(dx - avgdx) > 3 * ddx or abs(dy - avgdy)3 * ddy thencheckerr 二 trueend ifend functionprivate function checkerr2(mdx as double, mdy as double, lx() as double, ly () as double, rx() as double, ry

25、() as double, index () as long) as booleandim i as integer, id as integerdim dx as double, dy as doubledim ddx as double, ddy as doublefor i 二 1 to 4dx 二 dx + lx(index(i) 一 rx(index(i)dy 二 dy + ly (index(i) 一 ry (index(i)next dx 二 dx / 4dy 二 dy / 4if dx = 0 then dx = 1if dy 二 0 then dy 二 1if abs (md

26、x 一 dx) > abs (dx * 0. 2) or abs (mdy 一 dy) > abs (dy* 0. 2) thencheckerr2 二 trueend ifend functionprivate function checkerr3(mdx as double, mdy as double, lx() as double, ly () as double, rx() as double, ry () as double, index() as long) as booleandim i as integer, id as integerdim dx(l to 10

27、) as double, dy(1 to 10) as doubledim ddx as double, ddy as double, avgdx as double, avgdy as doublefor i = 1 to 10dx (i) = lx (index (i) 一 rx(index(i) dy(i) = ly(index(i) 一 ry(index(i)nextavgdx 二 0: avgdy 二 0for i = 1 to 10avgdx = avgdx + dx(i)avgdy = avgdy + dy(i)next iavgdx = avgdx / 10avgdy 二 av

28、gdy / 10dx (5) + dx (6) +dy(5) + dy (6) +' avgdx 二(dx(l) + dx(2) + dx(3) + dx(4) +dx(7) + dx (8) / 8' avgdy 二(dy(l) + dy(2) + dy(3) + dy(4) + dy(7) + dy(8) / 8for i = 1 to 10ddx 二 ddx + (dx(i) - avgdx) * (dx(i) - avgdx)ddy = ddy + (dy(i) - avgdy) * (dy(i) - avgdy)nextddx 二 sqr (ddx / 10)ddy

29、二 sqr (ddy / 10)if ddx 二 0 then ddx 二 1if ddy 二 0 then ddy 二 1if abs (mdx - avgdx) > 3 * ddx or abs (mdy - avgdy) > 3 * ddythencheckerr3 二 trueend ifend function5分析与结论(1)传统的粗差剔除方法都是基于数据本身的研究，而本文是从问题 的结果出发，采用反馈迭代的方法，利用当前的匹配结果來指导新一 轮的匹配过程。 粗差是一种错误，采用该方法可以有效地消除粗差。但是dem的质量除了与影像有关外，在很大程度上与算法有关，在整个迭

30、代和匹 配的过程中采用统一匹配方法，容易造成误差的积累，故不利于提高dem的精度。 从实验的结果来看，对于线阵推扫式卫星影像的匹配结果可以有 效的剔除粗差。该方法不能适用于任何地区，对于大面积纹理贫乏区域就不能有 效的使用。针对于该方法存在的不足之处，仍需进一步地努力探索，积极地 弥补。参考文献1 李志林，朱庆，数字高程模型m,武汉测绘科技大学出版社，20002 schenk.a.f.1989.determination of dem using iteratively rectifiedimage.photogrammetry teehnical report no.3,dept.of ge

31、odetic第14 jx （共16页）science and surveying ,the ohio state university,columbus,ohio3 lob on e.t. 1996.human.supervised tools for digital photogrammetricsystems,ph.d.thesis,purdue university4 novelle,f.r. 1996.using interative orthophoto refinement togen erate and correct digital elevati on models(dems

32、),asprs,digital potogrammetry.pp. 151-1555 王佩军，徐亚明，摄影测量学m,武汉大学出版社巩丹超，黄小波，邓雪清，数字高程模型的粗差探测与剔除信息工程大 学测绘学院abstract： the digital elevation model as an important digeital terrain surface expression, its quality issues more attention by people. various factors affected the accuracy of dem, including error

33、caused by the gross error, will cause serious distort!on in dem space, and even lead to serious distortion of dem and its products, so it is particularly important to detect and eliminate gross error.traditional gross error handling is based on the principle of the adjustment, but, if there is no ad

34、justment problem, also cannot adjustment problems in fixing the gross error. so check the 第15页(共16页)dem data error, can't simply use general adjustment method, and analysis of a single in depe ndent data only is not eno ugh, only as a whole or partial analysis of the data processing to make the problem be solved.by image matchi ng to build a digital elevati on model, gross error usually occurs whe n the automatic image correlation image matching error. on the dem generation, compared with manual op