Tidal flat hydrodynamic numerical modeling method and system based on nautical chart and underwater terrain inversion

Abstract

The invention provides a tidal flat hydrodynamic numerical modeling method and system based on a nautical chart and underwater terrain inversion, and relates to the field of coastal zone dynamic process modeling. The method aims at providing an efficient initial research means for tidal flat hydrodynamic process analysis and treatment. The method comprises the following steps: firstly, extracting sea chart topographic data as a model subtidal zone topography, and obtaining a model intertidal zone topography through water depth remote sensing inversion; then, two sets of topographic data are integrated, a sea chart topography is used in a subtidal zone, and a remote sensing inversion topography is used in an intertidal zone; secondly, only the sea chart topography is used, and intertidal zone topographic changes are ignored. Two sets of terrains are adopted to construct models respectively, and the finally selected model is determined through verification result comparison and error analysis. According to the method, the accuracy of tidal flat hydrodynamic simulation can be ensured, and meanwhile the economical efficiency and timeliness are improved.

Description

Technical Field

[0001] The present invention relates to the field of modeling of coastal zone dynamic processes, and particularly to a method and system for numerical modeling of tidal flat hydrodynamic based on chart and underwater terrain inversion. Background Art

[0002] As a widely distributed coastal wetland ecosystem, tidal flats have important ecological and economic values. Tidal flats usually have complex geomorphic forms, with the intertidal zone and subtidal zone intricately distributed. Coupled with various dynamic effects in the coastal zone (such as tides, runoff, wind, ocean currents, etc.), special hydrodynamic conditions are formed, which have important impacts on their geomorphic evolution and ecological functions. Mathematical models, with advantages such as the spatio-temporal continuity and convenience of obtaining information, have become an important means for studying and predicting tidal flat hydrodynamic processes.

[0003] However, in the process of establishing a tidal flat hydrodynamic model, the accuracy of terrain data and boundary information determines the prediction performance of the model, which often depends on high-precision measured data, increasing the modeling cost to a certain extent. Modern chart terrain data usually has high accuracy and resolution, and is relatively convenient to obtain, but lacks detailed information on the intertidal zone terrain. Tidal flats usually have a considerable intertidal area. The model established based on chart terrain ignores the influence of the intertidal zone terrain, which may reduce the accuracy of simulating tidal flat hydrodynamic processes. Summary of the Invention

[0004] Object of the Invention: To propose a method for numerical modeling of tidal flat hydrodynamic based on chart and underwater terrain inversion, and further propose a system that can execute the above method, aiming to improve the accuracy and timeliness of hydrodynamic simulation of tidal flats or similar coastal zone systems (such as lagoons, estuaries, and bays).

[0005] A method for numerical modeling of tidal flat hydrodynamic based on chart and underwater terrain inversion proposed by the present invention comprises the following steps:

[0006] Step 1: Plot the shoreline contour of the study area based on Geographic Information System (GIS), and extract chart data within a predetermined range around the study area;

[0007] Step 2: Geographically register the chart data, and extract the bathymetric contour and scattered point elevation information therein to obtain the chart terrain;

[0008] Step 3: Calculate the elevation of the intertidal zone area through a remote sensing inversion algorithm, including preliminary fitting, correction, and secondary fitting, to obtain the remotely sensed inversion terrain, and verify it through the chart data;

[0009] Step 4: Retain two sets of terrain data. For one set, use the chart terrain obtained in Step 2 in the subtidal zone and the remotely sensed inversion terrain obtained in Step 3 in the intertidal zone. For the other set, only use the chart terrain and ignore the terrain changes in the intertidal zone.

[0010] Step 5: Divide the computational domain grid and interpolate the two sets of terrain data obtained in Step 4 into the computational domain grid respectively.

[0011] Step 6: Based on the two sets of computational domain grids interpolated in Step 5, construct mathematical models respectively, set the model parameters, and generate boundary conditions.

[0012] Step 7: Conduct simulations using the two sets of mathematical models generated in Step 6 respectively, verify the mathematical models based on the measured data, evaluate the accuracy of the mathematical models, including water level and flow velocity verification, and determine the finally adopted mathematical model according to the verification results and error analysis.

[0013] In a further embodiment, the process of calculating the elevation of the intertidal zone area through the remotely sensed inversion algorithm specifically includes:

[0014] Identify the effective water area range of the study area, distinguish the foreground and background in the near-infrared band, and filter out the land and high-frequency noise areas.

