Real-time calculation method and system for geometric parameters and surface flow field of mountain torrent debris flow

Abstract

The invention discloses a mountain torrent debris flow geometric parameter and surface flow field real-time calculation method and system, and belongs to the technical field of geological disaster monitoring, and the method comprises the following three core steps: S1, precisely calculating the geometric parameters of a riverway, determining the main axis of the riverway through principal component analysis, and calculating the width through a multi-thread parallel convex hull algorithm; efficient depth inversion is realized through a KDTree spatial index, and an Alpha shape is utilized to extract a central area form so as to estimate flow; s2, dynamically analyzing a surface flow field, constructing a grid model by using Delaunay triangulation, rendering and enhancing micro-topographic features based on gradient topographic shadow, and extracting a surface flow velocity field from continuous topographic shadow by using a PIV technology; and S3, data stream processing and transmission architecture: constructing a directed acyclic graph topological structure and a dual-channel data transmission mechanism. According to the invention, the problems of real-time performance and accuracy of debris flow monitoring in a complex mountain area environment are solved, the system delay is less than 300 milliseconds, and reliable technical support can be provided for early warning of mountain torrents and debris flow disasters.

Description

Technical Field

[0001] The present invention belongs to the technical field of geological disaster monitoring, and particularly relates to a real-time calculation method and system for geometric parameters and surface flow fields of mountain torrents and debris flows. Background Art

[0002] As a highly frequent natural disaster, mountain torrents and debris flows are characterized by strong suddenness and great destructive power, seriously threatening the lives and property safety of people in mountainous areas. With the intensification of global climate change and the frequent occurrence of extreme precipitation events, mountain torrent and debris flow disasters show a multiple trend. Traditional monitoring methods for mountain torrents and debris flows mainly rely on single-point monitoring devices such as water level gauges and rain gauges. These methods can only obtain local and one-dimensional monitoring information and cannot comprehensively reflect the three-dimensional spatial characteristics and dynamic evolution process of mountain torrents and debris flows. For example, although traditional contact monitoring methods such as GNSS, total station, GPS, etc. can provide high-precision single-point monitoring data, limited by the number and distribution of monitoring devices, they can only monitor preset fixed points and it is difficult to obtain the overall geometric characteristics of the disaster body. At the same time, under extreme disaster conditions, the damage rate of traditional monitoring devices is extremely high, and the continuity and reliability of monitoring data are difficult to guarantee. In addition, traditional monitoring methods usually require manual on-site layout and maintenance, which not only has a large workload and high cost, but also brings safety hazards to monitoring personnel.

[0003] In recent years, with the rapid development of lidar and computer vision technologies, three-dimensional spatial monitoring methods based on point cloud data have gradually been applied to the field of mountain torrent and debris flow monitoring. Compared with traditional technologies, point cloud data can provide high-density and high-precision three-dimensional spatial information, and can more comprehensively describe the geometric shape and surface characteristics of mountain torrents and debris flows. However, mountain torrent and debris flow point cloud data has characteristics such as large data volume, complex spatio-temporal distribution, and rapid changes, which pose great challenges to real-time processing and analysis. Existing point cloud processing methods are mostly designed for static scenes, with high computational complexity and difficult to meet the real-time requirements in mountain torrent and debris flow monitoring; existing flow field analysis methods are mainly based on manual marker tracking or traditional image processing technologies, with limited accuracy and low automation, and cannot adapt to the continuous monitoring requirements under complex terrain conditions in mountainous areas.

[0004] Existing methods lack the ability to accurately calculate specific geometric parameters of mountain flood debris flows and extract surface flow fields. Existing point cloud segmentation and geometric feature extraction algorithms are mostly designed for regular objects and are difficult to cope with the irregular and dynamically changing characteristics of mountain flood debris flows. Traditional flow field extraction methods are mainly based on images or numerical simulations and lack effective integration with three-dimensional geometric features, resulting in the disconnection between flow field information and geometric information and making it difficult to comprehensively reflect the dynamic evolution process of mountain flood debris flows. Mountain flood debris flow monitoring systems are usually deployed in mountainous environments with complex terrains and poor network conditions, and the real-time transmission and processing of a large amount of high-precision point cloud data face severe challenges. Existing data transmission methods often use a single communication method and a fixed data compression strategy and are difficult to adapt to the complex mountain environment. The traditional batch processing mode is inefficient and cannot meet the requirements of millisecond-level response in mountain flood debris flow monitoring. In terms of data compression, general algorithms do not consider the spatial distribution characteristics of point cloud data. In terms of data transmission, there is a lack of differential transmission strategies for metadata and large data objects, resulting in low overall system efficiency.

[0005] In summary, there are obvious deficiencies in the existing technologies for monitoring the geometric parameters and surface flow fields of mountain flood debris flows based on point cloud data. It is urgent to develop new technical solutions to meet the needs of real-time monitoring and early warning of mountain flood debris flow disasters. Summary of the Invention

[0006] To solve the technical problems existing in the background art, the present invention aims to provide a real-time calculation method and system for geometric parameters and surface flow fields of mountain flood debris flows, which are used for real-time monitoring and analysis of the geometric morphological characteristics and dynamic change processes of mountain flood debris flows, and realize high-precision calculation of geometric parameters such as river channel width, depth, and central area and dynamic characteristics such as surface flow velocity, providing key technical support for mountain flood debris flow disaster monitoring and early warning.

[0007] To solve the technical problems, the technical solution of the present invention is as follows:

[0008] A real-time calculation method for geometric parameters and surface flow fields of mountain flood debris flows, the method comprising:

[0009] S1: Calculation of river channel geometric parameters and flow rate estimation. Collect the reference point cloud data of the riverbed and the dynamic point cloud of the water surface flow for data preprocessing. Use the principal component analysis PCA to determine the main axis of the river channel. Calculate the river channel width through a multi-threaded parallel convex hull algorithm. Use the KDTree spatial index for efficient depth inversion. Extract the central area morphology using the Alpha shape, and output the key geometric parameters such as the main axis, width, depth, and central area morphology of the river channel;

