CN113124834A - 一种结合多源数据的区域网平差方法、系统及存储介质 - Google Patents

一种结合多源数据的区域网平差方法、系统及存储介质 Download PDF

Info

Publication number
CN113124834A
CN113124834A CN202110668978.5A CN202110668978A CN113124834A CN 113124834 A CN113124834 A CN 113124834A CN 202110668978 A CN202110668978 A CN 202110668978A CN 113124834 A CN113124834 A CN 113124834A
Authority
CN
China
Prior art keywords
image
point
water surface
coordinate
elevation
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.)
Granted
Application number
CN202110668978.5A
Other languages
English (en)
Other versions
CN113124834B (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.)
China University of Geosciences
Original Assignee
China University of Geosciences
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 China University of Geosciences filed Critical China University of Geosciences
Priority to CN202110668978.5A priority Critical patent/CN113124834B/zh
Publication of CN113124834A publication Critical patent/CN113124834A/zh
Application granted granted Critical
Publication of CN113124834B publication Critical patent/CN113124834B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C11/00—Photogrammetry or videogrammetry, e.g. stereogrammetry; Photographic surveying
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C5/00—Measuring height; Measuring distances transverse to line of sight; Levelling between separated points; Surveyors' levels
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/48—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S17/00
    • G01S7/4802—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S17/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Multimedia (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Image Processing (AREA)

Abstract

本发明提供了一种结合多源数据的区域网平差方法、系统及存储介质,涉及测绘技术领域，该方法通过对目标区域的激光雷达点云数据进行滤波处理以获取目标区域高精度的水面高程数据；同时根据核线约束条件获取目标区域的立体像对的水边线上的同名像点的坐标，同名像点的坐标结合水面高程数据可获取控制点三维坐标，在获取控制点三维坐标后，可根据控制点三维坐标和RPC模型以两类未知数交替趋近法解算预设的像方空间补偿模型，直至满足精度要求，从而完成对目标区域影像的区域网平差。

Description

一种结合多源数据的区域网平差方法、系统及存储介质
技术领域
本发明涉及测绘技术领域，具体而言，涉及一种结合多源数据的区域网平差方法、系统及存储介质。
背景技术
遥感立体影像高精度几何定位需要可靠的控制信息，常规方法一般借助均匀分布的控制点进行联合平差计算，消除影像成像过程中的各种系统误差，从而达到相应比例尺地图测绘的精度要求。
对于难以进行实地测量的区域，例如境外的海岸带和海岛礁区域，难以采用人员实地测量的方法进行测量，而海岸带和岛礁浅海区域已有的数据如数字地面模型等，其精度往往不高且又因为各种原因难以获取。因此，难以进行实地测量的区域难以获取控制点，也无法获取控制点的相关数据，进而难以进行区域网平差。
发明内容
本发明解决的问题是难以进行实地测量的区域难以选取控制点，也无法获取控制点的相关数据。
为解决上述问题，本发明第一方面提供一种结合多源数据的区域网平差方法，包括：
对目标区域的激光雷达点云数据进行滤波处理，得到所述目标区域的水面高程数据；
根据核线约束条件获取所述目标区域的立体像对的水边线上的同名像点的坐标；
选取所述水边线上的同名像点对应的地面点作为控制点，并根据所述水面高程数据和所述同名像点的坐标获取控制点三维坐标；
根据所述控制点三维坐标和RPC模型以两类未知数交替趋近法解算预设的像方空间补偿模型，其中，两类未知数分别为控制点三维坐标的改正数和影像仿射变换参数的改正数。
进一步地，所述对目标区域的激光雷达点云数据进行滤波处理包括：
以预设的高程间隔将所述激光雷达点云数据在高程方向切分为K个切片区间；
获取各所述切片区间内的光子信号数量和各所述切片区间的中心高程；
根据各所述切片区间的光子信号数量和各所述切片区间的中心高程进行正态分布拟合，得到拟合函数；
对所述拟合函数进行最优解算，获取所述拟合函数的参数；
根据所述拟合函数的参数构建水面光子信号高程区间，并根据所述水面光子信号高程区间获取第一水面光子信号数据集。
进一步地，所述对目标区域的激光雷达点云数据进行滤波处理还包括：
对所述第一水面光子信号数据集进行曲线拟合，得到所述第一水面光子信号数据集的曲线方程；
对所述曲线方程求一阶导，获取导数为零的点的高程值；
将对应所述第一水面光子信号数据集的曲线的波谷的点的高程值作为分界值，根据所述分界值对所述第一水面光子信号数据集进行进一步滤波处理，得到第二水面光子信号数据集。
进一步地，所述对目标区域的激光雷达点云数据进行滤波处理还包括：
对所述第二水面光子信号数据集进行正态分布拟合，将得到的所述第二水面光子信号数据集的正态分布曲线的预设范围作为可靠性区域，并将所述可靠性区域的区域宽度作为圆形滤波窗口直径长度；
沿着轨道方向，依次以所述第二水面光子数据集中的每个光子信号作为中心，以所述圆形滤波窗口直径长度的一半为半径构建圆形滤波窗口，统计每个所述圆形滤波窗口内的光子信号数量，并计算光子信号数量均值；
获取内部光子信号数量大于所述光子信号数量均值的所述圆形滤波窗口，根据获取结果构建第三水面光子信号数据集。
进一步地，所述对目标区域的激光雷达点云数据进行滤波处理还包括：
对所述第三水面光子信号数据集中的光子信号进行潮汐改正，得到影像获取时刻的水面高程值。
进一步地，所述根据核线约束条件获取所述目标区域的立体像对的水边线上的同名像点像点坐标包括：
获取所述目标区域的立体像对，并根据核线约束条件，获取所述立体像对的核线影像；
根据归一化水体指数法分别将所述立体像对的核线影像分成水体区域和非水体区域；
根据所述水体区域和所述非水体区域分别提取所述立体像对的核线影像上的水边线区域；
选取所述立体像对的核线影像上的水边线区域上的相应点作为同名像点，获取所述同名像点的核线影像坐标；
根据所述立体像对的原始影像与核线影像之间的坐标映射关系获取所述同名像点的原始影像坐标。
进一步地，根据所述控制点三维坐标和两类未知数交替趋近法解算所述像点坐标误差方程包括：
将所述像方空间补偿模型转化为像点坐标误差方程，将控制点的三维坐标作为已知，以将所述像点坐标误差方程转化为影像仿射变换参数的误差方程，根据所述控制点三维坐标和所述RPC模型求解所述影像仿射变换参数的误差方程，得到影像仿射变换参数的改正数；
根据所述影像仿射变换参数的改正数对影像仿射变换参数进行改正；
根据改正后的影像仿射变换参数对所述控制点三维坐标进行改正，将所述影像仿射变换参数的改正数代入所述像点坐标误差方程得到控制点三维坐标的改正数的误差方程；
根据改正后的控制点三维坐标和所述RPC模型解算所述控制点三维坐标的改正数的误差方程，得到控制点三维坐标的改正数；
当所述影像仿射变换参数的改正数和所述控制点三维坐标的改正数均大于对应阈值时，则根据所述影像仿射变换参数的改正数和所述控制点三维坐标的改正数求解新的影像仿射变换参数的改正数和新的控制点三维坐标的改正数。
进一步地，所述根据改正后的控制点三维坐标和所述RPC模型解算所述控制点三维坐标的改正数的误差方程包括：
根据所述控制点的平面坐标精度和高程坐标精度确定所述控制点的平面坐标权重和高程坐标权重；
根据所述平面坐标权重和所述高程坐标权重以参数加权平差法解算所述控制点三维坐标的改正数的误差方程，得到所述控制点三维坐标的改正数。
本发明第二方面提供一种结合多源数据的区域网平差系统，包括：
滤波模块，用于对目标区域的激光雷达点云数据进行滤波处理，得到所述目标区域的水面高程数据；
第一获取模块，用于根据核线约束条件获取所述目标区域的立体像对的水边线上的同名像点的坐标；
第二获取模块，用于选取所述水边线上的同名像点对应的地面点作为控制点，并根据所述水面高程数据和所述同名像点的坐标获取控制点三维坐标；
计算模块，用于根据所述控制点三维坐标和RPC模型以两类未知数交替趋近法解算预设的像方空间补偿模型，其中，两类未知数为控制点三维坐标的改正数和影像仿射变换参数的改正数。
本发明第三方面提供一种存储介质，所述存储介质存储有计算机程序，所述计算机程序被处理器执行时，实现如上所述的结合多源数据的区域网平差方法。
本发明的有益效果：激光雷达点云的每一个点都包含了三维坐标信息，对目标区域的激光雷达点云数据进行滤波处理可以获取目标区域高精度的水面高程数据；同时根据核线约束条件可获取目标区域的立体像对的水边线上的同名像点坐标，同名像点坐标结合水面高程数据可获取控制点三维坐标，在获取控制点三维坐标后，可根据控制点三维坐标和RPC模型以两类未知数交替趋近法解算预设的像方空间补偿模型，直至满足精度要求，从而完成对目标区域影像的区域网平差。
附图说明
图1为本发明实施例的结合多源数据的区域网平差方法的流程图；
图2为本发明实施例的将激光雷达点云数据切分为K个切片区间的示意图；
图3为本发明实施例的第一水面光子信号数据集的拟合曲线示意图；
图4为本发明实施例的圆形滤波窗口滤波处理示意图；
图5为本发明实施例的核线影像提取同名点的示意图；
图6为本发明实施例的核线影像水边线上的同名像点获取示意图；
图7为本发明实施例的结合多源数据的区域网平差的结构示意图。
具体实施方式
为使本发明的上述目的、特征和优点能够更为明显易懂，下面结合附图对本发明的具体实施例做详细的说明。
本发明的说明书和权利要求书及上述附图中的术语“包括”以及它们任何变形，意图在于覆盖不排他的包含。例如包含一系列步骤或单元的过程、方法或系统、产品或设备没有限定于已列出的步骤或单元，而是可选地还包括没有列出的步骤或单元，或可选地还包括对于这些过程、方法、产品或设备固有的其它步骤或单元。
本说明书描述的“第一”、“第二”和“第三”等术语，仅用于区分装置/组件/子组件/部件等，不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量，由此，限定有如“第一”、“第二”和“第三”等的特征可以明示或者隐含地表示包括至少一个该特征，除非另有明确具体的限定，“多个”的含义是至少两个，例如两个，三个等，对于本领域技术人员而言，可以具体情况理解上述术语在本发明中的具体含义。
如图1所示，本发明实施例的一种结合多源数据的区域网平差方法包括：
S101：对目标区域的激光雷达点云数据进行滤波处理，得到所述目标区域的水面高程数据。
其中，通过单光子激光雷达卫星获取激光雷达点云数据。
可选地，所述对目标区域的激光雷达点云数据进行滤波处理包括：
以预设的高程间隔将所述激光雷达点云数据在高程方向切分为K个切片区间。
如图2所示，根据目标范围和光子数据的高程范围的实际情况，以一定的高程间隔Δh，沿高程方向从上至下将光子数据划分成K个切片区间，例如，假设光子高度的上限为U，则K个切片区间可依次为[U，U-Δh]、（U，U-2Δh]...（U，U-KΔh]。
获取各所述切片区间内的光子信号数量和各所述切片区间的中心高程。
根据各所述切片区间的光子信号数量和各所述切片区间的中心高程进行正态分布拟合，得到拟合函数。如式（1）：
Figure DEST_PATH_IMAGE002
（1）
n maxnum 表示包含光子信号数量最多的高程切片中的光子数量，h maxnum 表示光子信号数量最多的切片区间中的中心高程，n表示光子信号高程，d（n）表示光子信号数量。
对所述拟合函数进行最优解算，获取所述拟合函数的参数。
以切片区间内光子信号数量为坐标系的纵轴，以切片区间的中心高程为坐标系的横轴可得到K个点坐标，通过最小二乘法进行最优解算，可获取式（1）中的相关参数，即h maxnum 和σ。
根据所述拟合函数的参数构建水面光子信号高程区间，并根据所述水面光子信号高程区间获取第一水面光子信号数据集。
由水面高程预估区间的光子分布并结合实际可知，出现某些水下地形高程因接近水面高程，造成水下地形光子信号错分到水面光子信号的情况，因此，需要再进一步，进行精细滤波处理。在进行水面高程分界时，以该区间切片的中心高程h maxnum 为中心，取σ作为最小高程方向的阈值，取2σ作为最大高程方向的阈值，分别得到水面光子信号高程的下限H wmin 和上限H wmax ，如式（2）所示。
Figure DEST_PATH_IMAGE004
（2）
根据式（2）中的水面光子高程的上下限，获取水面光子信号，构建第一水面光子信号数据集SW和水下光子信号数据集UW。
Figure DEST_PATH_IMAGE006
（3）
可选地，所述对目标区域的激光雷达点云数据进行滤波处理还包括：
对所述第一水面光子信号数据集进行曲线拟合，得到所述第一水面光子信号数据集的曲线方程。
第一水面光子信号数据集中还存在一定数量非常靠近水面的水下地形光子信号和少量水下噪声光子信号（在后面，统一将水下地形光子信号和水下噪声光子信号称为水下光子信号），为获取高精度的水面光子信号，需对其进行更进一步的精细滤波处理。本实施例通过基于圆形搜索窗口的空间聚类密度滤波模型实现更进一步的精细滤波处理，具体如下：
首先，利用正态分布函数对第一水面光子信号数据集SW进行曲线拟合，由于存在一定数量的水下光子信号，因此在接近水面高度处，很可能会出现两个极为接近的波峰，分别是水面光子信号和水下光子信号，如图3所示。
对所述曲线方程求一阶导，获取导数为零的点的高程值；
将对应所述第一水面光子信号数据集的曲线的波谷的点的高程值作为分界值，根据所述分界值对所述第一水面光子信号数据集进行进一步滤波处理，得到第二水面光子信号数据集。
可对第一水面光子信号数据集SW的正态分布曲线求取一阶导，此时，波峰和波谷的导数均为零，但如图3所示，波谷上侧的曲线为递减曲线，波谷下侧的曲线为递增曲线，以及波谷位于两波峰之间，可确定波谷G点以及获取波谷G点所在的高程值H G ，并将H G 作为水面光子信号和水下光子信号所对应的两个波峰的分界线，则A区域表示水下光子信号，T区域表示水面光子信号，此时，水面光子信号数据集表示为SW B ，即第二水面光子信号数据集。
可选地，所述对目标区域的激光雷达点云数据进行滤波处理还包括：
对所述第二水面光子信号数据集进行正态分布拟合，将得到的所述第二水面光子信号数据集的正态分布曲线的预设范围作为可靠性区域，并将所述可靠性区域的区域宽度作为圆形滤波窗口直径长度。
具体地，对第二水面光子信号数据集SW B 进行正态分布拟合，优选将正态分布曲线的98%作为可靠性区域，并获取该区域范围的区域宽度作为圆形窗口直径长度D，其圆形窗口如式（4）所示，R为圆形窗口半径。
Figure DEST_PATH_IMAGE008
（4）
其中，S circle 表示圆形滤波窗口的面积。
沿着轨道方向，依次以所述第二水面光子信号数据集中的每个光子信号作为中心，以所述圆形滤波窗口直径长度的一半为半径构建圆形滤波窗口，统计每个所述圆形滤波窗口内的光子信号数量，并计算光子信号数量均值；
获取内部光子信号数量大于所述光子信号数量均值的所述圆形滤波窗口，根据获取结果构建第三水面光子信号数据集。
如图4所示，在第二水面光子信号数据集SW B 中，沿着单光子激光雷达轨道方向，依次以每个光子信号作为圆形滤波窗口中心，以R为半径构建圆形滤波窗口，统计每个滤波区域内的光子信号点数，记为P num ，
Figure DEST_PATH_IMAGE010
，其中，
Figure DEST_PATH_IMAGE012
为第二水面光子数据集SW B 中的总光子数，并计算点数光子信号均值
Figure DEST_PATH_IMAGE014
。
Figure DEST_PATH_IMAGE016
（5）
比较每个搜索区域内的光子数量与光子点数均值
Figure 350094DEST_PATH_IMAGE014
的大小。当搜索区域内的光子信号数量大于
Figure 329551DEST_PATH_IMAGE014
时，则认为搜索区域中心的光子为水面光子信号，否则为噪声光子信号。通过该方法，对SW B 中的每个光子信号进行滤波处理，提取水面光子信号，实现高精度水面光子信号精滤波，提取的水面光子信号数据集表示为
Figure DEST_PATH_IMAGE018
，即第三水面光子信号数据集。
可选地，所述对目标区域的激光雷达点云数据进行滤波处理还包括：
对所述第三水面光子信号数据集中的光子信号进行潮汐改正，得到影像获取时刻的水面高程值。
由于地球表面各处的海水受月球和太阳的吸引力作用进行周期性升降运动，而通常将开阔平静的海洋平面当作大地水准面。因此，为了得到影像获取时刻的水面高程值，需要对自然海面和作为大地水准面的平均海面之间海水层的吸引力与这段空间距离的重力正常变化要加以改正，称该改正过程为“潮汐改正”。
S102：根据核线约束条件获取所述目标区域的立体像对的水边线上的同名像点的坐标。
可选地，所述根据核线约束条件获取所述目标区域的立体像对的水边线上的同名像点像点坐标包括：
获取所述目标区域的立体像对，并根据核线约束条件，获取所述立体像对的核线影像；
根据归一化水体指数法分别将所述立体像对的核线影像分成水体区域和非水体区域；
根据所述水体区域和所述非水体区域分别提取所述立体像对的核线影像上的水边线区域；
选取所述立体像对的核线影像上的水边线区域上的相应点作为同名像点，获取所述同名像点的核线影像坐标；
根据所述立体像对的原始影像与核线影像之间的坐标映射关系获取所述同名像点的原始影像坐标。
首先，获取目标区域的立体像对，根据核线约束条件，获取立体像对的核线影像。立体像对的原始影像与核线影像之间存在坐标映射关系，如式（6）所示（由于不同影像的成像方式不同，因此核线影像与原始影像之间的坐标映射关系并不唯一，因此式（6）只用于表示两种影像之间存在函数关系）。同时，利用在核线对上的同名像点的左右影像纵坐标相等这一关系，可获取同名像点核线影像坐标
Figure DEST_PATH_IMAGE020
，
Figure DEST_PATH_IMAGE022
。
Figure DEST_PATH_IMAGE024
（6）
如图5所示，其中(a)，(b)表示核线影像上的一对同名像点，(c)，(d)表示取其中一对同名像点进行局部放大，便于清楚展示同名像点的核线影像纵坐标的关系。
由于立体像对影像获取的间隔时间较短，因此，可以认为在该时间间隔内立体像对的水边线基本保持不变，提取立体像对的水边线。首先，利用归一化水体指数法将核线影像分成2个区域，通过预设值（这里的预设值不是一个定值，需要考虑不同地区，不同地物类型等因素进行调整）判断水体区域和非水体区域，将水体部分设置为0，非水体部分设置为1。
以影像内任意一像素点t为例，对该像素点t进行八邻域判断，若所述像素点的八邻域均为0或1，则判定该像素点t不位于水边线上，该像素点t为岛礁内部点或水面点，遍历选取影像中的所有像素点，可提取得到一个像素宽度的水边线区域，并记录水边线区域上所有像素的像素中心坐标值。
在得到立体像对左右一个像素宽度的水边线区域之后，如图6所示，（图中白色部分表示水体区域，灰色表示岛礁内部区域，黑色表示水体与岛礁的分界的一个像素宽度的水边线区域，其中，图6只用于说明，图6中的水边线区域并不代表实际的水边线），使用黑色的环形圈分别表示左右核线影像中的一个像素宽度的水边线。以左核线影像上的水边线（水体部分与非水体部分的分界线）上的任意一个像素Q为例，在像素点Q左侧的边界上，任选一点W点，过W点做水平线与右核线影像的水边线上的像素Q’的左侧交会于W’点，根据核线影像的特性可知，W和W’ 的纵坐标相等，则可得W’的纵坐标，W’的横坐标可由W’所在的像素Q’的中心坐标
Figure DEST_PATH_IMAGE026
的横坐标
Figure DEST_PATH_IMAGE028
减去0.5个像素得到。按照此方法，可以得同名像点的核线影像坐标，然后根据坐标映射关系将同名像点的核线影像坐标转换为原始影像坐标，即通过公式（6）将同名像点的核线影像坐标转换为原始影像坐标，同名像点的原始影像坐标为上述的同名像点的坐标。
S103：选取所述水边线上的同名像点对应的地面点作为控制点，并根据所述水面高程数据和所述同名像点坐标获取控制点三维坐标。
具体地，水面高程数据是作为水边线上同名点的高程控制数据，而由水边线上同名像点进行前方交会获取的平面坐标则作为水边线上同名点的平面控制数据，则可得到同名像点对应的控制点的三维坐标，进而得到控制点的三维坐标。
S104：根据所述控制点三维坐标和RPC模型以两类未知数交替趋近法解算预设的像方空间补偿模型，其中，两类未知数分别为控制点三维坐标的改正数和影像仿射变换参数的改正数。
像方空间补偿模型如式（7）所示，考虑到是对立体像对的区域网平差，而立体像对包括左影像和右影像，因此在矩阵表达式中，用下标left，right区分左右影像的具体有关参数的表达（下同）。
Figure DEST_PATH_IMAGE030
（7）
式（7）中，
Figure DEST_PATH_IMAGE032
表示控制点对应的像点在左原始影像中的列坐标，
Figure DEST_PATH_IMAGE034
表示控制点对应的像点在左原始影像的行坐标，
Figure DEST_PATH_IMAGE036
表示控制点对应的像点在右原始影像的列坐标，
Figure DEST_PATH_IMAGE038
表示控制点对应的像点在右原始影像的行坐标，
Figure DEST_PATH_IMAGE040
，
Figure DEST_PATH_IMAGE042
，
Figure DEST_PATH_IMAGE044
，
Figure DEST_PATH_IMAGE046
表示将物方三维坐标带入RPC模型中计算得到的对应像点坐标近似值，
Figure DEST_PATH_IMAGE048
，
Figure DEST_PATH_IMAGE050
，
Figure DEST_PATH_IMAGE052
，
Figure DEST_PATH_IMAGE054
表示左右影像的像方补偿量，可进行如下表示：
Figure DEST_PATH_IMAGE056
（8）
将式（8）代入式（7）中，并进行整理可得式（9）：
Figure DEST_PATH_IMAGE058
（9）
式（9）中，
Figure DEST_PATH_IMAGE060
，
Figure DEST_PATH_IMAGE062
，
Figure DEST_PATH_IMAGE064
，
Figure DEST_PATH_IMAGE066
表示左右影像的像点近似值与像方补偿后的像点坐标的差值。
可选地，根据所述控制点三维坐标和两类未知数交替趋近法解算所述像点坐标误差方程包括：
将所述像方空间补偿模型转化为像点坐标误差方程，将控制点的三维坐标作为已知参数，以将所述像点坐标误差方程转化为影像仿射变换参数的误差方程，根据所述控制点三维坐标和所述RPC模型求解所述影像仿射变换参数的误差方程，得到影像仿射变换参数的改正数。
由于式（9）是非线性表达式，因此需要将
Figure 862951DEST_PATH_IMAGE060
，
Figure 320477DEST_PATH_IMAGE062
，
Figure 134849DEST_PATH_IMAGE064
，
Figure 19628DEST_PATH_IMAGE066
进行泰勒展开至一次项，才能使用最小二乘进行未知数解算。因此，像点坐标误差方程可以表示为：
V=M 1 X 1 + M 2 X 2 - L，P （10）
式（10）中X 1 表示控制点三维坐标的改正数，X 2 表示立体像对12个影像仿射变换参数的改正数（左右影像各6个影像仿射变换参数的改正数），L表示常数项阵，由像点量测坐标（
Figure 286662DEST_PATH_IMAGE032
、
Figure 547879DEST_PATH_IMAGE034
、
Figure 216757DEST_PATH_IMAGE036
、
Figure 210121DEST_PATH_IMAGE038
）与像点坐标近似值（
Figure 26768DEST_PATH_IMAGE040
，
Figure 763779DEST_PATH_IMAGE042
，
Figure 349482DEST_PATH_IMAGE044
，
Figure 779326DEST_PATH_IMAGE046
）之差计算得到；M 1表示控制点三维坐标的改正数的系数阵，由影像像点坐标误差方程对控制点三维坐标求偏导得到；M 2表示影像仿射变换参数的改正数的系数阵，由影像的像点坐标误差方程对影像仿射变换参数求偏导得到，P表示权阵。具体表达如下：
Figure DEST_PATH_IMAGE068
（11）
其中，
Figure DEST_PATH_IMAGE070
表示立体像对中的左影像列坐标的误差项，
Figure DEST_PATH_IMAGE072
表示立体像对中的左影像行坐标的误差项，
Figure DEST_PATH_IMAGE074
表示立体像对中的右影像列坐标的误差项，
Figure DEST_PATH_IMAGE076
表示立体像对中的右影像行坐标的误差项。
Figure DEST_PATH_IMAGE078
（12）
其中，dLat表示的是控制点的纬度的改正数，dLon表示的是控制点的经度的改正数，dH表示控制点高程的改正数。
Figure DEST_PATH_IMAGE080
（13）
其中，
Figure DEST_PATH_IMAGE082
，
Figure DEST_PATH_IMAGE084
，
Figure DEST_PATH_IMAGE086
表示的是立体像对中的左影像的列坐标方向上的仿射变换参数，
Figure DEST_PATH_IMAGE088
，
Figure DEST_PATH_IMAGE090
，
Figure DEST_PATH_IMAGE092
表示的是立体像对中的左影像的行坐标方向上的仿射变换参数，
Figure DEST_PATH_IMAGE094
，
Figure DEST_PATH_IMAGE096
，
Figure DEST_PATH_IMAGE098
表示的是立体像对中的右影像的列坐标方向上的仿射变换参数，
Figure DEST_PATH_IMAGE100
，
Figure DEST_PATH_IMAGE102
，
Figure DEST_PATH_IMAGE104
表示的是立体像对中的右影像的行坐标方向上的仿射变换参数。
Figure DEST_PATH_IMAGE106
（14）
Figure DEST_PATH_IMAGE108
（15）
Figure DEST_PATH_IMAGE110
（16）
根据目标区域的形状，大致在沿目标区域四周的水边线上选取分布均匀的同名像点作为控制点进行区域网平差。
由式（10）知，基于RPC模型的区域网平差中存在两类未知参数，分别是立体像对的影像仿射变换参数改正数和控制点三维坐标的改正数。使用交替趋近法在进行这两类未知数解算时，首先，假设控制点的三维坐标是已知的，即此时的误差方程就表示为：
Figure DEST_PATH_IMAGE112
（17）
根据所述影像仿射变换参数的改正数对影像仿射变换参数进行改正；
根据改正后的影像仿射变换参数对所述控制点三维坐标进行改正，将所述影像仿射变换参数的改正数代入所述像点坐标误差方程得到控制点三维坐标的改正数的误差方程。
由式（17）计算影像仿射变换参数的改正数，并对影像仿射变换参数进行改正。利用改正后的影像仿射变换参数对控制点三维坐标进行改正，则此时，控制点三维坐标改正的误差方程如式（18）所示：
Figure DEST_PATH_IMAGE114
（18）
根据改正后的控制点三维坐标和所述RPC模型解算所述控制点三维坐标的改正数的误差方程，得到控制点三维坐标的改正数。
可选地，所述根据所述控制点三维坐标和所述RPC模型求解所述影像仿射变换参数的误差方程包括：
根据所述控制点的平面坐标精度和高程坐标精度确定所述控制点的平面坐标权重和高程坐标权重；
根据所述平面坐标权重和所述高程坐标权重以参数加权平差法解算所述控制点三维坐标的改正数的误差方程，得到所述控制点三维坐标的改正数。
在式（18）中，当控制点三维坐标改正量作为未知参数进行解算时，由于控制点的平面坐标和高程坐标的精度不同，因此需要对未知参数赋予不同的权重进行解算。影像在无控制点的情况下，其平面定位精度为m h 米，而光子数据获取的水面高程的精度可表示为m v 米，则平面坐标和高程坐标的权重可分别用式（19）、式（20）表示。
Figure DEST_PATH_IMAGE116
（19）
Figure DEST_PATH_IMAGE118
（20）
Figure DEST_PATH_IMAGE120
（21）
Figure DEST_PATH_IMAGE122
（22）
若所述影像仿射变换参数的改正数和所述控制点三维坐标的改正数均大于对应阈值时，则根据所述影像仿射变换参数的改正数和所述控制点三维坐标的改正数求解新的影像仿射变换参数的改正数和新的控制点三维坐标的改正数。
具体地，根据控制点三维坐标的改正数对改正后的控制点三维坐标进行再次改正，根据再次改正后的控制点三维坐标和RPC模型求解所述影像仿射变换参数的误差方程，得到新的影像仿射变换参数的改正数；基于新的影像仿射变换参数的改正数再次以上述方法求解新的控制点三维坐标的改正数。
若所述影像仿射变换参数的改正数和所述控制点三维坐标的改正数均小于对应阈值时，则结束求解过程。
激光雷达点云的每一个点都包含了三维坐标信息，对目标区域的激光雷达点云数据进行滤波处理可以获取目标区域高精度的水面高程数据；同时根据核线约束条件可获取目标区域的立体像对的水边线上的同名像点坐标，同名像点坐标结合水面高程数据可获取控制点三维坐标，在获取控制点三维坐标后，可根据控制点三维坐标和RPC模型以两类未知数交替趋近法解算预设的像方空间补偿模型，直至满足精度要求，从而完成对目标区域影像的区域网平差。
如图7所示，本发明另一实施例提供一种结合多源数据的区域网平差系统，包括：
滤波模块71，用于对目标区域的激光雷达点云数据进行滤波处理，得到所述目标区域的水面高程数据；
第一获取模块72，用于根据核线约束条件获取所述目标区域的立体像对的水边线上的同名像点的坐标；
第二获取模块73，用于选取所述水边线上的同名像点对应的地面点作为控制点，并根据所述水面高程数据和所述同名像点坐标获取控制点三维坐标；
计算模块74，用于根据预设的像方空间补偿模型和RPC模型获取像点坐标误差方程，根据所述控制点三维坐标和两类未知数交替趋近法解算所述像点坐标误差方程。
所属领域的技术人员可以清楚地了解到，为了描述的方便和简洁，仅以上述各功能单元、模块的划分进行举例说明，实际应用中，可以根据需要而将上述功能分配由不同的功能单元、模块完成，即将所述装置的内部结构划分成不同的功能单元或模块，以完成以上描述的全部或者部分功能。实施例中的各功能单元、模块可以集成在一个处理单元中，也可以是各个单元单独物理存在，也可以两个或两个以上单元集成在一个单元中，上述集成的单元既可以采用硬件的形式实现，也可以采用软件功能单元的形式实现。另外，各功能单元、模块的具体名称也只是为了便于相互区分，并不用于限制本申请的保护范围。上述系统中单元、模块的具体工作过程，可以参考前述方法实施例中的对应过程，在此不再赘述。
在上述实施例中，对各个实施例的描述都各有侧重，某个实施例中没有详述或记载的部分，可以参见其它实施例的相关描述。
本发明实现上述实施例方法中的全部或部分流程，可以通过计算机程序来指令相关的硬件来完成，所述的计算机程序可存储于一存储介质中，该计算机程序在被处理器执行时，可实现上述各个方法实施例的步骤。其中，所述计算机程序包括计算机程序代码，所述计算机程序代码可以为源代码形式、对象代码形式、可执行文件或某些中间形式等。所述存储介质可以包括：能够携带所述计算机程序代码的任何实体或装置、记录介质、U盘、移动硬盘、磁碟、光盘、计算机存储器、只读存储器(ROM，Read-Only Memory)、随机存取存储器(RAM，Random Access Memory)、电载波信号、电信信号以及软件分发介质等。需要说明的是，所述存储介质包含的内容可以根据司法管辖区内立法和专利实践的要求进行适当的增减，例如在某些司法管辖区，根据立法和专利实践，存储介质不包括电载波信号和电信信号。
应理解，上述实施例中各步骤的序号的大小并不意味着执行顺序的先后，各过程的执行顺序应以其功能和内在逻辑确定，而不应对本发明实施例的实施过程构成任何限定。
虽然本公开披露如上，但本公开的保护范围并非仅限于此。本领域技术人员在不脱离本公开的精神和范围的前提下，可进行各种变更与修改，这些变更与修改均将落入本发明的保护范围。

Claims (10)

1.一种结合多源数据的区域网平差方法，其特征在于，包括：
对目标区域的激光雷达点云数据进行滤波处理，得到所述目标区域的水面高程数据；
根据核线约束条件获取所述目标区域的立体像对的水边线上的同名像点的坐标；
选取所述水边线上的同名像点对应的地面点作为控制点，并根据所述水面高程数据和所述同名像点的坐标获取控制点三维坐标；
根据所述控制点三维坐标和RPC模型以两类未知数交替趋近法解算预设的像方空间补偿模型；
将所述像方空间补偿模型转化为像点坐标误差方程，根据所述控制点三维坐标和两类未知数交替趋近法解算所述像点坐标误差方程;
其中，两类未知数分别为控制点三维坐标的改正数和影像仿射变换参数的改正数。
2.根据权利要求1所述的结合多源数据的区域网平差方法，其特征在于，所述对目标区域的激光雷达点云数据进行滤波处理包括：
以预设的高程间隔将所述激光雷达点云数据在高程方向切分为K个切片区间；
获取各所述切片区间内的光子信号数量和各所述切片区间的中心高程；
根据各所述切片区间的光子信号数量和各所述切片区间的中心高程进行正态分布拟合，得到拟合函数；
对所述拟合函数进行最优解算，获取所述拟合函数的参数；
根据所述拟合函数的参数构建水面光子信号高程区间，并根据所述水面光子信号高程区间获取第一水面光子信号数据集。
3.根据权利要求2所述的结合多源数据的区域网平差方法，其特征在于，所述对目标区域的激光雷达点云数据进行滤波处理还包括：
对所述第一水面光子信号数据集进行曲线拟合，得到所述第一水面光子信号数据集的曲线方程；
对所述曲线方程求一阶导，获取导数为零的点的高程值；
将对应所述第一水面光子信号数据集的曲线的波谷的点的高程值作为分界值，根据所述分界值对所述第一水面光子信号数据集进行进一步滤波处理，得到第二水面光子信号数据集。
4.根据权利要求3所述的结合多源数据的区域网平差方法，其特征在于，所述对目标区域的激光雷达点云数据进行滤波处理还包括：
对所述第二水面光子信号数据集进行正态分布拟合，将得到的所述第二水面光子信号数据集的正态分布曲线的预设范围作为可靠性区域，并将所述可靠性区域的区域宽度作为圆形滤波窗口直径长度；
沿着轨道方向，依次以所述第二水面光子数据集中的每个光子信号作为中心，以所述圆形滤波窗口直径长度的一半为半径构建圆形滤波窗口，统计每个所述圆形滤波窗口内的光子信号数量，并计算光子信号数量均值；
获取内部光子信号数量大于所述光子信号数量均值的所述圆形滤波窗口，根据获取结果构建第三水面光子信号数据集。
5.根据权利要求4所述的结合多源数据的区域网平差方法，其特征在于，所述对目标区域的激光雷达点云数据进行滤波处理还包括：
对所述第三水面光子信号数据集中的光子信号进行潮汐改正，得到影像获取时刻的水面高程值。
6.根据权利要求1所述的结合多源数据的区域网平差方法，其特征在于，所述根据核线约束条件获取所述目标区域的立体像对的水边线上的同名像点像点坐标包括：
获取所述目标区域的立体像对，并根据核线约束条件，获取所述立体像对的核线影像；
根据归一化水体指数法分别将所述立体像对的核线影像分成水体区域和非水体区域；
根据所述水体区域和所述非水体区域分别提取所述立体像对的核线影像上的水边线区域；
选取所述立体像对的核线影像上的水边线区域上的相应点作为同名像点，获取所述同名像点的核线影像坐标；
根据所述立体像对的原始影像与核线影像之间的坐标映射关系获取所述同名像点的原始影像坐标。
7.根据权利要求1所述的结合多源数据的区域网平差方法，其特征在于，根据所述控制点三维坐标和两类未知数交替趋近法解算所述像点坐标误差方程包括：
将控制点的三维坐标作为已知，以将所述像点坐标误差方程转化为影像仿射变换参数的误差方程，根据所述控制点三维坐标和所述RPC模型求解所述影像仿射变换参数的误差方程，得到影像仿射变换参数的改正数；
根据所述影像仿射变换参数的改正数对影像仿射变换参数进行改正；
根据改正后的影像仿射变换参数对所述控制点三维坐标进行改正，将所述影像仿射变换参数的改正数代入所述像点坐标误差方程得到控制点三维坐标的改正数的误差方程；
根据改正后的控制点三维坐标和所述RPC模型解算所述控制点三维坐标的改正数的误差方程，得到控制点三维坐标的改正数；
当所述影像仿射变换参数的改正数和所述控制点三维坐标的改正数均大于对应阈值时，则根据所述影像仿射变换参数的改正数和所述控制点三维坐标的改正数求解新的影像仿射变换参数的改正数和新的控制点三维坐标的改正数。
8.根据权利要求7所述的结合多源数据的区域网平差方法，其特征在于，所述根据改正后的控制点三维坐标和所述RPC模型解算所述控制点三维坐标的改正数的误差方程包括：
根据所述控制点的平面坐标精度和高程坐标精度确定所述控制点的平面坐标权重和高程坐标权重；
根据所述平面坐标权重和所述高程坐标权重以参数加权平差法解算所述控制点三维坐标的改正数的误差方程，得到所述控制点三维坐标的改正数。
9.一种结合多源数据的区域网平差系统，其特征在于，包括：
滤波模块，用于对目标区域的激光雷达点云数据进行滤波处理，得到所述目标区域的水面高程数据；
第一获取模块，用于根据核线约束条件获取所述目标区域的立体像对的水边线上的同名像点的坐标；
第二获取模块，用于选取所述水边线上的同名像点对应的地面点作为控制点，并根据所述水面高程数据和所述同名像点的坐标获取控制点三维坐标；
计算模块，用于根据所述控制点三维坐标和RPC模型以两类未知数交替趋近法解算预设的像方空间补偿模型，其中，两类未知数为控制点三维坐标的改正数和影像仿射变换参数的改正数。
10.一种存储介质，所述存储介质存储有计算机程序，其特征在于，所述计算机程序被处理器执行时，实现如权利要求1至8任一项所述的结合多源数据的区域网平差方法。
CN202110668978.5A 2021-06-17 2021-06-17 一种结合多源数据的区域网平差方法、系统及存储介质 Active CN113124834B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110668978.5A CN113124834B (zh) 2021-06-17 2021-06-17 一种结合多源数据的区域网平差方法、系统及存储介质

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110668978.5A CN113124834B (zh) 2021-06-17 2021-06-17 一种结合多源数据的区域网平差方法、系统及存储介质

Publications (2)

Publication Number Publication Date
CN113124834A true CN113124834A (zh) 2021-07-16
CN113124834B CN113124834B (zh) 2021-10-12

Family

ID=76783014

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110668978.5A Active CN113124834B (zh) 2021-06-17 2021-06-17 一种结合多源数据的区域网平差方法、系统及存储介质

Country Status (1)

Country Link
CN (1) CN113124834B (zh)

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108919319A (zh) * 2018-05-15 2018-11-30 中国人民解放军战略支援部队信息工程大学 海岛礁卫星影像无地面控制点定位方法及系统
CN111277243A (zh) * 2020-02-06 2020-06-12 中国地质大学（武汉） 一种面向水深提取的光子计数激光雷达自适应滤波算法
CN111522019A (zh) * 2020-07-02 2020-08-11 中国地质大学（武汉） 一种水下光子位置的误差修正方法及装置

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108919319A (zh) * 2018-05-15 2018-11-30 中国人民解放军战略支援部队信息工程大学 海岛礁卫星影像无地面控制点定位方法及系统
CN111277243A (zh) * 2020-02-06 2020-06-12 中国地质大学（武汉） 一种面向水深提取的光子计数激光雷达自适应滤波算法
CN111522019A (zh) * 2020-07-02 2020-08-11 中国地质大学（武汉） 一种水下光子位置的误差修正方法及装置

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
喻金桃等: "基于水边线等高的海岸带卫星影像定位方法 ", 《海洋测绘》 *
宋冰等: "线阵CCD卫星影像DSM自动生成技术 ", 《测绘通报》 *
王慧等: "潮汐数据支持的海岸带区域网平差 ", 《海洋测绘》 *
申家双等: "海岸带水边线等高条件控制下的航空影像外部定向 ", 《海洋测绘》 *

Also Published As

Publication number Publication date
CN113124834B (zh) 2021-10-12

Similar Documents

Publication Publication Date Title
Hodúl et al. Satellite derived photogrammetric bathymetry
Westaway et al. The development of an automated correction procedure for digital photogrammetry for the study of wide, shallow, gravel‐bed rivers
Stojic et al. The assessment of sediment transport rates by automated digital photogrammetry
Del Río et al. Error determination in the photogrammetric assessment of shoreline changes
Ellmann et al. The 5 mm geoid model for Estonia computed by the least squares modified Stokes’s formula
CN110335355B (zh) 一种大型浅水湖泊水面高自动计算方法
CN109100719B (zh) 基于星载sar影像与光学影像的地形图联合测图方法
Mason et al. Measurement of recent intertidal sediment transport in Morecambe Bay using the waterline method
CN112433227B (zh) 一种水容量变化监测方法、系统、终端设备及存储介质
Welch et al. A photogrammetric technique for measuring soil erosion
CN113358091A (zh) 一种利用三线阵立体卫星影像生产数字高程模型方法
CN111144350B (zh) 一种基于参考底图的遥感影像定位精度评价方法
Fryer et al. On the accuracy of heighting from aerial photographs and maps: implications to process modellers
CN110207676A (zh) 一种田沟塘参数的获取方法及装置
CN105571598B (zh) 一种卫星激光高度计足印相机姿态的测定方法
CN113238228B (zh) 基于水准约束的三维地表形变获取方法、系统及装置
CN110133655A (zh) 基于多源雷达遥感技术的河流径流监测反演方法
CN113124834B (zh) 一种结合多源数据的区域网平差方法、系统及存储介质
Partama et al. A simple and empirical refraction correction method for UAV-based shallow-water photogrammetry
CN116310901A (zh) 一种基于低空遥感的泥石流物源动态运移识别方法
CN113516764B (zh) 基于数字高程模型的湖库水下三维地形模拟方法及装置
CN114494625A (zh) 高精度地形图制作方法、装置和计算机设备
CN114255051A (zh) 基于立体测绘卫星的正射产品的真实性检验方法
CN113280789B (zh) 一种地表起伏区域激光测高点作为影像高程控制点方法
Gallant et al. Digital terrain analysis

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