[0015] Separate the red light layer partition L r 、green light layer partition L g and blue light layer partition L b respectively in the red, green, and blue bands;

[0016] Preliminary fitting: For each spectral layer partition, calculate the underwater terrain elevation H using the following formula:

[0017]

[0018] In the formula, λ r 、λ g and λ b represent the red, green, and blue bands of the multispectral pixel λ respectively; R is the reflectance of each band; l r 、l g and l b are the linear regression slopes of each layer respectively; b r 、b g and b b are the linear regression intercepts of each layer respectively; m is the adjustment coefficient.

[0019] Calibration: Compare the fitting results with the measured data extracted from the chart, count the deviation degree between the fitting results and the measured values of all measurement points, and eliminate the measurement points with the largest deviation degree (outliers) according to a certain proportion (generally not exceeding 20%).

[0020] Quadratic fitting: Fit the remaining measurement points again to obtain the linear regression slopes and intercepts of each spectral layer for calculating the underwater terrain elevation.

[0021] In a further embodiment, when distinguishing the foreground and background in the near-infrared band, the (Otsu algorithm) is adopted, that is, the optimal reflectance threshold δ for distinguishing the foreground and background is obtained through the traversal method R , so that the between-class variance reaches the maximum, and the between-class variance D R The calculation formula is as follows:

[0022] D R = ω 0 ω 1 (R 0 -R 1 ) 2

[0023] In the formula, ω 0 , ω 1 are the proportions of background and foreground pixel points in the total number of image pixels respectively; R 0 , R 1 are the average reflectances of the background and foreground respectively.

[0024] In a further embodiment, the boundary conditions described in step 6 include the tidal level process, ocean current process, water surface wind action, and runoff process at the river inlet of the sea area boundary.

[0025] In a further embodiment, the application of the ocean current process is achieved by interpolating the ocean current velocity components to the center points of the calculation domain boundary.

[0026] In a further embodiment, the two sets of terrain data are as follows:

[0027] The first set of terrain data: the integration of the subtidal terrain extracted from the nautical chart data and the remotely sensed inversion terrain of the intertidal zone;

[0028] The second set of terrain data: the terrain extracted from the nautical chart data is used for both the subtidal zone and the intertidal zone.

[0029] In a further embodiment, the control equation of the mathematical model described in step 6 is the depth-averaged Reynolds-averaged equation:

[0030] U t +E(U) x +H(U) y = S(U)+T(U)

[0031] In the formula, U t is the time derivative term; E(U) x , H(U) y are the horizontal convection terms; S(U) is the source term; T(U) is the turbulent viscous stress term.

[0032] In a further embodiment, the error analysis quantifies the accuracy of the mathematical model based on the model performance parameters. The calculation formula for the model performance parameter Skill is as follows:

[0033]

[0034] where x i represents the i-th simulated value; y i represents the i-th measured value; represents the arithmetic mean of the measured values; and n represents the number of measured values.

[0035] In addition, the present invention also provides a tidal flat hydrodynamic numerical modeling system, which automatically executes the above-mentioned tidal flat hydrodynamic numerical modeling method based on the inversion of the chart and underwater topography. The tidal flat hydrodynamic numerical modeling system includes:

[0036] A chart data extraction module, which is used to access the geographic information system, plot the shoreline contour of the study area, and extract the chart data within a predetermined range around the study area;

[0037] A chart topography generation module, which is used to perform geographic registration on the chart data, extract the isobath and scatter point elevation information therein, and obtain the chart topography;

[0038] A calculation module, which is used to calculate the elevation of the intertidal zone area, obtain the remotely sensed inverted topography, and verify it through the chart data;

[0039] A topographic data output module, which is used to output two sets of topographic data. One uses the chart topography in the subtidal zone area and the remotely sensed inverted topography in the intertidal zone area; the other only uses the chart topography;

[0040] A mesh generation module, which is used to divide the computational domain mesh and interpolate the two sets of topographic data into the computational domain mesh respectively;

[0041] A model construction module, which is used to construct mathematical models respectively based on the two sets of interpolated computational domain meshes, set the model parameters, and generate the boundary conditions;

[0042] A simulation verification module, which performs simulations respectively using the two generated mathematical models, verifies the mathematical models based on the measured data, evaluates the accuracy of the mathematical models, including water level and flow velocity verification, and determines the finally adopted mathematical model according to the verification results and error analysis.