[0010] S2: Dynamic analysis of the surface flow field. Based on the output result of step S1, a Delaunay triangulation is used to construct a surface flow field grid model to optimize the calculation accuracy. The micro-topography features are enhanced by terrain shading rendering based on the gradient to improve the visualization effect of the flow field. The particle image velocimetry (PIV) is adopted to calculate the surface velocity field in combination with continuous terrain shading, accurately capturing the dynamic changes of water flow, and outputting the dynamic parameters of the velocity field, flow direction distribution, average velocity, and maximum velocity.

[0011] S3: Data stream processing and transmission architecture. Based on the output of step S2, a directed acyclic graph (DAG) topological structure is constructed, and a dual-channel data transmission mechanism is adopted. By dynamically calculating the error value, key parameters are extracted and processed in the real-time data stream, and finally, the structured river channel dynamic characteristic data including the time series change rate, debris flow front position, volume estimation, and warning threshold judgment are output.

[0012] Furthermore, in the stage of benchmark point cloud data acquisition, first, high-precision point clouds of the riverbed debris flow channel are obtained through 3D laser scanning as the benchmark for subsequent geometric parameter calculations. In the real-time monitoring stage, dynamic point cloud data is extracted through a point cloud segmentation algorithm, and principal component analysis (PCA) is applied for principal axis orientation:

[0013]

[0014] Among them, p i and p j represent the three-dimensional coordinate point vectors in the point cloud, n is the total number of points in the point cloud, C is the covariance matrix, V is the eigenvector matrix, Λ is the eigenvalue diagonal matrix, and the eigenvalue decomposition equation CV = VΛ is used to solve the principal axis direction. Among the principal axis directions V = [v1, v2, v3] obtained by eigenvalue decomposition, v1 corresponds to the main flow direction of mountain flood debris flow, v2 corresponds to the river channel width direction, and v3 corresponds to the depth direction. The advantage of using PCA is that it can adaptively identify the main flow direction under complex terrain and provide a natural coordinate system for subsequent geometric parameters. The point cloud data is transformed into the principal axis coordinate system through coordinate transformation:

[0015]

[0016] Among them, p' i is the point coordinate after being transformed into the principal axis coordinate system, and V T is the rotation matrix composed of eigenvectors, that is, the transpose of the eigenvector matrix. This transformation aligns the river channel geometric features with the coordinate axes, significantly simplifying the subsequent width and depth calculations.

[0017] Furthermore, the multi-threaded parallel convex hull algorithm calculation includes:

[0018] The point cloud is evenly divided into m segments along the principal axis direction, and the segment boundaries and widths are defined as:

[0019]

[0020] where width j is the width of the j-th river channel segment, H j is the convex hull point set obtained by the Graham scan algorithm for the point cloud of this segment, p′ i [0] represents the coordinate component of the transformed point cloud in the main axis direction, min A is the minimum coordinate value in the main axis direction, range A is the coordinate range in the main axis direction, m is the number of river channel segments, ||h a -h b ||2 represents the Euclidean distance between two points on the convex hull;

[0021] To improve the calculation efficiency, ThreadPoolExecutor is used to implement multi-threaded parallel calculation, and the number of threads is adaptively adjusted according to the number of CPU cores:

[0022]

[0023] where W is the set of widths of all river channel segments, F parallel represents the parallel calculation function, S is the set of point cloud segments evenly divided along the main axis, s j is the main axis coordinate range of the j-th point cloud segment among them;

[0024] The reason for choosing the Graham scan algorithm to construct the convex hull lies in its O(n log n) time complexity and high numerical stability, which is especially suitable for the boundary extraction of irregular river channel cross-sections. Calculating the maximum distance between convex hull points instead of simply the coordinate range accurately represents the actual width of the irregular river channel and avoids measurement errors when the main flow direction is not parallel to the coordinate axes.

[0025] Furthermore, the use of KDTree spatial index for efficient depth inversion includes:

[0026] Constructing the KDTree spatial index of the riverbed reference point cloud to achieve efficient depth inversion from the dynamic point cloud of the water surface flow to the reference point cloud:

[0027]

[0028] where depth i is the depth value of the i-th water surface point, p i [2] represents the Z coordinate component of this point, represents the Z coordinate component of the point in the riverbed reference point cloud that is closest to the horizontal position of this water surface point, T KDIt is a KDTree spatial index constructed based on the horizontal coordinates (x, y) of the riverbed point cloud. pi[:2] represents taking the x and y coordinates of the point, and the query function returns the index of the nearest neighbor point. The KDTree spatial index has a query complexity of O(log n). Through depth statistical analysis, the point with the maximum depth and its position information are obtained:

[0029]

[0030] Among them, p max is the coordinate of the point with the maximum depth, and argmax i represents the point index i that makes the depth expression take the maximum value, and p base,max is the coordinate of the point in the riverbed reference point cloud that is closest to the horizontal position of the point with the maximum depth.

[0031] Furthermore, the extraction of the central region morphology using Alpha shape includes:

[0032] Using the method for constructing the central region contour based on Alpha shape:

[0033]

[0034] Among them, area center is the area of the central region, AlpohaShape represents the Alpha shape algorithm, which is used to construct a complex contour from a non-convex point set. p′ i [0] represents the coordinate of the point in the main axis direction, represents the average coordinate of the point cloud in the main axis direction, δ is the semi-width parameter of the central region, and its value is 1 / 3 of the maximum width. ɑ = 1.0 is the control parameter of the Alpha shape algorithm, which is used to balance the complexity and smoothness of the contour.

[0035] Furthermore, the construction of the surface flow field grid model using Delaunay triangulation includes:

[0036] Using the Delaunay triangulation algorithm to construct an irregular triangular network, which satisfies the circumcircle condition:

[0037]

