CN114998419B - 一种基于地形点云的沟壑体积变化切片方法 - Google Patents

一种基于地形点云的沟壑体积变化切片方法 Download PDF

Info

Publication number
CN114998419B
CN114998419B CN202210923370.7A CN202210923370A CN114998419B CN 114998419 B CN114998419 B CN 114998419B CN 202210923370 A CN202210923370 A CN 202210923370A CN 114998419 B CN114998419 B CN 114998419B
Authority
CN
China
Prior art keywords
point cloud
point
slice
area
section block
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
CN202210923370.7A
Other languages
English (en)
Other versions
CN114998419A (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.)
Xian University of Science and Technology
Original Assignee
Xian University of Science 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 Xian University of Science and Technology filed Critical Xian University of Science and Technology
Priority to CN202210923370.7A priority Critical patent/CN114998419B/zh
Publication of CN114998419A publication Critical patent/CN114998419A/zh
Application granted granted Critical
Publication of CN114998419B publication Critical patent/CN114998419B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/60—Analysis of geometric attributes
    • G06T7/62—Analysis of geometric attributes of area, perimeter, diameter or volume
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
    • G01B11/00—Measuring arrangements characterised by the use of optical techniques
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C7/00—Tracing profiles
    • G01C7/02—Tracing profiles of land surfaces
    • 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/05—Geographic models
    • 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/30—Polynomial surface description
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T3/00—Geometric image transformations in the plane of the image
    • G06T3/06—Topological mapping of higher dimensional structures onto lower dimensional surfaces
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T3/00—Geometric image transformations in the plane of the image
    • G06T3/40—Scaling of whole images or parts thereof, e.g. expanding or contracting
    • G06T3/4007—Scaling of whole images or parts thereof, e.g. expanding or contracting based on interpolation, e.g. bilinear interpolation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V10/00—Arrangements for image or video recognition or understanding
    • G06V10/70—Arrangements for image or video recognition or understanding using pattern recognition or machine learning
    • G06V10/762—Arrangements for image or video recognition or understanding using pattern recognition or machine learning using clustering, e.g. of similar faces in social networks
    • G06V10/763—Non-hierarchical techniques, e.g. based on statistics of modelling distributions
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10028—Range image; Depth image; 3D point clouds
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10032—Satellite or aerial image; Remote sensing
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30181—Earth observation
    • G06T2207/30184—Infrastructure
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Software Systems (AREA)
  • Geometry (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Remote Sensing (AREA)
  • Computer Graphics (AREA)
  • Multimedia (AREA)
  • Computing Systems (AREA)
  • Probability & Statistics with Applications (AREA)
  • Health & Medical Sciences (AREA)
  • Databases & Information Systems (AREA)
  • Evolutionary Computation (AREA)
  • General Health & Medical Sciences (AREA)
  • Medical Informatics (AREA)
  • Artificial Intelligence (AREA)
  • Algebra (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Physics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Length Measuring Devices With Unspecified Measuring Means (AREA)
  • Length Measuring Devices By Optical Means (AREA)

Abstract

本发明公开了一种基于地形点云的沟壑体积变化切片方法，该方法包括以下步骤：一、待测沟壑区域的点云获取；二、待测沟壑区域的点云的切片；三、待测沟壑区域的点云切片的聚类及插值处理；四、获取后一期点云切片相对前一期点云切片的切面变化面积；五、待测沟壑区域的体积变化量获取。本发明方法步骤简单，设计合理，获取待测沟壑区域沿高程方向各个点云切片的体积变化量，进而累加得到待测沟壑区域整体的体积变化量，提高了沟壑体积变化量计算的准确性，以适应沟壑地形。

Description

一种基于地形点云的沟壑体积变化切片方法
技术领域
本发明属于沟壑体积变化测量技术领域，尤其是涉及一种基于地形点云的沟壑体积变化切片方法。
背景技术
高精度地形变化监测为沟壑区域变化、地质灾害、地表沉降等过程研究提供了有效手段。随着遥感技术的快速发展，基于三维激光扫描等获取的三维点云已成为高精度地形变化监测的重要数据源。基于点云的地形变化监测方法可将获取的地形变化点云转换为体积变化量，因此，地形变化点云向体积变化量的转换是提升侵蚀量化精度的关键环节。传统二维方法将地形变化点云直接转换为规则的二维栅格，栅格属性值为其所包含变化点云值的函数，栅格面积与乘栅格属性值相乘得到栅格对应的体积变化量，将不同栅格的体积变化量进行加和可得研究区的体积变化。该方法在地形平缓区域适用性较好，然而对于地形陡峭区域，如沟壑区域，其适用性受到挑战。
因此，需要一种基于地形点云的沟壑体积变化切片方法，获取待测沟壑区域沿高程方向各个点云切片的体积变化量，进而累加得到待测沟壑区域整体的体积变化量，提高了沟壑体积变化量计算的准确性，以适应沟壑地形。
发明内容
本发明所要解决的技术问题在于针对上述现有技术中的不足，提供一种基于地形点云的沟壑体积变化切片方法，其方法步骤简单，设计合理且实现方便，获取待测沟壑区域沿高程方向各个点云切片的体积变化量，进而累加得到待测沟壑区域整体的体积变化量，提高了沟壑体积变化量计算的准确性，以适应沟壑地形。
为解决上述技术问题，本发明采用的技术方案是：一种基于地形点云的沟壑体积变化切片方法，其特征在于，该方法包括以下步骤：
步骤一、待测沟壑区域的点云获取：
步骤101、建立空间直角坐标系；其中，以待测沟壑区域外侧左下方的稳定区任意一点为原点O，以过原点O且沿架站区连线为X轴，过原点O且与X轴垂直的方向为Y轴，以过原点O垂直于由X轴和Y轴形成的平面OXY且沿高程方向为Z轴，建立空间直角坐标系OXYZ；
步骤102、采用地面三维激光扫描仪对待测沟壑区域进行扫描，获取第1期的激光点云，...，第
Figure 309312DEST_PATH_IMAGE001
期的激光点云，...，第
Figure 498985DEST_PATH_IMAGE002
期的激光点云；其中，
Figure 398808DEST_PATH_IMAGE003
和
Figure 964919DEST_PATH_IMAGE004
均为正整数，
Figure 204270DEST_PATH_IMAGE005
，
Figure 564844DEST_PATH_IMAGE006
；
步骤二、待测沟壑区域的点云的切片：
在OXYZ坐标系下，采用计算机按照相邻两期的点云切片厚度
Figure 951963DEST_PATH_IMAGE007
，将第
Figure 321765DEST_PATH_IMAGE003
期的激光点云沿Z轴方向做多个切平面并将多个切平面从下至上进行排序依次记作第1个切平面，...，第
Figure 212360DEST_PATH_IMAGE008
个切平面，...，第
Figure 697831DEST_PATH_IMAGE009
个切平面；其中，
Figure 306666DEST_PATH_IMAGE010
和
Figure 480159DEST_PATH_IMAGE009
均为正整数，且
Figure 490840DEST_PATH_IMAGE011
；任一个切平面平行平面OXY，第1个切平面和第
Figure 396479DEST_PATH_IMAGE012
个切平面之间相邻两个切平面的间距为
Figure 492611DEST_PATH_IMAGE013
，第1个切平面距离第
Figure 204215DEST_PATH_IMAGE014
期的激光点云的最小Z轴坐标处的间距为
Figure 69403DEST_PATH_IMAGE015
，第
Figure 208260DEST_PATH_IMAGE016
个切平面与第
Figure 978639DEST_PATH_IMAGE017
期的激光点云的最大Z轴坐标之间的间距记作
Figure 228355DEST_PATH_IMAGE018
,且
Figure 213628DEST_PATH_IMAGE019
不大于
Figure 523387DEST_PATH_IMAGE020
；
步骤三、待测沟壑区域的点云切片的聚类及插值处理：
步骤301、当
Figure 531794DEST_PATH_IMAGE021
处于1～
Figure 585201DEST_PATH_IMAGE022
时，采用计算机将第
Figure 424981DEST_PATH_IMAGE021
个切平面上方及下方Z轴方向
Figure 905641DEST_PATH_IMAGE023
范围内的点云投影至第
Figure 198082DEST_PATH_IMAGE021
个切平面，得到第
Figure 764103DEST_PATH_IMAGE017
期第
Figure 458389DEST_PATH_IMAGE021
个点云切片；
当
Figure 109951DEST_PATH_IMAGE021
取
Figure 155267DEST_PATH_IMAGE016
时，采用计算机将第
Figure 488159DEST_PATH_IMAGE016
个切平面上方Z轴方向的剩余点云及下方Z轴方向
Figure 36952DEST_PATH_IMAGE023
范围内的点云投影至第
Figure 124994DEST_PATH_IMAGE016
个切平面，得到第
Figure 392027DEST_PATH_IMAGE017
期第
Figure 777878DEST_PATH_IMAGE016
个点云切片；
步骤302、采用计算机利用FCM聚类算法对第
Figure 181178DEST_PATH_IMAGE017
期第
Figure 440121DEST_PATH_IMAGE024
个点云切片上的点进行聚类，将第
Figure 194450DEST_PATH_IMAGE017
期第
Figure 931462DEST_PATH_IMAGE024
个点云切片上各个点划分为
Figure 392530DEST_PATH_IMAGE025
个类簇，并得到
Figure 822375DEST_PATH_IMAGE026
个类簇的聚类中心点；其中，
Figure 798421DEST_PATH_IMAGE025
为正整数；
步骤303、采用计算机利用三次B样条插值函数对第
Figure 339124DEST_PATH_IMAGE017
期第
Figure 671010DEST_PATH_IMAGE024
个点云切片上
Figure 6177DEST_PATH_IMAGE025
个聚类中心点进行插值处理，得到第
Figure 735098DEST_PATH_IMAGE017
期第
Figure 813913DEST_PATH_IMAGE024
个点云切片上插值点集
Figure 311890DEST_PATH_IMAGE027
；
步骤四、获取后一期点云切片相对前一期点云切片的切面变化面积：
步骤401、按照步骤二至步骤三所述的方法，采用计算机对第
Figure 755641DEST_PATH_IMAGE028
期的激光点云进行处理，得到第
Figure 971859DEST_PATH_IMAGE028
期第
Figure 854364DEST_PATH_IMAGE029
个点云切片上插值点集
Figure 206848DEST_PATH_IMAGE030
；
步骤402、采用计算机将第
Figure 336347DEST_PATH_IMAGE031
期第
Figure 774281DEST_PATH_IMAGE032
个点云切片上插值点集
Figure 194898DEST_PATH_IMAGE033
和第
Figure 667468DEST_PATH_IMAGE034
期第
Figure 515338DEST_PATH_IMAGE035
个点云切片上插值点集
Figure 378252DEST_PATH_IMAGE036
同步绘制在平面OXY上，且
Figure 602560DEST_PATH_IMAGE037
和
Figure 929636DEST_PATH_IMAGE038
中相交的相邻两个交点将
Figure 213987DEST_PATH_IMAGE039
和
Figure 312000DEST_PATH_IMAGE040
之间围设的区域划分为第1个切面区块，..，第
Figure 74419DEST_PATH_IMAGE041
个切面区块，..，第
Figure 256002DEST_PATH_IMAGE042
个切面区块，并得到各个切面区块的面积；其中，第
Figure 711254DEST_PATH_IMAGE043
个点云切片上第
Figure 611077DEST_PATH_IMAGE044
个切面区块的面积记作
Figure 114871DEST_PATH_IMAGE045
；其中，
Figure 150960DEST_PATH_IMAGE046
和
Figure 777113DEST_PATH_IMAGE047
均为正整数，且
Figure 164232DEST_PATH_IMAGE048
，
Figure 455405DEST_PATH_IMAGE049
表示第
Figure 611580DEST_PATH_IMAGE050
个点云切片上切面区块的总数；
步骤五、待测沟壑区域的体积变化量获取：
步骤501、对第
Figure 408635DEST_PATH_IMAGE050
个点云切片上第
Figure 17470DEST_PATH_IMAGE046
个切面区块进行判断，如果第
Figure 190963DEST_PATH_IMAGE050
个点云切片上第
Figure 139327DEST_PATH_IMAGE046
个切面区块中第
Figure 107283DEST_PATH_IMAGE051
期插值点的Y轴坐标大于第
Figure 203415DEST_PATH_IMAGE052
期插值点的Y轴坐标，则第
Figure 649440DEST_PATH_IMAGE046
个切面区块为负地形区，则第
Figure 468623DEST_PATH_IMAGE050
个点云切片上第
Figure 607480DEST_PATH_IMAGE046
个切面区块的面积
Figure 190908DEST_PATH_IMAGE053
取正值；如果第
Figure 440624DEST_PATH_IMAGE050
个点云切片上第
Figure 425897DEST_PATH_IMAGE046
个切面区块中第
Figure 673339DEST_PATH_IMAGE054
期插值点的Y轴坐标小于等于第
Figure 744063DEST_PATH_IMAGE052
期插值点的Y轴坐标，则第
Figure 797470DEST_PATH_IMAGE050
个点云切片上第
Figure 637250DEST_PATH_IMAGE046
个切面区块为正地形区，则第
Figure 304860DEST_PATH_IMAGE050
个点云切片上第
Figure 597302DEST_PATH_IMAGE046
个切面区块的面积
Figure 454399DEST_PATH_IMAGE055
取负值；
步骤502、根据第
Figure 148686DEST_PATH_IMAGE050
个点云切片上第
Figure 800247DEST_PATH_IMAGE046
个切面区块的面积，获取第
Figure 783246DEST_PATH_IMAGE050
个点云切片上第
Figure 178456DEST_PATH_IMAGE046
个切面区块的变化体积，进而得到第
Figure 727249DEST_PATH_IMAGE050
个点云切片的变化体积；
步骤503、将各个点云切片的变化体积累加得到待测沟壑区域第
Figure 815290DEST_PATH_IMAGE056
期相对第
Figure 767809DEST_PATH_IMAGE057
期的体积变化量
Figure 701130DEST_PATH_IMAGE058
。
上述的一种基于地形点云的沟壑体积变化切片方法，其特征在于：步骤二中相邻两期的点云切片厚度
Figure 370009DEST_PATH_IMAGE059
的获取，具体过程如下：
步骤201、采用计算机从第
Figure 628952DEST_PATH_IMAGE060
期的激光点云中随机选择点云组成第
Figure 117702DEST_PATH_IMAGE060
个点云集合
Figure 57976DEST_PATH_IMAGE061
；其中，第
Figure 581362DEST_PATH_IMAGE062
个点云集合
Figure 745627DEST_PATH_IMAGE063
中第
Figure 987252DEST_PATH_IMAGE065
个点记作
Figure 714906DEST_PATH_IMAGE066
，
Figure 92797DEST_PATH_IMAGE067
，
Figure 693543DEST_PATH_IMAGE068
表示第
Figure 156885DEST_PATH_IMAGE069
个点云集合
Figure 501279DEST_PATH_IMAGE070
的总数；
步骤202、采用计算机获取第
Figure 936940DEST_PATH_IMAGE062
个点云集合
Figure 443007DEST_PATH_IMAGE070
中第
Figure 659225DEST_PATH_IMAGE065
个点
Figure 541730DEST_PATH_IMAGE071
与第
Figure 582630DEST_PATH_IMAGE072
期的激光点云中各个点的欧式距离，并将各个欧式距离按照从大到小的顺序排序，获取后
Figure 259599DEST_PATH_IMAGE073
个欧式距离；其中，与
Figure 697533DEST_PATH_IMAGE074
对应的后
Figure 383730DEST_PATH_IMAGE075
个欧式距离中第
Figure 528403DEST_PATH_IMAGE076
个欧式距离记作
Figure 641853DEST_PATH_IMAGE077
；
步骤203、根据公式
Figure 567083DEST_PATH_IMAGE078
，得到第
Figure 791391DEST_PATH_IMAGE079
期的激光点云的平均点间距
Figure 118468DEST_PATH_IMAGE080
；
步骤204、按照步骤201至步骤203所述的方法，得到第
Figure 589769DEST_PATH_IMAGE081
期的激光点云的平均点间距
Figure 736717DEST_PATH_IMAGE082
；并对
Figure 764715DEST_PATH_IMAGE083
和
Figure 946298DEST_PATH_IMAGE084
进行平均值处理，得到相邻两期的平均点间距
Figure 339233DEST_PATH_IMAGE085
；
步骤205、根据公式
Figure 239056DEST_PATH_IMAGE086
，得到相邻两期的点云切片厚度
Figure 805167DEST_PATH_IMAGE087
；其中，
Figure 841256DEST_PATH_IMAGE088
为常数且
Figure 467409DEST_PATH_IMAGE089
取值为0.5。
上述的一种基于地形点云的沟壑体积变化切片方法，其特征在于：步骤402中得到各个切面区块的面积的方法均相同，其中，得到第
Figure 274435DEST_PATH_IMAGE090
个点云切片上第
Figure 644236DEST_PATH_IMAGE091
个切面区块的面积
Figure 800411DEST_PATH_IMAGE092
，具体过程如下：
步骤4021、对第
Figure 597466DEST_PATH_IMAGE090
个点云切片上第
Figure 143985DEST_PATH_IMAGE091
个切面区块中包含的第
Figure 51898DEST_PATH_IMAGE093
期插值点拟合，得到第
Figure 62579DEST_PATH_IMAGE090
个点云切片上第
Figure 30535DEST_PATH_IMAGE091
个切面区块中第
Figure 126667DEST_PATH_IMAGE093
期的拟合函数段
Figure 25222DEST_PATH_IMAGE094
；其中，
Figure 155989DEST_PATH_IMAGE095
为自变量；
对第
Figure 294846DEST_PATH_IMAGE096
个点云切片上第
Figure 612695DEST_PATH_IMAGE097
个切面区块中包含的第
Figure 65673DEST_PATH_IMAGE098
期插值点拟合，得到第
Figure 50947DEST_PATH_IMAGE099
个点云切片上第
Figure 360705DEST_PATH_IMAGE100
个切面区块中第
Figure 165850DEST_PATH_IMAGE098
期的拟合函数段
Figure 219257DEST_PATH_IMAGE101
；
步骤4022、采用计算机根据公式
Figure 747453DEST_PATH_IMAGE102
，得到第
Figure 228113DEST_PATH_IMAGE103
个点云切片上第
Figure 786133DEST_PATH_IMAGE104
个切面区块的面积
Figure 377651DEST_PATH_IMAGE105
；其中，
Figure 275200DEST_PATH_IMAGE106
表示关于X轴坐标的积分，且
Figure 926761DEST_PATH_IMAGE107
的取值范围为
Figure 706498DEST_PATH_IMAGE108
，
Figure 101708DEST_PATH_IMAGE109
表示第
Figure 916080DEST_PATH_IMAGE110
个点云切片上第
Figure 925493DEST_PATH_IMAGE111
个切面区块的最小X轴坐标，
Figure 192526DEST_PATH_IMAGE112
表示第
Figure 391427DEST_PATH_IMAGE113
个点云切片上第
Figure 60305DEST_PATH_IMAGE111
个切面区块的最大X轴坐标；
Figure 991352DEST_PATH_IMAGE114
表示绝对值符合。
上述的一种基于地形点云的沟壑体积变化切片方法，其特征在于：步骤502中根据第
Figure 745682DEST_PATH_IMAGE113
个点云切片上第
Figure 748273DEST_PATH_IMAGE115
个切面区块的面积，获取第
Figure 271658DEST_PATH_IMAGE113
个点云切片上第
Figure 435923DEST_PATH_IMAGE115
个切面区块的变化体积，进而得到第
Figure 363034DEST_PATH_IMAGE113
个点云切片的变化体积，具体过程如下：
步骤5021、当
Figure 903737DEST_PATH_IMAGE113
处于1～
Figure 281629DEST_PATH_IMAGE116
时，则根据公式
Figure 882374DEST_PATH_IMAGE117
，得到第
Figure 283400DEST_PATH_IMAGE113
个点云切片上第
Figure 627794DEST_PATH_IMAGE118
个切面区块的变化体积
Figure 860192DEST_PATH_IMAGE119
；
当
Figure 631839DEST_PATH_IMAGE120
取
Figure 769428DEST_PATH_IMAGE121
时，则根据公式
Figure 651933DEST_PATH_IMAGE122
，得到第
Figure 4417DEST_PATH_IMAGE123
个点云切片上第
Figure 946965DEST_PATH_IMAGE124
个切面区块的变化体积
Figure 384900DEST_PATH_IMAGE125
；
步骤5022、当
Figure 8779DEST_PATH_IMAGE127
处于1～
Figure 215770DEST_PATH_IMAGE128
时，采用计算机根据
Figure 329219DEST_PATH_IMAGE129
，得到第
Figure 254450DEST_PATH_IMAGE130
个点云切片的变化体积
Figure 167173DEST_PATH_IMAGE131
；
当
Figure 494250DEST_PATH_IMAGE132
取
Figure 513021DEST_PATH_IMAGE133
时，采用计算机根据
Figure 925548DEST_PATH_IMAGE134
，得到第
Figure 953547DEST_PATH_IMAGE135
个点云切片的变化体积
Figure 72812DEST_PATH_IMAGE136
；
步骤503中将各个点云切片的变化体积累加得到待测沟壑区域第
Figure 528065DEST_PATH_IMAGE137
期相对第
Figure 162308DEST_PATH_IMAGE138
期的体积变化量
Figure 728419DEST_PATH_IMAGE139
，具体过程如下：
采用计算机根据
Figure 217038DEST_PATH_IMAGE140
，得到待测沟壑区域第
Figure 843191DEST_PATH_IMAGE141
期相对第
Figure 964731DEST_PATH_IMAGE142
期的体积变化量
Figure 334532DEST_PATH_IMAGE143
。
本发明与现有技术相比具有以下优点：
1、本发明方法步骤简单、设计合理且实现方便，精度高。
2、本发明采用地面三维激光扫描仪对待测沟壑区域进行扫描，获取多期激光点云，便于后续基于激光点云进行相邻两期沟壑体积变化量计算。
3、本发明对任一期的激光点云沿Z轴高程方向做切平面，并将切平面高程方向1/2点云切片厚度范围内的点云投影至切平面上，形成点云切片，从而将地面变形点云的三维问题简化为二维问题，可以更加精确的区分相同平面OXY位置但位于不同高程位置的点云之间的变化情况，进而有效获取和区分每个点云切片的正负地形变化量，实现复杂沟壑地形中各个点云切片体积变化的准确量化。
4、本发明采用FCM聚类算法对点云切片上的点进行聚类，是为了去z噪便于后续进行插值处理，得到各期点云切片的插值点集，从而便于将相邻两期之间围设的区域划分为各个切面区块。
5、本发明根据各个切面区块的面积，获取各个切面区块的变化体积，进而得到各个点云切片的变化体积，最后将各个点云切片的变化体积累加得到待测沟壑区域两期间的体积变化量，有效地适应沟壑地形。
6、本发明在获取各个切面区块的变化体积时，先对各个切面区块进行负地形区和正地形区的判断，从而根据负地形区的面积和正地形区的面积获取正地形区的体积变化量和负地形区的体积变化量，进而得到点云切片的变化体积，适应了沟壑地形的实际变化状态。
综上所述，本发明方法步骤简单，设计合理且实现方便，获取待测沟壑区域沿高程方向各个点云切片的体积变化量，进而累加得到待测沟壑区域整体的体积变化量，提高了沟壑体积变化量计算的准确性，以适应沟壑地形。
下面通过附图和实施例，对本发明的技术方案做进一步的详细描述。
附图说明
图1为本发明的方法流程图。
图2为本发明待测沟壑区域、基准点和架站区的结构示意图。
具体实施方式
如图1和图2所示的一种基于地形点云的沟壑体积变化切片方法，包括以下步骤：
该方法包括以下步骤：
步骤一、待测沟壑区域的点云获取：
步骤101、建立空间直角坐标系；其中，以待测沟壑区域外侧左下方的稳定区任意一点为原点O，以过原点O且沿架站区连线为X轴，过原点O且与X轴垂直的方向为Y轴，以过原点O垂直于由X轴和Y轴形成的平面OXY且沿高程方向为Z轴，建立空间直角坐标系OXYZ；
步骤102、采用地面三维激光扫描仪对待测沟壑区域进行扫描，获取第1期的激光点云，...，第
Figure 490707DEST_PATH_IMAGE001
期的激光点云，...，第
Figure 225445DEST_PATH_IMAGE002
期的激光点云；其中，
Figure 834281DEST_PATH_IMAGE003
和
Figure 742194DEST_PATH_IMAGE004
均为正整数，
Figure 752876DEST_PATH_IMAGE005
，
Figure 406317DEST_PATH_IMAGE006
；
步骤二、待测沟壑区域的点云的切片：
在OXYZ坐标系下，采用计算机按照相邻两期的点云切片厚度
Figure 502449DEST_PATH_IMAGE007
，将第
Figure 214053DEST_PATH_IMAGE003
期的激光点云沿Z轴方向做多个切平面并将多个切平面从下至上进行排序依次记作第1个切平面，...，第
Figure 79241DEST_PATH_IMAGE008
个切平面，...，第
Figure 218099DEST_PATH_IMAGE009
个切平面；其中，
Figure 739210DEST_PATH_IMAGE010
和
Figure 254505DEST_PATH_IMAGE009
均为正整数，且
Figure 974199DEST_PATH_IMAGE011
；任一个切平面平行平面OXY，第1个切平面和第
Figure 283958DEST_PATH_IMAGE012
个切平面之间相邻两个切平面的间距为
Figure 541632DEST_PATH_IMAGE013
，第1个切平面距离第
Figure 595039DEST_PATH_IMAGE014
期的激光点云的最小Z轴坐标处的间距为
Figure 434819DEST_PATH_IMAGE015
，第
Figure 915479DEST_PATH_IMAGE016
个切平面与第
Figure 473499DEST_PATH_IMAGE017
期的激光点云的最大Z轴坐标之间的间距记作
Figure 2701DEST_PATH_IMAGE018
,且
Figure 696987DEST_PATH_IMAGE019
不大于
Figure 614128DEST_PATH_IMAGE020
；
步骤三、待测沟壑区域的点云切片的聚类及插值处理：
步骤301、当
Figure 393865DEST_PATH_IMAGE021
处于1～
Figure 477490DEST_PATH_IMAGE022
时，采用计算机将第
Figure 291862DEST_PATH_IMAGE021
个切平面上方及下方Z轴方向
Figure 114324DEST_PATH_IMAGE023
范围内的点云投影至第
Figure 381358DEST_PATH_IMAGE021
个切平面，得到第
Figure 580258DEST_PATH_IMAGE017
期第
Figure 186820DEST_PATH_IMAGE021
个点云切片；
当
Figure 180184DEST_PATH_IMAGE021
取
Figure 934513DEST_PATH_IMAGE016
时，采用计算机将第
Figure 671525DEST_PATH_IMAGE016
个切平面上方Z轴方向的剩余点云及下方Z轴方向
Figure 381861DEST_PATH_IMAGE023
范围内的点云投影至第
Figure 811705DEST_PATH_IMAGE016
个切平面，得到第
Figure 53330DEST_PATH_IMAGE017
期第
Figure 328454DEST_PATH_IMAGE016
个点云切片；
步骤302、采用计算机利用FCM聚类算法对第
Figure 909608DEST_PATH_IMAGE017
期第
Figure 244774DEST_PATH_IMAGE024
个点云切片上的点进行聚类，将第
Figure 973696DEST_PATH_IMAGE017
期第
Figure 52510DEST_PATH_IMAGE024
个点云切片上各个点划分为
Figure 550488DEST_PATH_IMAGE025
个类簇，并得到
Figure 7621DEST_PATH_IMAGE026
个类簇的聚类中心点；其中，
Figure 958259DEST_PATH_IMAGE025
为正整数；
步骤303、采用计算机利用三次B样条插值函数对第
Figure 840765DEST_PATH_IMAGE017
期第
Figure 193248DEST_PATH_IMAGE024
个点云切片上
Figure 73480DEST_PATH_IMAGE025
个聚类中心点进行插值处理，得到第
Figure 511414DEST_PATH_IMAGE017
期第
Figure 932031DEST_PATH_IMAGE024
个点云切片上插值点集
Figure 404601DEST_PATH_IMAGE027
；
步骤四、获取后一期点云切片相对前一期点云切片的切面变化面积：
步骤401、按照步骤二至步骤三所述的方法，采用计算机对第
Figure 518051DEST_PATH_IMAGE028
期的激光点云进行处理，得到第
Figure 630232DEST_PATH_IMAGE028
期第
Figure 854540DEST_PATH_IMAGE029
个点云切片上插值点集
Figure 916037DEST_PATH_IMAGE030
；
步骤402、采用计算机将第
Figure 200388DEST_PATH_IMAGE031
期第
Figure 550597DEST_PATH_IMAGE032
个点云切片上插值点集
Figure 578596DEST_PATH_IMAGE033
和第
Figure 760179DEST_PATH_IMAGE034
期第
Figure 949852DEST_PATH_IMAGE035
个点云切片上插值点集
Figure 849675DEST_PATH_IMAGE036
同步绘制在平面OXY上，且
Figure 104201DEST_PATH_IMAGE037
和
Figure 405869DEST_PATH_IMAGE038
中相交的相邻两个交点将
Figure 766443DEST_PATH_IMAGE039
和
Figure 153562DEST_PATH_IMAGE040
之间围设的区域划分为第1个切面区块，..，第
Figure 461047DEST_PATH_IMAGE041
个切面区块，..，第
Figure 351643DEST_PATH_IMAGE042
个切面区块，并得到各个切面区块的面积；其中，第
Figure 148697DEST_PATH_IMAGE043
个点云切片上第
Figure 757533DEST_PATH_IMAGE044
个切面区块的面积记作
Figure 931026DEST_PATH_IMAGE045
；其中，
Figure 128658DEST_PATH_IMAGE046
和
Figure 96614DEST_PATH_IMAGE047
均为正整数，且
Figure 192746DEST_PATH_IMAGE048
，
Figure 904350DEST_PATH_IMAGE049
表示第
Figure 769537DEST_PATH_IMAGE050
个点云切片上切面区块的总数；
步骤五、待测沟壑区域的体积变化量获取：
步骤501、对第
Figure 846078DEST_PATH_IMAGE050
个点云切片上第
Figure 429506DEST_PATH_IMAGE046
个切面区块进行判断，如果第
Figure 679222DEST_PATH_IMAGE050
个点云切片上第
Figure 664495DEST_PATH_IMAGE046
个切面区块中第
Figure 659740DEST_PATH_IMAGE051
期插值点的Y轴坐标大于第
Figure 730464DEST_PATH_IMAGE052
期插值点的Y轴坐标，则第
Figure 783870DEST_PATH_IMAGE046
个切面区块为负地形区，则第
Figure 623650DEST_PATH_IMAGE050
个点云切片上第
Figure 104310DEST_PATH_IMAGE046
个切面区块的面积
Figure 334435DEST_PATH_IMAGE053
取正值；如果第
Figure 191532DEST_PATH_IMAGE050
个点云切片上第
Figure 885819DEST_PATH_IMAGE046
个切面区块中第
Figure 537380DEST_PATH_IMAGE054
期插值点的Y轴坐标小于等于第
Figure 769647DEST_PATH_IMAGE052
期插值点的Y轴坐标，则第
Figure 164856DEST_PATH_IMAGE050
个点云切片上第
Figure 713649DEST_PATH_IMAGE046
个切面区块为正地形区，则第
Figure 801691DEST_PATH_IMAGE050
个点云切片上第
Figure 6407DEST_PATH_IMAGE046
个切面区块的面积
Figure 205307DEST_PATH_IMAGE055
取负值；
步骤502、根据第
Figure 608607DEST_PATH_IMAGE050
个点云切片上第
Figure 867550DEST_PATH_IMAGE046
个切面区块的面积，获取第
Figure 621879DEST_PATH_IMAGE050
个点云切片上第
Figure 47307DEST_PATH_IMAGE046
个切面区块的变化体积，进而得到第
Figure 570692DEST_PATH_IMAGE050
个点云切片的变化体积；
步骤503、将各个点云切片的变化体积累加得到待测沟壑区域第
Figure 536DEST_PATH_IMAGE056
期相对第
Figure 976583DEST_PATH_IMAGE057
期的体积变化量
Figure 454968DEST_PATH_IMAGE058
。
本实施例中，步骤二中相邻两期的点云切片厚度
Figure 98439DEST_PATH_IMAGE059
的获取，具体过程如下：
步骤201、采用计算机从第
Figure 433606DEST_PATH_IMAGE060
期的激光点云中随机选择点云组成第
Figure 162527DEST_PATH_IMAGE060
个点云集合
Figure 241342DEST_PATH_IMAGE061
；其中，第
Figure 926270DEST_PATH_IMAGE062
个点云集合
Figure 432338DEST_PATH_IMAGE063
中第
Figure 648555DEST_PATH_IMAGE065
个点记作
Figure 531061DEST_PATH_IMAGE066
，
Figure 821228DEST_PATH_IMAGE067
，
Figure 763776DEST_PATH_IMAGE068
表示第
Figure 201711DEST_PATH_IMAGE069
个点云集合
Figure 622328DEST_PATH_IMAGE070
的总数；
步骤202、采用计算机获取第
Figure 803821DEST_PATH_IMAGE062
个点云集合
Figure 651691DEST_PATH_IMAGE070
中第
Figure 576922DEST_PATH_IMAGE065
个点
Figure 801229DEST_PATH_IMAGE071
与第
Figure 65989DEST_PATH_IMAGE072
期的激光点云中各个点的欧式距离，并将各个欧式距离按照从大到小的顺序排序，获取后
Figure 350340DEST_PATH_IMAGE073
个欧式距离；其中，与
Figure 762866DEST_PATH_IMAGE074
对应的后
Figure 525286DEST_PATH_IMAGE075
个欧式距离中第
Figure 706869DEST_PATH_IMAGE076
个欧式距离记作
Figure 349071DEST_PATH_IMAGE077
；
步骤203、根据公式
Figure 248894DEST_PATH_IMAGE078
，得到第
Figure 815005DEST_PATH_IMAGE079
期的激光点云的平均点间距
Figure 851094DEST_PATH_IMAGE080
；
步骤204、按照步骤201至步骤203所述的方法，得到第
Figure 414930DEST_PATH_IMAGE081
期的激光点云的平均点间距
Figure 802049DEST_PATH_IMAGE082
；并对
Figure 906272DEST_PATH_IMAGE083
和
Figure 62447DEST_PATH_IMAGE084
进行平均值处理，得到相邻两期的平均点间距
Figure 859501DEST_PATH_IMAGE085
；
步骤205、根据公式
Figure 156753DEST_PATH_IMAGE086
，得到相邻两期的点云切片厚度
Figure 330245DEST_PATH_IMAGE087
；其中，
Figure 340926DEST_PATH_IMAGE088
为常数且
Figure 308882DEST_PATH_IMAGE089
取值为0.5。
本实施例中，步骤402中得到各个切面区块的面积的方法均相同，其中，得到第
Figure 342697DEST_PATH_IMAGE090
个点云切片上第
Figure 788722DEST_PATH_IMAGE091
个切面区块的面积
Figure 919489DEST_PATH_IMAGE092
，具体过程如下：
步骤4021、对第
Figure 58347DEST_PATH_IMAGE090
个点云切片上第
Figure 641775DEST_PATH_IMAGE091
个切面区块中包含的第
Figure 78441DEST_PATH_IMAGE093
期插值点拟合，得到第
Figure 63715DEST_PATH_IMAGE090
个点云切片上第
Figure 373473DEST_PATH_IMAGE091
个切面区块中第
Figure 444197DEST_PATH_IMAGE093
期的拟合函数段
Figure 435287DEST_PATH_IMAGE094
；其中，
Figure 275067DEST_PATH_IMAGE095
为自变量；
对第
Figure 755727DEST_PATH_IMAGE096
个点云切片上第
Figure 48168DEST_PATH_IMAGE097
个切面区块中包含的第
Figure 590752DEST_PATH_IMAGE098
期插值点拟合，得到第
Figure 285038DEST_PATH_IMAGE099
个点云切片上第
Figure 936599DEST_PATH_IMAGE100
个切面区块中第
Figure 981916DEST_PATH_IMAGE098
期的拟合函数段
Figure 377125DEST_PATH_IMAGE101
；
步骤4022、采用计算机根据公式
Figure 863601DEST_PATH_IMAGE102
，得到第
Figure 951643DEST_PATH_IMAGE103
个点云切片上第
Figure 218676DEST_PATH_IMAGE104
个切面区块的面积
Figure 151997DEST_PATH_IMAGE105
；其中，
Figure 7827DEST_PATH_IMAGE106
表示关于X轴坐标的积分，且
Figure 266770DEST_PATH_IMAGE107
的取值范围为
Figure 755520DEST_PATH_IMAGE108
，
Figure 758111DEST_PATH_IMAGE109
表示第
Figure 219179DEST_PATH_IMAGE110
个点云切片上第
Figure 383444DEST_PATH_IMAGE111
个切面区块的最小X轴坐标，
Figure 625070DEST_PATH_IMAGE112
表示第
Figure 165772DEST_PATH_IMAGE113
个点云切片上第
Figure 543664DEST_PATH_IMAGE111
个切面区块的最大X轴坐标；
Figure 832825DEST_PATH_IMAGE114
表示绝对值符合。
本实施例中，步骤502中根据第
Figure 296168DEST_PATH_IMAGE113
个点云切片上第
Figure 640561DEST_PATH_IMAGE115
个切面区块的面积，获取第
Figure 138539DEST_PATH_IMAGE113
个点云切片上第
Figure 582290DEST_PATH_IMAGE115
个切面区块的变化体积，进而得到第
Figure 798507DEST_PATH_IMAGE113
个点云切片的变化体积，具体过程如下：
步骤5021、当
Figure 681013DEST_PATH_IMAGE113
处于1～
Figure 33497DEST_PATH_IMAGE116
时，则根据公式
Figure 710466DEST_PATH_IMAGE117
，得到第
Figure 335351DEST_PATH_IMAGE113
个点云切片上第
Figure 21547DEST_PATH_IMAGE118
个切面区块的变化体积
Figure 228537DEST_PATH_IMAGE119
；
当
Figure 341987DEST_PATH_IMAGE120
取
Figure 204901DEST_PATH_IMAGE121
时，则根据公式
Figure 429209DEST_PATH_IMAGE122
，得到第
Figure 756285DEST_PATH_IMAGE123
个点云切片上第
Figure 40636DEST_PATH_IMAGE124
个切面区块的变化体积
Figure 187583DEST_PATH_IMAGE125
；
步骤5022、当
Figure 901068DEST_PATH_IMAGE127
处于1～
Figure 82651DEST_PATH_IMAGE128
时，采用计算机根据
Figure 537903DEST_PATH_IMAGE129
，得到第
Figure 437726DEST_PATH_IMAGE130
个点云切片的变化体积
Figure 941519DEST_PATH_IMAGE131
；
当
Figure 977608DEST_PATH_IMAGE132
取
Figure 603762DEST_PATH_IMAGE133
时，采用计算机根据
Figure 725302DEST_PATH_IMAGE134
，得到第
Figure 95103DEST_PATH_IMAGE135
个点云切片的变化体积
Figure 438229DEST_PATH_IMAGE136
；
步骤503中将各个点云切片的变化体积累加得到待测沟壑区域第
Figure 235283DEST_PATH_IMAGE137
期相对第
Figure 844119DEST_PATH_IMAGE138
期的体积变化量
Figure 752032DEST_PATH_IMAGE139
，具体过程如下：
采用计算机根据
Figure 700397DEST_PATH_IMAGE140
，得到待测沟壑区域第
Figure 668353DEST_PATH_IMAGE141
期相对第
Figure 764485DEST_PATH_IMAGE142
期的体积变化量
Figure 476089DEST_PATH_IMAGE143
。
本实施例中，步骤一中待测沟壑区域的点云获取的方法，可参考申请日为2021年05月25日，申请号为CN202110569873.4的中国专利公开的一种基于地形点云的沟壑体积变化三维计算方法中步骤一至步骤四的方法，仅坐标系采用不同。
本实施例中，需要说明的是，架站区连线即为第一架站区S1、第二架站区S2和第三架站区S3的中心连线。
本实施例中，需要说明的是，一种基于地形点云的沟壑体积变化三维计算方法中步骤四中的二次滤波后第
Figure 606856DEST_PATH_IMAGE144
期激光点云即为本申请中第
Figure 434129DEST_PATH_IMAGE145
期的激光点云。
本实施例中，相邻两期的时间间隔为10days～20days，可以根据实际需要进行调整。
本实施例中，第
Figure 751978DEST_PATH_IMAGE147
期的激光点云中任一个点坐标数据记作
Figure 267272DEST_PATH_IMAGE148
，
Figure 252546DEST_PATH_IMAGE149
表示第
Figure 499988DEST_PATH_IMAGE150
期第
Figure 305133DEST_PATH_IMAGE151
个点在OXYZ坐标系下的X轴坐标，
Figure 358539DEST_PATH_IMAGE152
表示第
Figure 198319DEST_PATH_IMAGE153
期第
Figure 865930DEST_PATH_IMAGE154
个点在OXYZ坐标系下的Y轴坐标，
Figure 423950DEST_PATH_IMAGE155
表示第
Figure 15469DEST_PATH_IMAGE156
期第
Figure 975334DEST_PATH_IMAGE157
个点在OXYZ坐标系下的Z轴坐标；其中，
Figure 626896DEST_PATH_IMAGE158
，
Figure 344316DEST_PATH_IMAGE159
和
Figure 739525DEST_PATH_IMAGE160
均为正整数，且
Figure 553897DEST_PATH_IMAGE161
表示第
Figure 376360DEST_PATH_IMAGE162
期的激光点云总数。
本实施例中，步骤201中第
Figure 328879DEST_PATH_IMAGE164
个点云集合
Figure 527779DEST_PATH_IMAGE165
的总数
Figure 196658DEST_PATH_IMAGE166
的取值为
Figure 190022DEST_PATH_IMAGE167
。
本实施例中，实际使用时，
Figure 944351DEST_PATH_IMAGE168
的取值为25～30；
Figure 884625DEST_PATH_IMAGE168
个类簇的聚类中心点、插值点集均可用X轴坐标和Y轴坐标表示。
本实施例中，实际使用时，
Figure 408010DEST_PATH_IMAGE169
为正整数，
Figure 572276DEST_PATH_IMAGE169
小于
Figure 813901DEST_PATH_IMAGE170
，且
Figure 541554DEST_PATH_IMAGE169
的取值范围为3～6。
本实施例中，需要说明的是，第
Figure 919446DEST_PATH_IMAGE171
个点云切片上第
Figure 520192DEST_PATH_IMAGE172
个切面区块的面积
Figure 983534DEST_PATH_IMAGE173
取正值或者负值，仅代表地形变化的类型，不代表数值的正负。
本实施例中，负地形区是指沟壑区域侵蚀沉降，正地形区是指沟壑区域沉积抬升，沟壑区域侵蚀沉降量减去沟壑区域沉积抬升量，就是最终求得的沟壑区域体积改变量即产沙量。
本实施例中，需要说明的是，第
Figure 327928DEST_PATH_IMAGE174
期的激光点云和第
Figure 498009DEST_PATH_IMAGE175
期的激光点的切平面在Z轴方向的位置相同，以及总的切平面个数
Figure 269656DEST_PATH_IMAGE176
相同。
本实施例中，当
Figure 220295DEST_PATH_IMAGE177
取1时，
Figure 102800DEST_PATH_IMAGE178
表示第1个点云切片上切面区块的总数；当
Figure 143699DEST_PATH_IMAGE177
取2时，
Figure 86248DEST_PATH_IMAGE179
表示第2个点云切片上切面区块的总数；...；当
Figure 524182DEST_PATH_IMAGE177
取
Figure 210378DEST_PATH_IMAGE180
时，
Figure 355052DEST_PATH_IMAGE181
表示第
Figure 468501DEST_PATH_IMAGE180
个点云切片上切面区块的总数，且
Figure 393732DEST_PATH_IMAGE178
、
Figure 618040DEST_PATH_IMAGE179
和
Figure 945116DEST_PATH_IMAGE182
均为正整数。
综上所述，本发明方法步骤简单，设计合理且实现方便，获取待测沟壑区域沿高程方向各个点云切片的体积变化量，进而累加得到待测沟壑区域整体的体积变化量，提高了沟壑体积变化量计算的准确性，以适应沟壑地形。
以上所述，仅是本发明的较佳实施例，并非对本发明作任何限制，凡是根据本发明技术实质对以上实施例所作的任何简单修改、变更以及等效结构变化，均仍属于本发明技术方案的保护范围内。

Claims (3)

1.一种基于地形点云的沟壑体积变化切片方法，其特征在于，该方法包括以下步骤：
步骤一、待测沟壑区域的点云获取：
步骤101、建立空间直角坐标系；其中，以待测沟壑区域外侧左下方的稳定区任意一点为原点O，以过原点O且沿架站区连线为X轴，过原点O且与X轴垂直的方向为Y轴，以过原点O垂直于由X轴和Y轴形成的平面OXY且沿高程方向为Z轴，建立空间直角坐标系OXYZ；
步骤102、采用地面三维激光扫描仪对待测沟壑区域进行扫描，获取第1期的激光点云，...，第
Figure 20594DEST_PATH_IMAGE001
期的激光点云，...，第
Figure 347670DEST_PATH_IMAGE002
期的激光点云；其中，
Figure 569704DEST_PATH_IMAGE003
和
Figure 716651DEST_PATH_IMAGE004
均为正整数，
Figure 744650DEST_PATH_IMAGE005
，
Figure 926233DEST_PATH_IMAGE006
；
步骤二、待测沟壑区域的点云的切片：
在OXYZ坐标系下，采用计算机按照相邻两期的点云切片厚度
Figure 817703DEST_PATH_IMAGE007
，将第
Figure 717526DEST_PATH_IMAGE003
期的激光点云沿Z轴方向做多个切平面并将多个切平面从下至上进行排序依次记作第1个切平面，...，第
Figure 283637DEST_PATH_IMAGE008
个切平面，...，第
Figure 257409DEST_PATH_IMAGE009
个切平面；其中，
Figure 883562DEST_PATH_IMAGE010
和
Figure 5102DEST_PATH_IMAGE009
均为正整数，且
Figure 374903DEST_PATH_IMAGE011
；任一个切平面平行平面OXY，第1个切平面和第
Figure 468761DEST_PATH_IMAGE012
个切平面之间相邻两个切平面的间距为
Figure 265816DEST_PATH_IMAGE013
，第1个切平面距离第
Figure 874652DEST_PATH_IMAGE014
期的激光点云的最小Z轴坐标处的间距为
Figure 782565DEST_PATH_IMAGE015
，第
Figure 730930DEST_PATH_IMAGE016
个切平面与第
Figure 698886DEST_PATH_IMAGE017
期的激光点云的最大Z轴坐标之间的间距记作
Figure 795018DEST_PATH_IMAGE018
,且
Figure 942840DEST_PATH_IMAGE019
不大于
Figure 73607DEST_PATH_IMAGE020
；
步骤三、待测沟壑区域的点云切片的聚类及插值处理：
步骤301、当
Figure 212464DEST_PATH_IMAGE021
处于1～
Figure 530313DEST_PATH_IMAGE022
时，采用计算机将第
Figure 983291DEST_PATH_IMAGE021
个切平面上方及下方Z轴方向
Figure 968565DEST_PATH_IMAGE023
范围内的点云投影至第
Figure 278323DEST_PATH_IMAGE021
个切平面，得到第
Figure 83468DEST_PATH_IMAGE017
期第
Figure 74558DEST_PATH_IMAGE021
个点云切片；
当
Figure 914338DEST_PATH_IMAGE021
取
Figure 394998DEST_PATH_IMAGE016
时，采用计算机将第
Figure 953018DEST_PATH_IMAGE016
个切平面上方Z轴方向的剩余点云及下方Z轴方向
Figure 482220DEST_PATH_IMAGE023
范围内的点云投影至第
Figure 442085DEST_PATH_IMAGE016
个切平面，得到第
Figure 93647DEST_PATH_IMAGE017
期第
Figure 309602DEST_PATH_IMAGE016
个点云切片；
步骤302、采用计算机利用FCM聚类算法对第
Figure 704811DEST_PATH_IMAGE017
期第
Figure 519184DEST_PATH_IMAGE024
个点云切片上的点进行聚类，将第
Figure 341646DEST_PATH_IMAGE017
期第
Figure 546362DEST_PATH_IMAGE024
个点云切片上各个点划分为
Figure 745263DEST_PATH_IMAGE025
个类簇，并得到
Figure 414141DEST_PATH_IMAGE026
个类簇的聚类中心点；其中，
Figure 407505DEST_PATH_IMAGE026
为正整数；
步骤303、采用计算机利用三次B样条插值函数对第
Figure 99518DEST_PATH_IMAGE017
期第
Figure 102109DEST_PATH_IMAGE024
个点云切片上
Figure 625494DEST_PATH_IMAGE025
个聚类中心点进行插值处理，得到第
Figure 789759DEST_PATH_IMAGE017
期第
Figure 969068DEST_PATH_IMAGE024
个点云切片上插值点集
Figure 509770DEST_PATH_IMAGE027
；
步骤四、获取后一期点云切片相对前一期点云切片的切面变化面积：
步骤401、按照步骤二至步骤三所述的方法，采用计算机对第
Figure 887662DEST_PATH_IMAGE028
期的激光点云进行处理，得到第
Figure 488408DEST_PATH_IMAGE028
期第
Figure 387968DEST_PATH_IMAGE029
个点云切片上插值点集
Figure 732362DEST_PATH_IMAGE030
；
步骤402、采用计算机将第
Figure 964760DEST_PATH_IMAGE031
期第
Figure 736407DEST_PATH_IMAGE032
个点云切片上插值点集
Figure 624729DEST_PATH_IMAGE033
和第
Figure 507234DEST_PATH_IMAGE034
期第
Figure 859718DEST_PATH_IMAGE035
个点云切片上插值点集
Figure 739949DEST_PATH_IMAGE036
同步绘制在平面OXY上，且
Figure 177884DEST_PATH_IMAGE037
和
Figure 864080DEST_PATH_IMAGE038
中相交的相邻两个交点将
Figure 71071DEST_PATH_IMAGE039
和
Figure 122203DEST_PATH_IMAGE040
之间围设的区域划分为第1个切面区块，..，第
Figure 47434DEST_PATH_IMAGE041
个切面区块，..，第
Figure 271742DEST_PATH_IMAGE042
个切面区块，并得到各个切面区块的面积；其中，第
Figure 598818DEST_PATH_IMAGE043
个点云切片上第
Figure 77245DEST_PATH_IMAGE044
个切面区块的面积记作
Figure 489772DEST_PATH_IMAGE045
；其中，
Figure 517771DEST_PATH_IMAGE046
和
Figure 699354DEST_PATH_IMAGE047
均为正整数，且
Figure 92289DEST_PATH_IMAGE048
，
Figure 726533DEST_PATH_IMAGE049
表示第
Figure 292643DEST_PATH_IMAGE050
个点云切片上切面区块的总数；
步骤五、待测沟壑区域的体积变化量获取：
步骤501、对第
Figure 531995DEST_PATH_IMAGE050
个点云切片上第
Figure 158148DEST_PATH_IMAGE046
个切面区块进行判断，如果第
Figure 279688DEST_PATH_IMAGE050
个点云切片上第
Figure 587172DEST_PATH_IMAGE046
个切面区块中第
Figure 743347DEST_PATH_IMAGE051
期插值点的Y轴坐标大于第
Figure 540402DEST_PATH_IMAGE052
期插值点的Y轴坐标，则第
Figure 149238DEST_PATH_IMAGE046
个切面区块为负地形区，则第
Figure 493369DEST_PATH_IMAGE050
个点云切片上第
Figure 504050DEST_PATH_IMAGE046
个切面区块的面积
Figure 472006DEST_PATH_IMAGE053
取正值；如果第
Figure 568138DEST_PATH_IMAGE050
个点云切片上第
Figure 217426DEST_PATH_IMAGE046
个切面区块中第
Figure 82613DEST_PATH_IMAGE054
期插值点的Y轴坐标小于等于第
Figure 221471DEST_PATH_IMAGE052
期插值点的Y轴坐标，则第
Figure 742582DEST_PATH_IMAGE050
个点云切片上第
Figure 257877DEST_PATH_IMAGE046
个切面区块为正地形区，则第
Figure 977571DEST_PATH_IMAGE050
个点云切片上第
Figure 287330DEST_PATH_IMAGE046
个切面区块的面积
Figure 295737DEST_PATH_IMAGE055
取负值；
步骤502、根据第
Figure 349144DEST_PATH_IMAGE050
个点云切片上第
Figure 188924DEST_PATH_IMAGE046
个切面区块的面积，获取第
Figure 669584DEST_PATH_IMAGE050
个点云切片上第
Figure 663822DEST_PATH_IMAGE046
个切面区块的变化体积，进而得到第
Figure 255341DEST_PATH_IMAGE050
个点云切片的变化体积；
步骤503、将各个点云切片的变化体积累加得到待测沟壑区域第
Figure 949627DEST_PATH_IMAGE056
期相对第
Figure 804451DEST_PATH_IMAGE057
期的体积变化量
Figure 584188DEST_PATH_IMAGE058
；
步骤502中根据第
Figure 979397DEST_PATH_IMAGE059
个点云切片上第
Figure 793769DEST_PATH_IMAGE060
个切面区块的面积，获取第
Figure 553915DEST_PATH_IMAGE059
个点云切片上第
Figure 820948DEST_PATH_IMAGE060
个切面区块的变化体积，进而得到第
Figure 19848DEST_PATH_IMAGE059
个点云切片的变化体积，具体过程如下：
步骤5021、当
Figure 688727DEST_PATH_IMAGE059
处于1～
Figure 619774DEST_PATH_IMAGE061
时，则根据公式
Figure 374103DEST_PATH_IMAGE062
，得到第
Figure 111115DEST_PATH_IMAGE059
个点云切片上第
Figure 634500DEST_PATH_IMAGE063
个切面区块的变化体积
Figure 500563DEST_PATH_IMAGE064
；
当
Figure 742189DEST_PATH_IMAGE065
取
Figure 17312DEST_PATH_IMAGE066
时，则根据公式
Figure 598466DEST_PATH_IMAGE067
，得到第
Figure 933632DEST_PATH_IMAGE068
个点云切片上第
Figure 662554DEST_PATH_IMAGE069
个切面区块的变化体积
Figure 741369DEST_PATH_IMAGE070
；
步骤5022、当
Figure 177029DEST_PATH_IMAGE072
处于1～
Figure 948676DEST_PATH_IMAGE073
时，采用计算机根据
Figure 899314DEST_PATH_IMAGE074
，得到第
Figure 781820DEST_PATH_IMAGE075
个点云切片的变化体积
Figure 71987DEST_PATH_IMAGE076
；
当
Figure 14535DEST_PATH_IMAGE077
取
Figure 452470DEST_PATH_IMAGE078
时，采用计算机根据
Figure 873087DEST_PATH_IMAGE079
，得到第
Figure 781875DEST_PATH_IMAGE080
个点云切片的变化体积
Figure 895324DEST_PATH_IMAGE081
；
步骤503中将各个点云切片的变化体积累加得到待测沟壑区域第
Figure 820555DEST_PATH_IMAGE082
期相对第
Figure 982546DEST_PATH_IMAGE083
期的体积变化量
Figure 44043DEST_PATH_IMAGE084
，具体过程如下：
采用计算机根据
Figure 328394DEST_PATH_IMAGE085
，得到待测沟壑区域第
Figure 740920DEST_PATH_IMAGE086
期相对第
Figure 706602DEST_PATH_IMAGE087
期的体积变化量
Figure 888185DEST_PATH_IMAGE088
。
2.按照权利要求1所述的一种基于地形点云的沟壑体积变化切片方法，其特征在于：步骤二中相邻两期的点云切片厚度
Figure 77858DEST_PATH_IMAGE089
的获取，具体过程如下：
步骤201、采用计算机从第
Figure 977681DEST_PATH_IMAGE090
期的激光点云中随机选择点云组成第
Figure 481474DEST_PATH_IMAGE090
个点云集合
Figure 783143DEST_PATH_IMAGE091
；其中，第
Figure 143717DEST_PATH_IMAGE092
个点云集合
Figure 530836DEST_PATH_IMAGE093
中第
Figure 336856DEST_PATH_IMAGE095
个点记作
Figure 227451DEST_PATH_IMAGE096
，
Figure 24506DEST_PATH_IMAGE097
，
Figure 571025DEST_PATH_IMAGE098
表示第
Figure 744517DEST_PATH_IMAGE099
个点云集合
Figure 755199DEST_PATH_IMAGE100
的总数；
步骤202、采用计算机获取第
Figure 723155DEST_PATH_IMAGE092
个点云集合
Figure 756970DEST_PATH_IMAGE100
中第
Figure 468574DEST_PATH_IMAGE095
个点
Figure 333762DEST_PATH_IMAGE101
与第
Figure 472619DEST_PATH_IMAGE102
期的激光点云中各个点的欧式距离，并将各个欧式距离按照从大到小的顺序排序，获取后
Figure 993730DEST_PATH_IMAGE103
个欧式距离；其中，与
Figure 243446DEST_PATH_IMAGE104
对应的后
Figure 228719DEST_PATH_IMAGE105
个欧式距离中第
Figure 974696DEST_PATH_IMAGE106
个欧式距离记作
Figure 45420DEST_PATH_IMAGE107
；
步骤203、根据公式
Figure 98827DEST_PATH_IMAGE108
，得到第
Figure 938607DEST_PATH_IMAGE109
期的激光点云的平均点间距
Figure 356950DEST_PATH_IMAGE110
；
步骤204、按照步骤201至步骤203所述的方法，得到第
Figure 649391DEST_PATH_IMAGE111
期的激光点云的平均点间距
Figure 506489DEST_PATH_IMAGE112
；并对
Figure 200775DEST_PATH_IMAGE113
和
Figure 790020DEST_PATH_IMAGE114
进行平均值处理，得到相邻两期的平均点间距
Figure 835336DEST_PATH_IMAGE115
；
步骤205、根据公式
Figure 230545DEST_PATH_IMAGE116
，得到相邻两期的点云切片厚度
Figure 779338DEST_PATH_IMAGE117
；其中，
Figure 805063DEST_PATH_IMAGE118
为常数且
Figure 72096DEST_PATH_IMAGE119
取值为0.5。
3.按照权利要求1所述的一种基于地形点云的沟壑体积变化切片方法，其特征在于：步骤402中得到各个切面区块的面积的方法均相同，其中，得到第
Figure 270997DEST_PATH_IMAGE120
个点云切片上第
Figure 674296DEST_PATH_IMAGE121
个切面区块的面积
Figure 369457DEST_PATH_IMAGE122
，具体过程如下：
步骤4021、对第
Figure 123787DEST_PATH_IMAGE120
个点云切片上第
Figure 860799DEST_PATH_IMAGE121
个切面区块中包含的第
Figure 321867DEST_PATH_IMAGE123
期插值点拟合，得到第
Figure 751711DEST_PATH_IMAGE120
个点云切片上第
Figure 727758DEST_PATH_IMAGE121
个切面区块中第
Figure 268460DEST_PATH_IMAGE123
期的拟合函数段
Figure 849614DEST_PATH_IMAGE124
；其中，
Figure 184781DEST_PATH_IMAGE125
为自变量；
对第
Figure 913702DEST_PATH_IMAGE126
个点云切片上第
Figure 992517DEST_PATH_IMAGE127
个切面区块中包含的第
Figure 428177DEST_PATH_IMAGE128
期插值点拟合，得到第
Figure 934245DEST_PATH_IMAGE129
个点云切片上第
Figure 150463DEST_PATH_IMAGE130
个切面区块中第
Figure 32968DEST_PATH_IMAGE128
期的拟合函数段
Figure 821670DEST_PATH_IMAGE131
；
步骤4022、采用计算机根据公式
Figure 764218DEST_PATH_IMAGE132
，得到第
Figure 202153DEST_PATH_IMAGE133
个点云切片上第
Figure 560453DEST_PATH_IMAGE134
个切面区块的面积
Figure 33023DEST_PATH_IMAGE135
；其中，
Figure 880893DEST_PATH_IMAGE136
表示关于X轴坐标的积分，且
Figure 806124DEST_PATH_IMAGE137
的取值范围为
Figure 968115DEST_PATH_IMAGE138
，
Figure 295191DEST_PATH_IMAGE139
表示第
Figure 579542DEST_PATH_IMAGE140
个点云切片上第
Figure 992069DEST_PATH_IMAGE141
个切面区块的最小X轴坐标，
Figure 692171DEST_PATH_IMAGE142
表示第
Figure 873754DEST_PATH_IMAGE059
个点云切片上第
Figure 329006DEST_PATH_IMAGE141
个切面区块的最大X轴坐标；
Figure 228829DEST_PATH_IMAGE143
表示绝对值符号。
CN202210923370.7A 2022-08-02 2022-08-02 一种基于地形点云的沟壑体积变化切片方法 Active CN114998419B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210923370.7A CN114998419B (zh) 2022-08-02 2022-08-02 一种基于地形点云的沟壑体积变化切片方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210923370.7A CN114998419B (zh) 2022-08-02 2022-08-02 一种基于地形点云的沟壑体积变化切片方法

Publications (2)

Publication Number Publication Date
CN114998419A CN114998419A (zh) 2022-09-02
CN114998419B true CN114998419B (zh) 2022-11-04

Family

ID=83022226

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210923370.7A Active CN114998419B (zh) 2022-08-02 2022-08-02 一种基于地形点云的沟壑体积变化切片方法

Country Status (1)

Country Link
CN (1) CN114998419B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115824053B (zh) * 2023-02-16 2023-05-30 江苏南京地质工程勘察院 一种融合主成分分析与体素积分的滑坡体积量计算方法
CN117523124B (zh) * 2023-11-13 2024-04-26 西安科技大学 一种基于无人机激光雷达的大流域地形变化切片方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111595403A (zh) * 2020-05-15 2020-08-28 中交广州航道局有限公司 一种基于点云测量技术的工程土方计量方法
CN113160374A (zh) * 2021-05-25 2021-07-23 西安科技大学 一种基于地形点云的沟壑体积变化三维计算方法
CN114266987A (zh) * 2022-03-03 2022-04-01 水利部长江勘测技术研究所 一种无人机高边坡危岩体智能识别方法
CN114419130A (zh) * 2021-12-22 2022-04-29 中国水利水电第七工程局有限公司 一种基于图像特征和三维点云技术的散料体积测量方法
CN114429455A (zh) * 2022-01-05 2022-05-03 南京市测绘勘察研究院股份有限公司 一种基于机载激光点云辅助影像植被三维变化检测方法

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10151839B2 (en) * 2012-06-01 2018-12-11 Agerpoint, Inc. Systems and methods for determining crop yields with high resolution geo-referenced sensors
CN107464223B (zh) * 2017-07-19 2020-01-14 西安理工大学 一种基于切片的点云孔洞修补方法
CN110827233B (zh) * 2019-08-29 2022-06-14 杭州电子科技大学 一种牙齿三维点云数据表面窝沟区域的提取方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111595403A (zh) * 2020-05-15 2020-08-28 中交广州航道局有限公司 一种基于点云测量技术的工程土方计量方法
CN113160374A (zh) * 2021-05-25 2021-07-23 西安科技大学 一种基于地形点云的沟壑体积变化三维计算方法
CN114419130A (zh) * 2021-12-22 2022-04-29 中国水利水电第七工程局有限公司 一种基于图像特征和三维点云技术的散料体积测量方法
CN114429455A (zh) * 2022-01-05 2022-05-03 南京市测绘勘察研究院股份有限公司 一种基于机载激光点云辅助影像植被三维变化检测方法
CN114266987A (zh) * 2022-03-03 2022-04-01 水利部长江勘测技术研究所 一种无人机高边坡危岩体智能识别方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
Rapid clustering of colorized 3D point cloud data for reconstructing building interiors;Kuldeep K. Sareen et al;《 2010 International Symposium on Optomechatronic Technologies》;20110113;1-6 *
复杂地形中机载LiDAR点云构建DEM的插值算法对比;李朋飞等;《农业工程学报》;20210831;第37卷(第15期);146-153 *

Also Published As

Publication number Publication date
CN114998419A (zh) 2022-09-02

Similar Documents

Publication Publication Date Title
CN114998419B (zh) 一种基于地形点云的沟壑体积变化切片方法
CN108986048B (zh) 基于线激光扫描三维点云快速复合滤波处理方法
Robinson The accuracy of digital elevation models derived from digitised contour data
CN109816664B (zh) 一种三维点云分割方法及装置
CN107038717A (zh) 一种基于立体栅格自动分析3d点云配准误差的方法
CN110276732B (zh) 一种顾及地形特征线要素的山区点云空洞修复方法
CN111598780B (zh) 一种适用于机载LiDAR点云的地形自适应插值滤波方法
CN112215958B (zh) 一种基于分布式计算的激光雷达点云数据投影方法
CN111199549B (zh) 基于叶片型面测量点云的叶型提取方法
CN110726998B (zh) 一种激光雷达扫描测定矿区采煤塌陷盆地的方法
CN108765568A (zh) 一种基于激光雷达点云的多层次建筑物快速三维重建方法
CN107590347A (zh) 一种基于设计模型匹配孤立点识别与删除方法及系统
CN116486268A (zh) 一种激光点云地面滤波方法及装置
CN110458938B (zh) 一种散货料堆实时三维重建方法及系统
CN112666614B (zh) 基于电法勘探和数字高程模型泥石流物源静储量计算方法
CN112687005A (zh) 一种基于三维重建的煤堆体积测量方法
CN113160374B (zh) 一种基于地形点云的沟壑体积变化三维计算方法
CN107545601B (zh) 一种架空输电线路树高断面自动生成方法
CN113947628B (zh) 一种基于不同比尺河道地形图计算冲淤量的方法
CN112365534B (zh) 一种基于单目相机三维重建的大型煤堆体积测量方法
CN112381029B (zh) 一种基于欧氏距离的机载LiDAR数据建筑物提取方法
CN117523124B (zh) 一种基于无人机激光雷达的大流域地形变化切片方法
CN113496551B (zh) 一种基于地质露头三维模型的地形剖面线绘制方法
CN111563892A (zh) 基于点云数据立面采集精度优化的拟合算法
CN116645482B (zh) 一种大范围地面连续高程提取方法

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