[0043] In addition, the present invention also provides a computer-readable storage medium, which stores computer-executable instructions for causing a computer to execute the above-mentioned tidal flat hydrodynamic numerical modeling method based on the inversion of the chart and underwater topography.

[0044] Beneficial effects: The present invention proposes a tidal flat hydrodynamic numerical modeling method based on chart and underwater terrain inversion. The elevation information of the subtidal zone and intertidal zone in the study area is extracted and integrated through chart and optical remote sensing inversion. The data required for the model boundary conditions are obtained on an open access platform, and the spatio-temporal sequence of the model boundary conditions is obtained through interpolation, thereby establishing a tidal flat hydrodynamic model. This modeling method can improve its economy and timeliness while ensuring the accuracy of tidal flat hydrodynamic simulation, thus providing a scientific basis for tidal flat environmental governance. The method proposed by the present invention is not only applicable to the hydrodynamic simulation of tidal flats, but can also be extended to the research and governance of the dynamic processes of other similar coastal systems (such as lagoons, estuaries, and bays). Description of the Drawings

[0045] Figure 1 It is a flowchart for the implementation of the tidal flat hydrodynamic numerical modeling method in the embodiment.

[0046] Figure 2 It is a schematic diagram of the plot of the coastline and chart topography of the Xincun Lagoon.

[0047] Figure 3 It is a multi-spectral remote sensing image and spectral layer zoning map of the Xincun Lagoon.

[0048] Figure 4 It is a diagram of the underwater terrain inversion result and error analysis.

[0049] Figure 5 It is a schematic diagram of the grid division of the calculation domain of the Xincun Lagoon (85 elevation system).

[0050] Figure 6 It is a statistical chart of the water level verification results of the mathematical model at each measuring point at different times under two terrain conditions.

[0051] Figure 7 It is a statistical chart of the flow velocity verification results of the mathematical model at each measuring point at different times under two terrain conditions.

[0052] Figure 8 It is a schematic diagram of the composition of the tidal flat hydrodynamic numerical modeling in the embodiment. Detailed Embodiment

[0053] In the following description, a large number of specific details are given to provide a more thorough understanding of the present invention. However, it is obvious to those skilled in the art that the present invention can be implemented without one or more of these details. In other examples, in order to avoid confusion with the present invention, some technical features well known to the art are not described.

[0054] Taking the Xincun Lagoon as an example to illustrate the specific implementation steps of the present invention. The Xincun Lagoon is located in Lingshui County, southeastern Hainan Province, with a water area of 19.4 km 2。The internal terrain of the lagoon is complex, with an intricate distribution of intertidal zones and coral reef flats. There are two small rivers, the Qugang River and the Sancai River, flowing into the lagoon, and the average annual runoff is less than 10 m 3 / s. There are obvious seasonal variations in the wind speed and direction in the study area. The east wind and northeast wind are dominant in winter, while the south wind and southwest wind are dominant in summer. The wind speed is relatively high from July to November, and the maximum wind speed can reach 28 m / s in October.

[0055] Xincun Lagoon has a typical "big belly and small mouth" morphological feature. It is about 5.5 km long from east to west and about 4.5 km wide from north to south. It is connected to the ocean through a tidal channel about 1.7 km long and with an average width of 250 m (see Figure 2 ). There is a large area of subtidal plain in the center of the lagoon, with a depth of 5 - 11 m. The southern part of the lagoon is surrounded by intertidal flats, within a range of 500 - 600 m near the shoreline, and the beach slope is gentle; the northern part is a continuous distribution of coral reef flats with sediment deposition on them. The two constitute the main intertidal areas of the lagoon, accounting for 36.1% of its total area.

[0056] Step 1: Plotting the shoreline of the study area

[0057] In the GIS software, plot the shoreline of the study area through a polyline. First, it is necessary to obtain the satellite remote sensing images of the study area, and select the Sentinel - 2 (L2A) multispectral remote sensing images. The acquisition period of the selected remote sensing images should be as close as possible to the chart data described in Step 2. In addition, remote sensing images with lower cloud cover (no clouds covering the water area), higher tide levels (larger submerged range of the intertidal zone), and lower water turbidity should be selected as much as possible. The acquisition time of the selected remote sensing images in this example is 3:22 on December 28, 2018 (UTC), and the number is 49QDA.