[0038] Among them, T i represents a triangle, p k represents a point that does not belong to this triangle, c i represents the center point of the circumcircle of the triangle, and v j represents the vertex of the triangle. This formula expresses the circumcircle condition of the Delaunay triangulation: the distance from any point that does not belong to the triangle to the center of the circumcircle is greater than the radius of the circumcircle;

[0039]

[0040] Among them, z is the elevation value of the interpolation point, z1, z2, and z3 are the elevation values of the three vertices of the triangle, A1, A2, and a3 are the areas of the sub-triangles formed by the interpolation point and the other two vertices of the triangle respectively, and the whole expression is a linear interpolation formula based on barycentric coordinates;

[0041] The triangles generated by Delaunay triangulation have the property of maximizing the minimum angle. For the grid boundary or data void area, the inverse distance weighted interpolation method with an adaptive search radius is adopted:

[0042]

[0043] Among them, z is the elevation value of the interpolation point in the void area, z i is the elevation value of the i-th known point, d i is the distance from the interpolation point to the known point, p is the distance power parameter, defaulting to 2. This formula is the inverse distance weighted interpolation method, which is used to process the triangular mesh boundary or data void area.

[0044] Furthermore, the enhancement of micro-topographic features by terrain shading rendering based on gradients includes:

[0045] The terrain shading rendering algorithm based on gradients directly calculates the lighting model:

[0046]

[0047] Among them, hillshade is the terrain shading value, θ z is the altitude angle of the light source, defaulting to 45°, φ is the azimuth angle of the light source, defaulting to 315°, z factor is the elevation scaling factor, defaulting to 2.5, z(i,j) represents the elevation value at the grid point (i,j), res x and res y are the grid resolutions in the x and y directions, and are the gradients in the x and y directions calculated using the central difference method;

[0048] Calculating the gradient using the central difference method can more accurately capture the terrain change trend and reduce numerical errors; the introduction of the elevation scaling factor z_factor enables the slope calculation to have an adjustable vertical exaggeration degree and enhances the terrain expressiveness.

[0049] Furthermore, the adoption of the Particle Image Velocimetry (PIV) method specifically includes:

[0050] Applying the Particle Image Velocimetry (PIV) method to the extraction of the surface flow field of mountain floods and debris flows, processing the terrain shading image and calculating the relevant displacement:

[0051]

[0052] Among them, V is the calculated flow velocity magnitude, Δx and Δy are the displacements obtained from PIV cross-correlation analysis, C is the cross-correlation function, and C i，j represents the value of the cross-correlation function at the (i,j) position, and are sub-pixel displacement correction terms based on Gaussian fitting, scale is the spatial scale factor, dt is the time interval, and represent Fourier transform and inverse transform respectively, represents the conjugate Fourier transform, σ is the standard deviation of the Gaussian kernel, defaulting to 1, r is the radius of the Gaussian kernel, defaulting to 4, and I1 and I2 are two consecutive frames of terrain shadow images;

[0053] The selection of Gaussian filtering is based on the balance consideration of signal-to-noise ratio and detail retention. The setting of the PIV analysis window size and overlap rate achieves the best balance between spatial resolution and computational efficiency, while ensuring the continuity of the velocity field. The Fourier transform converts the convolution in the spatial domain into a product in the frequency domain, and the computational complexity is reduced from O(n 2 ) to O(n log n), significantly improving the real-time performance.

[0054] Furthermore, the step S3 includes:

[0055] Construct a streaming processing framework based on Streamz, decompose the point cloud processing task into a directed acyclic graph topology, and achieve efficient data stream processing:

[0056]

[0057] Among them, E is the event trigger flag, is the indicator function, which is 1 when the condition in the parentheses is satisfied and 0 otherwise. w is the sliding window size, defaulting to 2, and D k represents the point cloud data at time k, and f1, f2,..., f n represent the continuous processing functions in the data stream processing pipeline, represents the function composition operation, represents the derivative with respect to time, used to calculate the rate of change. τ is the trigger threshold, defaulting to 0.05. This formula expresses the event trigger mechanism based on the rate of change threshold, which is used to control the execution timing of complex analysis.

[0058] A real-time calculation system for the geometric parameters and surface flow field of mountain flood debris flow. The system belongs to the above method, and the system includes:

[0059] The river channel geometric parameter calculation and flow rate estimation module is used to collect the riverbed reference point cloud data and the water surface flow dynamic point cloud for data preprocessing. It uses the principal component analysis (PCA) to determine the main axis of the river channel, calculates the river channel width through the multi-threaded parallel convex hull algorithm, performs efficient depth inversion using the KDTree spatial index, and extracts the central region morphology using the Alpha shape, and outputs the key geometric parameters such as the main axis, width, depth, and central region morphology of the river channel;

[0060] The surface flow field dynamic analysis module is used to construct a surface flow field grid model using Delaunay triangulation to optimize the calculation accuracy, enhances the micro-topography features through gradient-based terrain shading rendering to improve the flow field visualization effect, uses the particle image velocimetry (PIV), combines with continuous terrain shading to calculate the surface velocity field, accurately captures the dynamic changes of water flow, and outputs the dynamic parameters such as the velocity field, flow direction distribution, average velocity, and maximum velocity;

[0061] The data stream processing and transmission architecture module is used to construct a directed acyclic graph (DAG) topological structure, adopts a dual-channel data transmission mechanism, extracts and processes key parameters in the real-time data stream by dynamically calculating the error value, and finally outputs the structured river channel dynamic feature data including the time series change rate, debris flow front position, volume estimation, and early warning threshold judgment.

[0062] A computer device includes a memory, a processor, and a computer program stored on the memory and executable on the processor. When the processor executes the program, it implements a method for real-time calculation of the geometric parameters and surface flow field of mountain torrents and debris flows as described in any one of the above.

[0063] A computer-readable storage medium stores a computer program thereon. When the program is executed by a processor, it implements a method for real-time calculation of the geometric parameters and surface flow field of mountain torrents and debris flows as described in any one of the above.

