Debris flow collapse object source bottom surface reconstruction and volume estimation method

Abstract

The invention provides a debris flow collapse sliding object source bottom surface reconstruction and volume estimation method, and relates to the technical field of geological disaster monitoring and prevention. The method comprises the following steps: acquiring a vector boundary and elevation data of a debris flow collapse source, and unifying coordinates to construct a basic data set; fitting a main trend constraint function in a sliding direction by utilizing elevation data points in a debris flow collapse and slip object source buffer area, and generating a trend fitting elevation vector; fitting the elevation vector based on the structure tensor and the trend, and performing interpolation propagation and expansion to obtain an initial bottom sliding surface elevation; constructing a sparse matrix system, and performing global optimization fitting on the elevation of the bottom sliding surface by using a generalized minimum residual method; and extending the bottom slip surface to the earth surface, and naturally intersecting with the top surface of the digital elevation model to form the volume of the closed slip mass. The problems that according to an existing debris flow collapse sliding object source bottom sliding surface interpolation reconstruction method, interpolation direction control is insufficient, boundary constraint is not rigorous, and a closed volume cannot be formed through natural intersection are solved.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of geological disaster monitoring and prevention, and particularly relates to a method for reconstructing a bottom surface of a debris flow and estimating a volume of the debris flow. BACKGROUND

[0002] Debris flow is a common and serious geological disaster in mountainous areas, and its activity process is often accompanied by dramatic changes in topography. Accurate acquisition of the volume of the debris flow source is of great significance for disaster assessment, prevention engineering design and management effect evaluation. Digital elevation model (DEM) is an important geographic information data for characterizing the topography of debris flow, and is the basic data source for volume estimation, risk assessment and disaster management. Building an accurate bottom sliding surface model is a key factor in determining the accuracy of the volume of the debris flow. At present, the methods for estimating the volume of the debris flow are generally divided into two categories: traditional field measurement methods and remote sensing data assisted methods. The traditional field measurement method usually determines the position of the sliding bottom surface through field reconnaissance, drilling or geological profile drawing. This method is limited by the complex field conditions, large workload, high cost and discrete data in mountainous areas, and it is difficult to meet the needs of large-scale, rapid disaster investigation and monitoring. In recent years, with the rapid development of remote sensing technology, high-precision digital elevation models (DEM) obtained by airborne LiDAR have gradually become an important means of debris flow investigation. The laser data of multiple echoes can effectively penetrate vegetation coverage and quickly and accurately obtain high-precision topographic data in debris flow areas, providing higher quality input information for three-dimensional reconstruction of the bottom sliding surface of the debris flow.

[0003] However, the existing remote sensing data assisted debris flow volume estimation method still has great limitations, especially in the reconstruction accuracy of the three-dimensional spatial structure of the bottom sliding surface. Conventional geological surface modeling methods such as inverse distance weighting (IDW), radial basis function (RBF) interpolation or Kriging interpolation method usually only consider the distance relationship of spatial data points, ignoring the obvious directionality caused by topographic slope, slope and geological structure in the formation process of the geological sliding body, resulting in local unreasonable concave or convex phenomenon of the fitted bottom sliding surface model, which seriously affects the accuracy of subsequent source volume calculation.

[0004] In recent years, some studies have attempted to use the discrete smooth interpolation (DSI) method with direction tensor constraints to construct geological surfaces in order to improve the spatial continuity and rationality of the bottom sliding surface. Although this method improves the fitting quality of the bottom sliding surface to some extent, it still has the following outstanding problems: (1) The single direction of interpolation propagation: existing discrete smooth interpolation (DSI) methods mostly use fixed direction tensors for interpolation, which makes it difficult to fully consider the complex topographic conditions and changes in the direction of geological structure in the sliding area, resulting in poor local fitting results of the interpolation results; (2) Lack of trend information of digital elevation model (DEM): the existing methods do not fully utilize the topographic trend information outside the sliding body area, and do not consider the effective constraint of the trend of the digital elevation model (DEM) surface on the interpolation of the bottom sliding surface, making it difficult for the interpolated surface to naturally connect and close with the digital elevation model (DEM); (3) Lack of effective global optimization strategy: the existing interpolation results are mostly local optimal solutions, and lack of global optimization strategy, making the overall and consistency of the interpolation results poor, and unable to directly form a closed body for volume calculation; (4) Unclear uncertainty range of interpolation results: due to the complex topographic conditions in mountainous areas, the slope and topographic relief are obvious, and the traditional method does not provide quantitative uncertainty of interpolation, which is difficult to use for risk assessment in engineering practice.

