CN117523130A - Method for generating city high-precision digital terrain and dividing sub-catchment area by LiDAR - Google Patents

Method for generating city high-precision digital terrain and dividing sub-catchment area by LiDAR Download PDF

Info

Publication number
CN117523130A
CN117523130A CN202311574228.7A CN202311574228A CN117523130A CN 117523130 A CN117523130 A CN 117523130A CN 202311574228 A CN202311574228 A CN 202311574228A CN 117523130 A CN117523130 A CN 117523130A
Authority
CN
China
Prior art keywords
sub
catchment area
point
points
dividing
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.)
Pending
Application number
CN202311574228.7A
Other languages
Chinese (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.)
Guangxi University
Original Assignee
Guangxi University
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 Guangxi University filed Critical Guangxi University
Priority to CN202311574228.7A priority Critical patent/CN117523130A/en
Publication of CN117523130A publication Critical patent/CN117523130A/en
Pending legal-status Critical Current

Links

Classifications

    • 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

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Geometry (AREA)
  • Software Systems (AREA)
  • Remote Sensing (AREA)
  • Computer Graphics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Image Processing (AREA)

Abstract

The method for generating the urban high-precision digital terrain and dividing the sub-catchment areas by using the LiDAR comprises the following steps: after data preprocessing is carried out on original LiDAR point cloud data and original image data of a certain city, a TIN model is built through a terraModel module, after TIN models of DEM and DSM are respectively generated, elevation values of target points are obtained through interpolation calculation, and high-precision DEM and DSM are generated; defining connection points by the preprocessed original LiDAR point cloud data and the original image data, and encrypting by using a connection point space triangle; defining Colorpoints for the original image data after space triangle encryption, and generating a high-precision DOM after correcting the spelling line; dividing a first-level sub-catchment area of the city by a hydrologic analysis method in combination with the high-precision DEM; on the basis of dividing the primary sub-catchment area of the city, combining a high-precision DSM, a road network and a pipe network to divide the secondary sub-catchment area of the city; and on the basis of the urban secondary sub-catchment area, combining a Thiessen polygon method and manually modifying and dividing the urban tertiary sub-catchment area by referring to a high-precision DOM. The method can improve the simulation accuracy of the urban rainfall flood process.

Description

