CN110322557B - 一种结合gps与srtm融合的游动平差方法 - Google Patents

一种结合gps与srtm融合的游动平差方法 Download PDF

Info

Publication number
CN110322557B
CN110322557B CN201910528737.3A CN201910528737A CN110322557B CN 110322557 B CN110322557 B CN 110322557B CN 201910528737 A CN201910528737 A CN 201910528737A CN 110322557 B CN110322557 B CN 110322557B
Authority
CN
China
Prior art keywords
grid
srtm
elevation
value
gps
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201910528737.3A
Other languages
English (en)
Other versions
CN110322557A (zh
Inventor
周璀
张贵
朱建军
李卫
刘志卫
张智填
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Central South University of Forestry and Technology
Original Assignee
Central South University of Forestry and Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Central South University of Forestry and Technology filed Critical Central South University of Forestry and Technology
Priority to CN201910528737.3A priority Critical patent/CN110322557B/zh
Publication of CN110322557A publication Critical patent/CN110322557A/zh
Application granted granted Critical
Publication of CN110322557B publication Critical patent/CN110322557B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00Pattern recognition
    • G06F18/20Analysing
    • G06F18/25Fusion techniques
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T17/00Three dimensional [3D] modelling, e.g. data description of 3D objects
    • G06T17/05Geographic models
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10032Satellite or aerial image; Remote sensing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10032Satellite or aerial image; Remote sensing
    • G06T2207/10044Radar image
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/30Subject of image; Context of image processing
    • G06T2207/30181Earth observation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Software Systems (AREA)
  • Geometry (AREA)
  • General Physics & Mathematics (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Evolutionary Computation (AREA)
  • Evolutionary Biology (AREA)
  • General Engineering & Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Artificial Intelligence (AREA)
  • Remote Sensing (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Computer Graphics (AREA)
  • Position Fixing By Use Of Radio Waves (AREA)

Abstract

一种结合GPS与SRTM融合的游动平差方法,包括以下步骤:1)在第一个已知GPS点建立已知第一个GPS点与对应SRTM格网九个相邻格网的平差模型;2)求解平差模型计算对应SRTM格网和四个角点的改进高程值;3)把四个角点作为下一轮改进SRTM的已知值,求每一层平差融合后的新值与原格网值的差值,若差值的均值<阀值(0.5米),则算法收敛;4)在第一个GPS点以外的区域继续搜索第二个已知GPS点,然后进行第二轮改进,方法和第一个GPS点一致,直至整个区域融合完成。本发明可以提高融合区域SRTM的高程精度。

Description