[0064] Compared with the prior art, the advantages of the present invention are as follows:

[0065] Compared with the traditional single-point monitoring technology, the present invention obtains the overall geometric shape of mountain torrents and debris flows through three-dimensional point cloud data, uses the principal component analysis and convex hull algorithm to calculate key parameters such as width and depth, avoids the deviation of single-point monitoring due to position limitations, and significantly improves the comprehensiveness and accuracy of monitoring, especially under complex terrain conditions.

[0066] Compared with the static point cloud processing method, the present invention innovatively combines the PIV technology to realize the dynamic extraction of the surface velocity field. Through terrain shading rendering and sub-pixel displacement correction, it overcomes the limitation that the traditional method can only monitor static geometric features, realizes the real-time monitoring of the dynamic process of mountain torrents and debris flows, and provides more comprehensive dynamic information for disaster early warning.

[0067] Compared with the traditional batch processing mode, the present invention constructs a streaming processing framework based on Streamz. Through a multi-level data stream topology structure and an event-triggered mechanism, combined with ZeroMQ dual-channel transmission and adaptive compression, it solves the problem of real-time processing and transmission of large-scale point cloud data in complex mountain environments, meeting the emergency monitoring requirements of millisecond-level response. BRIEF DESCRIPTION OF THE DRAWINGS

[0068] Figure 1 - The overall technical flow chart of the present invention;

[0069] Figure 2 - Schematic diagram of principal component analysis (PCA) and coordinate transformation;

[0070] Figure 3 - Schematic diagram of river channel cross-section analysis;

[0071] Figure 4 - Schematic diagram of surface flow field extraction. DETAILED DESCRIPTION OF THE INVENTION

[0072] The following describes the specific implementation manners of the present invention in conjunction with embodiments:

[0073] It should be noted that the structures, ratios, sizes, etc. shown in this specification are only used to cooperate with the content disclosed in the specification for those skilled in this technology to understand and read, and are not used to limit the limiting conditions under which the present invention can be implemented. Any modification of the structure, change of the proportional relationship, or adjustment of the size, without affecting the effects that the present invention can produce and the purposes that can be achieved, should still fall within the scope covered by the technical content disclosed in the present invention.

[0074] At the same time, the terms such as "upper", "lower", "left", "right", "middle", and "one" cited in this specification are only for the convenience of clear narration, and are not used to limit the scope under which the present invention can be implemented. The change or adjustment of their relative relationships, without substantial change in the technical content, should also be regarded as the scope within which the present invention can be implemented.

[0075] Example 1:

[0076] The present invention uses the open-source point cloud data of the Illgraben debris flow monitoring station in the Canton of Valais, Switzerland, for testing and verification. This monitoring station is located in the Swiss Alps and has typical mountain river terrain features. Debris flow activities occur periodically, making it an ideal place to study mountain floods and debris flows. In the present invention, "river channel" specifically refers to the narrow gully where debris flows occur in the mountains, including the gully bed (the basal part where the debris flow moves) and the dynamically flowing debris flow body. "Debris flow gully" has the same meaning as "river channel".

[0077] S1. Accurate calculation of river channel geometric parameters

[0078] S1.1. Principal Component Analysis and Coordinate Transformation

[0079] The riverbed reference point cloud is obtained by a 3D laser scanner with a scanning frequency of 10 Hz and a point density of no less than 1000 points per square meter. The static point cloud (riverbed) and the dynamic point cloud (water surface) are distinguished through a point cloud segmentation algorithm. And the principal component analysis (PCA) is applied to determine the main direction of change of the debris flow:

[0080]

[0081] Among the eigenvectors V = [v1, v2, v3] obtained after eigenvalue decomposition, v1 corresponds to the direction with the largest change (the main flow direction of the debris flow), v2 corresponds to the direction with the second largest change, and v3 corresponds to the direction with the smallest change. The calculation of the river channel width is carried out in the direction perpendicular to v1 and parallel to the horizontal plane, while the depth calculation directly uses the z-axis direction of the geodetic coordinate system to ensure consistency with the actual elevation. Through coordinate transformation:

[0082]

[0083] where R is the rotation matrix composed of eigenvectors. Standardization is not performed in the PCA calculation, and the geometric scale information of the original point cloud is retained to ensure that the main axis can accurately represent the river channel trend.

[0084] S1.2. Multi-threaded Parallel Width Calculation

[0085] The point cloud is divided into 20 uniform segments along the main flow direction, and the convex hull is constructed for each segment using the Graham scan algorithm. The Graham scan algorithm first finds the point with the smallest y coordinate as the starting point P0, then calculates the polar angles of the remaining points relative to P0 and sorts them by polar angle; then each point is checked in turn to ensure a counterclockwise turn and maintain the convex hull stack; finally, the convex hull point set CH(S i ) is generated. After the convex hull is constructed, the distances between all pairs of convex hull vertices are calculated, and the maximum distance is the river channel width:

[0086]

[0087] ThreadPoolExecutor is used to implement multi-threaded parallel calculation. The number of threads is adaptively set according to the number of CPU cores, with the core thread number set to 8, the maximum thread number set to 16, and the task queue size set to 100. Through parallel calculation, the processing time is reduced by about 75%, and a single width analysis only takes about 0.3 seconds, achieving centimeter-level width measurement accuracy.

[0088] S1.3. Depth Parameter Inversion and Maximum Depth Location

[0089] Construct a KDTree spatial index for the reference point cloud (riverbed) using the Euclidean distance metric, with a leaf node size of 16 and a search radius of 0.5 meters, to achieve efficient depth inversion from the dynamic point cloud (water surface) to the reference point cloud:

[0090]

[0091] By calculating the elevation difference between the water surface points and the nearest riverbed points, obtain the depth values for finding the maximum depth points:

[0092] d max = max i d i , p max = {p i | d i = d max} (1 - 5)