[0005] In fact, most current studies mainly use traditional interpolation methods to construct the debris flow and landslide bottom sliding surface, and do not effectively use the geological structure characteristics and digital elevation model (DEM) trend information to guide the interpolation of the bottom sliding surface, and lack appropriate structure tensor guidance strategies to significantly improve the interpolation accuracy and rationality. Accurate construction of the bottom sliding surface model of the landslide body requires consideration of complex and diverse factors: first, the slope and slope direction have a significant impact on the sliding process of the debris flow, and the direction and shape of the bottom sliding surface have significant spatial differences; second, the trend information of the digital elevation model (DEM) has an important constraint on the fitting of the boundary area of the sliding surface, and existing methods generally lack the use of digital elevation model (DEM) trend information to guide interpolation; in addition, the topography inside the landslide area changes dramatically, and the local propagation of the bottom sliding surface interpolation is prone to structural discontinuity, making it difficult for traditional methods to effectively control the interpolation error. Therefore, accurate closed sliding body volume calculation from the ideal landslide bottom sliding surface model requires consideration of many factors. SUMMARY

[0006] In view of the above problems in the prior art, the mudslide and landslide source bottom surface reconstruction and volume estimation method provided by the present application solves the problems of insufficient interpolation direction control, non-strict boundary constraint, and inability to naturally intersect to form a closed volume in the existing mudslide and landslide source bottom sliding surface interpolation reconstruction method.

[0007] To achieve the above objectives, the technical solution adopted by this invention is: a method for reconstructing the bottom surface and estimating the volume of debris flow landslide source, comprising the following steps: S1. Obtain the vector boundary and elevation data of the debris flow landslide source, and construct a basic dataset by unifying coordinates; S2. Based on the constructed basic dataset, using the elevation data points within the debris flow landslide source buffer zone, fit the main trend constraint function in the sliding direction to generate a trend-fitted elevation vector, where the elevation data points include boundary points. S3. Based on the structural tensor and the trend fitting elevation vector, perform direction-guided interpolation propagation to obtain the initial bottom sliding surface elevation; S4. Construct a sparse optimization system that integrates structural tensor direction constraints and trend function control terms. Use the generalized minimum residual method to globally optimize the initial bottom slip surface elevation to obtain the optimal bottom slip surface elevation vector and complete the reconstruction of the bottom surface of the debris flow landslide source. S5. Based on the optimal bottom slip surface elevation vector and the surface of the digital elevation model, construct vertical volume units, estimate the closed volume between them, and obtain the volume of the landslide body.

[0008] The beneficial effects of this invention are as follows: This invention proposes an interpolation method guided by structural tensors. Considering various geological structural factors such as slope, aspect, and elevation data trends, the interpolation direction and propagation process are guided by both boundary trends and structural tensors. This reconstructs a more accurate and reasonable sliding bottom surface, directly forming a closed volume that naturally intersects with the top surface of the digital elevation model (DEM) for volume estimation. This supports the risk analysis and decision-making needs of actual disaster investigation and mitigation projects. Further, S1 includes the following steps: S101. Extract the effective elevation point set from the high-resolution elevation data of the study area, remove null values ​​and outliers, and form a point set. , This represents the set of valid points for elevation data; S102. Obtain the vector boundary of the debris flow landslide source and ensure that the boundary polygon... Closed, where the boundary polygon All vertices satisfy closure and coordinate consistency, and the boundary polygon The set of vertices is: ; in, Represents the set of boundary control points. Indicates the first The x-coordinates of the boundary points Indicates the first The ordinates of the boundary points Indicates the original topographic elevation at the boundary. denotes the number of boundary points, and its value range is , M b denotes the number of boundary points; S103, constructing a control point set for controlling interpolation with the point set: in the buffer zone of the boundary polygon extended outward by a distance d extracted from the buffer zone and merging the boundary control point set and the buffer control point set obtaining the control point set for interpolation constraint : ; ; wherein, denotes the three-dimensional coordinates of the external constraint elevation data points extracted in the buffer zone of the boundary polygon extended outward by a distance d , denotes the number of points in the buffer zone, denotes the number of control points in the control point set , denotes the number of points in the buffer zone, denotes the number of control points in the control point set , denotes the boundary polygon and the buffer zone formed by extending the distance d outward, denotes the three-dimensional coordinates of the control points; S104, constructing a bottom sliding surface interpolation grid , and setting the unknown elevation data points screened from the bottom sliding surface interpolation grid as the interpolation target point set , and the interpolation target point set is the interpolation target point set after removing the known control points from the bottom sliding surface interpolation grid : ; ; wherein, denotes the two-dimensional plane coordinates of a point in the bottom sliding surface interpolation grid, belonging to the real number field , denotes the real number set, denotes the internal region of the boundary polygon , denotes the horizontal and vertical coordinates of the m th interpolation target point, corresponding to the position of unknown elevation in the bottom sliding surface interpolation grid, representing the undetermined elevation of the m interpolation target point, K representing the total number of the interpolation target point set , m representing the number in the interpolation target point set, m ∈{1,2,…,K}.