一种结合GPS与SRTM融合的游动平差方法
技术领域
本发明涉及一种图像处理方法,尤其涉及一种结合GPS与SRTM融合的游动平差方法。
背景技术
数字高程模型(Digital Elevation Model,简称DEM)是地表起伏形态的数字化表达,也是进行各类地学分析的重要基础数据,已被广泛用于地震、火山、滑坡、地面沉降、洪水灾害及军事等各个领域。伴随着对地观测技术的不断发展,获取DEM数据的能力日益增强,各研究领域对高质量DEM的需求也随之增长。目前,DEM数据获取方法主要包括野外测量、摄影测量和遥感、制图数字化、合成孔径雷达和机载激光扫描。如今,使用卫星图像生成DEM 有一个显著的优势因为它相对便宜而且生成DEM所需的时间更少。然而使用光谱范围的缺点是它需要高分辨率,良好的光照条件和无云视线才能得到较好的DEM精度。近年来,干涉合成孔径雷达(InSAR)因其不依赖于自然光照而作为一个主动系统在提取高程数据方面变成流行技术。2000年美国航天飞机雷达地形测绘任务(Shuttle radar topographymission,SRTM) 利用机载Insar技术,仅11d便获取了全球80%陆地1″分辨率的三维地形信息。2014年9 月,1″分辨率的SRTM数据逐步向全球用户免费开放,因其精度在平坦地区稳定而成为最常用的DEM数据源。然而受限于雷达侧视成像模式,仍有空隙和异常,直接影响数据的应用潜力。多源DEM融合方法通过综合不同数据之间的互补信息,能获取更准确、全面、可靠的DEM,实现现有数据集的质量提升。在我国,每个省市都建有CORS站,每天收集全省各地大量点的 GPS高程(或北斗高程),充分利用这些GPS高程数据具有重要的实际和实用价值。
发明内容
本发明要解决的技术问题是克服现有技术的不足,提供一种结合GPS与SRTM融合的游动平差方法。通过对研究区域已知GPS点的高程值与SRTM对应格网的高程值进行多项式系数加权游动平差处理,提高融合区域SRTM的高程精度。最后根据SRTM坡度选择地理特征明显的典型研究区域(包括平坦区、半平坦区和陡峭区)作为研究对象进行实验分析,证明该算法的可行性和有效性。
为解决上述技术问题,本发明提出的技术方案为:一种结合GPS与SRTM融合的游动平差方法,包括以下步骤:
1)在第一个已知GPS点建立已知第一个GPS点与对应SRTM格网九个相邻格网的平差模型; 2)求解平差模型计算对应SRTM格网和四个角点的改进高程值;
3)把四个角点作为下一轮改进SRTM的已知值,求每一层平差融合后的新值与原格网值的差值,若差值的均值<阀值(0.5米),则算法收敛;
4)在第一GPS点以外的区域继续搜索第二个已知GPS点,然后进行第二轮改进,方法和第一个GPS点一致,直至整个区域融合完成;
(5)精度评定。
以某区域一个已知GPS观测点M举例说明,首先把DEM中每一个格网的高程用一个多项式函数表达,然后为M点对应的SRTM格网建立平差模型,通过M点的GPS观测高程与对应SRTM格网周围的九个相邻格网联合解算出M点对应的SRTM格网改正高程值和该格网的四个角点高程改正值,接着把四个角点看做已知点,继续求解周围的格网改正值;待第一个已知GPS点修正完毕,在该已知点半径为9个像素范围以外的区域继续搜索第二个已知GPS点,然后进行第二轮改进,方法和第一个GPS点一致,直至整个区域融合完成。因为从已知GPS 点开始利用平差方法不断融合SRTM DEM,并且是逐个格网逐个点进行的运算,因此称为游动平差方法。具体的平差模型可按如下方法建立:
假设DEM中格网(i,j)的格网高程值
Figure 100002_1
用一个f(x,y)函数表示。f(x,y)表示任意一点(x,y) 的高程,即:
Figure BDA0002099040930000021
其中。DEM的各格网值一般是通过各种测绘手段获取的
Figure 100002_2
的观测值hi,j。考虑观测误差,有:
hi,j+v=f(x,y) (2)
M点的高程可表示为h(xM,yM),考虑到观测误差的存在,M点高程的观测值可表示为:
hM+v=f(xM,yM) (3)
因为实际地形一般都非常复杂,故难以用数学函数准确描述,所以很难求得f(x,y)的精确表述。这里采用多项式函数代替f(x,y)。
f(x,y)=a1+a2x+a3y+a4xy+a5x2+a6y2 (4)
然后根据式(1)利用M点的GPS观测值与对应SRTM格网周围相邻9个格网点的高程建立观测方程,再与式(2)联合建立平差模型,接下来联合求解f(x,y),由此可求解出融合后格网(i,j) 和四个角点的高程值:
L+V=AX (5)
其中
L=(hi-1,j-1 hi-1,j...hi+1,j hi+1,j+1 hM)T
V=(Vi-1,j-1 Vi-1,j...Vi+1,j Vi+1,j+1 VM)T
X=(a1 a2 a3 a4 a5 a6)T
Figure BDA0002099040930000031
根据最小二乘平差,有
Figure BDA0002099040930000032
其中P是权值,初始权值为900,采用反距离降权的方式定权。根据平差融合后,改进格网小区域内各点的高程值可表示为
Figure BDA0002099040930000033
格网(i,j)平差融合后的值分别为:
Figure BDA0002099040930000034
格网(i,j)四个角点的融合后的高程分别为
Figure BDA0002099040930000035
Figure BDA0002099040930000036
Figure BDA0002099040930000037
Figure BDA0002099040930000038
格网(i,j)平差融合后的值
Figure BDA00020990409300000311
即为该格网融合后的格网值。格网(i,j)四个角点融合后的高程值
Figure BDA0002099040930000039
对于相邻的格网就是新的观测值,如对于格网(i-1,j),
Figure BDA00020990409300000310
是2个新的观测,新观测值的权可按上式计算。按与上述一样的方法,取格网(i-1,j) 及周边9个格网值与这2个新观测一起进行平差,可以求得该格网融合后的高程值及格网4 个角点融合后的高程值。
与现有技术相比,本发明的优点在于:通过对研究区域已知GPS点的高程值与SRTM对应格网的高程值进行多项式系数加权游动平差处理,提高融合区域SRTM的高程精度。
附图说明
图1为实施例1中平坦区的SRTM影像图。
图2为实施例1中半平坦区的SRTM影像图。
图3为实施例1中陡峭区的SRTM影像图。
图4为平坦区SRTM高程与GPS高程直线拟合图。
图5为平坦区算法改正高程与GPS高程直线拟合图。
图6为半平坦区SRTM高程与GPS高程直线拟合图。
图7为半平坦区算法改正高程与GPS高程直线拟合图。
图8为陡峭区SRTM高程与GPS高程直线拟合图。
图9为陡峭区算法改正高程与GPS高程直线拟合图。
具体实施方式
为了便于理解本发明,下文将结合较佳的实施例对本发明作更全面、细致地描述,但本发明的保护范围并不限于以下具体的实施例。
需要特别说明的是,当某一元件被描述为“固定于、固接于、连接于或连通于”另一元件上时,它可以是直接固定、固接、连接或连通在另一元件上,也可以是通过其他中间连接件间接固定、固接、连接或连通在另一元件上。
除非另有定义,下文中所使用的所有专业术语与本领域技术人员通常理解的含义相同。本文中所使用的专业术语只是为了描述具体实施例的目的,并不是旨在限制本发明的保护范围。
实施例
选取中国某地区为研究区域,从网站http://earthexplorer.usgs.gov/上面下载该地区的SRTM数据,然后从该地区截取三个地理特征不同的区域作为实验数据,它们的大小均为 128*128,三个区域的坡度分别在0-10度,10-20度,20度以上,如图1-图3所示分别对应平坦区、有坡度区、陡峭区,地理特征极具代表性。该地区已知高程点的高程来自于CORS站记录的该区域大量的GPS流动站点的高程。每个研究区域的GPS数据集分为两部分:实验GPS 数据和校验GPS数据,如表1所示。然后对每个研究区域实验GPS高程数据运用游动平差算法改进对应研究区域的SRTM影像,从而达到提高研究区域DEM精度的目的。
表1研究区域的GPS数据
Figure BDA0002099040930000051
为了客观地评估改进的DEM精度,考虑到系统偏差和标准差等精度估计会受到异常值和误差非正态分布的影响,因此这里考虑采用稳健统计方法,采用误差均值(MeanError,ME)、均方根误差(Root Mean Square Error,RMSE)、归一化绝对偏差中位数(Normalized Median, NMAD)和误差的标准差(Standard deviation of error,SDE)来进行精度评定。
Figure BDA0002099040930000052
Figure BDA0002099040930000053
NMAD=1.4826×medianj(|Δhj-mΔh|)
Figure BDA0002099040930000054
其中,N是GPS精度评定样本的个数;Hi,Href分别是提出算法改进的GPS高程和GPS的高程;Δhi,mΔh分别是第i个精度评定点的高程差和高程差的中值;这些值都是精度的度量,越小表明改进算法效果越好。
对三个研究区域的实验GPS数据运用游动平差算法改进对应研究区域的SRTM影像,然后采用ME、RMSE、NMAD、SDE精度评定因子对研究区域校验GPS数据、SRTM DEM与游动平差算法改正DEM评估,得到精度评定如表2所示。
表2三个研究区域精度评定
Figure BDA0002099040930000055
Figure BDA0002099040930000061
由表2可得,在平坦区、过渡区和陡峭区采用所提算法得到的DEM具有明显优势,与SRTM DEM比较,改进的RMSE提高了1-2米,精度分别提高了25.9%、22.4%和21.5%,表明所提出算法能够有效改善SRTM DEM,三个区域的精度提高程度相差不大,但在陡峭区域,因为受到地势的影响,融合的精度受到了影响。为了更好地展示和分析所提出算法的优势,这里对校验GPS数据的测量高程和游动平差融合之后的高程之间的测量关系进行详细分析,分别列出 GPS测量高程与相应SRTM的高程的直线拟合图,和GPS测量高程与游动平差算法改正高程的直线拟合图,如图4-图9所示。
图4、图6和图8是三个研究区域表明精度评定点的SRTM高程值与GPS观测高程之间的相关性;图5、图7和图9反映的是算法改正高程与GPS观测高程之间的相关性。这六个图显示的均是线性与正斜率的关系,表明两两变量之间朝着同一个方向发展。与图4、图6和图8相比,图5、图7和图9拥有更高的正斜率,而且散点的分布更加朝拟合直线聚拢,说明提出的游动平差算法改进DEM效果明显,拥有相对较高的置信限分析。
本实施例主要结合GPS实测高程数据,提出一种结合GPS与SRTM融合的游动平差方法,该方法通过对研究区域已知GPS点的高程值与SRTM对应格网的高程值进行多项式系数加权游动平差处理,已达到提高研究区域SRTM DEM精度的目的。该方法已通过三个典型的研究区域 (平坦区、半平坦区和陡峭区)实验验证,实验结果表明所提出游动平差算法能够有效改善三个实验区域的SRTM DEM,精度大约能提高21%以上。