[0093] The KDTree adopts a balanced tree structure with a search complexity of O(log n), significantly superior to the O(n) complexity of brute-force search. In the processing of point clouds in the millions, the query speed is increased by about 200 times, achieving millimeter-level depth measurement accuracy.

[0094] S1.4, Central Region Morphology Analysis and Discharge Estimation

[0095] To accurately estimate the debris flow discharge, construct a central region interface based on the Alpha shape algorithm:

[0096]

[0097] where δ is the central region half-width parameter, taking a value of 1 / 3 of the maximum width, and the α parameter is set to 1.0. This method effectively filters noise points while retaining the detailed features of the debris flow. By calculating the interface area and average velocity, the discharge can be estimated:

[0098]

[0099]

[0100] where N is the number of valid points in the velocity field; V i is the velocity value of the i-th point. Among them, Q is the estimated discharge, A center is the central interface area, is the average velocity. The voxel downsampling parameter is set to 0.1 meters, which significantly reduces the computational complexity while maintaining the morphological features, and the number of nodes is reduced by about 85%.

[0101] S2. Dynamic Analysis of the Surface Flow Field

[0102] S2.1, Triangulation and Interpolation from Point Cloud to Mesh

[0103] The Delaunay triangulation algorithm is used to construct an irregular triangular network, which satisfies the circumcircle condition:

[0104]

[0105] where c ijl and r ijl are the center and radius of the circumcircle of the triangle respectively. Based on the triangular network, linear interpolation using barycentric coordinates is performed to generate a regular grid with a resolution of 0.2 meters:

[0106]

[0107] where λ i is the barycentric coordinate, A i is the area of the sub-triangle formed by a point and the other two points of the triangle, and A is the area of the whole triangle. For the data void area, inverse distance weighted interpolation is used, and the search radius is set to 3 times the local point density, and the power parameter p = 2 to ensure the continuity and accuracy of the interpolation result.

[0108] S2.2. Terrain Shadow Rendering Algorithm

[0109] Design a terrain shadow rendering algorithm based on gradients, and use the central difference method to calculate the gradients:

[0110]

[0111] Calculate the slope and aspect based on the gradients, and apply the lighting model:

[0112] S = 255·[sin(A alt )·sin(S slope ) + cos(A alt )·cos(S slope )·cos(A az - S aspect )](2 - 4)

[0113] The altitude angle of the light source is set to 45°, the azimuth angle is 315°, and the elevation scaling factor z_factor is set to 2.5. These parameters are determined through experiments and can enhance the terrain details to the greatest extent, facilitating the identification of minor changes and texture features of debris flows.

[0114] S2.3. Extraction of Surface Flow Field Based on PIV

[0115] Apply particle image velocimetry (PIV) to the analysis of continuous terrain shadow images. First, perform Gaussian filtering preprocessing (σ = 1, kernel radius 4), and then perform cross-correlation analysis:

[0116]

[0117] The PIV analysis window size is set to 32×32 pixels, the overlap rate is 50%, and sub-pixel Gaussian fitting is used to determine the peak position of the displacement field, achieving an accuracy of 0.1 pixel. For post-processing, global thresholding (-50, 50) is used to remove abnormal vectors, and the local mean method (4×4 neighborhood) is used to repair missing values.

[0118] It can be understood that by enhancing micro-topographic features through terrain shadow rendering, regarding the micro-topographic changes between consecutive frames as "particles" in the flow field, using Gaussian filtering (σ = 1) for noise reduction, a 32×32 pixel window and a 50% overlap rate are used to obtain the best balance between spatial resolution and continuity, enabling the PIV technology to be successfully adapted to the analysis of debris flow surface flow fields.

[0119] S3. Streaming Processing Architecture and Data Transmission

[0120] A streaming processing framework is constructed based on the Streamz library, decomposing the point cloud processing task into multiple mutually dependent processing nodes. Each node represents a processing operation, and directed edges are formed between nodes through input-output relationships, forming an acyclic graph structure as a whole to ensure the unidirectional flow of data and efficient parallel processing.

[0121] Based on Streamz, a streaming processing framework is constructed to form a directed acyclic graph topological structure. A sliding window (size = 2) is used to detect plane changes, and a complex analysis is triggered when the threshold exceeds 0.05. The data transmission adopts a dual-channel mechanism (JSON metadata + binary data), and zlib compression is enabled for data over 100KB. The overall system latency is <300ms, meeting the requirements of real-time monitoring.

[0122] Example 2:

[0123] This embodiment provides a terminal device, which includes a processor and a memory. The memory is used to store a computer program, and the computer program includes program instructions. The processor is used to execute the program instructions stored in the computer storage medium. The processor may be a Central Processing Unit (CPU), or may also be other general-purpose processors, Digital Signal Processors (DSPs), Application Specific Integrated Circuits (ASICs), Field-Programmable Gate Arrays (FPGAs), or other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components, etc. It is the computing core and control core of the terminal, and is suitable for implementing one or more instructions. Specifically, it is suitable for loading and executing one or more instructions to implement the corresponding method flow or corresponding function. The processor described in the embodiments of the present invention can be used for the operation of a method for real-time calculation of geometric parameters and surface flow fields of mountain torrents and debris flows based on point cloud data, including the following steps:

[0124] S1: Accurate calculation of river channel geometric parameters

[0125] S1.1. Apply principal component analysis to determine the main flow direction of the debris flow, and obtain the main direction through eigenvector decomposition, where v1 represents the main flow direction, the width calculation direction is perpendicular to v1 and parallel to the horizontal plane, and the depth direction uses the z-axis of the geodetic coordinate system;

[0126] S1.2. Divide the point cloud segments along the main flow direction, construct convex hulls for each segment using the Graham scan, calculate the maximum distance of the vertices as the width, and use multi-threaded parallel processing (core thread number 8, maximum thread number 16);

[0127] S1.3. Construct a KDTree index of the reference point cloud (leaf node size 16), calculate the depth through the elevation difference between the dynamic point cloud and the reference point, and determine the position of the maximum depth;