[0058] Import the remote sensing images into the GIS software (taking ArcGIS10.5 as an example), and plot the shoreline through polyline vectorization. The shoreline should fit the non - water boundary of the study area according to the actual situation in the figure, and the retained morphological features should meet the accuracy requirements of the simulation (such as Figure 2 (a)).

[0059] Step 2: Extracting the chart topography

[0060] Taking the chart drawn by a certain party as an example, for this study area, select the harbor chart (number C1516311, updated on March 21, 2016). If the chart is in paper form, it needs to be scanned into an electronic picture format (such as PNG format). Import the chart into the GIS software (such as ArcGIS10.5). The chart picture itself has no coordinate system, and it is necessary to make the chart shoreline overlap with the lagoon shoreline plotted in Step 1 through georeferencing and scaling operations.

[0061] After registration is completed, the bathymetric contours in the nautical chart are plotted by polyline vectors. For easy distinction, different colored lines can be used for plotting. Subsequently, the elevation points in the nautical chart are plotted by point vectors (see Figure 2 (b)). The vector plotting must include geographical location information, which can be longitude and latitude or projected coordinates. For elevation points, the position and elevation information (X, Y, Z coordinates) need to be saved in ASCII text format. In ArcGIS 10.5, the X and Y coordinates can be added through the field calculator in the layer attribute table. In addition, the elevation information (Z coordinate) needs to be added manually.

[0062] Step 3: Terrain remote sensing inversion

[0063] Filter the multispectral remote sensing images. Preprocessing operations such as atmospheric correction and geometric correction can be completed through remote sensing image processing software (such as ENVI 5.6) to obtain the reflectance data of the study area. The Sentinel-2 (L2A) remote sensing images have completed the atmospheric correction and geometric correction steps at the time of release, and the pixel values of each band are the reflectance.

[0064] The Stumpf ratio method based on spectral stratification is used for underwater terrain inversion. This method utilizes the differences in the water penetration of different spectral bands. For different water depth ranges, different bands are used for underwater terrain inversion, which can effectively improve the inversion accuracy. First, the effective water area range of the study area must be identified, which requires distinguishing the foreground and background of the near-infrared band through the Otsu algorithm and filtering out the land and high-frequency noise areas (i.e., the near-infrared layer partition L nr ). The Otsu algorithm obtains the optimal reflectance threshold δ R for distinguishing the foreground and background through the traversal method, so that the between-class variance reaches the maximum. The calculation formula of the between-class variance D R is as follows:

[0065] D R = ω 0 ω 1 (R 0 - R 1 ) 2 (1)

[0066] In the formula: ω 0 , ω 1 are the proportions of background and foreground pixel points in the total number of pixels in the image respectively; R 0 , R 1 are the average reflectances of the background and foreground respectively.