Method for generating city high-precision digital terrain and dividing sub-catchment area by LiDAR
Technical Field
The invention belongs to the technical field of high-precision digital terrain models and urban sub-catchment areas, and particularly relates to a method for generating urban high-precision digital terrain and dividing the sub-catchment areas by using LiDAR (Chinese is a laser radar).
Background
The former person earlier carries out urban sub-catchment area division according to the distribution of terrains, blocks and catch basins by a manual division method. And e.g. Barco divides the catchment area according to the pipe network, the catchment area data and the water flow direction. With the development of GIS technology, the division of urban sub-catchment areas is gradually changed into a semi-automatic and automatic mode. Zhang Shuliang and the like establish an automatic dividing method of urban rainwater drainage basin catchment areas by using drainage basin topography segmentation technology and professional GIS topography processing software, so that the efficiency and quality of urban catchment area division are greatly improved, but only topography factors are considered. Xue Fengchang A technical method for classifying catchments for urban plain areas is proposed, the method divides cities into central urban areas and suburban areas by urban land classification, divides primary catchments according to urban drainage main water systems, combines different runoff influencing factors with DEM, divides secondary catchments according to refined DEM, and combines Voronoi diagrams according to actual confluence conditions on the basis that the method for dividing the central urban three-level catchments does not consider influencing factors such as buildings. Wang Yimei A method for dividing urban rainwater catchment areas with multiple factors is provided, and is based on a hydrologic analysis method and a GIS technology, and the urban subcollection areas are divided in three stages by considering influence factors such as topography, roads, buildings, drainage pipe networks and the like, so that the urban subcollection areas are further thinned, but the accuracy of DEM does not reach higher level.
In the past, sub-meter level cannot be realized due to the accuracy of the DEM, the low-accuracy DEM is generally used for dividing the urban sub-catchment area as a dividing basis, and factors affecting the urban sub-catchment area relate to multiple aspects and multiple factors, so that all influencing factors cannot be fully considered to divide the sub-catchment area. Under the background that the occurrence frequency and the intensity of urban inland inundation are gradually increased, the simulation of the urban rainfall flood process is particularly important, so that the urban sub-catchment area division is refined, and the simulation accuracy of the urban rainfall flood process is improved.
In summary, based on the dividing requirement of the urban sub-catchment area, a sub-meter grade even higher-precision digital terrain model needs to be manufactured, and the urban sub-catchment area is finely divided by combining with the high-precision digital terrain model results and considering various influencing factors, so that a key scientific foundation is laid for improving the simulation result of the urban rainfall flood process and strengthening the urban rainfall warning, management and forecasting system.
Disclosure of Invention
Aiming at the problems existing in the prior art, the invention provides a method for generating urban high-precision digital topography and dividing a sub-catchment area by using LiDAR, which aims to generate a high-precision digital topography model (DEM, DSM, DOM) based on LiDAR (Chinese laser radar) point cloud data by using terrasoid software under a MicroStation platform, thereby further refining the urban sub-catchment area by considering various factors such as topography, road network, pipe network, building and the like.
In order to achieve the above object, the present invention is specifically as follows:
a method for generating city high-precision digital terrains and dividing sub-catchment areas by using LiDAR comprises the following steps:
s1, carrying out data preprocessing on original LiDAR point cloud data and original image data of a certain city;
s2, constructing a TIN model by using the original LiDAR point cloud data preprocessed in the step S1 through a terra model module in terra solid software, and respectively generating a TIN model of a DEM and a DSM;
s3, respectively carrying out interpolation calculation on the DEM and the DSM TIN model generated in the step S2 to obtain the elevation value of the target point, and generating a high-precision DEM and DSM;
s4, defining connection points for the original LiDAR point cloud data and the original image data which are preprocessed in the step S1, and performing space triangular encryption through the connection points to match the original LiDAR point cloud data and the original image data;
s5, defining Color points for the original image data subjected to space triangle encryption in the step S4, and correcting the spelling line to generate a high-precision DOM;
s6, dividing a first-level sub-catchment area of the city by a hydrologic analysis method in combination with the high-precision DEM in the step S2;
s7, dividing a city secondary sub-catchment area by combining a high-precision DSM, a road network and a pipe network on the basis of dividing the city primary sub-catchment area in the step S6;
and S8, dividing the urban three-level sub-catchment area on the basis of the urban two-level sub-catchment area in the step S7 by combining a Thiessen polygon method, and manually modifying and dividing the urban three-level sub-catchment area by referring to the high-precision DOM.
Further, the original LiDAR point cloud data in step S1 is a data set containing a large number of discrete points generated by LiDAR collecting the surface features of the ground, the building, the vegetation and the river of the urban area, wherein the data set comprises coordinates of the points, attributes of the points, density and a space reference system; the properties of the points include color, reflection intensity, and classification information; the original image data is unprocessed urban area satellite remote sensing image data and is stored in a TIFF format.
Further, the step S1 of preprocessing the original LiDAR point cloud data includes the following steps:
s11, data import: the method comprises the steps of importing original LiDAR point cloud data through a terraScan module of TerraSolid software, performing thinning treatment, and importing original image data through a terraPhoto module;
s12, density detection: displaying a density map of original LiDAR point cloud data through a visualization tool of terraSolid software;
s13, dividing the dense point cloud data into blocks: original LiDAR point cloud data imported in the Terra Scan module are referred to the Terra photo module, and dense original LiDAR point cloud data are divided into blocks through block polygons;
s14, filtering and classifying irregular dense point cloud data: the point cloud filtering and classifying function of the Terra Scan module loads corresponding filtering and classifying algorithm according to the type of the point cloud data to carry out progressive encryption filtering processing on the original LiDAR point cloud data, remove isolated points andrespectively calculating average distances d 'of each point and k' nearest neighbors of each point after removing isolated pointsAnd standard deviation s, according to the set threshold valuePerforming filtering operation when->In the case of ground points, when +.>And is a non-ground point, where k' and k s All are preset parameters, the ground points are separated from the non-ground points, so that the ground points are extracted, and the non-ground points are reclassified through a ground clearance algorithm;
s15, homogenizing and color-homogenizing and splicing original image data: the terraPhoto module is used for adjusting the numerical value of the image color band and the brightness, the brightness and the contrast of the image to the imported original image data; and automatically splicing the plurality of image pictures into an image in a splicing line mode through a terraPhoto module, and completing the splicing of the images.
Further, in step S14, the calculation formulas of the average distances d 'between each point and k' nearest neighbors of each point after removing the isolated points are as follows:
in the formula (2), k' is a preset parameter; x is x i ’、y i ’、z i ' represents the coordinates of each point after removal of the isolated point; x is x ij ’、y ij ’、z ij 'represents the coordinates of each point corresponding to the respective k' nearest neighbors after removal of the outlier;
k' of each pointAverage value of average distance d' of nearest neighborThe calculation formula of (2) is as follows:
in the formula (3), m represents the total number of remaining points after removing isolated points;
the calculation formula of the standard deviation s is as follows:
further, the step of constructing the TIN model in step S2 is as follows:
s21, performing triangulation algorithm on the original LiDAR point cloud data set to obtain a ground point set DT (P) ground ) The ground point set DT (P ground ) Connecting to generate a TIN model;
s22, adding new ground points by using a random sampling consistency algorithm, and updating a TIN model;
s23, checking triangles in the TIN model by using an elevation difference method;
s24, repeating the steps S22 and S23, continuously adding new ground points, and checking triangles in the TIN model until all the ground points are processed.
Further, the ground point set DT (P ground ) The formula of (2) is as follows:
TIN=DT(P ground ),(10)
wherein TIN represents a TIN model; DT (P) ground ) Representing a set of ground points.
Further, the calculation formula of the elevation value of the target point in step S3 is as follows:
wherein z represents the elevation value of the target point; z i Representing the elevation value of the interpolation point i; w (w) i The weight between the target point and the interpolation point i is represented.
Further, the specific steps of dividing the first-level sub-catchment area of the city by the high-precision DEM in step S6 through the hydrologic analysis method are as follows:
s61, filling the depressions, namely correcting the low-lying areas in the high-precision DEM data through a depression filling tool so as to simulate the water flow path, the water flow direction and the water flow rate;
s62, calculating the flow direction and flow rate of water flow: simulating the water flow direction by a flow direction tool, so that the water flow flows from the area of Gao Chengda to the area of Gao Chengxiao, simulating the water flow by a flow direction tool, and accumulating the flow from upstream to downstream to obtain the water flow direction and flow data of the whole area;
s63, generating a grid river network: setting a corresponding flow threshold value through a grid calculator tool, screening and reserving rivers with flow larger than the threshold value, and generating a grid river network so as to acquire grid-form river network data;
s64, river link: according to the water flow direction and flow data of the whole area, discrete water flow paths in the grid river network are connected into a continuous river line through a river link tool;
s65, grid river network vectorization: converting river network data in a grid form into vector element representation to obtain water system distribution conditions of the whole urban area;
s66, extracting watershed: firstly, selecting the intersection point of tributaries and main flows in a river network as a water outlet, then using a pouring point capturing tool to set the water outlet as a pouring point, and finally extracting a watershed by a watershed tool;
s67, dividing a catchment area: based on the steps S61-S66, collecting area distribution of the whole urban area is obtained through a collecting area tool, a plurality of adjacent collecting areas are combined to generate sub-collecting areas according to the water system distribution condition of the whole urban area, each collecting area is guaranteed to be included in a unique sub-collecting area until all the collecting areas are processed, and finally the first-level sub-collecting area of the city is obtained.
Further, in step S7, the step of dividing the urban secondary sub-catchment is as follows:
s71, respectively extracting and integrating a trunk road and a trunk pipeline according to road network and pipe network data to generate a trunk drainage network;
s72, dividing the urban neighborhood sub-catchment area based on the high-precision DEM by a hydrologic analysis method according to the step S6;
s73, dividing a main road surface according to the main drainage network generated in the step S71, dividing the main road surface into a plurality of road ranges by the width of the road, and if an overlapping phenomenon occurs at the road intersection, distributing the road with the overlapping road intersection into the road ranges according to the main drainage network, wherein the road is used as a drainage system of a secondary sub-catchment area of the city, and each road range is a sub-catchment area of the road;
s74, determining the flow accumulation condition of each raster data in the city neighborhood sub-catchment area in the step S72 when dividing, determining the water outlet of each city neighborhood sub-catchment area according to the principle that the maximum flow is the water outlet, wherein the neighborhood sub-catchment area with the water outlet at the junction of the neighborhood and the road is called a direct catchment area, the neighborhood sub-catchment area with the water outlet at the junction of the same neighborhood sub-catchment area is called an indirect catchment area, merging the indirect catchment area with the direct catchment area according to the flow direction, and merging the merged sub-catchment area with the road sub-catchment area in the step S73 again according to the position of the water outlet to form a city secondary sub-catchment area.
Further, in step S8, the calculation formula for dividing the urban three-level sub-catchment area based on the urban two-level sub-catchment area in step S7 by combining the Thiessen polygon method is as follows:
wherein L is i The coordinates representing the i point are (x i ,y i );L j Representing an arbitrary point; d (L) j ,L i ) Representing an arbitrary point L j With other points L i Distance between them.
THE ADVANTAGES OF THE PRESENT INVENTION
Compared with the traditional method for manufacturing the digital terrain model based on the satellite image, the method for dividing the LiDAR generating city high-precision digital terrain and the sub-catchment area is used for obtaining the digital terrain model reaching sub-meter level and even millimeter level high precision in an automatic and semi-automatic mode based on LiDAR point cloud data through terra solid software, so that the hardware processing difficulty is reduced, the real-time processing time is shortened, and the digital terrain model generating precision is improved to a great extent. In the aspect of urban sub-catchment area division, the urban sub-catchment area is divided by a hydrologic analysis method based on the high-precision DEM, the urban secondary sub-catchment area is divided by combining the high-precision DSM, the road network and the pipe network, the urban tertiary sub-catchment area is corrected according to the water diversion ridge generated by the high-precision DEM and the building outline in the DOM, the urban sub-catchment area is divided by taking various complex factors such as urban topography, buildings, road network and pipe network into consideration, the urban sub-catchment area division is refined, and the simulation accuracy of the urban rainfall flood process is improved.
Drawings
FIG. 1 is a flow chart of a method for generating urban high-precision digital terrain and dividing sub-catchment areas by using LiDAR according to the invention.
Fig. 2 is a diagram of original point cloud data of fig. 1.
Fig. 3 is a ground point map extracted from the point cloud data of fig. 1.
FIG. 4 is a diagram of the result of a high-precision digital elevation model of a region based on LiDAR point clouds.
Detailed Description
The invention is further illustrated in the following drawings and description of specific embodiments, it being noted that the specific embodiments are not intended to limit the scope of the invention.
As shown in fig. 1 to 4, the method for generating urban high-precision digital terrain and dividing sub-catchment areas by using the LiDAR according to the present embodiment includes the following steps:
s1, carrying out data preprocessing on original LiDAR point cloud data and original image data of a certain city;
the original LiDAR point cloud data is a data set which is generated by collecting the surface characteristics of the ground, buildings, vegetation and rivers of a city area through LiDAR (Chinese name: laser radar), wherein the data set comprises coordinates of points, attributes of the points, density and a space reference system; the properties of the points include color, reflection intensity, and classification information; the original image data is unprocessed urban area satellite remote sensing image data, is stored in a TIFF (Chinese name: bitmap) format, has the characteristic of high resolution, and meets the requirement of manufacturing high-precision digital topography. The data preprocessing comprises data importing, density detecting, dense point cloud data dividing into blocks, irregular dense point cloud data filtering, classification processing, original image data dodging and splicing, wherein the data preprocessing steps of the original LiDAR point cloud data and the original image data are as follows:
s11, data import: the TerraScan module of terraslid software is used for importing original LiDAR point cloud data and performing thinning treatment, the thinning purpose is to reduce the processing time and not influence the processing result, and the data of the specific embodiment uses original LiDAR point cloud data in a foreign place, as shown in fig. 2; setting a coordinate system and coordinate system conversion through a terraPhoto module (Chinese name: orthophoto generation module), and importing original image data containing information such as northeast coordinates, rotation angle coordinates and the like after camera parameters;
s12, density detection: displaying a density map of original LiDAR point cloud data through a visualization tool of terraSolid software; preferably, the density of the original LiDAR point cloud data is controlled within a reasonable range according to the requirement and the existing condition;
s13, dividing the dense point cloud data into blocks: original LiDAR point cloud data imported in the Terra Scan module are referred to the Terra photo module, and dense original LiDAR point cloud data are divided into blocks through block polygons, so that complexity of processing the whole point cloud data set is avoided, and therefore the speed and efficiency of data processing are improved;
s14, filtering and classifying irregular dense point cloud data: the point cloud filtering and classifying function of the Terra Scan module loads corresponding filtering and classifying algorithm according to the type of the point cloud data to carry out progressive encryption filtering processing on the original LiDAR point cloud data so as to make the original LiDAR point cloud data smooth in the mode thatTraversing each point, removing isolated points, and respectively calculating average distances d 'of each point and k' nearest neighbors of each point after removing isolated pointsAnd standard deviation s, according to the set threshold +.>Performing filtering operation when->In the case of ground points, when When it is a non-ground point, where k s The method is characterized in that the method is a preset parameter, the ground points and the non-ground points are separated, so that the ground points are extracted, the non-ground points are reclassified through a ground clearance algorithm, the ground point extraction result is shown in fig. 3, and the specific steps are as follows:
assume that the original LiDAR point cloud data is P= { P 1 ,p 2 ,...,p n The coordinates of each point are p i =(x i ,y i ,z i ) Where i=1, 2, n, calculating the average distance d between each point and the k nearest neighbors of each point i Where k is a preset parameter, the formula is as follows:
in the formula (1), k represents a preset parameter; n represents the total number of point cloud data; x is x i 、y i 、z i Representing coordinates of each point of the point cloud data; x is x ij 、y ij 、z ij Representing coordinates of each point of the point cloud data corresponding to respective k nearest neighbors;
Setting a distance threshold value u according to actual conditions, when d i Marking the points as isolated points and removing the points when the points are more than u;
the coordinates of each point remaining after the isolated point is removed are p i ’,=(x i ’,y i ’,z i '), wherein i=1, 2,..m, calculating the average distance d' between each point and k 'nearest neighbors of each point after removing isolated points' i And d' i Average value of (2)And standard deviation s, specifically as follows:
in the formula (2), k' is a preset parameter; m represents the total number of remaining points after removing isolated points; x is x i ’、y i ’、z i ' represents the coordinates of each point after removal of the isolated point; x is x ij ’、y ij ’、z ij 'represents the coordinates of each point corresponding to the respective k' nearest neighbors after removal of the outlier;
the average distance d ' of the k ' nearest neighbors of each point ' i Average value of (2)The calculation formula of (2) is as follows:
the calculation formula of the standard deviation s is as follows:
d 'for any point after removal of the outlier' i Greater than a threshold valueIt is classified as non-ground points including building points, vegetation points, etc., otherwise it is classified as ground points, where k s Is a preset parameter, the result is as follows:
classification formula of non-ground points:
classification formula of ground points:
for non-ground points, the non-ground points are divided into building points, vegetation points and the like, and the non-ground points are classified by a ground clearance algorithm, and the processing flow is as follows: firstly, all non-ground points are classified into low vegetation types and a judgment threshold value is set; generating a triangular model according to the non-ground points; the points in the low-level vegetation class are carried into a triangular model for calculation, the calculation result is compared with a threshold value, the points are low vegetation points when the calculation result is smaller than the threshold value, and the points are classified as medium-level vegetation class when the calculation result is higher than the threshold value; repeating the previous step until all points in the low-grade altitude value class are traversed; results were obtained for medium-height vegetation. The treatment process of separating the high vegetation from the medium-height vegetation is the same as that of the high vegetation, and finally, building points are separated from the high vegetation, wherein the threshold value of the low vegetation is set to be 0.5m, the medium-height vegetation is set to be 0.5-1 m, the high vegetation is set to be 1-2 m, and the height of the high vegetation is set to be higher than 2m to be a building;
because the filtering operation can change the distribution and density of the point cloud data, the filtering operation is needed to be performed again on the filtered data until the preset requirement is met;
s15, homogenizing and color-homogenizing and splicing original image data: the terraPhoto module is used for adjusting the numerical value of the image color band and the brightness, the brightness and the contrast of the image to the imported original image data, and a plurality of image images are automatically spliced into an image in the form of splicing lines to finish the splicing of the images, so that the later-stage working processing is convenient.
S2, constructing a TIN (Chinese name: irregular triangle net) model by using the original LiDAR point cloud data preprocessed in the step S1 through a terra model module (Chinese name: terrain model generating module) in terra solid software, and respectively generating a TIN model of a DEM (Chinese name: digital elevation model) and a DSM (Chinese name: digital surface model), wherein the TIN model of the DEM is a TIN model of a ground point set, and the TIN model of the DSM is a TIN model of a ground and non-ground point set;
specifically, the method for constructing the TIN model of the DEM is consistent with the method for constructing the TIN model of the DSM, and the method comprises the following steps of:
s21, triangulating the original LiDAR point cloud data set to obtain a ground point set DT (P) ground ) The formula is as follows:
TIN=DT(P ground )(5)
set of ground points DT (P) ground ) Connecting the triangle to generate a TIN model, and ensuring that each triangle in the TIN meets the condition by a triangulation algorithm, namely that each triangle circumscribed circle does not contain any other points;
s22, adding new ground points by using a random sampling consistency algorithm, and updating a TIN model;
the random sampling consistency algorithm selects some ground points in a random sampling mode so as to generate a new triangle;
assume that the set of ground points that have been currently processed is P processed The ground point set to be processed is P unprocessed
Randomly selecting k unprocessed ground points, denoted as P k ={p k1 ,p k2 ,...,p kk The coordinates of each point are p ki =(x i ,y i ,z i ) K is a preset parameter, representing the number of ground points sampled randomly each time, and a least square method is used for fitting a plane equation:
z=ax+by+c(6)
wherein a, b, c are the x, y, z components of the planar normal vector, respectively;
calculating a set of all unprocessed ground points P having a distance from the plane less than a threshold d add The formula is as follows:
P add ={p ki ∈P unprocessed :|ax i +by i +c-z i |<d}(7)
In the formula (7), x i 、y i 、z i Representing coordinates of each unprocessed ground point;
will P add Adding all points in (1) to TIN and adding P add Marking as processed;
repeating the above process until all the ground points have been processed;
s23, checking the ground points, and checking the triangles in the TIN by using an elevation difference method;
calculating the height difference between all points in the triangle and the triangle, and if a certain triangle is found to contain non-ground points, dividing the triangle into a plurality of small triangles;
for triangle t in each TIN, calculate all points p i Distance to t, find the maximum distance d max If d max Greater than a set threshold d', the triangle is considered to contain non-ground points, in which case t is split into three sub-triangles;
specifically, if three vertices of triangle t are p m 、p n 、p j Then calculate their distance to plane t, set as d m 、d n 、d j Take the maximum value d max The formula is as follows:
d max =max{d m ,d n ,d j }(8)
if d max If the triangle t is larger than the set threshold d', splitting the triangle t into three sub-triangles t 1 、t 2 、t 3 The method comprises the following steps:
t 1 from p m 、p n And p m 、p n One point of connection;
t 2 from p m 、p j And p m 、p j One point of connection;
t 3 from p n 、p j And p n 、p j One point of connection;
s24, repeating the steps S22 and S23, continuously adding new ground points, and checking triangles in the TIN until all the ground points are processed.
S3, respectively carrying out interpolation calculation on the DEM and the DSM TIN model generated in the step S2 to obtain the elevation value of the target point, and generating a high-precision DEM and DSM;
specifically, the method for generating the high-precision DEM is consistent with the method for generating the DSM, taking the high-precision DEM as an example, and the specific steps are as follows:
s31, dividing each triangle on the TIN model into a plurality of small grids, and selecting a proper interpolation algorithm according to the characteristics of the TIN model to determine elevation values of four vertex points of all grids so as to generate the DEM. Taking an inverse distance weighted interpolation algorithm as an example, calculating the weight between each point and the target point, calculating the elevation value of the target point according to the weight, and w i Representing the weight between the target point and the interpolation point i, a i Representing the distance between the target point and the interpolation point i, p representing the index value, typically taking 2 or 3, z representing the elevation value of the target point, z i The elevation value representing the interpolation point i is given by:
performing interpolation calculation on the target points on each small grid, verifying and adjusting the results, and splicing elevation values of the target points on all the small grids to obtain a DEM of the whole area, wherein the high-precision DEM result is shown in fig. 4;
and S32, performing post-processing on the generated DEM and DSM, such as removing isolated points, smoothing, filling holes and the like, namely repeating the operations of filtering, classifying and generating a TIN model until the high-precision DEM and DSM are obtained.
S4, defining connection points for the original LiDAR point cloud data and the original image data which are preprocessed in the step S1, and performing space triangular encryption through the connection points to match the original LiDAR point cloud data and the original image data;
the method comprises the following specific steps:
s41, defining connection points, manually adding the connection points by a terraPhoto module, automatically adding the function of the connection points to define the connection points, and manually adjusting the connection points with larger deviation to reduce the mismatching degree of the connection points so as to achieve the required effect, wherein the connection points are distributed in four corner areas of an image as much as possible, and are uniformly distributed to maximize the control range;
s42, performing space triangulation through the connection points, optimizing camera parameters, and resolving point coordinates of the external azimuth element and the encryption point. The matching degree between the point cloud data and the image can be improved by the space triangle encryption, the accuracy of DOM (Chinese name: digital orthophoto image) can be greatly influenced by the result, and in actual work, control points can be obtained through field operation to correct the result of the space triangle encryption, so that the accuracy of DOM is improved.
S5, defining Color points (Chinese names: color points) on the original image data subjected to space triangle encryption in the step S4, and correcting the splice lines to generate a high-precision DOM; the step of generating the high-precision DOM is as follows:
s51, defining Color points, manually adding the Color points, automatically adding the Color points through a terraPhoto module, and finally manually adjusting the Color points which do not meet the requirements to realize the purposes of modifying and balancing the Color and brightness difference of the images in the overlapping area, wherein when the Color points are manually added, the Color points are properly added on two sides of a spliced line, and places with obvious ground feature changes, obvious Color and inconsistent image Color tone are also avoided, so that noise points such as Color cast, late, bright spots and the like are avoided, the Color tone is uniform and the Color tone is excessively natural;
s52, editing the spliced line, wherein the spliced line in the step S15 is automatically generated by a terraPhoto module, and the problems of color mismatch and the like can occur at the joint of the image, so that the spliced line is manually edited, the integrity of a ground object unit is ensured when the spliced line is edited, the errors of misplacement of the ground object, obvious loss of the ground object and the like do not occur, the spliced line does not pass through a building, and the spliced line is along the road, the water surface, the field and other flat terrains as much as possible;
s53, parameters such as required resolution, generation format and the like are adjusted to generate the high-precision DOM.
S6, dividing a first-level sub-catchment area of the city by a hydrologic analysis method in combination with the high-precision DEM in the step S2;
because the urban primary sub-catchment area mainly considers the terrain and water system distribution factors, the urban primary sub-catchment area is divided by an ArcGIS software (Chinese name: geographic information system software) hydrologic analysis method based on the high-precision DEM, and the specific steps comprise filling, calculating flow direction and flow, generating grid river network, river link, grid river network vectorization, extracting watershed, and dividing catchment area, wherein tools used in the steps are all in the ArcGIS software tool box.
The method for dividing the urban primary sub-catchment area by the high-precision DEM combined with the step S2 through the hydrologic analysis method specifically comprises the following steps:
and S61, filling the depression, namely correcting the depression area in the high-precision DEM data through a depression filling tool to simulate the water flow path, the water flow direction and the water flow, so that the water flow path is more accurately simulated and the simulated water flow direction and water flow are more practical.
S62, calculating the flow direction and flow rate of water flow: simulating the water flow direction by a flow direction tool, so that the water flow flows from the area of Gao Chengda to the area of Gao Chengxiao, simulating the water flow by a flow direction tool, and accumulating the flow from upstream to downstream to obtain the water flow direction and flow data of the whole area;
s63, generating a grid river network: setting a corresponding flow threshold value through a grid calculator tool, screening and reserving rivers with flow larger than the threshold value, and generating a grid river network so as to acquire grid-form river network data;
s64, river link: according to the water flow direction and flow data of the whole area, discrete water flow paths in the grid river network are connected into a continuous river line through a river link tool;
s65, grid river network vectorization: converting river network data in a grid form into vector element representation to obtain water system distribution conditions of the whole urban area;
s66, extracting watershed: firstly, selecting the intersection point of tributaries and main flows in a river network as a water outlet, then using a pouring point capturing tool to set the water outlet as a pouring point, and finally extracting a watershed by a watershed tool;
s67, dividing a catchment area: based on the steps S61-S66, collecting area distribution of the whole urban area is obtained through a collecting area tool, a plurality of adjacent collecting areas are combined to generate sub-collecting areas according to the water system distribution condition of the whole urban area, each collecting area is guaranteed to be included in a unique sub-collecting area until all the collecting areas are processed, and finally the first-level sub-collecting area of the city is obtained.
S7, dividing a city secondary sub-catchment area by combining a high-precision DSM, a road network and a pipe network on the basis of the city primary sub-catchment area divided in the step S6; the urban secondary sub-catchment area is a refinement of the primary sub-catchment area, and compared with the urban primary sub-catchment area, the main factors considered by the urban secondary sub-catchment area are terrain, road network, pipe network, buildings and the like, and the dividing steps comprise: respectively extracting and integrating a trunk road and a trunk pipeline according to road network and pipe network data to generate a trunk drainage network; dividing a city neighborhood sub-catchment area by utilizing an ArcGIS hydrologic analysis method according to the high-precision DSM; dividing a main road surface according to a main drainage network; the street block and the road sub-catchment area are combined to generate a city secondary sub-catchment area; the method specifically comprises the following steps:
s71, respectively extracting and integrating a trunk road and a trunk pipeline according to road network and pipe network data to generate a trunk drainage network;
the method comprises the following specific steps: collecting road network and pipe network data stored by a node-arc model, connecting the road arc segments according to attribute information of roads of the road network data, deflection angles of the road arc segments and overall trends of the arc segments, wherein the connection mode of the pipeline arc segments is basically the same as that of the road arc segments, the difference between the road arc segments is that the attribute information of the road is a road name and a road width, and the attribute information of the pipeline is pipe diameter, pipe materials and a road to which the pipeline belongs;
the main road and the main pipeline are extracted according to the length and the width of the connecting road, the length and the pipe diameter characteristics of the connecting pipeline, and the main road and the main pipeline are integrated as the ground surface of the road surface on which the pipelines are arranged is communicated with underground confluence; the main drainage network is formed according to a road center line.
S72, dividing a city neighborhood sub-catchment area by using an ArcGIS hydrologic analysis method according to the high-precision DSM;
specifically, the city block refers to an area surrounded by a road in the city, which is generally rectangular, and the division of the city block sub-catchment area refers to a hydrological analysis method of the first-level sub-catchment area, and the steps of filling a depression, calculating flow direction and flow, generating grid river networks, river links, vectorizing the grid river networks, extracting watershed, dividing the catchment area and the like are performed on a high-precision DSM containing building and vegetation data, and after the division of the city block sub-catchment area is completed, the situation that the city block sub-catchment area is converged to the road is discussed;
s73, dividing a main road surface according to a main drainage network, dividing the main road into a plurality of road ranges by the width of the road, and reasonably distributing the road at the intersection into the road ranges according to the main drainage network when the phenomenon of overlapping occurs at the intersection of the roads;
the method comprises the following specific steps:
for the junction of four road sections or three road sections, connecting the intersection point of a trunk drainage network with the top point of a city block, and distributing the road range according to the connection line; when no urban block top point exists at the junction of two road sections or between two trunk drainage lines, the two trunk drainage lines are distributed by taking the angular bisector as a dividing line, so that the road surface of the trunk road is completely divided into a plurality of road ranges, the road is used as a drainage system of a secondary sub-catchment area of the city, and each road range is a road sub-catchment area;
s74, merging the neighborhood sub-catchment area and the road sub-catchment area to generate a city secondary sub-catchment area, wherein the neighborhood sub-catchment area determines the flow accumulation condition of each grid data when dividing, and determines the water outlet of each neighborhood sub-catchment area according to the principle that the maximum flow is the water outlet, wherein one water outlet is positioned at the neighborhood boundary, namely the junction of the neighborhood and the road, and the other water outlet is positioned at the junction of the same neighborhood sub-catchment area;
the sub-catchment areas of the neighborhood with the water outlets at the junction of the neighborhood and the road are called direct catchment areas, the sub-catchment areas of the neighborhood with the water outlets at the junction of the same neighborhood are called indirect catchment areas, the indirect catchment areas and the direct catchment areas are combined according to the flow direction, and then the combined sub-catchment areas and the road sub-catchment areas are combined again according to the positions of the water outlets to form the urban secondary sub-catchment areas.
And S8, dividing the three-level sub-catchment area of the city by combining a Thiessen polygon method on the basis of the two-level sub-catchment area of the city in the step S7, and manually modifying and dividing the three-level sub-catchment area of the city by referring to the high-precision DOM.
Specifically, the division of the three-level sub-catchment area of the city mainly considers the distribution situation of drainage points, the Thiessen polygon method is a space data analysis method and is used for dividing a point set on a plane into a plurality of areas, wherein each area surrounds a specific point, the point is the nearest point to any point in the area, the three-level sub-catchment area of the city is generated by using the Thiessen polygon method, and the final three-level sub-catchment area of the city is formed by referring to the water diversion ridge generated by the DEM and the manual modification of the building outline in the DOM.
The specific steps of dividing the urban three-level sub-catchment area based on the urban two-level sub-catchment area in the step S7 by combining the Thiessen polygon method are as follows:
s81, constructing a three-level sub-catchment area by using a Thiessen polygon method, wherein the principle is as follows:
a set of discrete points L= { L is arranged 1 ,L 2 ,...,L n }, wherein L i The coordinates representing the i point are (x i ,y i ) For any point L j Calculate it and other points L i Distance d (L) j ,L i ) The calculation formula is as follows:
for each point L j Find the nearest point L to it i (i.noteq.j), this is then appliedThe two points are connected to form an edge, the edge becomes a Thiessen edge, and polygons formed by all adjacent Thiessen edges are Thiessen polygons of the point and represent the nearest neighbor relation between the point and surrounding points;
selecting drainage points as an initial discrete point set, generating Thiessen polygons, and ensuring that the distance between the point in each sub-catchment area and the corresponding drainage point is nearest;
s82, manually modifying Thiessen polygons by referring to the watershed generated by the DEM and the building outline in the DOM, and equally dividing the part which does not contain the water outlet into surrounding sub-catchments for the sub-catchments divided by the watershed and the building outline, so that all the sub-catchments are not divided any more, and forming the final three-level sub-catchments of the city.

Claims (10)

1. A method for generating city high-precision digital terrains and dividing sub-catchment areas by using LiDAR is characterized by comprising the following steps:
s1, carrying out data preprocessing on original LiDAR point cloud data and original image data of a certain city;
s2, constructing a TIN model by using the original LiDAR point cloud data preprocessed in the step S1 through a terra model module in terra solid software, and respectively generating a TIN model of a DEM and a DSM;
s3, respectively carrying out interpolation calculation on the DEM and the DSM TIN model generated in the step S2 to obtain the elevation value of the target point, and generating a high-precision DEM and DSM;
s4, defining connection points for the original LiDAR point cloud data and the original image data which are preprocessed in the step S1, and performing space triangular encryption through the connection points to match the original LiDAR point cloud data and the original image data;
s5, defining Color points for the original image data subjected to space triangle encryption in the step S4, and correcting the spelling line to generate a high-precision DOM;
s6, dividing a first-level sub-catchment area of the city by a hydrologic analysis method in combination with the high-precision DEM in the step S3;
s7, dividing a city secondary sub-catchment area by combining a high-precision DSM, a road network and a pipe network on the basis of dividing the city primary sub-catchment area in the step S6;
and S8, dividing the three-level sub-catchment area of the city by combining a Thiessen polygon method on the basis of the two-level sub-catchment area of the city in the step S7, and manually modifying and dividing the three-level sub-catchment area of the city by referring to the high-precision DOM in the step S5.
2. The method for generating urban high-precision digital terrain and sub-catchment area division according to claim 1, wherein the original LiDAR point cloud data in step S1 is a data set comprising a plurality of discrete points generated by the surface features of the ground, buildings, vegetation and rivers of the urban area collected by the LiDAR, the data set comprising coordinates of the points, attributes of the points, density and a spatial reference system; the properties of the points include color, reflection intensity, and classification information; the original image data is unprocessed urban area satellite remote sensing image data and is stored in a TIFF format.
3. The method for generating urban high-precision digital terrain and sub-catchment area division by using LiDAR according to claim 1, wherein the step S1 of preprocessing the original LiDAR point cloud data and the original image data comprises the following steps:
s11, data import: the method comprises the steps of importing original LiDAR point cloud data through a terraScan module of TerraSolid software, performing thinning treatment, and importing original image data through a terraPhoto module;
s12, density detection: displaying a density map of original LiDAR point cloud data through a visualization tool of terraSolid software;
s13, dividing the dense point cloud data into blocks: original LiDAR point cloud data imported in the Terra Scan module are referred to the Terra photo module, and dense original LiDAR point cloud data are divided into blocks through block polygons;
s14, filtering and classifying irregular dense point cloud data: the method comprises the steps of loading corresponding filtering and classifying algorithms according to the types of point cloud data through the point cloud filtering and classifying functions of a Terra Scan module, carrying out progressive encryption filtering processing on original LiDAR point cloud data, removing isolated points, and respectively calculating each point and k' of each point after the isolated points are removedAverage distance d 'of nearest neighbors' i Average distance d ' of k ' nearest neighbors of each point ' i Average value of (2)And standard deviation s, according to the set threshold valuePerforming filtering operation when->In the case of ground points, when +.>And is a non-ground point, where k' and k s All are preset parameters, the ground points are separated from the non-ground points, so that the ground points are extracted, and the non-ground points are reclassified through a ground clearance algorithm;
s15, homogenizing and color-homogenizing and splicing original image data: the terraPhoto module is used for adjusting the numerical value of the image color band and the brightness, the brightness and the contrast of the image to the imported original image data; and automatically splicing the plurality of image pictures into an image in a splicing line mode through a terraPhoto module, and completing the splicing of the images.
4. The method for generating urban high-precision digital terrain and sub-catchment area division according to claim 3, wherein in step S14, the average distance d ' between each point and k ' nearest neighbors of each point after removing isolated points ' i The calculation formulas of (a) are respectively as follows:
in the formula (2), k' is a preset parameter; x is x i ’、y i ’、z i ' represents the coordinates of each point after removal of the isolated point; x is x ij ’、y ij ’、z ij 'represents the coordinates of each point corresponding to the respective k' nearest neighbors after removal of the outlier;
the average distance d ' of the k ' nearest neighbors of each point ' i Average value of (2)The calculation formula of (2) is as follows:
in the formula (3), m represents the total number of remaining points after removing isolated points;
the calculation formula of the standard deviation s is as follows:
5. the method for generating urban high-precision digital terrain and sub-catchment area division according to claim 1, wherein the step of constructing the TIN model in the step S2 is as follows:
s21, performing triangulation algorithm on the original LiDAR point cloud data set to obtain a ground point set DT (P) ground ) The ground point set DT (P ground ) Connecting to generate a TIN model;
s22, adding new ground points by using a random sampling consistency algorithm, and updating a TIN model;
s23, checking triangles in the TIN model by using an elevation difference method;
s24, repeating the steps S22 and S23, continuously adding new ground points, and checking triangles in the TIN model until all the ground points are processed.
6. The method for generating urban high-precision digital terrain and sub-catchment area division according to claim 4, wherein the steps ofGround point set DT (P) ground ) The formula of (2) is as follows:
TIN=DT(P ground ) (5)
in the formula (5), TIN represents a TIN model; DT (P) ground ) Representing a set of ground points.
7. The method for generating urban high-precision digital terrain and sub-catchment area division by using LiDAR according to claim 1, wherein the calculation formula of the elevation value of the target point in the step S3 is as follows:
wherein z represents the elevation value of the target point; z i Representing the elevation value of the interpolation point i; w (w) i The weight between the target point and the interpolation point i is represented.
8. The method for generating urban high-precision digital terrains and dividing sub-catchment areas by using LiDAR according to claim 1, wherein the specific steps of dividing the first-level sub-catchment areas of the city by using the high-precision DEM in the step S6 through a hydrologic analysis method are as follows:
s61, filling: correcting the low-lying areas in the high-precision DEM data through a filling tool so as to simulate the water flow path, the water flow direction and the water flow rate;
s62, calculating the flow direction and flow rate of water flow: simulating the water flow direction by a flow direction tool, so that the water flow flows from the area of Gao Chengda to the area of Gao Chengxiao, simulating the water flow by a flow direction tool, and accumulating the flow from upstream to downstream to obtain the water flow direction and flow data of the whole area;
s63, generating a grid river network: setting a corresponding flow threshold value through a grid calculator tool, screening and reserving rivers with flow larger than the threshold value, and generating a grid river network so as to acquire grid-form river network data;
s64, river link: according to the water flow direction and flow data of the whole area, discrete water flow paths in the grid river network are connected into a continuous river line through a river link tool;
s65, grid river network vectorization: converting river network data in a grid form into vector element representation to obtain water system distribution conditions of the whole urban area;
s66, extracting watershed: firstly, selecting the intersection point of tributaries and main flows in a river network as a water outlet, then using a pouring point capturing tool to set the water outlet as a pouring point, and finally extracting a watershed by a watershed tool;
s67, dividing a catchment area: based on the steps S61-S66, collecting area distribution of the whole urban area is obtained through a collecting area tool, a plurality of adjacent collecting areas are combined to generate sub-collecting areas according to the water system distribution condition of the whole urban area, each collecting area is guaranteed to be included in a unique sub-collecting area until all the collecting areas are processed, and finally the first-level sub-collecting area of the city is obtained.
9. The method for generating urban high-precision digital terrains and sub-catchment area division by using LiDAR according to claim 1, wherein the step of dividing the urban secondary sub-catchment area in step S7 is as follows:
s71, respectively extracting and integrating a trunk road and a trunk pipeline according to road network and pipe network data to generate a trunk drainage network;
s72, dividing the urban neighborhood sub-catchment area based on the high-precision DEM by a hydrologic analysis method according to the step S6;
s73, dividing a main road surface according to the main drainage network generated in the step S71, dividing the main road surface into a plurality of road ranges by the width of the road, and if an overlapping phenomenon occurs at the road intersection, distributing the road with the overlapping road intersection into the road ranges according to the main drainage network, wherein the road is used as a drainage system of a secondary sub-catchment area of the city, and each road range is a sub-catchment area of the road;
s74, determining the flow accumulation condition of each raster data in the city neighborhood sub-catchment area in the step S72 when dividing, determining the water outlet of each city neighborhood sub-catchment area according to the principle that the maximum flow is the water outlet, wherein the neighborhood sub-catchment area with the water outlet at the junction of the neighborhood and the road is called a direct catchment area, the neighborhood sub-catchment area with the water outlet at the junction of the same neighborhood sub-catchment area is called an indirect catchment area, merging the indirect catchment area with the direct catchment area according to the flow direction, and merging the merged sub-catchment area with the road sub-catchment area in the step S73 again according to the position of the water outlet to form a city secondary sub-catchment area.
10. The method for generating urban high-precision digital terrains and dividing sub-catchment areas by using LiDAR according to claim 1, wherein the calculation formula for dividing the urban three-level sub-catchment areas on the basis of the urban two-level sub-catchment areas in step S7 by combining the Thiessen polygon method in step S8 is as follows:
wherein L is i The coordinates representing the i point are (x i ,y i );L j Representing an arbitrary point; d (L) j ,L i ) Representing an arbitrary point L j With other points L i Distance between them.
CN202311574228.7A 2023-11-23 2023-11-23 Method for generating city high-precision digital terrain and dividing sub-catchment area by LiDAR Pending CN117523130A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202311574228.7A CN117523130A (en) 2023-11-23 2023-11-23 Method for generating city high-precision digital terrain and dividing sub-catchment area by LiDAR

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202311574228.7A CN117523130A (en) 2023-11-23 2023-11-23 Method for generating city high-precision digital terrain and dividing sub-catchment area by LiDAR

Publications (1)

Publication Number Publication Date
CN117523130A true CN117523130A (en) 2024-02-06

Family

ID=89749150

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202311574228.7A Pending CN117523130A (en) 2023-11-23 2023-11-23 Method for generating city high-precision digital terrain and dividing sub-catchment area by LiDAR

Country Status (1)

Country Link
CN (1) CN117523130A (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN118071118A (en) * 2024-04-18 2024-05-24 深圳市中科云驰环境科技有限公司 Multi-node sewage discharge management method, device and system

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN118071118A (en) * 2024-04-18 2024-05-24 深圳市中科云驰环境科技有限公司 Multi-node sewage discharge management method, device and system

Similar Documents

Publication Publication Date Title
CN111553963B (en) Meta-grid generation method and device based on geographic information
CN116310192A (en) Urban building three-dimensional model monomer reconstruction method based on point cloud
CN107679498A (en) A kind of airborne laser point cloud downtown roads recognition methods
CN104574512A (en) Multi-scale DEM (digital elevation model) construction method considering topographical semantic information
CN117523130A (en) Method for generating city high-precision digital terrain and dividing sub-catchment area by LiDAR
CN111652241B (en) Building contour extraction method integrating image features and densely matched point cloud features
CN111858810B (en) Modeling elevation point screening method for road DEM construction
CN105956542B (en) High-resolution remote sensing image road extraction method based on statistical matching of structural wire harnesses
CN114580126A (en) Construction method and construction system of urban drainage waterlogging prevention system
CN114648617A (en) Water system extraction method based on digital elevation model DEM
Borisov et al. An automated process of creating 3D city model for monitoring urban infrastructures
CN112700411B (en) Karst tunnel water collection area calculation method based on remote sensing technology
CN113987969A (en) Land-water dynamic simulation model of multi-elevation scale water flow network
CN117012004B (en) Urban waterlogging risk early warning method and system
CN105894553A (en) Street spatial form layout method based on grid selection
CN108062438A (en) A kind of wire selection system and selection method for wind power plant operation road
CN113032877A (en) BIM technology-based optimization method for construction site temporary construction scheme
CN114969944B (en) High-precision road DEM construction method
Xu et al. Methods for the construction of DEMs of artificial slopes considering morphological features and semantic information
Gargallo-Peiró et al. Representing urban geometries for unstructured mesh generation
CN110737931A (en) ArcGIS-based railway bridge crossing hydrological key parameter extraction method
Ganendra et al. The role of airborne LiDAR survey technology in digital transformation
Mandlburger et al. Derivation of a countrywide river network based on Airborne Laser Scanning DEMs-results of a pilot study
Kaubukowski et al. Classification of Large and Heterogeneous LiDAR Data Sets
Zhao et al. Rapid construction algorithm of 3D urban road network from raster maps

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