Method for determining riverbed morphological feature parameters based on river digital terrain
By constructing an orthogonal curvilinear coordinate system and combining remote sensing images and digital elevation models, the riverbed morphology feature parameters are automatically extracted, solving the problem of discontinuity in the extraction of existing technologies and achieving efficient and accurate acquisition of riverbed morphology feature parameters.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- 中国雅江集团有限公司
- Filing Date
- 2026-03-25
- Publication Date
- 2026-06-26
AI Technical Summary
Existing technologies struggle to efficiently and automatically extract continuous riverbed morphology parameters from natural river channels, resulting in issues such as strong human subjectivity, low efficiency, and discontinuous results.
By integrating multi-source remote sensing images with a digital elevation model of the river channel, an orthogonal curvilinear coordinate system is constructed to automatically extract cross-sections perpendicular to the direction of water flow and calculate riverbed morphological characteristic parameters.
It enables continuous extraction of riverbed morphology parameters across the entire river section, improving the accuracy and automation of cross-section extraction, reducing data acquisition and processing costs, and is suitable for precise analysis of meandering river sections.
Smart Images
Figure CN122289353A_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of river channel management, and more specifically, to a method for determining riverbed morphological characteristic parameters based on digital topography of rivers. Background Technology
[0002] Riverbed morphology parameters are core indicators describing the geometric shape of rivers, encompassing three main categories: planar morphology, cross-sectional morphology, and thoroughfare longitudinal profile morphology. These parameters play a fundamental role in engineering practice and scientific research, including river classification, analysis of riverbed evolution patterns, flood capacity assessment, and numerical simulation of water and sediment transport. Traditional methods for obtaining riverbed morphology parameters mainly rely on the distance between the starting points and elevation data of fixed cross-sections. However, natural river topography is complex and varied, and fixed cross-sections are often spaced far apart, resulting in a discrete distribution of parameters that fail to accurately reflect the continuous changes in riverbed morphology between adjacent cross-sections, leading to significant spatial representativeness issues.
[0003] Existing technologies attempt to extract cross-sectional data based on digital river topography, but they generally suffer from drawbacks such as the need for manual processing of river types and sections, the fact that manually determined cross-sectional lines are not perpendicular to the water flow direction, low work efficiency, and the inability to directly calculate various morphological characteristic parameters. Other technical solutions propose extracting cross-sections based on point cloud data, but for natural rivers, the cost of simultaneously acquiring high-precision point cloud data both above and below water is extremely high, and data processing is difficult, making it difficult to promote and apply in large-scale engineering projects. Therefore, there is an urgent need for a method that can work with digital river topography to automatically extract cross-sections perpendicular to the water flow direction from large volumes of irregularly arranged river topographic data, and to calculate the continuous changes in riverbed morphological characteristic parameters along the river in real time, overcoming the shortcomings of existing technologies such as strong human subjectivity, low efficiency, and discontinuous results. Summary of the Invention
[0004] This invention provides a method for determining riverbed morphology feature parameters based on digital river topography. Addressing the problem that existing hydrological measurement data cannot directly and efficiently extract continuous riverbed morphology feature parameters, this invention integrates multi-source remote sensing images with a digital elevation model of the river channel to construct an orthogonal curvilinear coordinate system, thereby achieving automated, continuous, and high-precision extraction of riverbed morphology feature parameters.
[0005] The technical solution adopted in this invention is: A method for determining riverbed morphology characteristic parameters based on digital river topography includes the following steps: Acquire remote sensing image data, determine the corresponding time period of the study river section's left and right bank boundary lines, and obtain a planar scatter set of the bank boundary lines; Obtain the river topographic map, convert it into a regular rectangular digital elevation model file, and extract the elevation data of the riverbank boundary lines based on the planar scatter point set to construct a three-dimensional coordinate set of the riverbank boundary lines; The three-dimensional coordinate set of the two riverbank boundary lines is fitted to the riverbank line. The scatter point set of the river center is determined by searching the closest scatter point pair on the riverbank line. The river center line is fitted to generate the river center line and the cross-sectional line perpendicular to the river center line is extracted. An orthogonal curve coordinate system is established along the direction of the river center line, perpendicular to the direction of the river center line and the direction of water depth. Under the orthogonal curvilinear coordinate system, the river channel planar morphology characteristic parameters are determined; wherein, the river channel planar morphology characteristic parameters include the width of the flat riverbed at each cross-section, the curvature coefficient along the river centerline, the relative curvature along the river centerline, and the tortuosity coefficient of the bend. Based on the digital elevation model file, the elevation point matrix is obtained, and the morphological characteristic parameters of the river channel cross section and the deep channel longitudinal section are determined in the orthogonal curvilinear coordinate system. Among them, the morphological characteristic parameters of the river channel cross section include the flat area, flat water depth, width-to-depth ratio and wetted perimeter of each section; the morphological characteristic parameters of the deep channel longitudinal section include the location and elevation of the deep channel points of each section and the deep channel longitudinal ratio.
[0006] Furthermore, the acquisition of remote sensing image data, determination of the boundary lines on both banks of the river section under study at the corresponding time, and acquisition of a planar scatter set of the boundary lines on both banks are specifically as follows: Download remote sensing images of the river section during the high water period of the flood season in the corresponding year, and synthesize multiple different bands based on remote sensing image processing software; By selecting the green and mid-infrared bands, an improved normalized differential water index was calculated based on the water index calculation formula; the calculation formula is as follows: MNDWI = (Green - MIR) / (Green + MIR) In the formula, MNDWI is the improved normalized differential water index, with a value of [-1,1]. This index threshold is used to identify land and water areas and extract water areas; Green and MIR represent the reflectance of the green light band and the mid-infrared band, respectively; in Landsat 5, 7, and 8 satellites, Green corresponds to the B2, B2, and B3 bands, respectively, and MIR corresponds to the B5, B6, and B6 bands, respectively. A binary river topographic map is generated based on the improved normalized differential water index calculation results, and the boundary lines of the left and right banks are extracted using vectorization tools. The extracted left and right boundary lines are divided into scattered points at equal intervals. After being assigned planar coordinates, planar scattered point sets are obtained, namely the scattered point set of the left boundary line and the scattered point set of the right boundary line.
[0007] Furthermore, the process of acquiring the river topographic map, converting it into a regular rectangular digital elevation model file, and extracting the elevation data of the riverbank boundaries based on the planar scatter set to construct a three-dimensional coordinate set of the riverbank boundaries specifically involves: Obtain the river topographic map, determine the minimum and maximum values of the horizontal and vertical coordinates of the elevation points on the river topographic map, and construct a regular rectangular range containing four virtual corner points; The elevation of the virtual corner point is assigned the maximum elevation value on the topographic map; Spatial interpolation is performed between the measured elevation points and the virtual corner points to generate a regular rectangular digital elevation model file in Geo-TIFF format; The elevation values of the corresponding locations of the plane scatter point set are extracted from the digital elevation model file and combined to form a three-dimensional coordinate set of the two bank boundary lines containing plane coordinates and elevation information, namely the three-dimensional coordinate set of the left bank boundary line and the three-dimensional coordinate set of the right bank boundary line.
[0008] Furthermore, the three-dimensional coordinate set of the two riverbank boundaries is fitted to the riverbank line. By searching for the nearest scatter point pair on the riverbank line, the scatter point set of the river center is determined. The river centerline is then fitted and extracted, and the cross-sectional line perpendicular to the river centerline is extracted. An orthogonal curve coordinate system is established along the direction of the river centerline, perpendicular to the river centerline, and in the direction of water depth. Specifically: Using the spline curve tool in Matlab, the original scatter set of the boundary lines on both sides of the river was fitted into two riverbank lines. Search for the nearest scatter point pair on the two fitted riverbanks to obtain the scatter point set at the center of the river; The spline curve tool is used to fit the scatter set of points at the center of the river channel to the center line of the river channel. The scatter points are then evenly divided along the fitted center line to obtain a new scatter set of points for the center line of the river channel, which is denoted as the coordinate set of the center line of the river channel. Draw perpendicular lines to the river centerline from the new scatter set of points passing through the river centerline to obtain the corresponding multiple cross-sectional lines; Determine the intersection points of each cross-section line with the corresponding left and right bank boundary lines, ensuring that the scattered points on the left bank boundary line, right bank boundary line and river center line correspond one-to-one and that the line connecting the scattered points on the left and right bank boundary lines is basically perpendicular to the river center line. Based on the direction of each cross-section line, the coordinate sets of the left and right bank boundary lines and the river centerline are transformed into an orthogonal curvilinear coordinate system.
[0009] Furthermore, the process of searching for the nearest scatter pair is specifically as follows: The points are evenly distributed along the two riverbanks after fitting, forming a new set of scattered points for the two riverbanks. For each scatter point on one riverbank, search for the scatter point closest to that scatter point in the scatter point set on the other riverbank; for bends, search only the scatter point set on the concave bank for the scatter point closest to each scatter point on the convex bank. Pair the searched scattered points with scattered points along one of the riverbanks, determine the midpoint of each pair of scattered points, delete outlier midpoints, and obtain the scattered point set at the center of the river.
[0010] Furthermore, under the orthogonal curvilinear coordinate system, the river channel planar morphological characteristic parameters are determined; wherein, the river channel planar morphological characteristic parameters include the width of the flat riverbed at each cross-section, the curvature coefficient along the river centerline, the relative curvature along the river centerline, and the tortuosity coefficient of the bends, specifically: Based on the intersection of each cross-section line with the boundary lines on both banks, the length of each cross-section line is calculated to obtain the width of the flat river channel for each cross-section line. In an orthogonal curvilinear coordinate system, the curvature coefficient along the river centerline is calculated based on scattered points along the centerline; the calculation formula is as follows:
[0011] In the formula, The curvature coefficient along the centerline of the river channel; and These are scattered points along the center line of the river channel. X c (s), Y c (s)} for distance along the route s The first derivative, and These are scattered points along the center line of the river channel. X c (s), Y c (s)} for distance along the route s The second derivative; Based on the width of the flat river channel and the curvature along the river centerline at each cross-section, the relative curvature along the river centerline is calculated; the calculation formula is as follows:
[0012] In the formula, The width of the flat riverbed at each cross-section line; The curvature coefficient along the centerline of the river channel; and These are scattered points along the center line of the river channel. X c (s), Y c (s)} for distance along the route s The first derivative, and These are scattered points along the center line of the river channel. X c (s), Yc (s)} for distance along the route s The second derivative of is defined as follows: the curvature bending to the left bank is negative, and the curvature bending to the right bank is positive. Determine the location of the bend crest and inflection point of the curved section; Based on the locations of the bend crest and inflection point of the curved section, the length of the river channel centerline and the straight-line distance of the curved section are determined, and the tortuosity coefficient of the curved section is calculated; the calculation formula is as follows:
[0013] In the formula, S The tortuosity coefficient of the curved section; L c This is the length of the centerline of the bend in the river channel, which usually starts at the center point of the upstream transition section and ends at the center point of the downstream transition section. L s The straight-line distance between the starting point and the ending point.
[0014] Furthermore, determining the location of the bend apex and inflection point of the curved segment specifically involves: 1) Traverse the relative curvature data along the river centerline and determine the location of the section with the largest absolute value of relative curvature as the bend apex; 2) The cross-sectional positions where the relative curvature value is zero or the relative curvature sign changes are determined as inflection points.
[0015] Furthermore, based on the digital elevation model file, an elevation point matrix is obtained, and the morphological characteristic parameters of the river channel cross section and the thalweg longitudinal profile are determined in an orthogonal curvilinear coordinate system. The morphological characteristic parameters of the river channel cross section include the flat area, flat water depth, width-to-depth ratio, and wetted perimeter of each section; the morphological characteristic parameters of the thalweg longitudinal profile include the location and elevation of the thalweg points at each section, and the thalweg longitudinal gradient, specifically: In orthogonal curvilinear coordinate system ( s , n , z Generate evenly spaced planar meshes in ) According to each n The orientation angle of the directional grid lines determines the orientation of each grid node. snz coordinates and xyz The correspondence between coordinates; Use Matlab to extract elevation data from DEM files, generate an elevation point matrix, and draw a contour map. The extracted elevation point matrix xyz Elevation interpolation in coordinate system to snz On the grid nodes of the coordinate system, obtain the elevation points on each grid node; The area of the flat beach, the water depth of the flat beach, the width-to-depth ratio of the flat beach channel and the wetted perimeter are calculated on each cross section to determine the characteristic parameters of the river channel cross section. Determine the location and elevation of the thalweg on each cross-section, in xyz Plot the plan view and longitudinal profile of the thalweg in the coordinate system, and calculate the longitudinal gradient of the thalweg. J In order to determine the morphological characteristics of the deepwater longitudinal profile; Furthermore, the characteristic parameters of the river channel cross-section morphology are determined, specifically: Assume that there are a total of There are several grid nodes, which are numbered sequentially from one bank to the other, with serial numbers ranging from 1 to... ; The flat beach elevation of each cross section is set as the arithmetic mean of the elevation points on both banks. Based on the elevation of each grid node, the flat beach elevation of the cross section is calculated using the following formula:
[0016] In the formula, The elevation of the flat beach at the cross-section; , Elevations of grid nodes 1 and m respectively; Based on the flat elevation of the cross-section and the coordinate set of grid nodes on each cross-section, the flat area of the cross-section is calculated; the calculation formula is as follows:
[0017] In the formula, The flat area of the cross-section; The elevation of the flat beach at the cross-section; , Grid nodes and Elevation; , Grid nodes i and In Coordinates of direction; Based on the width of the flat river channel and the flat area of each cross-section, the flat water depth of the cross-section is calculated using the following formula:
[0018] In the formula, The flat area of the cross-section; The width of the flat riverbed at each cross-section line; The depth of the flat beach at the cross-section; Based on the width of the river channel at each cross-section and the water depth at the cross-section, the width-to-depth ratio is calculated using the following formula:
[0019] In the formula, Aspect ratio; The width of the flat riverbed at each cross-section line; The depth of the flat beach at the cross-section; Based on the elevation and coordinate point set of the grid nodes on each cross-section, the wetted perimeter is calculated using the following formula:
[0020] In the formula, For wet period; , Grid nodes and Elevation; , Grid nodes i and In Coordinates of direction.
[0021] Furthermore, the morphological characteristic parameters of the deepwater longitudinal profile were determined, specifically as follows: Set the thalweg point as the point with the lowest elevation on the cross section, and find the elevation set {z} of the grid nodes on each cross section. 1, z2,……,z m The minimum value z in} k The coordinates of the thalweg point are obtained. s k , n k , z k ,) and their corresponding xyz Coordinates in a coordinate system x k , y k , z k ,), that is, to determine the location and elevation of the thalweg on each cross section; Based on determining the location and elevation of the thalweg at each cross-section, in xyz Draw the plan view and longitudinal section view of the thalweg in the coordinate system; Assuming grid nodes on the cross-section i and i The thalweg elevations of +1 are respectively z i,min , z i+1,min Calculate the grid nodes on the cross section i and i The distance Δs between +1 is used to calculate the longitudinal gradient of the cross-section, and the calculation formula is as follows: J = ( z i,min — zi+1,min ) / △s In the formula, J The longitudinal ratio of the deep channel in the cross-section decreases; z i,min , z i+1,min These are the grid nodes on the cross-section. i and i +1 represents the depth elevation; △s represents the grid nodes on the cross-section. i and i The distance between +1.
[0022] Compared with the prior art, the present invention has the following advantages: 1) This invention overcomes the limitation of the traditional fixed section method, which can only obtain discrete points, by constructing an orthogonal curve coordinate system along the center line of the river channel. It can fully reflect the continuous change law of the riverbed morphology along the course and realize the continuous extraction of the riverbed morphology characteristic parameters of the entire river section. 2) The method of the present invention automatically searches for the nearest point pair on both banks to generate the center line of the river channel, and generates a cross-sectional line perpendicular to the direction of water flow in strict accordance with geometric relationships, eliminating the subjective error caused by manually delineating the cross-section. It is particularly suitable for the precise analysis of curved river sections and improves the accuracy and automation of cross-section extraction. 3) This invention only requires conventional remote sensing images and digital elevation models of rivers to achieve the entire process calculation, without the need for expensive underwater point cloud scanning equipment. Moreover, the standardized algorithm process greatly improves work efficiency, has wide engineering applicability, and reduces data acquisition and processing costs. 4) This invention covers a complete set of characteristic parameters from planar morphology to cross-section to longitudinal profile of the thalweg. The parameters are associated through a unified coordinate system to ensure the consistency and reliability of the data, providing a solid data foundation for riverbed evolution analysis. Attached Figure Description
[0023] The present invention will be further described below with reference to the accompanying drawings and specific embodiments: Figure 1 This is a flowchart illustrating a method for determining riverbed morphological characteristic parameters based on digital river topography according to the present invention. Figure 2 These are remote sensing images of the river section studied in this invention and schematic diagrams of the riverbank boundaries. Figure 3 This is a schematic diagram of the Geo-TIFF format DEM file of the river section studied in this invention; Figure 4 This is a schematic diagram comparing the fitted lines of the riverbank boundaries on both sides of the river section studied in this invention with the original measured points; Figure 5 This is a schematic diagram of the riverbank boundary lines, the river centerline, and the cross-sectional line perpendicular to the river centerline of the river section studied in this invention. Figure 6 This is a schematic diagram showing the calculation results of the river channel planar morphology characteristic parameters involved in this invention; Figure 7 It is a contour map of the digital topography of the river channel involved in this invention; Figure 8 This is a schematic diagram illustrating the definition of cross-sectional morphological feature parameters in an embodiment of the present invention; Figure 9 This is a diagram showing the calculation results of the cross-sectional morphological characteristic parameters involved in the embodiments of the present invention; Figure 10 This is a diagram showing the calculated results of the deep-water plane position, longitudinal section, and longitudinal gradient involved in the embodiments of the present invention. Detailed Implementation
[0024] Example 1
[0025] like Figure 1 As shown, a method for determining riverbed morphology characteristic parameters based on digital river topography includes: S1. Acquire remote sensing image data, determine the boundary lines of the left and right banks of the river section under study at the corresponding time, and obtain the planar scatter set of the boundary lines of the two banks. Specifically, S11, download remote sensing images of the corresponding year's flood season high water period for the research river section, and synthesize multiple different bands based on remote sensing image processing software; S12. Selecting the green light band and mid-infrared band, the improved normalized differential water index is calculated based on the water index calculation formula; the calculation formula is as follows: MNDWI = (Green - MIR) / (Green + MIR) In the formula, MNDWI is the improved normalized differential water index, with a value of [-1,1]. This index threshold is used to identify land and water areas and extract water areas; Green and MIR represent the reflectance of the green light band and the mid-infrared band, respectively; in Landsat 5, 7, and 8 satellites, Green corresponds to the B2, B2, and B3 bands, respectively, and MIR corresponds to the B5, B6, and B6 bands, respectively. S13. Generate a binary river topographic map based on the improved normalized differential water index calculation results, and extract the boundary lines of the left and right banks using vectorization tools. S14. Divide the extracted left and right boundary lines into scattered points at equal intervals, assign plane coordinates to obtain plane scattered point sets, namely the scattered point set of the left bank boundary and the scattered point set of the right bank boundary.
[0026] In this example, Landsat 8 satellite remote sensing images of a certain section of a river in the middle reaches of a certain basin during the high water level of the flood season in August 2018 were downloaded. The Green band corresponds to B3, and the MIR band to B6. Then, 11-band composite images were created using remote sensing image processing software. Figure 2 This document presents a remote sensing image of a section of a river in the middle reaches of a certain basin, taken on July 20, 2018. Using QGIS software, the following steps are performed: First, select Raster Analysis - Raster Calculator in the toolbox. Then, select the appropriate band according to the water index formula in the Raster Calculator to output the water distribution structure. Next, select Raster - Transform - Raster Vectorization in the toolbar, and select the newly generated raster file as the input feature. Then, select Vector Geometry - Add Geometric Attributes in the toolbox, select the face containing the sandbar as the input feature, and extract its area and perimeter. Open the attribute table, delete the area corresponding to the sandbar, and obtain the river surface features. Finally, select Vector Geometry - Polygon to Line in the toolbox to convert the river surface features into line features. Delete redundant line features and break the inlet and outlet lines to obtain the left and right bank boundary lines, as shown below. Figure 2 As shown by the red line; select Vector Geometry - Points along the Geometry, select the left bank shoreline in the input features, and set the spacing of the scatter points to 100 m as needed to ensure the accuracy of the description of the boundary lines on both banks. The generated features are the scatter point set of the left bank boundary line, as shown below. Figure 2 As shown by the yellow dot, select the "Add Coordinates to Points" plugin in SAGA, choose the scatter plot layer of the left bank shoreline in the input features, select an appropriate projected coordinate system, and obtain the scatter plot of the left bank boundary line with projected coordinates. X L , Y L Similarly, the scatter set of the right bank boundary line is obtained. X R , Y R}
[0027] S2. Obtain the river topographic map, convert it into a regular rectangular digital elevation model file, and extract the elevation data of the riverbank boundary lines based on the planar scatter point set to construct a three-dimensional coordinate set of the riverbank boundary lines. Specifically, S21, determine the minimum and maximum values of the horizontal and vertical coordinates of the elevation points on the river topographic map, and construct a regular rectangular range containing four virtual corner points; S22. Assign an elevation value to the virtual corner point, taking the maximum elevation value from the river topographic map; S23. Spatial interpolation is performed between the measured elevation points and the virtual corner points to generate a regular rectangular digital elevation model file (DEM file) in Geo-TIFF format; S24. Extract the elevation values of the two bank boundary lines corresponding to the planar scatter set from the DEM file, and combine them to form a three-dimensional coordinate set of the two bank boundary lines containing planar coordinates and elevation information, namely the three-dimensional coordinate set of the left bank boundary line and the three-dimensional coordinate set of the right bank boundary line.
[0028] In this example, we collected a river topographic map of a section of a river in the middle reaches of a certain basin, measured after the flood season in 2018. In AutoCAD, we selected Tools - Data Extraction, selected the layer with the corresponding elevation points, and extracted the elevation points from the river topographic map. At the same time, we determined the minimum x-coordinate of these elevation points. X min Maximum value of x-coordinate X max Minimum value of the ordinate Y min and the maximum value of the ordinate Y max Combine them into four virtual corner points, namely ( X min , Y min ), ( X min , Y max ), ( X max , Y min ), ( X max , Y max The elevations of the four virtual corner points are assigned the maximum elevation value of 40m from the entire river channel topographic map. The measured elevation points and the four virtual corner points are imported into QGIS and named TG-LJP. In the toolbox, select Interpolation - Irregular Triangular Network Interpolation or Inverse Distance Weighted Interpolation. In the input features, select the "TG-LJP" layer, set the relevant interpolation parameters, and obtain a regular rectangular DEM layer. Extract the layer as a Geo-TIFF format DEM file, i.e., TG-LJP.tif. Figure 3 As shown; select Vector Geometry - Extract Z-value, select the left and right bank boundary line layers and DEM layer respectively in the input features, extract the elevation of the left and right bank boundary lines, and obtain the three-dimensional coordinate set of the bank boundary lines { X L , Y L , Z L}、{ X R , Y R , Z R} is denoted as the original three-dimensional coordinate set of the boundary lines between the two sides.
[0029] S3. Fit the three-dimensional coordinate set of the two riverbank boundary lines to the riverbank line. Determine the scatter point set of the river center by searching for the closest scatter point pair on the riverbank line. Fit and generate the river center line and extract the cross-sectional line perpendicular to the river center line. Establish an orthogonal curve coordinate system along the direction of the river center line, perpendicular to the direction of the river center line and the direction of water depth. Specifically, S31, using the spline curve tool in Matlab, fits the original scatter set of the boundary lines on both sides into two riverbank lines respectively. S32. Search for the closest scatter point pair on the two fitted riverbanks to obtain the scatter point set at the center of the river. The process of searching for the nearest scatter point pair is as follows: The points are evenly distributed along the two riverbanks after fitting, forming a new set of scattered points for the two riverbanks. For each scatter point on one riverbank, search for the scatter point closest to that scatter point in the scatter point set on the other riverbank; for bends, search only the scatter point set on the concave bank for the scatter point closest to each scatter point on the convex bank. Pair the searched scattered points with scattered points along one of the riverbanks, determine the midpoint of each pair of scattered points, delete outlier midpoints, and obtain the scattered point set at the center of the river. S33. Use the spline curve tool to fit the scatter set of the river center to the river centerline, and divide the scatter points evenly along the fitted river centerline to obtain a new scatter set of the river centerline, which is denoted as the coordinate set of the river centerline. S34. Draw perpendicular lines to the river centerline from the new scatter set of points passing through the river centerline to obtain the corresponding multiple cross-sectional lines. S35. Determine the intersection point of each cross-section line with the corresponding left and right bank boundary lines, and ensure that the scattered points on the left bank boundary line, right bank boundary line and river center line correspond one-to-one and that the line connecting the scattered points on the left and right bank boundary lines is basically perpendicular to the river center line. S36. Based on the direction of each cross-section line, transform the coordinate system of the left and right bank boundary lines and the river centerline to the orthogonal curve coordinate system. s , n , z ) in; among them, s Let the direction be along the centerline of the river channel, and let the direction of water flow be positive; n The direction perpendicular to the center line of the river channel is defined as positive, with the left bank as the positive direction. z Let the direction of water depth be denoted by "upwards" as positive.
[0030] In this example, the `spaps` function in Matlab is used to fit the three-dimensional coordinate sets of the riverbank boundaries into lines, resulting in two riverbank lines. Furthermore, to verify the reliability of the fitting, the fitted two riverbank lines and the original scattered points in the three-dimensional coordinate sets of the riverbank boundaries are plotted on the same graph, as shown below. Figure 4 As shown, this is a schematic diagram comparing the fitted boundary lines and the original measured points on both banks. The fnval function is used to evenly divide the points along the fitted riverbanks, forming new scatter sets for the two riverbanks, denoted as the target scatter sets for the boundary lines, respectively {...} X L ', Y L ', Z L '}、{ X R ', Y R ', Z R '}; For straight sections, search the scattered points on the right bank for the point with the smallest Euclidean distance to each scattered point on the left bank; for curved sections, search only the scattered points on the concave bank for the point with the smallest Euclidean distance to each scattered point on the convex bank, grouping them in pairs, calculating the straight-line distance between each group of scattered points, and determining the midpoint; plot the obtained midpoints and the boundary lines of the left and right banks on the same map, delete outlier midpoints, and the coordinates corresponding to the remaining midpoints form the scattered point set at the center of the river channel. X C , Y C , Z C}; The spaps function is used again to fit the scatter set of points along the river centerline to obtain the river centerline, and the fnval function is used to uniformly divide the scatter points along the river centerline to obtain a new scatter set of points along the river centerline. X C ', Y C ', Z C In this example, the curve length between two adjacent scattered points is set to 500m, thus determining the cross-sectional line perpendicular to the river centerline corresponding to each scattered point; the intersection points of each cross-sectional line with the boundary lines on both banks are solved, so that the scattered points on the left bank boundary, right bank boundary, and river centerline correspond one-to-one, and the line connecting the scattered points on the left and right bank boundary lines is basically perpendicular to the river centerline, such as... Figure 5 As shown, the red and blue lines represent the riverbank boundaries; the pink line represents the river centerline; and the black line represents the cross-sectional lines. Based on the direction of each cross-sectional line, the coordinate sets of the left and right bank boundaries and the river centerline are transformed to an orthogonal curvilinear coordinate system. s , n , zIn this context, scattered points along the centerline of the river channel can all be represented as { X c (s), Y c (s)}.
[0031] S4. Under the orthogonal curvilinear coordinate system, determine the river channel planar morphology characteristic parameters; wherein, the river channel planar morphology characteristic parameters include the width of the flat riverbed at each cross-section, the curvature coefficient along the river centerline, the relative curvature along the river centerline, and the tortuosity coefficient of the bend. Specifically, S41, based on the intersection of each cross-section line with the boundary lines on both banks, calculate the length of each cross-section line to obtain the width of the flat river channel for each cross-section line. ; S42. In an orthogonal curvilinear coordinate system, the curvature coefficient along the river centerline is calculated based on scattered points along the river centerline. The calculation formula is as follows:
[0032] In the formula, The curvature coefficient along the centerline of the river channel; and These are scattered points along the center line of the river channel. X c (s), Y c (s)} for distance along the route s The first derivative, and These are scattered points along the center line of the river channel. X c (s), Y c (s)} for distance along the route s The second derivative; S43. Based on the width of the flat river channel and the curvature along the centerline of the river at each cross-section, the relative curvature along the centerline of the river is calculated. The calculation formula is as follows:
[0033] In the formula, The width of the flat riverbed at each cross-section line; The curvature coefficient along the centerline of the river channel; and These are scattered points along the center line of the river channel. X c (s), Y c (s)} for distance along the route s The first derivative, and These are scattered points along the center line of the river channel. X c (s), Y c (s)} for distance along the route s The second derivative of is defined as follows: the curvature bending to the left bank is negative, and the curvature bending to the right bank is positive. S44. Determine the location of the bend apex and inflection point of the curved section; specifically: 1) Traverse the relative curvature data along the river centerline and determine the location of the section with the largest absolute value of relative curvature as the bend apex; 2) The cross-sectional locations where the relative curvature value is zero or the sign of the relative curvature changes are determined as inflection points; S45. Based on the location of the bend crest and inflection point of the curved section, determine the length of the river channel centerline and the straight-line distance of the curved section, and calculate the tortuosity coefficient of the curved section; the calculation formula is as follows:
[0034] In the formula, S The tortuosity coefficient of the curved section; L c This is the length of the centerline of the bend in the river channel, which usually starts at the center point of the upstream transition section and ends at the center point of the downstream transition section. L s The straight-line distance between the starting point and the ending point.
[0035] In this embodiment, the width of the flat river channel at each cross-section is calculated. B bf The curvature coefficient along the centerline of the river channel 1 / R C and relative curvature B bf / R C Determine the location of the bend apex and inflection point of the curved section and calculate the tortuosity coefficient of the curved section. S The result is as follows Figure 6 As shown.
[0036] S5. Based on the digital elevation model file, obtain the elevation point matrix and determine the morphological characteristic parameters of the river channel cross section and the deep channel longitudinal section in the orthogonal curvilinear coordinate system. Among them, the morphological characteristic parameters of the river channel cross section include the flat area, flat water depth, width-to-depth ratio and wetted perimeter of each section; the morphological characteristic parameters of the deep channel longitudinal section include the location and elevation of the deep channel points of each section and the deep channel longitudinal ratio. Specifically, S51, in the orthogonal curvilinear coordinate system ( s , n , z Generate an evenly spaced planar mesh in ); where, sThe grid spacing in the direction is given based on requirements, and is generally taken as the spacing between section lines. n The spacing in each direction is given based on the resolution of the DEM file; then, according to each... n The orientation angle of the directional grid lines determines the orientation of each grid node. snz coordinates and xyz The correspondence between coordinates; S52. Use Matlab to extract elevation data from DEM files, generate an elevation point matrix, and draw a contour map. S53. Extract the elevation point matrix. xyz Elevation interpolation in coordinate system to snz On the grid nodes of the coordinate system, obtain the elevation points on each grid node; S54. Calculate the flat beach area on each cross section. A bf The water is deep on the flat beach. H bf The width-to-depth ratio of the flat river channel B bf / H bf And the wetted perimeter χ, and then determine the characteristic parameters of the river channel cross-section morphology; specifically: Assume that there are a total of There are several grid nodes, which are numbered sequentially from one bank to the other, with serial numbers ranging from 1 to... ; The flat beach elevation of each cross section is set as the arithmetic mean of the elevation points on both banks. Based on the elevation of each grid node, the flat beach elevation of the cross section is calculated using the following formula:
[0037] In the formula, The elevation of the flat beach at the cross-section; , Elevations of grid nodes 1 and m respectively; Based on the flat elevation of the cross-section and the coordinate set of grid nodes on each cross-section, the flat area of the cross-section is calculated; the calculation formula is as follows:
[0038] In the formula, The flat area of the cross-section; The elevation of the flat beach at the cross-section; , Grid nodes and Elevation; , Grid nodes i and In Coordinates of direction; Based on the width of the flat river channel and the flat area of each cross-section, the flat water depth of the cross-section is calculated using the following formula:
[0039] In the formula, The flat area of the cross-section; The width of the flat riverbed at each cross-section line; The depth of the flat beach at the cross-section; Based on the width of the river channel at each cross-section and the water depth at the cross-section, the width-to-depth ratio is calculated using the following formula:
[0040] In the formula, Aspect ratio; The width of the flat riverbed at each cross-section line; The depth of the flat beach at the cross-section; Based on the elevation and coordinate point set of the grid nodes on each cross-section, the wetted perimeter is calculated using the following formula:
[0041] In the formula, For wet period; , Grid nodes and Elevation; , Grid nodes i and In Coordinates of direction.
[0042] S55. Determine the location of the thalweg on each cross-section ( s t , n t ) and its elevation z t ,exist xyz Plot the plan view and longitudinal profile of the thalweg in the coordinate system, and calculate the longitudinal gradient of the thalweg. J This allows for the determination of the morphological characteristics of the deepwater longitudinal profile; specifically: Set the thalweg point as the point with the lowest elevation on the cross section, and find the elevation set {z} of the grid nodes on each cross section. 1, z2,……,z m The minimum value z in} k The coordinates of the thalweg point are obtained. s k , n k ,z k ,) and their corresponding xyz Coordinates in a coordinate system x k , y k , z k ,), that is, to determine the location and elevation of the thalweg on each cross section; Based on determining the location and elevation of the thalweg at each cross-section, in xyz Draw the plan view and longitudinal section view of the thalweg in the coordinate system; Assuming grid nodes on the cross-section i and i The thalweg elevations of +1 are respectively z i,min , z i+1,min Calculate the grid nodes on the cross section i and i The distance Δs between +1 is used to calculate the longitudinal gradient of the cross-section, and the calculation formula is as follows: J = ( z i,min — z i+1,min ) / △s In the formula, J The longitudinal ratio of the deep channel in the cross-section decreases; z i,min , z i+1,min These are the grid nodes on the cross-section. i and i +1 represents the depth elevation; △s represents the grid nodes on the cross-section. i and i The distance between +1.
[0043] In this example, we first generate a planar mesh with equal spacing. S , N The grid spacing along the river centerline S The spacing is set to 500m, perpendicular to the centerline of the river channel. N Set to 20m, depending on each n The orientation angle of the directional grid lines determines the orientation of each grid node. snz coordinates and xyz The coordinate correspondence; using the `readgeoraster` function in Matlab to read the regular rectangular digital elevation model file (TG-LJP.tif file) in Geo-TIFF format, and automatically generating the elevation point matrix { X , Y , ZUse the `contourf` function in Matlab to draw contour plots, such as... Figure 7 As shown, after the terrain is drawn, a preliminary check is performed to see if it is reasonable: if it is reasonable, proceed to the next step; if it is not reasonable, repeat step S24, adding multiple scattered points with an elevation of 40 on the four sides of the rectangle formed by the four virtual corner points, and interpolating again; then use the interp2 function in Matlab to obtain the elevation point matrix { X , Y , Z Elevation values in} are interpolated to the planar grid. S , N On}, a new elevation point matrix is formed { S , N , Z};exist snz Calculate each line in the coordinate system n Cross-sectional morphological characteristic parameters of the grid lines in the direction, including the flat area A bf Pingtan River is wide B bf The water is deep on the flat beach. H bf The width-to-depth ratio of the river channel in the Pingtan area B bf / H bf And the wetted period χ, the definitions of each parameter are as follows Figure 8 As shown, the calculation results are as follows Figure 9 As shown; determine the location of the thalweg on each cross-section ( s t , n t ) and its elevation z t ,exist xyz Plot the plan view and longitudinal profile of the thalweg in the coordinate system, and calculate the longitudinal gradient of the thalweg. J This allows for the determination of the characteristic parameters of the longitudinal profile of the deepwater channel. The calculation results are as follows: Figure 10 As shown.
[0044] This invention first utilizes satellite remote sensing and geographic information system (GIS) technology to acquire the boundary lines of the left and right banks of the studied river section and their planar scatter points. Next, the river topographic map is converted into a regular rectangular DEM file in Geo-TIFF format, and the elevation points and three-dimensional coordinate sets of the left and right bank boundary lines are obtained. Then, Matlab programming is used to calculate and determine the river's planar morphological characteristic parameters. Finally, elevation data is extracted from the river topographic map to obtain the elevation point matrix, and the morphological characteristic parameters of the river cross-section and the longitudinal section of the thalweg are calculated and determined. This invention, based on remote sensing imagery and river topographic data, solves the problems of previous methods that only allowed extraction by river type and river section, excessively large fixed cross-section spacing, unreasonable cross-section layout, and the inability to directly calculate riverbed morphological characteristic parameters. It overcomes the technical bottleneck of automatically and efficiently extracting riverbed morphological characteristic parameters from digital river topography. This invention has a strong theoretical and professional foundation, and is highly operable and applicable to engineering projects, providing an effective tool for extracting riverbed morphological characteristic parameters for various types of rivers.
[0045] The embodiments described above are merely preferred embodiments of the present invention and are not intended to limit the scope of the present invention. Various modifications and improvements made by those skilled in the art to the technical solutions of the present invention without departing from the principles and essence of the present invention should fall within the protection scope defined by the claims of the present invention.
Claims
1. A method for determining riverbed morphology characteristic parameters based on digital river topography, characterized in that, Includes the following steps: Acquire remote sensing image data, determine the corresponding time period of the study river section's left and right bank boundary lines, and obtain a planar scatter set of the bank boundary lines; Obtain the river topographic map, convert it into a regular rectangular digital elevation model file, and extract the elevation data of the riverbank boundary lines based on the planar scatter point set to construct a three-dimensional coordinate set of the riverbank boundary lines; The three-dimensional coordinate set of the two riverbank boundary lines is fitted to the riverbank line. The scatter point set of the river center is determined by searching the closest scatter point pair on the riverbank line. The river center line is fitted to generate the river center line and the cross-sectional line perpendicular to the river center line is extracted. An orthogonal curve coordinate system is established along the direction of the river center line, perpendicular to the direction of the river center line and the direction of water depth. Under the orthogonal curvilinear coordinate system, the river channel planar morphology characteristic parameters are determined; wherein, the river channel planar morphology characteristic parameters include the width of the flat riverbed at each cross-section, the curvature coefficient along the river centerline, the relative curvature along the river centerline, and the tortuosity coefficient of the bend. Based on the digital elevation model file, the elevation point matrix is obtained, and the morphological characteristic parameters of the river channel cross section and the deep channel longitudinal section are determined in the orthogonal curvilinear coordinate system. Among them, the morphological characteristic parameters of the river channel cross section include the flat area, flat water depth, width-to-depth ratio and wetted perimeter of each section; the morphological characteristic parameters of the deep channel longitudinal section include the location and elevation of the deep channel points of each section and the deep channel longitudinal ratio.
2. The method for determining riverbed morphology characteristic parameters based on digital river topography according to claim 1, characterized in that, The process of acquiring remote sensing image data, determining the boundary lines of the left and right banks of the river section under study at the corresponding time, and obtaining a planar scatter set of the boundary lines is as follows: Download remote sensing images of the river section during the high water period of the flood season in the corresponding year, and synthesize multiple different bands based on remote sensing image processing software; By selecting the green and mid-infrared bands, an improved normalized differential water index was calculated based on the water index calculation formula; the calculation formula is as follows: MNDWI = (Green - MIR) / (Green + MIR) In the formula, MNDWI is the improved normalized differential water index, with a value of [-1,1]. This index threshold is used to identify land and water areas and extract water areas; Green and MIR represent the reflectance of the green light band and the mid-infrared band, respectively; in Landsat 5, 7, and 8 satellites, Green corresponds to the B2, B2, and B3 bands, respectively, and MIR corresponds to the B5, B6, and B6 bands, respectively. A binary river topographic map is generated based on the improved normalized differential water index calculation results, and the boundary lines of the left and right banks are extracted using vectorization tools. The extracted left and right boundary lines are divided into scattered points at equal intervals. After being assigned planar coordinates, planar scattered point sets are obtained, namely the scattered point set of the left boundary line and the scattered point set of the right boundary line.
3. The method for determining riverbed morphology characteristic parameters based on digital river topography according to claim 1, characterized in that, The process of acquiring the river topographic map, converting it into a regular rectangular digital elevation model file, and extracting the elevation data of the riverbank boundaries based on the planar scatter point set to construct a three-dimensional coordinate set for the riverbank boundaries is as follows: Obtain the river topographic map, determine the minimum and maximum values of the horizontal and vertical coordinates of the elevation points on the river topographic map, and construct a regular rectangular range containing four virtual corner points; The elevation of the virtual corner point is assigned the maximum elevation value on the topographic map; Spatial interpolation is performed between the measured elevation points and the virtual corner points to generate a regular rectangular digital elevation model file in Geo-TIFF format; The elevation values of the corresponding locations of the plane scatter point set are extracted from the digital elevation model file and combined to form a three-dimensional coordinate set of the two bank boundary lines containing plane coordinates and elevation information, namely the three-dimensional coordinate set of the left bank boundary line and the three-dimensional coordinate set of the right bank boundary line.
4. The method for determining riverbed morphology characteristic parameters based on digital river topography according to claim 1, characterized in that, The three-dimensional coordinate set of the riverbank boundary lines is fitted to the riverbank line. The scatter point set of the river center is determined by searching for the nearest scatter point pair on the riverbank line. The river centerline is then fitted to generate the river centerline, and the cross-sectional line perpendicular to the river centerline is extracted. An orthogonal curve coordinate system is established along the direction of the river centerline, perpendicular to the river centerline, and in the direction of water depth. Specifically: Using the spline curve tool in Matlab, the original scatter set of the boundary lines on both sides of the river was fitted into two riverbank lines. Search for the nearest scatter point pair on the two fitted riverbanks to obtain the scatter point set at the center of the river; The spline curve tool is used to fit the scatter set of points at the center of the river channel to the center line of the river channel. The scatter points are then evenly divided along the fitted center line to obtain a new scatter set of points for the center line of the river channel, which is denoted as the coordinate set of the center line of the river channel. Draw perpendicular lines to the river centerline from the new scatter set of points passing through the river centerline to obtain the corresponding multiple cross-sectional lines; Determine the intersection points of each cross-section line with the corresponding left and right bank boundary lines, ensuring that the scattered points on the left bank boundary line, right bank boundary line and river center line correspond one-to-one and that the line connecting the scattered points on the left and right bank boundary lines is basically perpendicular to the river center line. Based on the direction of each cross-section line, the coordinate sets of the left and right bank boundary lines and the river centerline are transformed into an orthogonal curvilinear coordinate system.
5. The method for determining riverbed morphology characteristic parameters based on digital river topography according to claim 3, characterized in that, The process of searching for the nearest scatter pair is as follows: The points are evenly distributed along the two riverbanks after fitting, forming a new set of scattered points for the two riverbanks. For each scatter point on one riverbank, search for the scatter point closest to that scatter point in the scatter point set on the other riverbank; for bends, search only the scatter point set on the concave bank for the scatter point closest to each scatter point on the convex bank. Pair the searched scattered points with scattered points along one of the riverbanks, determine the midpoint of each pair of scattered points, delete outlier midpoints, and obtain the scattered point set at the center of the river.
6. The method for determining riverbed morphology characteristic parameters based on digital river topography according to claim 1, characterized in that, Under the orthogonal curvilinear coordinate system, the river channel planar morphological characteristic parameters are determined; wherein, the river channel planar morphological characteristic parameters include the width of the flat riverbed at each cross-section, the curvature coefficient along the river centerline, the relative curvature along the river centerline, and the tortuosity coefficient of the bends, specifically: Based on the intersection of each cross-section line with the boundary lines on both banks, the length of each cross-section line is calculated to obtain the width of the flat river channel for each cross-section line. In an orthogonal curvilinear coordinate system, the curvature coefficient along the river centerline is calculated based on scattered points along the centerline; the calculation formula is as follows: In the formula, The curvature coefficient along the centerline of the river channel; and These are scattered points along the center line of the river channel. X c (s), Y c (s)} for distance along the route s The first derivative, and These are scattered points along the center line of the river channel. X c (s), Y c (s)} for distance along the route s The second derivative; Based on the width of the flat river channel and the curvature along the river centerline at each cross-section, the relative curvature along the river centerline is calculated; the calculation formula is as follows: In the formula, The width of the flat riverbed at each cross-section line; The curvature coefficient along the centerline of the river channel; and These are scattered points along the center line of the river channel. X c (s), Y c (s)} for distance along the route s The first derivative, and These are scattered points along the center line of the river channel. X c (s), Y c (s)} for distance along the route s The second derivative of is defined as follows: the curvature bending to the left bank is negative, and the curvature bending to the right bank is positive. Determine the location of the bend crest and inflection point of the curved section; Based on the locations of the bend crest and inflection point of the curved section, the length of the river channel centerline and the straight-line distance of the curved section are determined, and the tortuosity coefficient of the curved section is calculated; the calculation formula is as follows: In the formula, S The tortuosity coefficient of the curved section; L c This is the length of the centerline of the bend in the river channel, which usually starts at the center point of the upstream transition section and ends at the center point of the downstream transition section. L s The straight-line distance between the starting point and the ending point.
7. The method for determining riverbed morphology characteristic parameters based on digital river topography according to claim 6, characterized in that, The determination of the bend apex and inflection point positions of the curved section specifically involves: 1) Traverse the relative curvature data along the river centerline and determine the location of the section with the largest absolute value of relative curvature as the bend apex; 2) The cross-sectional positions where the relative curvature value is zero or the relative curvature sign changes are determined as inflection points.
8. The method for determining riverbed morphology characteristic parameters based on digital river topography according to claim 6, characterized in that, Based on the digital elevation model file, an elevation point matrix is obtained, and the morphological characteristic parameters of the river channel cross section and the thalweg longitudinal section are determined in an orthogonal curvilinear coordinate system. The morphological characteristic parameters of the river channel cross section include the flat area, flat water depth, width-to-depth ratio, and wetted perimeter of each section. The morphological characteristic parameters of the thalweg longitudinal section include the location and elevation of the thalweg points at each section, and the thalweg longitudinal gradient, specifically: In orthogonal curvilinear coordinate system ( s , n , z Generate evenly spaced planar meshes in ) According to each n The orientation angle of the directional grid lines determines the orientation of each grid node. snz coordinates and xyz The correspondence between coordinates; Use Matlab to extract elevation data from DEM files, generate an elevation point matrix, and draw a contour map. The extracted elevation point matrix xyz Elevation interpolation in coordinate system to snz On the grid nodes of the coordinate system, obtain the elevation points on each grid node; The area of the flat beach, the water depth of the flat beach, the width-to-depth ratio of the flat beach channel and the wetted perimeter are calculated on each cross section to determine the characteristic parameters of the river channel cross section. Determine the location and elevation of the thalweg on each cross-section, in xyz Plot the plan view and longitudinal profile of the thalweg in the coordinate system, and calculate the longitudinal gradient of the thalweg. J This allows us to determine the morphological characteristics of the deep thalweg longitudinal profile.
9. The method for determining riverbed morphology characteristic parameters based on digital river topography according to claim 8, characterized in that, The specific parameters for determining the cross-sectional morphological characteristics of the river channel are as follows: Assume that there are a total of There are several grid nodes, which are numbered sequentially from one bank to the other, with serial numbers ranging from 1 to... ; The flat beach elevation of each cross section is set as the arithmetic mean of the elevation points on both banks. Based on the elevation of each grid node, the flat beach elevation of the cross section is calculated using the following formula: In the formula, The elevation of the flat beach at the cross-section; , Elevations of grid nodes 1 and m respectively; Based on the flat elevation of the cross-section and the coordinate set of grid nodes on each cross-section, the flat area of the cross-section is calculated; the calculation formula is as follows: In the formula, The flat area of the cross-section; The elevation of the flat beach at the cross-section; , Grid nodes and Elevation; , Grid nodes i and In Coordinates of direction; Based on the width of the flat river channel and the flat area of each cross-section, the flat water depth of the cross-section is calculated using the following formula: In the formula, The flat area of the cross-section; The width of the flat riverbed at each cross-section line; The depth of the flat beach at the cross-section; Based on the width of the river channel at each cross-section and the water depth at the cross-section, the width-to-depth ratio is calculated using the following formula: In the formula, Aspect ratio; The width of the flat riverbed at each cross-section line; The depth of the flat beach at the cross-section; Based on the elevation and coordinate point set of the grid nodes on each cross-section, the wetted perimeter is calculated using the following formula: In the formula, For wet period; , Grid nodes and Elevation; , Grid nodes i and In Coordinates of direction.
10. The method for determining riverbed morphology characteristic parameters based on digital river topography according to claim 8, characterized in that, The morphological characteristics of the deepwater longitudinal profile are determined as follows: Set the thalweg point as the point with the lowest elevation on the cross section, and find the elevation set {z} of the grid nodes on each cross section. 1, z2,……,z m The minimum value z in} k The coordinates of the thalweg point are obtained. s k , n k , z k ,) and their corresponding xyz Coordinates in a coordinate system x k , y k , z k ,), that is, to determine the location and elevation of the thalweg on each cross section; Based on determining the location and elevation of the thalweg at each cross-section, in xyz Draw the plan view and longitudinal section view of the thalweg in the coordinate system; Assuming grid nodes on the cross-section i and i The thalweg elevations of +1 are respectively z i,min , z i+1,min Calculate the grid nodes on the cross section i and i The distance Δs between +1 is used to calculate the longitudinal gradient of the cross-section, and the calculation formula is as follows: J =（ z i,min — z i+1,min ） / △s In the formula, J The longitudinal ratio of the deep channel in the cross-section decreases; z i,min , z i+1,min These are the grid nodes on the cross-section. i and i +1 depth level; △s represents the grid nodes on the cross-section. i and i The distance between +1.