[0067] Subsequently, the Otsu algorithm is continued to screen out the red light layer, green light layer, and blue light layer partitions L r , L g and L bThese three partitions respectively correspond to water depth regions of 0 - 5m, 5 - 15m, and >15m. The spectral stratification of Xincun Lagoon and its surrounding waters is shown in Figure 3 ( Figure 3 In (a) is the multispectral image; (b) is the spectral layer partition. Since the water depth of the lagoon and adjacent waters is less than 15m, the blue light layer partition does not appear in the figure.

[0068] For each spectral layer partition, the following method is used to calculate the underwater topography:

[0069]

[0070] In the formula: H is the elevation of the underwater topography to be obtained; λ r , λ g and λ b respectively represent the red, green, and blue light bands of the multispectral pixel λ; R is the reflectance of each band; l r , l g and l b are the linear regression slopes of each layer respectively; b r , b g and b b are the linear regression intercepts of each layer; m is the adjustment coefficient, and in this example, m = 1.0. In the formula, the linear fitting slopes l and intercepts b of different partitions are obtained by linearly fitting the reflectance ratios at each measuring point on the nautical chart (such as for the red light layer, it is ln[mR(λ r )] / ln[mR(λ g )]) with the elevation extracted from the nautical chart.

[0071] To exclude the interference of obstacles or artificial light sources on the reflectance of some measuring points, first, a preliminary fitting is carried out, and the deviation degree of the fitting results of each measuring point from the measured values is statistically analyzed. The measuring points with the largest deviation degree (outliers) are excluded at a ratio of 10%, that is, correction. Subsequently, a secondary fitting is carried out on the remaining measuring points to obtain the linear regression slopes and intercepts of each layer, which are used to calculate the elevation of the underwater topography. The following root mean square error RMSE and mean absolute error MAE are used for error analysis:

[0072]

[0073] In the formula: x and y respectively represent the true value and the predicted value. Table 1 shows the linear fitting parameters and error parameters of each spectral layer. Figure 4 is the inversion result and error distribution of the underwater topography of Xincun Lagoon and adjacent waters, Figure 4 (a) is the elevation inversion result; Figure 4 (b) is the deviation of the inversion results of each spectral layer partition from the measured results on the nautical chart; Figure 4 (c) is the error frequency distribution.

[0074] Table 1 Linear fitting parameters and error parameters of each spectral layer

[0075]

[0076] Step 4: Integration of Topographic Data

[0077] Two sets of topographic data are retained. One is that the subtidal zone adopts the chart topography and the intertidal zone adopts the remote sensing inversion results, denoted as Topography A; the other is that only the chart topography is adopted and the water depth below the lowest tide level in the intertidal zone is counted as 0 m, denoted as Topography B.

[0078] Step 5: Mesh Generation

[0079] For the two topographies (Topography A and Topography B) in Step 4, the same mesh generation is used to calculate the domain. It covers the entire lagoon and a 10-km range around the lagoon mouth. The triangular mesh is used to divide the calculation domain, including 35,855 units and 54,301 unit edges (such as Figure 5 ). The side length of the unit gradually increases from 3 - 5 m near the lagoon mouth to the interior of the lagoon and the open boundary of the sea area, and the maximum side length of the open boundary is 823 m. In areas where the topographic shoreline changes violently and in narrow areas, the mesh is locally refined. The elevations of the mesh units and edges are obtained by linear interpolation from the above topographic data, and inverse distance weighted interpolation is used at the mesh boundary.

[0080] Step 6: Generation of Boundary Conditions

[0081] The boundary conditions include the tidal level and ocean current processes at the open boundary of the sea area, the flow rate process at the river inlet boundary, and the wind vector process applied to the entire calculation domain.

[0082] (1) Tidal level process: The tidal level process at the open boundary of the sea area is obtained by using the Chinatide tidal prediction model for the Chinese sea area. The model considers a total of eight astronomical partial tides, namely Q1, P1, O1, K1, N2, M2, S2, and K2, and the long-period meteorological partial tide Sa, and has a high prediction accuracy for the tidal level in the Chinese sea area.

[0083] (2) Ocean current process: The model considers the ocean current process in the adjacent sea area and realizes it by applying the velocity vector U Sea =(u, v) on the open boundary edge of the sea area. The velocity data on the sea area boundary comes from the hindcast results of the 3.1 version of the US Global Ocean Forecast System (GOFS). It covers the global spatial range with a latitude of 80S - 90N. The spatial resolution within the range of 40S - 40N is 0.08° lon × 0.08° lat, and it is divided into 41 layers vertically. Using the 3-hour ocean current field analysis data of GOFS 3.1, the velocity vector at the center of the open boundary edge of the outer sea is obtained through vertical averaging and linear interpolation and applied to the boundary edge.

[0084] (3) Runoff: The runoff injection of the river is achieved by specifying the flow process at the river inlet boundary. The flow process is specified at a daily time resolution and is obtained by the following method:

[0085] For the two main rivers (Qugang River, Sancai River), the runoff coefficient is calculated by the following formula:

[0086]

[0087] In the formula, R ave,i is the average annual runoff depth (mm) of each river, which is the ratio of the average annual runoff of the river to the basin area, and I ave,i is the average annual rainfall (mm) of the basin. Thus, the calculation formula for the runoff of each river is as follows:

[0088]

[0089] In the formula: A i is the catchment area (km 2 ) of the river basin, and I i is the rainfall intensity (mm / h) averaged monthly or daily in the basin, obtained from the measured data of nearby rain gauges. The hydrological characteristic parameters of each inflowing river of Xincun Lagoon are shown in Table 2.

[0090] Table 2 Hydrological Characteristic Parameters of Xincun Lagoon

[0091]

[0092] Step 7: Construct a mathematical model

[0093] A two-dimensional hydrodynamic mathematical model is used to simulate the hydrodynamic and mass transport processes in the lagoon under the action of runoff, astronomical tide, ocean current and wind. The governing equations of the model are the depth-averaged Reynolds-averaged equations:

[0094] U t +E(U) x +H(U) y =S(U)+T(U)(7)

[0095] The meanings of each term in the above formula are as follows:

[0096]

[0097] In the formula: u and v are the flow velocities in the x and y directions respectively, h is the water depth, g is the acceleration due to gravity, f is the Coriolis force coefficient, and the subscripts t, x, and y represent the partial derivatives in the time and x, y directions respectively. The effects of bed friction, bed slope and wind stress are represented by source terms (S fx , S fy ), (S 0x , S 0y ) and (Swx , S wy ) is applied to the model in the form of:

[0098]

[0099] where: C D = gn 2 H -1 / 3 is the bed surface resistance coefficient, where n is the Manning roughness coefficient. According to the sediment particle size distribution of the Xin Village Lagoon and previous research results, n = 0.019 is taken; z b is the bed surface elevation; C DW is the wind drag coefficient. For weak winds (U 10 ≤5 m / s), C DW = 0.0013, while for moderate and strong winds, C DW is taken as 0.0025, (u 10 , v 10 ) is the wind speed vector at 10 m above the water surface; T xx , T xy , T yy is the viscous stress term, where the molecular viscosity coefficient ν is replaced by the horizontal turbulent viscosity coefficient ν H , that is:

[0100]

[0101] In the above formula, the turbulent viscosity coefficient ν H is parameterized by the Smagorinsky method to close the control equations. Its calculation formula is as follows:

[0102]

[0103] where: A ν is a constant, with a value of 0.2; Ω e is the control cell area.

[0104] The model discretizes the computational domain using the finite volume method based on the Godunov format. The solution of the control equations takes the unit as the control volume, the conserved variables are stored at the center of each unit, and the edges of each unit are defined as the faces of the control volume. The calculation of the fluxes on the faces uses Roe's approximate Riemann solution. The time integration uses the TVD Runge - Kutta method with second - order accuracy. The topographic source term uses the upwind difference format, and other source terms are discretized by the semi - implicit format of the algorithm. The two sets of grids divided in step 5 are used to establish the model respectively. The tidal level and ocean current processes are given at the open boundary of the sea area, the flow rate process is given at the river inlet, and the wind stress source term is applied to the entire computational domain and calculated through (Equation 9).

[0105] Step 8: Model verification and selection

[0106] The simulated results are compared with the measured tidal levels and flow velocity data to calibrate and verify the model parameters. The verification is carried out in three time periods, and the information of the verification time periods and points is shown in Table 3. To evaluate the accuracy of the model, the following root mean square error RMSE (see Equation 3) and model performance parameter (Skill) are used to quantify the accuracy of the model. The calculation formula of the latter is as follows:

[0107]

[0108] In the formula: x and y represent the measured and simulated data respectively. Skill is a dimensionless parameter, and its values in the ranges of 0.65 - 1.0, 0.5 - 0.65, 0.2 - 0.5, and 0 - 0.2 represent excellent, good, qualified, and poor model performances respectively.

[0109] Table 3 Information of verification time periods and points of Xincun Lagoon

[0110]

[0111]

[0112] Figure 6 、 Figure 7 The comparison between the simulated results and the measured data at each measurement point of Xincun Lagoon under two topographic conditions (Topography A and Topography B) is given. It can be seen that the simulated results of water level and flow velocity under the two topographic conditions are basically consistent with the measured data. Among them, the Skill of the water level simulation results > 0.9, and the Skill of the flow velocity simulation results is between 0.56 and 0.89, meeting the simulation accuracy requirements. If the combination of intertidal zone topographic inversion and subtidal zone chart topography (Topography A) is adopted, the Skill of the simulated results as a whole increases and the RMSE as a whole decreases at each measurement point and each time period compared with only using the chart topography. Therefore, the model built based on Topography A is finally selected.

[0113] In addition, this embodiment also discloses a tidal flat hydrodynamic numerical modeling system 900, as shown in Figure 8 . This system consists of a chart data extraction module 901, a chart topography generation module 902, a calculation module 903, a topographic data output module 904, a grid division module 905, a model construction module 906, and a simulation verification module 907.

[0114] The chart data extraction module 901 accesses the geographic information system, which is used to plot the shoreline contour of the research area and extract the chart data within a predetermined range around the research area. The chart terrain generation module 902 is used to perform geographic registration on the chart data, extract the isobath and scatter point elevation information therein, and obtain the chart terrain. The calculation module 903 is used to calculate the elevation of the intertidal zone area to obtain the remotely sensed inversion terrain, and verify it through the chart data. The terrain data output module 904 outputs two sets of terrain data. One uses the chart terrain in the subtidal zone area and the remotely sensed inversion terrain in the intertidal zone area; the other only uses the chart terrain. The grid division module 905 is used to divide the computational domain grid and interpolate the two sets of terrain data into the computational domain grid respectively. The model construction module 906 respectively constructs mathematical models based on the two sets of computational domain grids obtained by interpolation, sets model parameters, and generates boundary conditions. The simulation verification module 907 performs simulations respectively using the two sets of generated mathematical models, verifies the mathematical models based on the measured data, evaluates the accuracy of the mathematical models, including water level and flow velocity verification, and determines the finally adopted mathematical model according to the verification results and error analysis.

[0115] When the above tidal flat hydrodynamic numerical modeling system 900 is running, it can automatically execute all or part of the process flow of the tidal flat hydrodynamic numerical modeling method based on chart and underwater terrain inversion disclosed in the above embodiments, which will not be elaborated here.

[0116] Moreover, the tidal flat hydrodynamic numerical modeling method based on chart and underwater terrain inversion disclosed in the above embodiments can be implemented in whole or in part by software, hardware, firmware or any other combination. When implemented using software, the above embodiments can be implemented in whole or in part in the form of a computer program product. The computer program product includes one or more computer instructions or computer programs.

[0117] When the computer instructions or computer programs are loaded or executed on a computer, the processes or functions described in the embodiments of the present application are generated in whole or in part. The computer can be a general-purpose computer, a special-purpose computer, a computer network, or other programmable devices. The computer instructions can be stored in a computer-readable storage medium or transmitted from one computer-readable storage medium to another. For example, the computer instructions can be transmitted from one website, computer, server or data center to another website, computer, server or data center by wired (such as infrared, wireless, microwave, etc.) means. The computer-readable storage medium can be any available medium that the computer can access or a data storage device such as a server or data center containing one or more sets of available media. The available medium can be a magnetic medium (such as a floppy disk, hard disk, magnetic tape), an optical medium (such as a DVD), or a semiconductor medium. The semiconductor medium can be a solid-state drive.

[0118] It should be understood that in various embodiments of the present application, the magnitude of the sequence numbers of the above processes does not mean the order of execution. The order of execution of each process should be determined by its function and internal logic, and should not constitute any limitation on the implementation process of the embodiments of the present application.

[0119] Although the present invention has been shown and described with reference to specific preferred embodiments, it should not be construed as a limitation of the present invention itself. Various changes can be made in its form and details without departing from the spirit and scope of the present invention defined by the appended claims.

Claims

1. A tidal flat hydrodynamic numerical modeling method based on nautical chart and underwater topography inversion, characterized in that: The steps include: Step 1: Plot the coastline contour of the study area based on the Geographic Information System (GIS) and extract the nautical chart data of the predetermined range around the study area; Step 2: geo-reference the nautical chart data, extract the isobaths and scattered elevation information therein, and obtain the nautical chart topography; Step 3, calculating the intertidal zone elevation in the study area by remote sensing inversion algorithm, including preliminary fitting, correction and secondary fitting, to obtain remote sensing inversion terrain, and verifying it by the nautical chart data; Step 4, retaining two sets of terrain data, one using the chart terrain obtained in step 2 in the subtidal zone, and the remote sensing inversion terrain obtained in step 3 in the intertidal zone; the other using only the chart terrain, ignoring the intertidal zone terrain changes; Step 5, dividing the computational domain into grids, and interpolating the two sets of terrain data obtained in step 4 into the computational domain grids respectively; Step 6: Based on the two sets of computational domain grids interpolated in step 5, mathematical models are constructed respectively, model parameters are set, and boundary conditions are generated; Step 7: Use the two sets of mathematical models generated in step 6 to perform simulations respectively, verify the mathematical models based on measured data, evaluate the accuracy of the mathematical models, including water level and flow rate verification, and determine the final mathematical model based on the verification results and error analysis.

2. The method for numerical modeling of tidal flat hydrodynamics based on nautical chart and underwater topography inversion according to claim 1 is characterized in that: The calculation of the intertidal zone elevation by remote sensing inversion algorithm specifically includes: Identify the effective water area of ​​the target tidal flat, distinguish the foreground and background of the near-infrared band, and filter out the land and high-frequency noise areas; The red light layer partition L is screened out in the red, green and blue bands respectively r , Green light layer partition L g and Blu-ray Layer Partition L b ; Preliminary fitting: For each spectral layer partition, the underwater terrain elevation H is calculated using the following formula: In the formula, λ r , h and λ b Represent the red, green and blue light bands of the multispectral pixel λ respectively; R is the reflectivity of each band; l r , l g and l b are the linear regression slopes of each stratum; b r , b g and b b is the linear regression intercept of each stratum; m is the adjustment coefficient; Correction: Compare the fitting results with the measured data extracted from the nautical chart, count the deviations between the fitting results and the measured values ​​of all measuring points, and remove the measuring points with the largest deviations according to a predetermined ratio; Secondary fitting: Fit the remaining measurement points again to obtain the linear regression slope and intercept of each spectral layer, which are used to calculate the underwater terrain elevation.

3. The tidal flat hydrodynamic numerical modeling method based on nautical chart and underwater terrain inversion according to claim 1 is characterized in that: The boundary conditions described in step 6 include the tidal process at the sea area boundary, the ocean current process, the surface wind effect, and the runoff process at the river entrance.

4. The method for numerical modeling of tidal flat hydrodynamics based on nautical chart and underwater terrain inversion according to claim 4 is characterized in that: The application of the ocean current process is achieved by interpolating the ocean current velocity component to the center point of the boundary edge of the calculation domain.

5. The method for numerical modeling of tidal flat hydrodynamics based on nautical chart and underwater topography inversion according to claim 1 is characterized in that: The two sets of terrain data include: The first set of topographic data: integration of subtidal topography extracted from chart data and intertidal topography inverted from remote sensing; The second set of terrain data: The subtidal zone and intertidal zone are both extracted using the terrain extracted from the nautical chart data.

6. The method for numerical modeling of tidal flat hydrodynamics based on nautical chart and underwater topography inversion according to claim 1 is characterized in that: The governing equation of the mathematical model described in step 6 is the depth-averaged Reynolds time-averaged equation: IN t +E(U) x +H(U) y =S(U)+T(U) Where U t is the time derivative term; E(U) x 、H(U) y is the horizontal convection term; S(U) is the source term; T(U) is the turbulent viscous stress term.

7. The method for numerical modeling of tidal flat hydrodynamics based on nautical chart and underwater topography inversion according to claim 1 is characterized in that: The error analysis is based on the accuracy of the model quantified by the model performance parameter Skill, which is calculated as follows: In the formula, x i represents the i-th simulation value; y i represents the i-th measured value; represents the arithmetic mean of the measured values; n represents the number of measured values.

8. A tidal flat hydrodynamic numerical modeling system, capable of executing the tidal flat hydrodynamic numerical modeling method based on nautical chart and underwater terrain inversion as claimed in any one of claims 1 to 7, characterized in that: include: Chart data extraction module; The chart data extraction module is connected to the geographic information system and is used to plot the coastline contour of the study area and extract chart data of a predetermined range around the study area; A nautical chart terrain generation module; the nautical chart terrain generation module is used to perform geo-registration on the nautical chart data, extract the isobaths and scattered elevation information therein, and obtain the nautical chart terrain; Calculation module: The calculation module is used to calculate the elevation of the intertidal zone area, obtain the remote sensing inversion terrain, and verify it through the nautical chart data; A terrain data output module; the terrain data output module outputs two sets of terrain data, one of which uses the chart terrain in the subtidal zone and the remote sensing inversion terrain in the intertidal zone; the other uses only the chart terrain; Grid division module; the grid division module is used to divide the computational domain grid and interpolate the two sets of terrain data into the computational domain grid respectively; Model building module; the model building module builds mathematical models based on two sets of computational domain grids obtained by interpolation, sets model parameters, and generates boundary conditions; Simulation verification module: The simulation verification module uses the two sets of generated mathematical models to perform simulations respectively, verifies the mathematical models based on measured data, evaluates the accuracy of the mathematical models, including water level and flow rate verification, and determines the final mathematical model based on the verification results and error analysis.

9. A computer-readable storage medium, characterized in that: The storage medium stores computer-executable instructions, and the computer-executable instructions are used to enable a computer to execute the tidal flat hydrodynamic numerical modeling method based on nautical chart and underwater terrain inversion as described in any one of claims 1 to 7.

Images