Claims (3)

1.一种结合GPS与SRTM融合的游动平差方法,其特征在于:包括以下步骤:
1)在第一个已知GPS点建立已知第一个GPS点与对应SRTM格网九个相邻格网的平差模型;
2)求解平差模型计算对应SRTM格网和四个角点的改进高程值;
3)把四个角点作为下一轮改进SRTM的已知值,求每一层平差融合后的新值与原格网值的差值,若差值的均值<阀值(0.5米),则算法收敛;
4)在第一GPS点以外的区域继续搜索第二个已知GPS点,然后进行第二轮改进,方法和第一个GPS点一致,直至整个区域融合完成。
2.根据权利要求1所述的结合GPS与SRTM融合的游动平差方法,其特征在于:所述步骤1)中首先将待融合区域的DEM中每一格网的高程用一个多项式表达;然后为GPS观测点M对应的SRTM格网建立平差模型;
DEM中格网(i,j)的格网高程值
Figure 1
用一个f(x,y)函数表示;f(x,y)表示任意一点(x,y)的高程,即:
Figure FDA0002099040920000012
其中;DEM的各格网值一般是通过各种测绘手段获取的
Figure 2
的观测值hi,j;考虑观测误差,有:
hi,j+v=f(x,y) (2)
M点的高程可表示为h(xM,yM),考虑到观测误差的存在,M点高程的观测值可表示为:
hM+v=f(xM,yM) (3)
因为实际地形一般都非常复杂,故难以用数学函数准确描述,所以很难求得f(x,y)的精确表述;这里采用多项式函数代替f(x,y);
f(xy)=a1+a2x+a3y+a4xy+a5x2+a6y2 (4)
然后根据式(1)利用M点的GPS观测值与对应SRTM格网周围相邻9个格网点的高程建立观测方程,再与式(2)联合建立平差模型,
L+V=AX (5)
其中
L=(hi-1,j-1 hi-1,j...hi+1,j hi+1,j+1 hM)T
V=(Vi-1,j-1 Vi-1,j...Vi+1,j Vi+1,j+1 VM)T
X=(a1 a2 a3 a4 a5 a6)T
Figure FDA0002099040920000021
根据最小二乘平差,有
Figure FDA0002099040920000028
其中P是权值,初始权值为900,采用反距离降权的方式定权;根据平差融合后,改进格网小区域内各点的高程值可表示为
Figure FDA0002099040920000022
3.根据权利要求1所述的结合GPS与SRTM融合的游动平差方法,其特征在于:所述步骤2)中格网(i,j)平差融合后的值分别为:
Figure FDA0002099040920000023
格网(i,j)四个角点的融合后的高程分别为
Figure FDA0002099040920000024
Figure FDA0002099040920000029
Figure FDA00020990409200000210
Figure FDA00020990409200000211
格网(i,j)平差融合后的值
Figure FDA0002099040920000025
即为该格网融合后的格网值;
格网(i,j)四个角点融合后的高程值
Figure FDA0002099040920000026
对于相邻的格网就是新的观测值,如对于格网(i-1,j),
Figure FDA0002099040920000027
是2个新的观测,新观测值的权可按式(6)计算;按与上述一样的方法,取格网(i-1,j)及周边9个格网值与这2个新观测一起进行平差,可以求得该格网融合后的高程值及格网4个角点融合后的高程值。
CN201910528737.3A 2019-06-18 2019-06-18 一种结合gps与srtm融合的游动平差方法 Active CN110322557B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910528737.3A CN110322557B (zh) 2019-06-18 2019-06-18 一种结合gps与srtm融合的游动平差方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910528737.3A CN110322557B (zh) 2019-06-18 2019-06-18 一种结合gps与srtm融合的游动平差方法