[0009] The beneficial effects of the above further scheme are: the application realizes the closed control of the interpolation region and the structured organization of the elevation data by constructing the boundary and buffer control point set and extracting the interpolation target point, and provides high-quality data support for the subsequent precise fitting of the bottom sliding surface and trend constraint.

[0010] Further, the S2 comprises the following steps: S201, using the boundary control point set as the sample set of the main trend constraint function; S202, based on the sample set of the main trend constraint function, adopting a second-order polynomial function to perform regression fitting processing on the elevation data points, to form the main trend constraint function : ; wherein, represents the trend elevation of any two-dimensional point in the bottom sliding surface region, represents the two-dimensional coordinates of any to-be-predicted position in the bottom sliding surface interpolation grid, represents the polynomial coefficient to be regressed; S203, solving the fitting coefficient by the least square method, so that the main trend constraint function can fit the continuous form of the terrain outside the boundary: ; wherein, j represents the number of the points in the buffer zone, represents the number of the points in the buffer zone, represents the three-dimensional coordinates of the external constraint elevation data points extracted in the buffer zone outside the boundary polygon by the outward expansion distance d , represents the elevation numerical prediction value at the point ; S204, applying the main trend constraint function to the interpolation target point set , to obtain the trend fitting elevation vector : ; wherein, Indicates the first m The x and y coordinates of each interpolation target point correspond to the location of the unknown elevation in the bottom sliding surface interpolation grid. m Represents the interpolation target point set The number in m ∈{1,2,…,K}, K Represents the interpolation target point set Total points.

[0011] The beneficial effect of the above-mentioned further solution is that the present invention effectively avoids unreasonable fluctuations or depressions at the boundary of the bottom slip surface interpolation by introducing the main trend function constraint.

[0012] Furthermore, step S3 includes the following steps: S301, Based on interpolation target point set At each target point to be interpolated, extract its spatial neighborhood. k neighborhood points And calculate the structure tensor It is used to characterize the main extension direction of the neighborhood point set. This represents the three-dimensional spatial coordinates of the neighborhood points of the m-th interpolation target point. i This represents the first interpolation target point. i Neighboring points, This represents the set of control points used for interpolation constraints. S302. The structure tensor constructed at each interpolation target point Perform eigenvalue decomposition to extract the structure tensor. Maximum eigenvalue corresponding feature vector , to feature vector The main propagation direction for the current interpolation target point: ; in, Represents the structure tensor matrix. Represents the principal eigenvector of unit length; S303, in the main direction vector Under propagation control, a direction-guided interpolation weighting function is constructed by combining Gaussian distance weights for each neighborhood point. Assign weights : ; ; in, Indicates the first m The nth interpolation target point is related to its nth i The structure tensor guides the weights assigned to each neighboring point. Indicates the distance from the current interpolation target point to the th i A vector of neighborhood points, This represents the Gaussian kernel width control parameter. express With the principal direction vector The included angle, Indicates the first i Two-dimensional coordinates of a neighboring point Indicates the first m Two-dimensional coordinates of the interpolation target point K Represents the interpolation target point set Total points m ∈{1,2,…,K}; S304. Combining the weighted elevations of each neighboring point, calculate the elevation value of the interpolated target point to obtain the initial bottom slip surface elevation. .

[0013] The beneficial effects of the above-mentioned further solutions are: the present invention realizes automatic guided interpolation propagation based on the principal direction of the structural tensor, which effectively improves the directional continuity and terrain adaptability of the bottom sliding surface fitting.

[0014] Furthermore, the structural tensor matrix The expression is as follows: ; in, k Indicates the number of neighboring points. i Indicates the first i 1 neighboring point; The expression for the initial bottom surface is as follows: ; in, This represents the initial interpolated elevation value of the bottom slip surface. This represents the predicted elevation value using the trend function. This indicates the trend fusion weight.

[0015] Furthermore, step S4 includes the following steps: S401. Construct a sparse Laplacian matrix weighted by the structure tensor. , used to describe the propagation constraints of elevation; S402. Construct the trend constraint term matrix C and the trend elevation vector. Apply the fitted value of the trend function as a soft constraint to each interpolation target point; S403, Based on Sparse Laplacian Matrix and trend elevation vector A joint optimization system is constructed that integrates structural tensor direction constraints and trend function control terms, and the final solution form is defined. The trend constraint term matrix C is used to construct the weight coefficient matrix of the constraint terms in the joint optimization system equation. S404, Based on Sparse Laplacian Matrix The joint optimization system is transformed into a sparse linear system of equations. Based on the initial bottom slip surface elevation, a generalized minimum residual method is used to perform global optimization fitting of the bottom slip surface, and the optimal bottom slip surface elevation vector is iteratively solved. .

[0016] Furthermore, the expression for the joint optimization system is as follows: ; ; in, This represents the column vector of initial elevation predictions after interpolation propagation. This represents the column vector of the bottom slip surface elevation solution obtained by joint optimization. Indicates the initial bottom slip surface elevation. This represents the trend elevation of any two-dimensional point in the bottom slip surface region, and T represents the transpose operation; For each target point to be interpolated sparse Laplace matrix Each line in The row is constructed as follows: diagonal elements off-diagonal elements ,in, Indicates the first m The nth interpolation target point is related to its nth i The structure tensor guides the weights assigned to each neighboring point. k Indicates the number of neighboring points. Indicates the first i The global index of each neighboring point, except for diagonal elements and off-diagonal elements, is 0; The expression for the trend constraint term matrix C is as follows: ; in, This represents the penalty factor for the trend term, and I represents the identity matrix. K Represents the interpolation target point set Total points.

[0017] The beneficial effects of the above-mentioned further solutions are: the present invention proposes an interpolation strategy jointly guided by structural tensor and DEM trend, which improves the local continuity and global accuracy of bottom slip surface interpolation.

[0018] Furthermore, step S5 includes the following steps: S501, Calculate the optimal bottom sliding surface elevation vector Grid interpolation is performed within the boundary to construct a continuous bottom surface; S502. Within the fitted bottom slip surface and the top surface of the collapse material source, construct vertical unit columns and use the column volume method to estimate the closed volume between the two surfaces. S503. Based on the estimation results, the top surface of the landslide source and the fitted bottom sliding surface are plotted together in a three-dimensional coordinate system, and the landslide boundary line is marked to show the three-dimensional closed structure of the entire landslide body, forming the landslide volume of the closed landslide body, and completing the volume estimation of the landslide volume.

[0019] The beneficial effects of the above-mentioned further solutions are: the present invention provides a solution for the bottom slip surface to naturally intersect with the top surface of the digital elevation model (DEM) and form a closed volume, which improves the reliability and practicality of the volume estimation of landslide material sources.

[0020] Furthermore, the volume of the landslide body in S5 is calculated as follows: ; in, Indicates the volume of the landslide. Indicates the first Area of ​​each grid cell This represents the elevation data of the earth's surface. This indicates the fitted bottom sliding surface elevation. K Represents the interpolation target point set Total points m Indicates the number in the set of interpolation target points. m ∈{1,2,…,K}. Attached Figure Description

[0021] Figure 1 This is a flowchart of the method of the present invention. Detailed Implementation

[0022] The specific embodiments of the present invention are described below to enable those skilled in the art to understand the present invention. However, it should be understood that the present invention is not limited to the scope of the specific embodiments. For those skilled in the art, various changes are obvious as long as they are within the spirit and scope of the present invention as defined and determined by the appended claims. All inventions utilizing the concept of the present invention are protected.

[0023] Example like Figure 1 As shown, this invention provides a method for reconstructing the bottom surface and estimating the volume of debris flow landslide source, the implementation method of which is as follows: S1. Obtain the vector boundary and elevation data of the debris flow landslide source, and construct a basic dataset using unified coordinates. The implementation method is as follows: S101. Extract the effective elevation point set from the high-resolution elevation data of the study area, remove null values ​​and outliers, and form a point set. , This represents the set of valid points for elevation data. S102. Obtain the vector boundary of the debris flow landslide source and ensure that the boundary polygon... Closed, where the boundary polygon All vertices satisfy closure and coordinate consistency, and the boundary polygon The set of vertices is: ; in, Represents the set of boundary control points. Indicates the first The x-coordinates of the boundary points Indicates the first The ordinates of the boundary points Indicates the original topographic elevation at the boundary. Indicates the number of the boundary point, with a value range of 1. , M b Indicates the number of boundary points; S103. Constructing a point set for control and interpolation: in the boundary polygon Expand outward distance d The set of buffer control points extracted within the buffer and merged boundary control point set and buffer control point set The resulting set of control points for interpolation constraints : ; ; in, Indicates the boundary polygon Expand outward distance d The three-dimensional coordinates of the external constraint elevation data points extracted within the buffer. Indicates the number of the point within the buffer. Represents the control point set Numbering of all control points in the middle, This indicates the number of points within the buffer. Represents the control point set The number of all control points in the system Represents the boundary polygon and distance d The buffer zone formed by outward expansion Represents the three-dimensional coordinates of the control points; S104, construct the bottom sliding surface interpolation grid , and set the unknown elevation data points screened from the bottom sliding surface interpolation grid as the interpolation target point set , the interpolation target point set after excluding the known control points from the bottom sliding surface interpolation grid : ; ; wherein, represents the two-dimensional plane coordinates of a point in the bottom sliding surface interpolation grid, belonging to the real number field , represents a real set, represents the internal area of the boundary polygon , represents the horizontal and vertical coordinates of the Kth interpolation target point, corresponding to the position of the unknown elevation in the bottom sliding surface interpolation grid, m represents the pending elevation of the Kth interpolation target point, represents the total number of interpolation target points m , K represents the number in the interpolation target point set, , m ∈{1,2,…,K}. m

[0024] In this embodiment, a digital elevation model (DEM) of the debris flow watershed is constructed based on airborne laser radar data, a vector boundary of the debris flow collapse and landslide source is extracted, and coordinate unification processing is performed to obtain a basic data set.

[0025] S2, based on the constructed basic data set, using the elevation data points in the debris flow collapse and landslide source buffer area, fitting the main trend constraint function in the sliding direction, generating a trend fitting elevation vector, wherein the elevation data points include boundary points, and the implementation method is as follows: S201, the boundary control point set is used as a sample set of the main trend constraint function to fit the trend of elevation change outside the boundary, participates in the construction of the constraint term, the sample set is applied in subsequent S202 and S203, using the sample points as regression input data, fitting the main trend constraint function by least square method. The main trend constraint function is used to: guide the sliding surface interpolation propagation (improve the terrain continuity); construct the trend term in the generalized minimum residual method; control the interpolation surface to avoid excessive depression or arching; S202, based on the sample set of the main trend constraint function, a second order polynomial function is used to regress and fit the elevation data points to form the main trend constraint function : ​ ; in, Represents any two-dimensional point in the bottom sliding surface region. The trend of elevation, This represents the two-dimensional coordinates of any position to be predicted within the bottom sliding surface interpolation grid. Represents the coefficients of the polynomial to be regressed; S203. Solve for the fitting coefficients using the least squares method. Make the main trend constraint function Capable of fitting continuous terrain morphology outside the boundary: ; in, j Indicates the number of the point within the buffer. This indicates the number of points within the buffer. Indicates the boundary polygon Expand outward distance d The three-dimensional coordinates of the external constraint elevation data points extracted within the buffer. Indicates at point Predicted elevation values ​​at the location; S204, Apply the main trend constraint function Applied to interpolation target point set The trend-fitted elevation vector is obtained. : ; in, Indicates the first m The x and y coordinates of each interpolation target point correspond to the location of the unknown elevation in the bottom sliding surface interpolation grid. m Represents the interpolation target point set The number in m ∈{1,2,…,K}, K Represents the interpolation target point set Total points.

[0026] In this embodiment, a trend-fitted elevation vector is constructed by utilizing the boundary of the landslide source vector and the elevation data trend information.

[0027] S3. Based on the structural tensor and the trend-fitted elevation vector, direction-guided interpolation propagation is performed to obtain the initial bottom slip surface elevation. The implementation method is as follows: S301, Based on interpolation target point set At each target point to be interpolated, extract its spatial neighborhood. k neighborhood points And calculate the structure tensor It is used to characterize the main extension direction of the neighborhood point set. Indicates the firstm The three-dimensional spatial coordinates of the neighborhood points of the interpolation target point i This represents the first interpolation target point. i Neighboring points, This represents the set of control points used for interpolation constraints. S302. The structure tensor constructed at each interpolation target point Perform eigenvalue decomposition to extract the structure tensor. Maximum eigenvalue Corresponding feature vector , to feature vector The main propagation direction for the current interpolation target point: ; in, Represents the structure tensor matrix. Represents the principal eigenvector of unit length; S303, in the main direction vector Under propagation control, a direction-guided interpolation weighting function is constructed by combining Gaussian distance weights for each neighborhood point. Assign weights : ; ; in, Indicates the first m The nth interpolation target point is related to its nth i The structure tensor guides the weights assigned to each neighboring point. Indicates the distance from the current interpolation target point to the th i A vector of neighborhood points, This represents the Gaussian kernel width control parameter. express With the principal direction vector The included angle, Indicates the first i Two-dimensional coordinates of neighboring points Indicates the first m Two-dimensional coordinates of the interpolation target point K Represents the interpolation target point set Total points m ∈{1,2,…,K}; S304. Combining the weighted elevations of each neighboring point, calculate the elevation value of the interpolated target point to obtain the initial bottom slip surface elevation. .

[0028] In this embodiment, the KD-Tree spatial indexing method is used in S301 at each target point to be interpolated. At the place, using KD-Tree index structure constructed by known point set, calling query function query() to search its nearest k neighbor point in space quickly , as its interpolation support point set.

[0029] In this embodiment, the expression of structure tensor matrix is as follows: ; Wherein, k N represents the number of neighbor points, i and i represents the i-th neighbor point.

[0030] In this embodiment, a round-by-round propagation interpolation strategy is adopted, and the interpolation target point calculated in the last round is taken as the neighbor control point in the new round to participate in the subsequent interpolation process, so as to gradually promote the interpolation range and effectively solve the problem of interpolation interruption caused by the sparsity of the middle control point. The expression of the interpolation target point is as follows: ; Wherein, represents the initial bottom sliding surface interpolation elevation value, represents the trend function elevation prediction value, and represents the trend fusion weight.

[0031] S4, a sparse optimization system fusing the structure tensor direction constraint and the trend function control term is constructed, and the generalized minimum residual method is adopted to globally optimize the initial bottom sliding surface elevation, so as to obtain the optimal bottom sliding surface elevation vector and complete the reconstruction of the debris flow collapse source bottom surface. The implementation method is as follows: S401, a structure tensor weighted sparse Laplace matrix is constructed, which is used to describe the propagation constraint relationship of the elevation; S402, a trend constraint term matrix C and a trend elevation vector are constructed, and the trend function fitting value is applied to each interpolation target point as a soft constraint; S403, based on the sparse Laplace matrix and the trend elevation vector , a joint optimization system fusing the structure tensor direction constraint and the trend function control term is constructed, and the final solving form is defined, wherein the trend constraint term matrix C is used to construct the weight coefficient matrix of the constraint term in the joint optimization system equation; S404, based on the sparse Laplace matrix , the joint optimization system is converted into a sparse linear equation set, and based on the initial bottom sliding surface elevation, the generalized minimum residual method is used to globally optimize and fit the bottom sliding surface, and the optimal bottom sliding surface elevation vector is iteratively solved.

[0032] In this embodiment, the expression for the joint optimization system is as follows: ; ; in, This represents the column vector of initial elevation predictions after interpolation propagation. This represents the column vector of the bottom slip surface elevation solution obtained by joint optimization. Indicates the initial bottom slip surface elevation. This represents the trend elevation of any two-dimensional point in the bottom slip surface region, and T represents the transpose operation; For each target point to be interpolated sparse Laplace matrix Each line in The row is constructed as follows: diagonal elements off-diagonal elements ,in, Indicates the first m The nth interpolation target point is related to its nth i The structure tensor guides the weights assigned to each neighboring point. k Indicates the number of neighboring points. Indicates the first i The global index of each neighboring point, except for diagonal elements and off-diagonal elements, is 0; The expression for the trend constraint term matrix C is as follows: ; in, This represents the penalty factor for the trend term, and I represents the identity matrix. K Represents the interpolation target point set Total points; S5. Based on the optimal bottom slip surface elevation vector and the surface of the digital elevation model, vertical volume units are constructed, and the closed volumes between them are estimated to obtain the volume of the landslide body. The implementation method is as follows: S501, Calculate the optimal bottom sliding surface elevation vector Grid interpolation is performed within the boundary to construct a continuous bottom surface; S502. Within the fitted bottom slip surface and the top surface of the collapse material source, construct vertical unit columns and use the column volume method to estimate the closed volume between the two surfaces. S503. Based on the estimation results, the top surface of the landslide source and the fitted bottom sliding surface are plotted together in a three-dimensional coordinate system, and the landslide boundary line is marked to show the three-dimensional closed structure of the entire landslide body, forming the landslide volume of the closed landslide body, and completing the volume estimation of the landslide volume.

[0033] In this embodiment, the volume of the landslide body is calculated as follows: wherein, represents the volume of the landslide, represents the area of the first grid unit, represents the ground surface elevation of the elevation data, represents the fitted bottom sliding surface elevation, K represents the total number of the interpolation target point set, represents the number in the interpolation target point set, m ∈{1,2,…,K}. m

[0034] In summary, the present application solves the problems of insufficient control of interpolation direction, non-strict boundary constraint and inability to naturally intersect to form a closed volume in the existing landslide source bottom sliding surface interpolation reconstruction method through the above design.​​

Claims

1. A method for reconstructing the source bottom surface of a debris flow and estimating the volume of the debris flow, characterized in that, The method comprises the following steps: S1, obtaining vector boundary and elevation data of the debris flow collapse source, and constructing a basic data set by unifying coordinates; S2, based on the constructed basic data set, fitting a main trend constraint function in the sliding direction using elevation data points in the buffer zone of the debris flow collapse source, to generate a trend-fitted elevation vector, wherein the elevation data points include boundary points; S3, based on the structure tensor and the trend-fitted elevation vector, performing direction-guided interpolation propagation to obtain an initial bottom sliding surface elevation; S4, constructing a sparse optimization system that fuses the direction constraint of the structure tensor and the control term of the trend function, and performing global optimization on the initial bottom sliding surface elevation using a generalized minimum residual method to obtain an optimal bottom sliding surface elevation vector, thereby completing reconstruction of the bottom surface of the debris flow collapse source; S5, based on the optimal bottom sliding surface elevation vector and a digital elevation model ground surface, constructing a vertical volume unit to estimate the closed volume therebetween, thereby obtaining a collapse body volume.

2. The method according to claim 1, characterized in that, The S1 comprises the following steps: S101、Extract the effective elevation point set in the high-resolution elevation data of the study area, eliminate null values and abnormal points, and form a point set , represent the effective point set of the elevation data; S102, acquire the vector boundary of the debris flow collapse source, and ensure the boundary polygon is closed, wherein all vertices of the boundary polygon satisfy the closure and the coordinate consistency, the vertex set of the boundary polygon is: ； wherein, denotes the set of boundary control points, denotes the horizontal coordinate of the th boundary point, denotes the vertical coordinate of the th boundary point, denotes the original terrain elevation at the boundary, denotes the number of the boundary point, which takes the value range of , M b denotes the number of boundary points; S103, constructing a control point set for controlling the interpolation with the point set: in the boundary polygon outward expansion distance d buffer control point set extracted in the buffer zone and merging the boundary control point set and the buffer control point set the resulting control point set for interpolation constraints : ； ； in, Indicates the boundary polygon Expand outward distance d The three-dimensional coordinates of the external constraint elevation data points extracted within the buffer. Indicates the number of the point within the buffer. Represents the control point set Numbering of all control points in the middle, This indicates the number of points within the buffer. Represents the control point set The number of all control points in the system Represents the boundary polygon and distance d The buffer zone formed by outward expansion Represents the three-dimensional coordinates of the control points; S104, constructing a bottom sliding surface interpolation grid , and setting unknown elevation data points screened from the bottom sliding surface interpolation grid as an interpolation target point set , the interpolation target point set after excluding known control points from the bottom sliding surface interpolation grid ​ ； ； wherein, represents the two-dimensional plane coordinates of a certain point in the bottom surface interpolation grid, belonging to the real number field , represents the real set, represents the internal area of the boundary polygon , represents the horizontal and vertical coordinates of the first m interpolation target point, corresponding to the position of the unknown elevation in the bottom surface interpolation grid, represents the pending elevation of the first m interpolation target point, K represents the total number of points in the interpolation target point set , m represents the number in the interpolation target point set, m ∈{1,2,…,K}.

3. The method according to claim 2, characterized in that, The S2 comprises the following steps: S201、obtain a set of boundary control points a sample set used as a main trend constraint function; S202, based on the sample set of the main trend constraint function, a second-order polynomial function is used to regress and fit the elevation data points, and a main trend constraint function is formed : ； wherein, denotes the trend elevation of an arbitrary two-dimensional point in the base-slip surface region, denotes the two-dimensional coordinates of an arbitrary location to be predicted in the base-slip surface interpolation grid, denotes the polynomial coefficients to be regressed; S203, solving the fitting coefficients by least square method making the main trend constraint function enable to fit the continuous form of the terrain outside the boundary ； wherein, j denotes the number of points within the buffer, denotes the number of points within the buffer, denotes the elevation value at the point the outward expansion distance d the three-dimensional coordinates of the external constraint elevation data points extracted within the buffer, denotes the elevation number prediction value at the point denotes the elevation number prediction value at the point S204, the main trend constraint function applied to the interpolation target point set , to obtain the trend fitting elevation vector : ； wherein, denotes the horizontal and vertical coordinates of the m interpolation target point, corresponding to the position of the unknown elevation in the base-slip surface interpolation grid, m denotes the number of the interpolation target point in the set of interpolation target points , m ∈{1,2,…,K}, K denotes the total number of the set of interpolation target points .

4. The method of claim 1, wherein, The S3 comprises the following steps: S301, Based on interpolation target point set At each target point to be interpolated, extract its spatial neighborhood. k neighborhood points And calculate the structure tensor It is used to characterize the main extension direction of the neighborhood point set. This represents the three-dimensional spatial coordinates of the neighborhood points of the m-th interpolation target point. i This represents the first interpolation target point. i Neighboring points, This represents the set of control points used for interpolation constraints. S302, construct the structure tensor at each interpolation target point Perform eigenvalue decomposition to extract the structure tensor The maximum eigenvalue The corresponding eigenvector The eigenvector is taken as the main propagation direction of the current interpolation target point The maximum eigenvalue ； wherein denotes the structure tensor matrix, denotes the principal eigenvector per unit length; S303、In the main direction vector Under the control of the propagation, combined with the Gaussian distance weight, the interpolation weighted function of direction guidance is constructed, and the interpolation weighted function of direction guidance is constructed Assign weights : ； ； wherein, denotes the structure tensor steering weight assigned by the m th interpolation target point to its i th neighborhood point, denotes the vector from the current interpolation target point to the i th neighborhood point, denotes the Gaussian kernel width control parameter, denotes the angle between the principal direction vector and the main direction vector denotes the two-dimensional coordinate of the i th neighborhood point, denotes the two-dimensional coordinate of the m th interpolation target point, K denotes the total number of interpolation target points in the set m ∈ {1,2,..., K}; S304、Combining the weighted elevations of each neighborhood point, calculate the elevation value of the interpolation target point, get the initial bottom sliding surface elevation .

5. The method according to claim 4, characterized in that, The structure tensor matrix The expression is as follows: ； wherein, k represents the number of neighborhood points, i represents the first i neighborhood point; The expression of the initial bottom sliding surface is as follows: ； wherein, represents an initial base-surface interpolation elevation value, represents a trend function elevation prediction value, represents a trend blending weight.

6. The method of Figure 1, wherein the method further comprises: The S4 comprises the following steps: S401、constructing a sparse Laplacian matrix weighted by structure tensor , for describing the propagation constraint relationship of the elevation; S402、constructing a trend constraint term matrix C and a trend elevation vector applying the trend function fitting value as a soft constraint to each interpolation target point; S403、based on sparse Laplacian matrix and trend elevation vector , construct a joint optimization system which combines the structural tensor direction constraint and the trend function control term, and define the final solution form, wherein the trend constraint term matrix C is used to construct the weight coefficient matrix of the constraint term in the equation of the joint optimization system; S404、based on the sparse Laplacian matrix The joint optimization system is converted into a sparse linear equation set, and based on the initial bottom sliding surface elevation, a global optimization fitting of the bottom sliding surface is performed using a generalized minimum residual method to iteratively solve an optimal bottom sliding surface elevation vector .

7. The method of claim 6, wherein, The expression of the joint optimization system is as follows: ； ； wherein, represents the initial elevation prediction vector after interpolation propagation, represents the joint optimization obtained bottom-surface elevation solution vector, represents the initial bottom-surface elevation, represents the trend elevation of any two-dimensional point in the bottom-surface region, and T represents the transpose operation; For each target point to be interpolated sparse Laplace matrix Each line in The row is constructed as follows: diagonal elements off-diagonal elements ,in, Indicates the first m The nth interpolation target point is related to its nth i The structure tensor guides the weights assigned to each neighboring point. k Indicates the number of neighboring points. Indicates the first i The global index of each neighboring point, except for diagonal elements and off-diagonal elements, is 0; The expression of the trend constraint term matrix C is as follows: ； wherein, is a penalty factor for the trend term, I is an identity matrix, K is the total number of interpolation target points. is the total number of interpolation target points.

8. The method of claim 1, wherein, The S5 comprises the following steps: S501、constructing the optimal base-slip surface elevation vector Grid interpolation is performed within the boundary to construct a continuous base-slip surface. S502, constructing a vertical unit column within the fitted bottom sliding surface and the top surface of the collapse source, and using a column volume method to estimate the closed volume between the two surfaces; S503, based on the estimation result, drawing the top surface of the collapse source and the fitted bottom sliding surface in a three-dimensional coordinate system, marking a collapse boundary line, displaying the entire sliding body three-dimensional closed structure, forming a collapse body volume of the closed sliding body, and completing volume estimation of the collapse volume.

9. The method of claim 8, wherein, The calculation of the collapse body volume in the S5 is as follows: ； in, Indicates the volume of the landslide. Indicates the first Area of ​​each grid cell This represents the elevation data of the earth's surface. This indicates the fitted bottom sliding surface elevation. K Represents the interpolation target point set Total points m Indicates the number in the set of interpolation target points. m ∈{1,2,…,K}.

Images