[0128] S1.4. Apply the Alpha shape algorithm to construct the central region interface (δ is 1 / 3 of the maximum width, α = 1.0), and estimate the flow rate in combination with the average flow velocity;

[0129] S2: Dynamic analysis of the surface flow field

[0130] S2.1. Use Delaunay triangulation to construct an irregular triangular network, and interpolate to generate a grid with a resolution of 0.2 meters;

[0131] S2.2. Calculate the terrain shadow based on the gradient, with the light source altitude angle of 45°, azimuth angle of 315°, and elevation scaling factor of 2.5;

[0132] S2.3. Apply particle image velocimetry (window 32×32 pixels, overlap rate 50%), and use sub-pixel fitting to extract the flow velocity field;

[0133] S3: Flow processing architecture and data transmission

[0134] Construct a flow processing architecture (directed acyclic graph topology, sliding window detection), perform dual-channel data transmission (compression is enabled for data above 100KB), and the system latency is less than 300ms.

[0135] Example 3:

[0136] This example provides a storage medium, specifically a computer-readable storage medium (Memory). The computer-readable storage medium is a memory device in a terminal device, used to store programs and data. It can be understood that the computer-readable storage medium here can include both the built-in storage medium in the terminal device and, of course, the extended storage medium supported by the terminal device. The computer-readable storage medium provides storage space, and this storage space stores the operating system of the terminal. And, in this storage space, there is also stored one or more instructions suitable for being loaded and executed by a processor. These instructions can be one or more computer programs (including program codes). It should be noted that the computer-readable storage medium here can be a high-speed RAM memory or a non-volatile memory, such as at least one disk memory.

[0137] One or more instructions stored in the computer-readable storage medium can be loaded and executed by a processor to implement the corresponding steps in the above example regarding a real-time calculation method for geometric parameters and surface flow fields of mountain torrents and debris flows based on point cloud data; one or more instructions in the computer-readable storage medium are loaded and executed by a processor to perform the following steps:

[0138] S1: Accurately calculate river channel geometric parameters

[0139] S1.1. Apply principal component analysis to determine the main flow direction of the debris flow, and obtain the main direction through eigenvector decomposition, where v1 represents the main flow direction. The width calculation direction is perpendicular to v1 and parallel to the horizontal plane, and the depth direction uses the z-axis of the geodetic coordinate system;

[0140] S1.2. Divide the point cloud segments along the main flow direction, construct a convex hull for each segment using Graham scan, calculate the maximum distance of the vertices as the width, and use multi-threaded parallel processing (8 core threads, 16 maximum threads);

[0141] S1.3. Construct a KDTree index for the reference point cloud (leaf node size 16), calculate the depth based on the elevation difference between the dynamic point cloud and the reference point, and determine the maximum depth position.

[0142] S1.4. Apply the Alpha shape algorithm to construct the interface of the central region (δ is 1 / 3 of the maximum width, α = 1.0), and estimate the flow rate in combination with the average flow velocity.

[0143] S2: Dynamic analysis of the surface flow field

[0144] S2.1. Use Delaunay triangulation to construct an irregular triangular network and interpolate to generate a grid with a resolution of 0.2 meters.

[0145] S2.2. Calculate the terrain shadow based on the gradient, with the light source elevation angle of 45°, azimuth angle of 315°, and elevation scaling factor of 2.5.

[0146] S2.3. Apply particle image velocimetry (window 32×32 pixels, overlap rate 50%), and extract the flow velocity field by sub-pixel fitting.

[0147] S3: Flow processing architecture and data transmission

[0148] Construct a flow processing architecture (directed acyclic graph topology, sliding window detection), perform dual-channel data transmission (compression is enabled for data above 100KB), and the system latency is less than 300 ms.

[0149] Those skilled in the art should understand that the embodiments of the present invention can be provided as a method, a system, or a computer program product. Therefore, the present invention can take the form of a complete hardware embodiment, a complete software embodiment, or an embodiment combining software and hardware aspects. Moreover, the present invention can take the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk memory, CD-ROM, optical memory, etc.) containing computer-usable program code.

[0150] The present invention is described with reference to the flowcharts and / or block diagrams of methods, devices (systems), and computer program products according to the embodiments of the present invention. It should be understood that each flow and / or block in the flowcharts and / or block diagrams, and the combination of flows and / or blocks in the flowcharts and / or block diagrams, can be implemented by computer program instructions. These computer program instructions can be provided to the processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing devices to generate a machine, so that the instructions executed by the processor of the computer or other programmable data processing devices generate means for realizing the functions specified in Figure 1 one or more of the flows or multiple flows and / or blocks Figure 1 one or more of the blocks or multiple blocks.

[0151] These computer program instructions can also be stored in a computer-readable memory that can direct a computer or other programmable data processing apparatus to work in a particular manner, such that the instructions stored in the computer-readable memory produce a manufacture including an instruction device that implements the functions specified in one process Figure 1 one or more processes and / or blocks Figure 1 specified in a block or blocks.

[0152] These computer program instructions can also be loaded onto a computer or other programmable data processing apparatus, such that a series of operational steps are performed on the computer or other programmable apparatus to produce a computer-implemented process, whereby the instructions executed on the computer or other programmable apparatus provide steps for implementing the functions specified in one process Figure 1 one or more processes and / or blocks Figure 1 specified in a block or blocks.

[0153] The preferred embodiments of the present invention have been described in detail above. However, the present invention is not limited to the above embodiments, and various changes can be made without departing from the spirit of the present invention within the knowledge of those of ordinary skill in the art.

[0154] Many other changes and modifications can be made without departing from the concept and scope of the present invention. It should be understood that the present invention is not limited to the specific embodiments, and the scope of the present invention is defined by the appended claims.

Claims

1. A real-time calculation method for the geometric parameters and surface flow field of mountain flood debris flow, characterized in that The method includes: S1: River channel geometric parameter calculation and flow estimation. Collect the riverbed reference point cloud data and the dynamic point cloud of the water surface flow for data preprocessing. Use the principal component analysis (PCA) to determine the main axis of the river channel. Calculate the river channel width through the multi-threaded parallel convex hull algorithm. Use the KDTree spatial index for efficient depth inversion. Extract the central region morphology using the Alpha shape. Output the key geometric parameters of the river channel main axis, width, depth, and central region morphology. S2: Dynamic analysis of the surface flow field. Based on the output result of step S1, use the Delaunay triangulation to construct the surface flow field grid model to optimize the calculation accuracy. Enhance the micro-topography features through terrain shading rendering based on the gradient to improve the visualization effect of the flow field. Use the particle image velocimetry (PIV) method, combined with continuous terrain shading to calculate the surface velocity field, accurately capture the dynamic changes of water flow, and output the dynamic parameters of the velocity field, flow direction distribution, average velocity, and maximum velocity. S3: Data stream processing and transmission architecture. Based on the output of step S2, construct a directed acyclic graph (DAG) topological structure. Use a dual-channel data transmission mechanism. Extract and process the key parameters in the real-time data stream by dynamically calculating the error value. Finally, output the structured river channel dynamic feature data including the time series change rate, debris flow front position, volume estimation, and warning threshold judgment.

2. The real-time calculation method for the geometric parameters and surface flow field of mountain flood and debris flow according to claim 1, characterized in that, In the stage of collecting the reference point cloud data, first obtain the high-precision point cloud of the riverbed debris flow channel through 3D laser scanning as the reference for subsequent geometric parameter calculation. In the real-time monitoring stage, extract the dynamic point cloud data through the point cloud segmentation algorithm and apply the principal component analysis (PCA) for main axis orientation: where p i and p j represent the three-dimensional coordinate point vectors in the point cloud, n is the total number of points in the point cloud, C is the covariance matrix, V is the eigenvector matrix, Λ is the eigenvalue diagonal matrix, and the eigenvalue decomposition equation CV = VΛ is used to solve the main axis direction; in the main axis direction V = [v1, v2, v3] obtained by eigenvalue decomposition, v1 corresponds to the main flow direction of mountain torrents and debris flows, v2 corresponds to the river width direction, and v3 corresponds to the depth direction; the advantage of using PCA is that it can adaptively identify the main flow direction under complex terrains and provide a natural coordinate system for subsequent geometric parameters; the point cloud data is transformed to the main axis coordinate system through coordinate transformation: where p i ′ is the point coordinate after being transformed into the principal axis coordinate system, and V T is the rotation matrix composed of eigenvectors, that is, the transpose of the eigenvector matrix. This transformation aligns the river channel geometric features with the coordinate axes, significantly simplifying the subsequent width and depth calculations.

3. A real-time calculation method for geometric parameters and surface flow fields of mountain torrents and debris flows according to claim 1, characterized in that The multi-threaded parallel convex hull algorithm calculation includes: Evenly divide the point cloud into m segments along the main axis direction, and define the segment boundary and width as: where width j is the width of the j-th river channel segment, H j is the convex hull point set obtained by the Graham scan algorithm for the point cloud of this segment, p i ′[0] represents the coordinate component of the transformed point cloud in the main axis direction, min A is the minimum coordinate value in the main axis direction, range A is the coordinate range in the main axis direction, m is the number of river channel segments, ‖h a -h b ‖2 represents the Euclidean distance between two points on the convex hull; To improve the calculation efficiency, use ThreadPoolExecutor to implement multi-threaded parallel calculation, and the number of threads is adaptively adjusted according to the number of CPU cores: Among them, W is the set of widths of all river sections, and F parallel represents the parallel computing function, S is the set of point cloud segments evenly divided along the main axis, and s j is the main axis coordinate range of the j-th point cloud segment among them; The reason for choosing the Graham scan algorithm to construct the convex hull is its O(n log n) time complexity and high numerical stability, which is especially suitable for the boundary extraction of irregular river channel cross-sections. Calculate the maximum distance between the convex hull points instead of simply the coordinate range difference, accurately representing the actual width of the irregular river channel and avoiding measurement errors when the main flow direction is not parallel to the coordinate axes.

4. A real-time calculation method for the geometric parameters and surface flow field of mountain torrents and debris flows according to claim 1, characterized in that The efficient depth inversion using the KDTree spatial index includes: Construct the KDTree spatial index of the riverbed reference point cloud to achieve efficient depth inversion of the dynamic point cloud of the water surface flow to the reference point cloud: where depth i is the depth value of the i-th water surface point, and p i [2] represents the Z coordinate component of this point, represents the Z coordinate component of the point in the riverbed reference point cloud that is closest to the horizontal position of this water surface point, and T KD is the KDTree spatial index constructed based on the horizontal coordinates (x, y) of the riverbed point cloud. p i [:2] represents taking the x and y coordinates of this point. The query function returns the index of the nearest neighbor point. The KDTree spatial index has a query complexity of O(log n). Through depth statistical analysis, the maximum depth point and its position information are obtained: where p max is the coordinate of the maximum depth point, and argmax i represents the point index i that maximizes the depth expression, and p base,max is the coordinate of the point in the riverbed reference point cloud that is closest to the horizontal position of the maximum depth point.

5. A real-time calculation method for the geometric parameters and surface flow field of mountain torrents and debris flows according to claim 1, characterized in that The extraction of the central region morphology using the Alpha shape includes: Use the method of constructing the central region contour based on the Alpha shape: Among them, area center is the area of the central region. AlphaShape represents the Alpha shape algorithm, which is used to construct a complex contour from a non-convex point set. p i ′[0] represents the coordinate of the point in the main axis direction, represents the average coordinate of the point cloud in the main axis direction. δ is the semi-width parameter of the central region, and its value is 1 / 3 of the maximum width. α = 1.0 is the control parameter of the Alpha shape algorithm, which is used to balance the complexity and smoothness of the contour.

6. A real-time calculation method for geometric parameters and surface flow fields of mountain torrents and debris flows according to claim 1, characterized in that, The construction of the surface flow field grid model using the Delaunay triangulation includes: Use the Delaunay triangulation algorithm to construct an irregular triangular network that meets the circumcircle condition: Among them, T i represents a triangle, p k represents a point not belonging to the triangle, c i represents the center point of the circumcircle of the triangle, v j represents the vertex of the triangle. This formula expresses the circumcircle condition of Delaunay triangulation: the distance from any point not belonging to the triangle to the center of the circumcircle is greater than the radius of the circumcircle; Among them, z is the elevation value of the interpolation point, z1, z2, and z3 are the elevation values of the three vertices of the triangle, A1, A2, and A3 are the areas of the sub-triangles formed by the interpolation point and the other two vertices of the triangle respectively, and the whole expression is a linear interpolation formula based on barycentric coordinates; The triangles generated by Delaunay triangulation have the property of maximizing the minimum angle. For the grid boundary or data hole area, an inverse distance weighted interpolation method with an adaptive search radius is adopted: Among them, z is the elevation value of the interpolation point in the void area, z i is the elevation value of the i-th known point, d i is the distance from the interpolation point to the known point, p is the distance power parameter, defaulting to 2. This formula is the inverse distance weighted interpolation method, which is used to handle the boundary of the triangular mesh or the data void area.

7. A real-time calculation method for the geometric parameters and surface flow field of mountain torrents and debris flows according to claim 1, characterized in that The enhancement of micro-topographic features by terrain shading rendering based on gradients includes: The terrain shading rendering algorithm based on gradients directly calculates the lighting model: where hillshade is the terrain shadow value, θ z is the light source altitude angle, default 45°, φ is the light source azimuth angle, default 315°, z factor is the elevation scaling factor, default 2.5, z(i,j) represents the elevation value at grid point (i,j), res x and res y are the grid resolutions in the x and y directions, and are the gradients in the x and y directions calculated using the central difference method; Using the central difference method to calculate the gradient can more accurately capture the terrain change trend and reduce numerical errors; the introduction of the elevation scaling factor z_factor enables the slope calculation to have an adjustable vertical exaggeration degree and enhances the terrain expressiveness.

8. A real-time calculation method for the geometric parameters and surface flow field of mountain torrents and debris flows according to claim 1, characterized in that, The adoption of the particle image velocimetry method PIV specifically includes: Applying the particle image velocimetry method PIV to the extraction of the surface flow field of mountain torrents and debris flows, processing the terrain shading image and calculating the relevant displacement: Among them, V is the calculated flow velocity magnitude, Δx and Δy are the displacements obtained from PIV cross-correlation analysis, C is the cross-correlation function, and C i,j represents the value of the cross-correlation function at the (i,j) position, and are the sub-pixel displacement correction terms based on Gaussian fitting, scale is the spatial scale factor, dt is the time interval, and represent the Fourier transform and the inverse transform respectively, represents the conjugate Fourier transform, σ is the standard deviation of the Gaussian kernel, defaulting to 1, r is the radius of the Gaussian kernel, defaulting to 4, and I1 and I2 are two consecutive frames of terrain shadow images; The selection of Gaussian filtering is based on the balance consideration of signal-to-noise ratio and detail preservation. The setting of the PIV analysis window size and overlap rate achieves the best balance between spatial resolution and computational efficiency, while ensuring the continuity of the velocity field. The Fourier transform converts the convolution in the spatial domain into a product in the frequency domain, and the computational complexity is reduced from O(n 2 ) to O(n log n), significantly improving the real-time performance.

9. A real-time calculation method for the geometric parameters and surface flow field of mountain torrents and debris flows according to claim 1, characterized in that, The step S3 includes: Constructing a streaming processing framework based on Streamz, decomposing the point cloud processing task into a directed acyclic graph topological structure, and realizing efficient data stream processing: Among them, E is the event trigger flag, is an indicator function that has a value of 1 when the condition in the parentheses is satisfied and 0 otherwise. w is the size of the sliding window, defaulting to 2. D k represents the point cloud data at time k, f1, f2,..., f n represent the continuous processing functions in the data stream processing pipeline, represents the function composition operation, represents the derivative with respect to time, used to calculate the rate of change. τ is the trigger threshold, defaulting to 0.

05. This formula expresses the event trigger mechanism based on the rate of change threshold, which is used to control the execution timing of complex analysis.

10. A real-time calculation system for the geometric parameters and surface flow field of mountain torrents and debris flows, characterized in that, The system is applied to the method described in any one of claims 1-9. The system includes: A river channel geometric parameter calculation and flow rate estimation module, which is used to collect the riverbed reference point cloud data and the water surface flow dynamic point cloud for data preprocessing, determine the main axis of the river channel by using the principal component analysis PCA, calculate the river channel width through a multi-threaded parallel convex hull algorithm, perform efficient depth inversion by using the KDTree spatial index, and extract the central area morphology by using the Alpha shape, and output the key geometric parameters such as the main axis, width, depth, and central area morphology of the river channel; A surface flow field dynamic analysis module, which is used to construct a surface flow field grid model by using Delaunay triangulation, optimize the calculation accuracy, enhance the micro-topographic features by terrain shading rendering based on gradients, improve the visualization effect of the flow field, adopt the particle image velocimetry method PIV, and combine the continuous terrain shading to calculate the surface velocity field, accurately capture the dynamic changes of water flow, and output the kinetic parameters such as the velocity field, flow direction distribution, average velocity, and maximum velocity; A data stream processing and transmission architecture module, which is used to construct a directed acyclic graph DAG topological structure, adopt a dual-channel data transmission mechanism, extract and process the key parameters in the real-time data stream by dynamically calculating the error value, and finally output the structured river channel dynamic feature data including the time series change rate, debris flow front position, volume estimation, and warning threshold judgment.

Images