Publications (2)

Publication Number Publication Date
CN110322557A CN110322557A (zh) 2019-10-11
CN110322557B true CN110322557B (zh) 2022-11-08

Family

ID=68119795

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910528737.3A Active CN110322557B (zh) 2019-06-18 2019-06-18 一种结合gps与srtm融合的游动平差方法

Country Status (1)

Country Link
CN (1) CN110322557B (zh)

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102213762B (zh) * 2011-04-12 2012-11-28 中交第二公路勘察设计研究院有限公司 基于rfm模型的多源星载sar影像自动匹配方法
WO2013121340A1 (en) * 2012-02-13 2013-08-22 Stellenbosch University Digital elevation model
CN108305237B (zh) * 2018-01-23 2021-09-21 中国科学院遥感与数字地球研究所 考虑不同光照成像条件的多立体影像融合制图方法

Also Published As

Publication number Publication date
CN110322557A (zh) 2019-10-11

Similar Documents

Publication Publication Date Title
Pieczonka et al. Generation and evaluation of multitemporal digital terrain models of the Mt. Everest area from different optical sensors
Tang et al. Verification of ZY-3 satellite imagery geometric accuracy without ground control points
Capaldo et al. DSM generation from high resolution imagery: applications with WorldView-1 and Geoeye-1.
Favalli et al. LIDAR strip adjustment: Application to volcanic areas
CN109100719B (zh) 基于星载sar影像与光学影像的地形图联合测图方法
Conte et al. Structure from Motion for aerial thermal imagery at city scale: Pre-processing, camera calibration, accuracy assessment
Martha et al. Effect of sun elevation angle on DSMs derived from Cartosat-1 data
CN109696152B (zh) 一种低相干性区域地面沉降量估算方法
Chrysoulakis et al. Validation of ASTER GDEM for the Area of Greece
Cenci et al. Describing the quality assessment workflow designed for DEM products distributed via the Copernicus Programme. Case study: The absolute vertical accuracy of the Copernicus DEM dataset in Spain
CN110310370B (zh) 一种gps与srtm点面融合的方法
Feng et al. Automatic selection of permanent scatterers-based GCPs for refinement and reflattening in InSAR DEM generation
Baltsavias et al. Geometric and radiometric investigations of Cartosat-1 data
Weifeng et al. Multi-source DEM accuracy evaluation based on ICESat-2 in Qinghai-Tibet Plateau, China
CN116299466B (zh) 一种输电通道地质形变监测方法和装置
Kirk et al. The effect of illumination on stereo dtm quality: simulations in support of europa exploration
CN110322557B (zh) 一种结合gps与srtm融合的游动平差方法
Spore et al. Collection, processing, and accuracy of mobile terrestrial lidar survey data in the coastal environment
Isioye et al. Comparison and validation of ASTER-GDEM and SRTM elevation models over parts of Kaduna State, Nigeria
Sefercik et al. Area-based quality control of airborne laser scanning 3D models for different land classes using terrestrial laser scanning: sample survey in Houston, USA
Yastikli et al. Quantitative assessment of remotely sensed global surface models using various land classes produced from Landsat data in Istanbul
CN115346128A (zh) 一种光学立体卫星dem高程改正和融合方法
Hnila et al. Quality Assessment of Digital Elevation Models in a Treeless High-Mountainous Landscape: A Case Study from Mount Aragats, Armenia
Papasaika et al. Fusion of LIDAR and photogrammetric generated Digital Elevation Models
Abileah et al. Bathymetry from fusion of multi-temporal Landsat and radar altimetery

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant