CN114463523B - 一种最小豪斯多夫距离面特征约束的点云航带平差方法 - Google Patents
一种最小豪斯多夫距离面特征约束的点云航带平差方法 Download PDFInfo
- Publication number
- CN114463523B CN114463523B CN202210014763.6A CN202210014763A CN114463523B CN 114463523 B CN114463523 B CN 114463523B CN 202210014763 A CN202210014763 A CN 202210014763A CN 114463523 B CN114463523 B CN 114463523B
- Authority
- CN
- China
- Prior art keywords
- building
- laser radar
- point
- adjacent
- airborne laser
- 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
Links
- 238000000034 method Methods 0.000 title claims abstract description 36
- 239000013598 vector Substances 0.000 claims abstract description 49
- 238000006243 chemical reaction Methods 0.000 claims abstract description 16
- 238000012545 processing Methods 0.000 claims abstract description 12
- 238000001914 filtration Methods 0.000 claims abstract description 11
- 230000000750 progressive effect Effects 0.000 claims abstract description 6
- 238000004364 calculation method Methods 0.000 claims description 22
- 230000011218 segmentation Effects 0.000 claims description 10
- 238000012216 screening Methods 0.000 claims description 6
- 238000004891 communication Methods 0.000 claims description 3
- 239000011159 matrix material Substances 0.000 description 8
- 238000001514 detection method Methods 0.000 description 4
- 238000010586 diagram Methods 0.000 description 4
- 238000012937 correction Methods 0.000 description 2
- 238000012946 outsourcing Methods 0.000 description 2
- 230000000452 restraining effect Effects 0.000 description 2
- ATJFFYVFTNAWJD-UHFFFAOYSA-N Tin Chemical compound [Sn] ATJFFYVFTNAWJD-UHFFFAOYSA-N 0.000 description 1
- 238000007792 addition Methods 0.000 description 1
- 239000007795 chemical reaction product Substances 0.000 description 1
- 238000000354 decomposition reaction Methods 0.000 description 1
- 238000003384 imaging method Methods 0.000 description 1
- 230000010354 integration Effects 0.000 description 1
- 238000002372 labelling Methods 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 230000009897 systematic effect Effects 0.000 description 1
- 230000009466 transformation Effects 0.000 description 1
- 238000013519 translation Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T19/00—Manipulating 3D models or images for computer graphics
- G06T19/20—Editing of 3D images, e.g. changing shapes or colours, aligning objects or positioning parts
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/23—Clustering techniques
- G06F18/232—Non-hierarchical techniques
- G06F18/2321—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions
- G06F18/23213—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions with fixed number of clusters, e.g. K-means clustering
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
- G06T17/20—Finite element generation, e.g. wire-frame surface description, tesselation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Data Mining & Analysis (AREA)
- General Physics & Mathematics (AREA)
- Software Systems (AREA)
- Computer Graphics (AREA)
- General Engineering & Computer Science (AREA)
- Evolutionary Computation (AREA)
- Evolutionary Biology (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Bioinformatics & Computational Biology (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Artificial Intelligence (AREA)
- Life Sciences & Earth Sciences (AREA)
- Geometry (AREA)
- Probability & Statistics with Applications (AREA)
- Architecture (AREA)
- Computer Hardware Design (AREA)
- Image Processing (AREA)
Abstract
本发明公开了一种最小豪斯多夫距离面特征约束的点云航带平差方法，主要包括：使用渐进三角网加密滤波算法对每条航带点云进行处理，得到地面点与非地面点；通过半抑制模糊C均值聚类与区域增长从非地面点中探测到每条航带中建筑物点；根据每条航带最小外包矩形确定相邻航带重叠区域内建筑物点，并对重叠区域内建筑物点进行处理得到单栋建筑物；使用由粗到精的方法从单栋建筑物中分割出面片；使用最小豪斯多夫距离与法向量约束从相邻航带重叠区域内探测出同名面片；使用同名面片构建误差方程并解算相邻航带转换参数。本发明以数据驱动的方式对航带点云进行平差处理，可为后续各类应用提供可靠数据。
Description
技术领域
本发明属于机载激光雷达数据处理技术领域,尤其涉及一种最小豪斯多夫距离面特征约束的点云航带平差方法。
背景技术
机载激光雷达(Airborne Light Detection And Ranging，LiDAR)是一种主动式航空遥感对地观测系统，其不受天气环境影响，可以获取高精度、高密度的点云数据，其为高时空高分辨率的地球空间信息提供了一种全新的技术手段。
在实际应用中，LiDAR以单条航带的方式搜集测区内地形数据，其中单条航带的范围远远小于整个测区范围，因此一般需要多条航带才能完整覆盖整个测区。在实际应用中，LiDAR系统需要在搜集数据前由系统提供商和最终用户进行校准，以消除系统各部件之间集成而引起的系统误差。由于缺少精确的系统校准或者不规范的操作，在实际获取的相邻数据中仍存在偏差。因此要获取高精度的最终产品，如DEM，则必须进行航带拼接，以消除航带间偏差。然而由于机载点云密度低，很难从点云中直接获取同名点，因此如何获取高精度的同名特征是航带平差中的一个难点。另外，传统的基于人工手动选取同名特征的方式耗时不准确，其难以适应当今数据自动处理的需要。
发明内容
本发明针对基于点云数据本身的航带平差方法中难以自动搜索同名特征的问题，在此提出一种最小豪斯多夫距离面特征约束的点云航带平差方法，实现机载LiDAR航带间自动平差计算，消除由系统检校不完全而引起的误差。
实现本发明所说的一种最小豪斯多夫距离面特征约束的点云航带平差方法的技术方案具体如下：
步骤1：引入多条机载激光雷达航带点云数据；
步骤2：将每条机载激光雷达航带点云数据使用渐进三角网加密滤波算法进行滤波处理，获取每条机载激光雷达航带点云数据中非地面点云以及每条机载激光雷达航带点云数据中地面点云；
步骤3：将每条机载激光雷达航带点云数据中非地面点云通过半抑制模糊C均值聚类方法得到每条机载激光雷达航带点云数据中多个初始建筑物种子点，将每条机载激光雷达航带点云数据中多个初始建筑物种子点通过区域增长聚类方法，得到每条机载激光雷达航带点云数据中非地面点云中的多个建筑物点；
步骤4：搜索每条机载激光雷达航带点云数据中非地面点云中的多个建筑物点的X轴坐标最大值、X轴坐标最小值、Y轴坐标最大值、Y轴坐标最小值，根据X轴与Y轴最大值与最小值确定两边平行于X轴与Y轴的包含每条机载激光雷达航带的最小外包矩形，将相邻两条机载激光雷达航带的最小外包矩形进行相交处理得到相邻两条机载激光雷达航带的重叠区域，将相邻两条机载激光雷达航带的重叠区域进行格网化处理得到相邻两条机载激光雷达航带的重叠区域二值图像，将相邻两条机载激光雷达航带的重叠区域二值图像通过四连通区域标签法进行处理得到相邻两条机载激光雷达航带的重叠区域多个图像块，在相邻两条机载激光雷达航带的重叠区域每个图像块中搜索步骤3所述的每条机载激光雷达航带点云数据中非地面点云中的多个建筑物点，得到相邻两条机载激光雷达航带的重叠区域每个图像块中多个建筑物点；
步骤5：使用主成分分析法计算相邻两条机载激光雷达航带的重叠区域每个图像块中每个建筑物点的法向量、每个建筑物点的平面、每个建筑物点的曲率，将相邻两条机载激光雷达航带的重叠区域每个图像块中建筑物点的曲率最小的建筑物点作为相邻两条机载激光雷达航带的重叠区域对应图像块的种子点，将相邻两条机载激光雷达航带的重叠区域每个图像块的除种子点外的建筑物点构建相邻两条机载激光雷达航带的重叠区域每个图像块的非种子点集，使用点到平面计算方法与向量间夹角计算方法分别计算相邻两条机载激光雷达航带的重叠区域中每个图像块的非种子点集中每个建筑物点与种子点的法向量夹角、相邻两条机载激光雷达航带的重叠区域中每个图像块的非种子点集中每个建筑物点到种子点的平面的距离，将相邻两条机载激光雷达航带的重叠区域每个图像块的非种子点集中建筑物点与种子点的法向量夹角大于等于夹角阈值或到种子点的平面的距离大于等于距离阈值的建筑物点筛选出来得到相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集，将相邻两条机载激光雷达航带的重叠区域每个图像块的非种子点集中建筑物点与种子点的法向量夹角小于夹角阈值且到种子点的平面的距离小于距离阈值的建筑物点筛选出来得到相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物平面点集，进一步与相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物平面点集通过最小二乘方法构建初始平面，使用向量间夹角计算方法计算相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集每个建筑物点到种子点的法向量夹角，使用点到平面距离计算方法计算相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集每个建筑物点到初始平面的距离，将相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集每个建筑物点与种子点的法向量夹角小于距离阈值且到初始平面的距离小于距离阈值的建筑物点筛选出来与建筑物平面点集、种子点共同构成相邻两条机载激光雷达航带的重叠区域每个图像块的最终建筑物平面点集，将相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集每个建筑物点到种子点的法向量夹角大于等于距离阈值或者到初始平面的距离大于等于距离阈值的建筑物点筛选出来构建相邻两条机载激光雷达航带的重叠区域每个图像块的最终非建筑物平面点集；
步骤5所述相邻两条机载激光雷达航带的重叠区域每个图像块的最终建筑物平面点集中每个建筑物点均属于初始平面；
步骤5所述相邻两条机载激光雷达航带的重叠区域每个图像块的最终非建筑物平面点集每个建筑物点均有对应的平面；
步骤6：对于一条机载激光雷达航带的重叠区域每个图像块的建筑物点集，使用豪斯多夫距离计算方法计算其到相邻机载激光雷达航带的重叠区域每个图像块的建筑物点集的豪斯多夫距离，若一条机载激光雷达航带的重叠区域的一个图像块的建筑物点集到相邻机载激光雷达航带的重叠区域中每个图像块的建筑物点集的最小豪斯多夫距离小于阈值，将一条机载激光雷达航带的重叠区域建筑物点集与最小豪斯多夫距离对应的相邻机载激光雷达航带的重叠区域图像块的建筑物点集作为一对同名建筑物点集。使用豪斯多夫距离计算方法与向量间夹角计算方法计算一对同名建筑物点集中一栋建筑物点集分割对应的平面到一对同名建筑物点集中另一栋建筑物点集中每个面片的豪斯多夫距离以及法向量夹角。若一对同名建筑物点集中一栋建筑物点集分割对应的平面到最小豪斯多夫距离对应的一对同名建筑物点集中另一栋建筑物点集分割的面片对应的法向量夹角小于夹角阈值，则将一对同名建筑物点集中一栋建筑物点集分割对应的平面与一对同名建筑物点集中另一栋建筑物点集分割的面片作为同名面片。
步骤7：以相邻航带探测得到的同名面片方程以及同名面片对应的点云数据，根据同名面片中一个面片中点云数据满足同名面片中另一面片方程构造误差方程，根据最小二乘原理解算转换参数。以相邻航带中一条航带为基准，使用解算的转换参数将相邻航带中另一条航带中点云数据进行旋转与平移实现将相邻航带点云数据转换到同一坐标系下。
本发明与现有技术相比，具有以下优点：
现有人工手动获取同名特征面的方法，难以选取满足平差条件的同名平面，且自动化程度低，选取的同名面片数量有限，使用最小豪斯多夫距离和法向量一致性进行约束，能够从建筑物屋面中探测出满足航带平差的同名特征面。由于使用了严格的约束限制条件，从而加强了该方法在探测获取准确同名特征面的鲁棒性。
采用由粗到精的建筑物屋面分割方法，可以精确分割出常见建筑物结构的屋顶面片，方便后续同名面片的探测。
附图说明
图1：一种最小豪斯多夫距离面特征约束的点云航带平差方法的流程图；
图2：PTD中初始种子点选取并构建TIN示意图；
图3：计算角度与距离示意图；
图4：单栋建筑物分割示意图；
图5：由粗到精的建筑物屋面分割示意图；
图6：同名特征面探测结果示意图；
图7：航带平差结果示意图。
具体实施方式
为了便于本领域普通技术人员理解和实施本发明，下面结合附图及实施例对本发明作进一步的详细描述，应当理解，此处所描述的实施示例仅用于说明和解释本发明，并不用于限定本发明。
下面通过实施例，并结合附图，对本发明的技术方案作进一步具体的说明本发明第一实施例以及第二实施例。
如图1所示，本发明第一实施例为一种最小豪斯多夫距离面特征约束的云航带的流程图，包括以下步骤：
步骤1：引入多条机载激光雷达航带点云数据；
步骤2：将每条机载激光雷达航带点云数据使用渐进三角网加密滤波算法进行滤波处理，获取每条机载激光雷达航带点云数据中非地面点云以及每条机载激光雷达航带点云数据中地面点云；
步骤3：将每条机载激光雷达航带点云数据中非地面点云通过半抑制模糊C均值聚类方法得到每条机载激光雷达航带点云数据中多个初始建筑物种子点，将每条机载激光雷达航带点云数据中多个初始建筑物种子点通过区域增长聚类方法，得到每条机载激光雷达航带点云数据中非地面点云中的多个建筑物点；
步骤4：搜索每条机载激光雷达航带点云数据中非地面点云中的多个建筑物点的X轴坐标最大值、X轴坐标最小值、Y轴坐标最大值、Y轴坐标最小值，根据X轴与Y轴最大值与最小值确定两边平行于X轴与Y轴的包含每条机载激光雷达航带的最小外包矩形，将相邻两条机载激光雷达航带的最小外包矩形进行相交处理得到相邻两条机载激光雷达航带的重叠区域，将相邻两条机载激光雷达航带的重叠区域进行格网化处理得到相邻两条机载激光雷达航带的重叠区域二值图像，将相邻两条机载激光雷达航带的重叠区域二值图像通过四连通区域标签法进行处理得到相邻两条机载激光雷达航带的重叠区域多个图像块，在相邻两条机载激光雷达航带的重叠区域每个图像块中搜索步骤3所述的每条机载激光雷达航带点云数据中非地面点云中的多个建筑物点，得到相邻两条机载激光雷达航带的重叠区域每个图像块中多个建筑物点；
步骤5：使用主成分分析法计算相邻两条机载激光雷达航带的重叠区域每个图像块中每个建筑物点的法向量、每个建筑物点的平面、每个建筑物点的曲率，将相邻两条机载激光雷达航带的重叠区域每个图像块中建筑物点的曲率最小的建筑物点作为相邻两条机载激光雷达航带的重叠区域对应图像块的种子点，将相邻两条机载激光雷达航带的重叠区域每个图像块的除种子点外的建筑物点构建相邻两条机载激光雷达航带的重叠区域每个图像块的非种子点集，使用点到平面计算方法与向量间夹角计算方法分别计算相邻两条机载激光雷达航带的重叠区域中每个图像块的非种子点集中每个建筑物点与种子点的法向量夹角、相邻两条机载激光雷达航带的重叠区域中每个图像块的非种子点集中每个建筑物点到种子点的平面的距离，将相邻两条机载激光雷达航带的重叠区域每个图像块的非种子点集中建筑物点与种子点的法向量夹角大于等于夹角阈值或到种子点的平面的距离大于等于距离阈值的建筑物点筛选出来得到相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集，将相邻两条机载激光雷达航带的重叠区域每个图像块的非种子点集中建筑物点与种子点的法向量夹角小于夹角阈值且到种子点的平面的距离小于距离阈值的建筑物点筛选出来得到相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物平面点集，进一步与相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物平面点集通过最小二乘方法构建初始平面，使用向量间夹角计算方法计算相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集每个建筑物点到种子点的法向量夹角，使用点到平面距离计算方法计算相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集每个建筑物点到初始平面的距离，将相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集每个建筑物点与种子点的法向量夹角小于距离阈值且到初始平面的距离小于距离阈值的建筑物点筛选出来与建筑物平面点集、种子点共同构成相邻两条机载激光雷达航带的重叠区域每个图像块的最终建筑物平面点集，将相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集每个建筑物点到种子点的法向量夹角大于等于距离阈值或者到初始平面的距离大于等于距离阈值的建筑物点筛选出来构建相邻两条机载激光雷达航带的重叠区域每个图像块的最终非建筑物平面点集；
步骤5所述相邻两条机载激光雷达航带的重叠区域每个图像块的最终建筑物平面点集中每个建筑物点均属于初始平面；
步骤5所述相邻两条机载激光雷达航带的重叠区域每个图像块的最终非建筑物平面点集每个建筑物点均有对应的平面；
步骤6：对于一条机载激光雷达航带的重叠区域每个图像块的建筑物点集，使用豪斯多夫距离计算方法计算其到相邻机载激光雷达航带的重叠区域每个图像块的建筑物点集的豪斯多夫距离，若一条机载激光雷达航带的重叠区域的一个图像块的建筑物点集到相邻机载激光雷达航带的重叠区域中每个图像块的建筑物点集的最小豪斯多夫距离小于阈值，将一条机载激光雷达航带的重叠区域建筑物点集与最小豪斯多夫距离对应的相邻机载激光雷达航带的重叠区域图像块的建筑物点集作为一对同名建筑物点集。使用豪斯多夫距离计算方法与向量间夹角计算方法计算一对同名建筑物点集中一栋建筑物点集分割对应的平面到一对同名建筑物点集中另一栋建筑物点集中每个面片的豪斯多夫距离以及法向量夹角。若一对同名建筑物点集中一栋建筑物点集分割对应的平面到最小豪斯多夫距离对应的一对同名建筑物点集中另一栋建筑物点集分割的面片对应的法向量夹角小于夹角阈值，则将一对同名建筑物点集中一栋建筑物点集分割对应的平面与一对同名建筑物点集中另一栋建筑物点集分割的面片作为同名面片。
步骤7：以相邻航带探测得到的同名面片方程以及同名面片对应的点云数据，根据同名面片中一个面片中点云数据满足同名面片中另一面片方程构造误差方程，根据最小二乘原理解算转换参数。以相邻航带中一条航带为基准，使用解算的转换参数将相邻航带中另一条航带中点云数据进行旋转与平移实现将相邻航带点云数据转换到同一坐标系下。
本发明第二实施例具体如下：
步骤1，在航飞载体上搭载一套LiDAR系统，包括惯性测量单元(IMU)、差分GPS、激光扫描测距系统和成像装置，采集得到建筑物场景的机载激光雷达点云数据；
步骤2，对步骤1中的激光雷达数据利用渐进三角网加密滤波算法对点云进行滤波操作，获取地面点与非地面点；
步骤3，对步骤2中的非地面点，使用半抑制模糊C均值和限制区域增长算法从非地面点中提取建筑物；
步骤4，将包含航带最小外包矩形划分成规则格网，并将步骤3中获取的建筑物点依据其坐标划分到相应的格网内，构建二值图，并使用连通区域标签法对二值图进行处理，以此对建筑物进行分割，得到单栋建筑物；
步骤5，使用由粗到精的方法对每一栋建筑物进行屋面分割得到建筑物面片；
步骤6，使用最小豪斯多夫距离以及法向量一致性进行约束，从同名建筑物中探测获取同名特征面；
步骤7，使用步骤6中获取得到的同名特征面构建航带平差模型，并依据一条航带中点经转换到另一航带中对应面的欧氏距离最小作为限制条件，解算转换参数，实现航带平差。
如上所述的步骤2中，使用渐进三角网加密滤波算法对点云进行滤波操作，得到非地面点的方法为：
步骤2.1，首先将整个数据划分成规则格网，并将每个格网中最低点指定为初始地面点(如图2所示)；
步骤2.2，将步骤2.1中获取的地面点构建成三角网作为参考地面(如图2所示)；
步骤2.3，对于每个三角形，对位于该三角形中未标记的点，计算该点到三角网中最近的距离d，以及该点和该三角形顶点的连线与三角形面之间的夹角α(如图3所示)；
其中为未标记点与三角网顶点构成的向量，为该三角网面对应的法向量，三角网平面方程ax+by+cz+d＝0，arccos()为反余弦函数；
步骤2.4，若步骤2.3中计算的角度与距离均小于设定的阈值，则将该点标记为地面点；
步骤2.5，将当前迭代中标记的所有地面点添加到先前的地面点集中并重建三角网，执行步骤2.3-2.4，直至所有点标记为地面点或者非地面点时，停止迭代。
如上所述的步骤3中，利用半抑制模糊C均值和限制区域增长算法从非地面点中提取建筑物的方法为：
步骤3.1，特征提取，每个点依据其圆柱体近邻点，计算获取5个特征，其中包括最高点与最低点的高差Δh、平均高程hmean、平面粗糙度σ、平面性planity以及点密度density：
Δh＝hmax-hmin (3)
其中hmax与hmin为圆柱体近邻点中最高点与最低点高程；hi为近邻点中某一点的高程，di为点到面的距离，n为圆柱体近邻点的个数，s为圆柱体的底面面积，λ1、λ2、λ3为近邻点构成协方差矩阵分解对应的3个特征值，其中λ1≤λ2≤λ3；
步骤3.2，将所有点的特征集合作为模糊C均值分类器的输入部分，进行无监督聚类，根据下式原则，将点集划分成两类，得到建筑物初始种子点：
其中dij为采样点xj到第i个聚类中心的距离，uij为模糊成员矩阵；
步骤3.3，以步骤3.2中探测得到的建筑物种子点，进行区域增长。计算其他点到种子点拟合的平面距离d、点间高差Δh以及该点与平面法向量之间夹角θ，以此提取出建筑物点；
Δh＝hmax-hmin (10)
其中，ax+by+cz+d＝0为种子点拟合的平面方程；hmax与hmin为最高点与最低点高程；分别为待判断点与种子点的法向量；
步骤3.4，根据点间连通性对所有屋顶点进行聚类，将所占面积过小的点集进行剔除，完成建筑物检测。
如上所述的步骤4中，利用构建的规则格网，将建筑物点投影到XOY面上提取出单栋建筑物的方法为：
步骤4.1，统计出建筑物点坐标，获取x、y坐标的最大值与最小值，记为xmin、xmax、ymin、ymax，并依据坐标值计算其最小外包矩形的长L与宽W；
步骤4.2，依据建筑物点每一点的坐标xi、yi，将其放入到相应的格网内构建二值图，其中该点对应的像素位置为(ui，vi)，该格网像素设置成0，否则设置成255，得到二值化图像(如图4(b)所示)。
其中l为格网的尺寸；
步骤4.3，选取像素为0的一个格网作为种子点Pseed，若有其他像素为0的格网与Pseed相邻，则将其归为种子点集，更新种子点集；
步骤4.4，对更新的种子点集重复步骤4.3，直至无种子点加入则结束，将该格网全部标记为相同的标签号；
步骤4.5，选取新的未标记像素为0的格网为像素种子点，重复步骤4.3-4.4，直至所有像素为0的格网均进行标记(如图4(c)所示)；
步骤4.6，将标签号相同的格网内点提取出来，即对应为同一栋建筑物点，实现建筑物单栋分割(如图4(d)所示)。
如上所述的步骤5中，基于由粗到精的建筑物屋面分割的方法为：
步骤5.1，点云预处理，依据每一点p的近邻点，通过协方差矩阵计算该点的法向量，以及该点拟合平面方程；
其中qi∈N＝(q1，q1，...，qk)为p点的近邻点；
zi+vzi＝a(x+vxi)+b(y+vyi)+c (15)
其中(a，b，c)为拟合平面参数，vxi、vyi、vzi为该点在系统上的观测误差；
步骤5.2，选取平面拟合误差最小的一点作为种子点，计算其他非种子点法向量与种子点法向量间夹角θ，以及该点到种子点所在面的距离d。若均小于设定阈值，则将其归为一类，提取出所有满足上述条件的点；
其中，ax+by+cz+d＝0为种子点拟合的平面方程；分别为待判断点与种子点的法向量；
步骤5.3，从剩余的未分割的点集中重新选择拟合误差最小的点作为种子点，重复步骤5.2，直至将所有点归为一平面点(如图5(a)所示)；
步骤5.4，从粗分割平面中选择包含点数最多的点集拟合成面，再计算所有点到该拟合面的距离，若小于阈值则认为属于该面；
步骤5.5，步骤5.4中提取的点集拟合成面，使用坐标转换公式(18)将三维点集转换到拟合平面上，并进行格网化，将包含点的格网像素设置成255，否则为0(如图5(c)(d))；
其中[x′ y′ z′]T为投影到拟合平面上的坐标，[x y z]T为转换投影前坐标，R为坐标转换矩阵；
步骤5.6，根据步骤5.5得到的投影坐标，将其进行格网化得到二值图，采用步骤4.3-4.6进行连通区域标签法得到多块点集，选取包含点数最多的点集，得到精确的属于同一面片的点集(如图5(e)(f)所示)；
步骤5.7，对于剩余的点集，重复步骤5.3-5.6，提取出屋面，直至所有屋面上点集提取完毕结束。
如上所述的步骤6，使用最小豪斯多夫距离以及法向量一致性进行约束，从同名建筑物中探测获取同名特征面的方法为：
步骤6.1，对于一航带中所有的单栋建筑物，分别进行步骤6.2-6.4：
步骤6.2，使用豪斯多夫距离公式，计算该单栋建筑物与另一航带所有建筑物的豪斯多夫距离；
H(A，B)＝max(h(A，B)，h(B，A)) (19)
其中
步骤6.3，统计步骤6.2中计算得到的所有豪斯多夫距离最小值，若该值小于设定阈值则认为最小豪斯多夫距离所对应的建筑物为该建筑物的同名建筑物；
步骤6.4，经过步骤6.1-6.3在相邻航带中获取多个同名建筑物，对每一个同名建筑物，分别进行步骤6.5-6.7，获取得到同名建筑物；
步骤6.5，使用由粗到精的建筑物屋面分割的方法提取出同名建筑物的屋面，剔除点数少于一定数量的屋面；
步骤6.6，对于建筑物中一屋面片，计算该屋面与同名建筑物中所有屋面的豪斯多夫距离，若最小豪斯多夫距离对应的面片距离以及两个平面法向量之间的夹角，若均小于设定的阈值则认为该面片为同名面片；
步骤6.7，对剩下的面片重复步骤6.6，直至该同名建筑物中所有同名平面被探测出(如图6所示)。
如上所述的步骤7，利用探测得到的同名屋面进行航带平差计算的方法为：
步骤7.1，依据航带平差采用的刚体变换模型，列出转换模型，
X1＝T+RX2 (20)
其中X1为一同名面上点集，X2为对应面片上点集，为旋转矩阵，为平移矩阵；
步骤7.2，对模型进行线性化；
其中Δx＝dx0+dxΔy＝dy0+dyΔz＝dz0+dz
步骤7.3，构建误差方程；
V＝BX-L (23)
其中B为系数矩阵，L为近似值；
步骤7.4，迭代求解该次改正数，当改正数小于设定阈值或者迭代次数达到一定限制后迭代终止：
X＝(BTPB)-1BTPL (24)
其中P为权重矩阵；
步骤7.5，根据公式计算精度指标，单位权中误差σ0；
其中n为所有的观测方程，t为必要观测数；
步骤7.6，根据求得的转换参数，将一航带点集转换到同一坐标系下，即完成航带平差。
本发明中所描述的具体实施例仅仅是对本发明精神作举例说明。本发明所属技术领域的技术人员可以对所描述的具体实施例做各种各样的修改或补充或采用类似的方式替代，但并不会偏离本发明的精神或者超越所附权利要求书所定义的范围。
Claims (1)
1.一种最小豪斯多夫距离面特征约束的点云航带平差方法，通过引入最小豪斯多夫距离实现点云航带平差，消除航带间偏差，其特征在于包括以下步骤：
步骤1：引入多条机载激光雷达航带点云数据；
步骤2：将每条机载激光雷达航带点云数据使用渐进三角网加密滤波算法进行滤波处理，获取每条机载激光雷达航带点云数据中非地面点云以及每条机载激光雷达航带点云数据中地面点云；
步骤3：将每条机载激光雷达航带点云数据中非地面点云通过半抑制模糊C均值聚类方法得到每条机载激光雷达航带点云数据中多个初始建筑物种子点，将每条机载激光雷达航带点云数据中多个初始建筑物种子点通过区域增长聚类方法，得到每条机载激光雷达航带点云数据中非地面点云中的多个建筑物点；
步骤4：搜索每条机载激光雷达航带点云数据中非地面点云中的多个建筑物点的X轴坐标最大值、X轴坐标最小值、Y轴坐标最大值、Y轴坐标最小值，根据X轴与Y轴最大值与最小值确定两边平行于X轴与Y轴的包含每条机载激光雷达航带的最小外包矩形，将相邻两条机载激光雷达航带的最小外包矩形进行相交处理得到相邻两条机载激光雷达航带的重叠区域，将相邻两条机载激光雷达航带的重叠区域进行格网化处理得到相邻两条机载激光雷达航带的重叠区域二值图像，将相邻两条机载激光雷达航带的重叠区域二值图像通过四连通区域标签法进行处理得到相邻两条机载激光雷达航带的重叠区域多个图像块，在相邻两条机载激光雷达航带的重叠区域每个图像块中搜索步骤3所述的每条机载激光雷达航带点云数据中非地面点云中的多个建筑物点，得到相邻两条机载激光雷达航带的重叠区域每个图像块中多个建筑物点；
步骤5：使用主成分分析法计算相邻两条机载激光雷达航带的重叠区域每个图像块中每个建筑物点的法向量、每个建筑物点的平面、每个建筑物点的曲率，将相邻两条机载激光雷达航带的重叠区域每个图像块中建筑物点的曲率最小的建筑物点作为相邻两条机载激光雷达航带的重叠区域对应图像块的种子点，将相邻两条机载激光雷达航带的重叠区域每个图像块的除种子点外的建筑物点构建相邻两条机载激光雷达航带的重叠区域每个图像块的非种子点集，使用点到平面计算方法与向量间夹角计算方法分别计算相邻两条机载激光雷达航带的重叠区域中每个图像块的非种子点集中每个建筑物点与种子点的法向量夹角、相邻两条机载激光雷达航带的重叠区域中每个图像块的非种子点集中每个建筑物点到种子点的平面的距离，将相邻两条机载激光雷达航带的重叠区域每个图像块的非种子点集中建筑物点与种子点的法向量夹角大于等于夹角阈值或到种子点的平面的距离大于等于距离阈值的建筑物点筛选出来得到相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集，将相邻两条机载激光雷达航带的重叠区域每个图像块的非种子点集中建筑物点与种子点的法向量夹角小于夹角阈值且到种子点的平面的距离小于距离阈值的建筑物点筛选出来得到相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物平面点集，进一步与相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物平面点集通过最小二乘方法构建初始平面，使用向量间夹角计算方法计算相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集每个建筑物点到种子点的法向量夹角，使用点到平面距离计算方法计算相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集每个建筑物点到初始平面的距离，将相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集每个建筑物点与种子点的法向量夹角小于距离阈值且到初始平面的距离小于距离阈值的建筑物点筛选出来与建筑物平面点集、种子点共同构成相邻两条机载激光雷达航带的重叠区域每个图像块的最终建筑物平面点集，将相邻两条机载激光雷达航带的重叠区域每个图像块的建筑物非平面点集每个建筑物点到种子点的法向量夹角大于等于距离阈值或者到初始平面的距离大于等于距离阈值的建筑物点筛选出来构建相邻两条机载激光雷达航带的重叠区域每个图像块的最终非建筑物平面点集；
步骤5所述相邻两条机载激光雷达航带的重叠区域每个图像块的最终建筑物平面点集中每个建筑物点均属于初始平面；
步骤5所述相邻两条机载激光雷达航带的重叠区域每个图像块的最终非建筑物平面点集每个建筑物点均有对应的平面；
步骤6：对于一条机载激光雷达航带的重叠区域每个图像块的建筑物点集，使用豪斯多夫距离计算方法计算其到相邻机载激光雷达航带的重叠区域每个图像块的建筑物点集的豪斯多夫距离，若一条机载激光雷达航带的重叠区域的一个图像块的建筑物点集到相邻机载激光雷达航带的重叠区域中每个图像块的建筑物点集的最小豪斯多夫距离小于阈值，将一条机载激光雷达航带的重叠区域建筑物点集与最小豪斯多夫距离对应的相邻机载激光雷达航带的重叠区域图像块的建筑物点集作为一对同名建筑物点集，使用豪斯多夫距离计算方法与向量间夹角计算方法计算一对同名建筑物点集中一栋建筑物点集分割对应的平面到一对同名建筑物点集中另一栋建筑物点集中每个面片的豪斯多夫距离以及法向量夹角，若一对同名建筑物点集中一栋建筑物点集分割对应的平面到最小豪斯多夫距离对应的一对同名建筑物点集中另一栋建筑物点集分割的面片对应的法向量夹角小于夹角阈值，则将一对同名建筑物点集中一栋建筑物点集分割对应的平面与一对同名建筑物点集中另一栋建筑物点集分割的面片作为同名面片；
步骤7：以相邻航带探测得到的同名面片方程以及同名面片对应的点云数据，根据同名面片中一个面片中点云数据满足同名面片中另一面片方程构造误差方程，根据最小二乘原理解算转换参数，以相邻航带中一条航带为基准，使用解算的转换参数将相邻航带中另一条航带中点云数据进行旋转与平移实现将相邻航带点云数据转换到同一坐标系下。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210014763.6A CN114463523B (zh) | 2022-01-07 | 2022-01-07 | 一种最小豪斯多夫距离面特征约束的点云航带平差方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210014763.6A CN114463523B (zh) | 2022-01-07 | 2022-01-07 | 一种最小豪斯多夫距离面特征约束的点云航带平差方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN114463523A CN114463523A (zh) | 2022-05-10 |
| CN114463523B true CN114463523B (zh) | 2024-08-20 |
Family
ID=81409435
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202210014763.6A Active CN114463523B (zh) | 2022-01-07 | 2022-01-07 | 一种最小豪斯多夫距离面特征约束的点云航带平差方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114463523B (zh) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115752448B (zh) * | 2022-11-03 | 2024-05-10 | 武汉际上导航科技有限公司 | 结合点云匹配和传感器数据的激光点云航带平差方法 |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110992341A (zh) * | 2019-12-04 | 2020-04-10 | 沈阳建筑大学 | 一种基于分割的机载LiDAR点云建筑物提取方法 |
| CN113658190A (zh) * | 2021-06-29 | 2021-11-16 | 桂林理工大学 | 张量投票的面特征航带平差方法 |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US8224097B2 (en) * | 2008-06-12 | 2012-07-17 | Sri International | Building segmentation for densely built urban regions using aerial LIDAR data |
-
2022
- 2022-01-07 CN CN202210014763.6A patent/CN114463523B/zh active Active
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110992341A (zh) * | 2019-12-04 | 2020-04-10 | 沈阳建筑大学 | 一种基于分割的机载LiDAR点云建筑物提取方法 |
| CN113658190A (zh) * | 2021-06-29 | 2021-11-16 | 桂林理工大学 | 张量投票的面特征航带平差方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN114463523A (zh) | 2022-05-10 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN111527467B (zh) | 使用机器学习、图像分析和/或计算机视觉自动限定计算机辅助设计文件的方法和装置 | |
| US8427505B2 (en) | Geospatial modeling system for images and related methods | |
| CN104123730B (zh) | 基于道路特征的遥感影像与激光点云配准方法及系统 | |
| CN110006408B (zh) | LiDAR数据“云控制”航空影像摄影测量方法 | |
| Zeybek | Accuracy assessment of direct georeferencing UAV images with onboard global navigation satellite system and comparison of CORS/RTK surveying methods | |
| CN108428220B (zh) | 静止轨道卫星序列遥感影像海岛礁区域自动几何校正方法 | |
| CN113593017A (zh) | 露天矿地表三维模型构建方法、装置、设备及存储介质 | |
| US20160321818A1 (en) | Geospatial imaging system providing segmentation and classification features and related methods | |
| CN114755661A (zh) | 一种移动激光扫描系统参数标定方法及标定装置 | |
| Salah et al. | Evaluation of the self‐organizing map classifier for building detection from lidar data and multispectral aerial images | |
| CN117387580B (zh) | 一种基于倾斜摄影大比例尺地形图的测绘方法及系统 | |
| CN114463523B (zh) | 一种最小豪斯多夫距离面特征约束的点云航带平差方法 | |
| CN117078870A (zh) | 融合高精地图与激光稀疏点云的道路环境三维重建方法 | |
| Yu et al. | Automatic extrinsic self-calibration of mobile LiDAR systems based on planar and spherical features | |
| Lu et al. | A lightweight real-time 3D LiDAR SLAM for autonomous vehicles in large-scale urban environment | |
| Belfiore et al. | Orthorectification and pan-sharpening of worldview-2 satellite imagery to produce high resolution coloured ortho-photos | |
| CN117029870A (zh) | 一种基于路面点云的激光里程计 | |
| CN112561981A (zh) | 一种融合影像信息的摄影测量点云滤波方法 | |
| CN110148218B (zh) | 一种大批量机载LiDAR点云数据整体优化的方法 | |
| CN113516059B (zh) | 固体废弃物的识别方法、装置、电子设备及存储介质 | |
| Sefercik et al. | Point-based and model-based geolocation analysis of airborne laser scanning data | |
| Klapa et al. | Charting topographic maps based on UAV data using the image classification method | |
| Ye et al. | Gaussian mixture model of ground filtering based on hierarchical curvature constraints for airborne lidar point clouds | |
| Zhang | Photogrammetric point clouds: quality assessment, filtering, and change detection | |
| Huang et al. | Ground filtering algorithm for mobile LIDAR using order and neighborhood point information |
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 |