Data processing method based on ocean tide quasi-harmonic analysis

Abstract

The invention provides a data processing method based on ocean tide quasi-harmonic analysis, and relates to the technical field of physical oceanography. The data processing method based on ocean power flow quasi-harmonic analysis comprises the following steps: S1, based on satellite altimeter surface flow velocity data, coastal radar horizontal flow field data and buoy vertical profile data, aligning a space-time reference by adopting a space-time registration algorithm, dynamically fusing three types of data sources through Bayesian weighted Kalman filtering, eliminating instrument system errors, and obtaining a time-space reference; and generating a three-dimensional power flow reference data set. By fusing multi-source ocean observation data and applying space-time registration and dynamic filtering technologies, high-precision three-dimensional power flow reference data is constructed, the data coverage range and reliability are remarkably expanded, in a partial tide parameter analysis stage, a least square fitting and complex vector superposition method is adopted, elliptic motion elements of a key partial tide are accurately solved, and the maximum power flow is obtained. And the integrity of tide motion rule description is ensured.

Description

Technical Field

[0001] This application relates to the field of physical oceanography, specifically a data processing method based on ocean current quasi-harmonic analysis. Background Technology

[0002] Physical oceanography is a core branch of marine science, focusing on the study of the physical properties of seawater (including temperature, salinity, and density distribution) and its dynamic patterns (such as tides, currents, waves, and turbulence), while also exploring the interaction mechanisms between seawater and the atmosphere, lithosphere, and biosphere. This field utilizes principles of fluid mechanics and thermodynamics, combined with mathematical modeling, field observations (buoys, moorings, shipboard instruments, satellite remote sensing), and numerical simulation techniques to quantitatively analyze ocean dynamic processes and their impacts on the climate system, ecosystems, and human activities. One data processing method based on ocean tidal current harmonic analysis refers to a technical solution that uses mathematical methods to decompose the periodic tidal components dominated by astronomical tides from ocean observation data. The core application of this method is to accurately extract the elliptical motion parameters (major axis, minor axis, and inclination angle) of key tidal constituents (such as O1, K1, M2, and S2), thereby serving scenarios such as the formulation of marine engineering safety design standards (e.g., calculation of scour for offshore wind farm foundations), prediction of channel sediment transport, pollution diffusion simulation, and tidal power generation site selection, providing a scientific basis for decision-making in marine activities.

[0003] The lack of a unified standard for determining marine dynamic conditions in existing technologies leads to chaotic selection of dynamic parameters in scour calculations, resulting in inconsistent results and engineering risks. Significant differences exist in key parameters calculated using different methods, causing inaccuracies in the strength design of engineering structures and increasing safety hazards. Current operational models do not fully consider the nonlinear effects of shallow water topography and cannot correct for shallow water influences, leading to systematic biases in nearshore tidal current predictions and affecting the reliability of sediment transport and pollution diffusion simulations. Extreme event handling relies on fixed tidal constituent combinations without dynamically linking meteorological event characteristics, resulting in inaccurate predictions under typhoon or storm surge scenarios. Real-world examples show that some marine structures have suffered damage due to underestimating extreme conditions. The interaction mechanism between tidal currents and waves is not integrated, and the contribution of wave radiation stress to the disturbance of the tidal field is ignored, causing predicted particle transport distances to deviate significantly from measured results, reducing the effectiveness of coastal protection engineering design.

[0004] Application content To address the shortcomings of existing technologies, this application provides a data processing method based on ocean current harmonic analysis. This method solves the problem of the lack of a unified standard for determining ocean dynamic conditions in existing technologies, which leads to chaotic selection of dynamic parameters in scour calculations, resulting in inconsistent results and engineering risks. Significant differences exist in the key parameters calculated by different methods, causing inaccuracies in the strength design of engineering structures and increasing safety hazards.

[0005] To achieve the above objectives, this application provides the following technical solution: a data processing method based on ocean current harmonic analysis, comprising the following steps: S1: Based on surface flow velocity data from satellite altimeters, horizontal flow field data from coastal radar, and vertical profile data from buoys, a spatiotemporal registration algorithm is used to align the spatiotemporal reference, and Bayesian weighted Kalman filtering is used to dynamically fuse the three types of data sources to eliminate instrument system errors and generate a three-dimensional tidal current reference dataset. S2: Based on the three-dimensional tidal reference dataset, the least squares quasi-harmonic analysis algorithm is used to solve the elliptical elements of the six tidal constituents O1, K1, M2, S2, M4, and MS4. A complex vector superposition model is introduced to calculate the amplitude and lag angle of each tidal constituent and generate the elliptical element matrix of the tidal constituents. S3: Based on the M2, M4, and MS4 tidal parameters in the tidal ellipse element matrix, the shallow water effect coefficients are constructed using the Legendre polynomial fitting algorithm. δ=α·(M4 / M2) 2 +β·(MS4 / M2) 2 α and β are generated by inversion through the topographic gradient dependence model, which corrects the tidal amplitude and generates a shallow water corrected tidal parameter set; S4: Based on the shallow water corrected tidal parameter set and the historical typhoon event database, the random forest contribution evaluation algorithm is used to dynamically identify highly correlated tidal constituents. Monte Carlo vector synthesis is used to simulate and generate tidal constituent combinations under extreme scenarios, and to generate optimized combinations of tidal constituents for extreme events. S5: Based on the M2 / S2 semi-major axis vector of the extreme event tidal optimization combination, the boundary conditions of the third-generation wave model are input, and the radiation stress-tidal shear coupling equation is adopted: Solve for the wave disturbance intensity on the tidal field to generate a coupled shear stress field; S6: Based on the coupled shear stress field and the shallow water corrected tidal parameter set, the complex vector integral algorithm is used to calculate the trajectory of tidal water particles: Output the migration distance vectors under three scenarios: spring tide, mid tide, and neap tide, and generate multi-scenario particle migration distance vectors; S7: Based on multi-scenario particle transport distance vectors, a GPU-accelerated vector synthesis algorithm is used to perform parallel calculations on millions of grids for: the maximum possible flow velocity of the tidal current, the maximum possible transport distance of particles, and the 3D tidal current elliptical visualization map, generating an engineering-grade tidal current feature parameter package.

[0006] Preferably, generating a three-dimensional power flow benchmark dataset based on S1 includes the following steps: S101: Based on the sea surface height anomaly data from the satellite altimeter, a radial basis function interpolation algorithm is used for spatial resampling to generate a gridded sea surface velocity field and generate a satellite surface velocity grid. S102: Based on satellite surface velocity grid and buoy ADCP vertical profile data, a dynamic optimal interpolation method is used for vertical correction to generate a vertically corrected velocity field. S103: Based on the vertically corrected velocity field and coastal radar horizontal flow data, a three-dimensional tidal current benchmark dataset is generated by fusing spatiotemporal outliers using Bayesian weighted Kalman filtering and removing outliers.

[0007] Preferably, generating the tidal constituent ellipse element matrix based on S2 includes the following steps: S201: Based on the three-dimensional tidal reference dataset, the fundamental frequency components are extracted using the fast Fourier transform preprocessing algorithm to generate the tidal frequency characteristic spectrum. S202: Based on the tidal frequency characteristic spectrum, the constrained least squares fitting method is used to solve for the initial tidal values ​​of O1 / K1 / M2 / S2 / M4 / MS6, and the initial solution of the tidal parameters is generated. S203: Based on the initial solution of the tidal parameters, the major semi-axis, minor semi-axis, inclination angle and rotation direction are calculated by the complex ellipse vector iteration algorithm to generate the tidal ellipse element matrix.

[0008] Preferably, the generation of shallow water corrected tidal parameter set based on S3 includes the following steps: S301: Based on the M2 and M4 tidal constituents in the tidal constituent ellipse element matrix, the shallow water intensity factor is calculated using a nonlinear ratio model: γ=(M4 amp / M2 amp ) 2 , generating a scalar field of shallow water intensity factor; S302: Based on the shallow water intensity factor scalar field and seabed topographic gradient data, the Legendre multinomial regression algorithm is used to fit the topographic dependence parameters α and β to generate the topographic coupling correction coefficient. S303: Based on the terrain coupling correction coefficient and the tidal ellipse element matrix, through the amplitude correction equation: A corr =A ori ·(1+αγ+βγ 2 Update the tidal parameters and generate a shallow water corrected tidal parameter set.

[0009] Preferably, the S4-based extreme event tidal optimization combination includes the following steps: S401: Based on the historical typhoon event database and the shallow water corrected tidal parameter set, the Pearson correlation matrix algorithm is used to calculate the correlation between tidal constituents and extreme current velocities, and generate the tidal constituent event correlation matrix. S402: Based on the correlation matrix of tidal events, a random forest feature importance assessment is used to screen high-contribution tidal periods and generate a ranking of key tidal contributions; S403: Based on the ranking of key tidal contributions, generate extreme scenario tidal combinations through Monte Carlo vector synthesis simulation, and generate optimized combinations of extreme event tidal combinations.

[0010] Preferably, the generation of the coupled shear stress field based on S5 includes the following steps: S501: Based on the M2 / S2 tidal vector in the extreme event tidal optimization combination, a spectral boundary condition generation algorithm is used to construct the wave model input field and generate a set of wave spectral boundary conditions. S502: Based on the wave spectrum boundary condition set, run the third-generation wave numerical model to solve the radiation stress tensor and generate the wave radiation stress field. S503: Solving the coupled shear stress equation based on the wave radiation stress field and the shallow water modified tidal parameter set: Generate a coupled shear stress field.

[0011] Preferably, the generation of multi-scenario particle transport distance vectors based on S6 includes the following steps: S601: Based on the coupled shear stress field, the velocity correction model is used to update the flow velocity vector of the split current. κ is an empirical coefficient used to generate a stress-corrected velocity field; S602: Calculate particle displacement using a complex harmonic integral algorithm based on a stress-corrected velocity field. Generate the original set of particle trajectories; S603: Based on the original set of particle trajectories, a spring tide / mid tide / neap tide scenario separation algorithm is used to classify and statistically analyze the transport distance, generating multi-scenario particle transport distance vectors.

[0012] Preferably, the generation of engineering-grade power flow characteristic parameters based on S7 includes the following steps: S701: Based on the multi-scenario particle transport distance vector, the extreme value statistical analysis algorithm is used to calculate the quantile transport distance and generate the maximum possible transport distance of the particle. S702: Based on the shallow water corrected tidal parameter set, the maximum possible tidal velocity is solved by the vector superposition extreme value search algorithm to generate the engineering extreme value velocity vector. S703: Based on the engineering extreme velocity vector, the three-dimensional power flow ellipse is reconstructed using the elliptic parameterization algorithm to generate a three-dimensional elliptic dynamic topology map; S704: Integrates the maximum possible transport distance of mass points, engineering extreme velocity vectors, and three-dimensional elliptical dynamic topology diagrams, and generates a standardized output package through the CUDA parallel compression encapsulation algorithm to generate an engineering-grade power flow characteristic parameter package.

[0013] In summary, this application includes at least one of the following beneficial technical effects: This application constructs high-precision three-dimensional tidal current benchmark data by fusing multi-source ocean observation data and applying spatiotemporal registration and dynamic filtering techniques, significantly expanding the data coverage and reliability. In the tidal parameter analysis stage, least squares fitting and complex vector superposition methods are employed to accurately solve the elliptical motion elements of key tidal constituents, ensuring the completeness of the tidal current motion description. A topographically dependent shallow-water effect correction model is introduced, and nonlinear correction coefficients are retrieved based on seafloor gradient characteristics, effectively improving the accuracy of tidal current characterization in estuaries and continental shelf areas. Through machine learning-driven tidal contribution assessment and stochastic simulation techniques, tidal combination under extreme events is dynamically optimized. This system enhances predictive adaptability by integrating the physical mechanisms of wave-tidal coupling, quantifying the disturbance effect of radiated stress on the tidal field, generating a stress field that better reflects the actual marine environment, applying harmonic integral algorithms to calculate particle migration trajectories, and combining multi-tidal scenario separation outputs to achieve refined prediction of migration distances. It employs a parallel computing architecture to accelerate engineering parameter generation, efficiently processes large-scale grid data, and comprehensively improves computational efficiency and visualization capabilities. The overall process forms a closed-loop optimization from data to application, achieving substantial breakthroughs in tidal parameter accuracy, extreme event adaptability, and engineering output efficiency, providing reliable technical support for marine engineering safety design and marine resource development. Attached Figure Description

[0014] Figure 1 This is a schematic diagram of the main steps of this application; Figure 2 This is a detailed schematic diagram of S1 in this application; Figure 3 This is a detailed schematic diagram of S2 in this application; Figure 4 This is a detailed schematic diagram of S3 in this application; Figure 5 This is a detailed schematic diagram of S4 in this application; Figure 6 This is a detailed schematic diagram of S5 in this application; Figure 7 This is a detailed schematic diagram of S6 in this application; Figure 8 This is a detailed schematic diagram of S7 in this application. Detailed Implementation

[0015] The following is in conjunction with the appendix Figure 1-8 This application will be described in further detail.

[0016] like Figure 1 As shown in the figure, this application provides a data processing method based on ocean current harmonic analysis, including the following steps: S1: Based on surface flow velocity data from satellite altimeters, horizontal flow field data from coastal radar, and vertical profile data from buoys, a spatiotemporal registration algorithm is used to align the spatiotemporal reference, and Bayesian weighted Kalman filtering is used to dynamically fuse the three types of data sources to eliminate instrument system errors and generate a three-dimensional tidal current reference dataset. S2: Based on the three-dimensional tidal reference dataset, the least squares quasi-harmonic analysis algorithm is used to solve the elliptical elements of the six tidal constituents O1, K1, M2, S2, M4, and MS4. A complex vector superposition model is introduced to calculate the amplitude and lag angle of each tidal constituent and generate the elliptical element matrix of the tidal constituents. S3: Based on the M2, M4, and MS4 tidal parameters in the tidal ellipse element matrix, the shallow water effect coefficients are constructed using the Legendre polynomial fitting algorithm. δ=α·(M4 / M2) 2 +β·(MS4 / M2) 2 α and β are generated by inversion through the topographic gradient dependence model, which corrects the tidal amplitude and generates a shallow water corrected tidal parameter set; S4: Based on the shallow water corrected tidal parameter set and the historical typhoon event database, the random forest contribution evaluation algorithm is used to dynamically identify highly correlated tidal constituents. Monte Carlo vector synthesis is used to simulate and generate tidal constituent combinations under extreme scenarios, and to generate optimized combinations of tidal constituents for extreme events. S5: Based on the M2 / S2 semi-major axis vector of the extreme event tidal optimization combination, the boundary conditions of the third-generation wave model are input, and the radiation stress-tidal shear coupling equation is adopted: Solve for the wave disturbance intensity on the tidal field to generate a coupled shear stress field; S6: Based on the coupled shear stress field and the shallow water corrected tidal parameter set, the complex vector integral algorithm is used to calculate the trajectory of tidal water particles: Output the migration distance vectors under three scenarios: spring tide, mid tide, and neap tide, and generate multi-scenario particle migration distance vectors; S7: Based on multi-scenario particle transport distance vectors, a GPU-accelerated vector synthesis algorithm is used to perform parallel calculations on millions of grids for: the maximum possible flow velocity of the tidal current, the maximum possible transport distance of particles, and the 3D tidal current elliptical visualization map, generating an engineering-grade tidal current feature parameter package.

[0017] See Figure 2 The generation of a three-dimensional power flow benchmark dataset based on S1 includes the following steps: S101: Based on the sea surface height anomaly data from the satellite altimeter, a radial basis function interpolation algorithm is used for spatial resampling to generate a gridded sea surface velocity field and generate a satellite surface velocity grid. Based on sea surface height anomaly data from satellite altimeters, a radial basis function interpolation algorithm is used for spatial resampling to generate a gridded sea surface velocity field. Specifically, satellite altimeter data from the northern South China Sea over a 24-hour period, including longitude, latitude, and sea surface height anomalies, is acquired. A Gaussian kernel function is selected as the radial basis function, with a kernel radius of 50 kilometers, and interpolation calculations are performed on 1°×1° grid cells. Taking a grid point (20°N, 115°E) as an example, sea surface height anomalies from its 30 neighboring points are extracted, and the weight of each point is calculated: weight = exp(-distance² / kernel radius²), where distance is in kilometers. For a neighboring point (20.1°N, 115.2°E), 25 kilometers from the target grid point, the weight = exp(-625 / 2500) ≈ 0.78. The weighted average of all neighboring point height anomalies is then calculated, combined with the geostrophic flow formula. (g-force acceleration 9.8 m / s²) 2 f Coriolis parameter 5×10 - 5 s -1 The density of seawater is 1025 kg / m³. 3 , The surface velocity at a given grid point is calculated using the height anomaly gradient. In this example, a height anomaly gradient of 0.1 m / 100 km yields a velocity of approximately 0.4 m / s. After traversing all grid points, a regularized surface velocity field covering the study area is output.

[0018] S102: Based on satellite surface velocity grid and buoy ADCP vertical profile data, a dynamic optimal interpolation method is used for vertical correction to generate a vertically corrected velocity field. Based on sea surface height anomaly data from satellite altimeters, a radial basis function interpolation algorithm is used for spatial resampling to generate a gridded sea surface velocity field. Specifically, satellite altimeter data from the northern South China Sea over a 24-hour period, including longitude, latitude, and sea surface height anomalies, is acquired. A Gaussian kernel function is selected as the radial basis function, with a kernel radius of 50 kilometers, and interpolation calculations are performed on 1°×1° grid cells. Taking a grid point (20°N, 115°E) as an example, sea surface height anomalies from its 30 neighboring points are extracted, and the weight of each point is calculated: weight = exp(-distance² / kernel radius²), where distance is in kilometers. For a neighboring point (20.1°N, 115.2°E), 25 kilometers from the target grid point, the weight = exp(-625 / 2500) ≈ 0.78. The weighted average of all neighboring point height anomalies is then calculated, combined with the geostrophic flow formula. (g-force acceleration 9.8 m / s²) 2 f Coriolis parameter 5×10 - 5 s -1 The density of seawater is 1025 kg / m³. 3 , The surface velocity at a given grid point is calculated using the height anomaly gradient. In this example, a height anomaly gradient of 0.1 m / 100 km yields a velocity of approximately 0.4 m / s. After traversing all grid points, a regularized surface velocity field covering the study area is output.

[0019] S103: Based on the vertically corrected velocity field and coastal radar horizontal flow data, a three-dimensional tidal current benchmark dataset is generated by fusing spatiotemporal outliers using Bayesian weighted Kalman filtering and removing outliers.

[0020] Based on the vertically corrected current velocity field and the horizontal current data from coastal radar, a Bayesian weighted Kalman filter was used to fuse and remove spatiotemporal outliers. The implementation process used a 24-hour dataset from the Pearl River Estuary as an example, with a time step of 1 hour. At a certain time t, the vertically corrected current velocity at the target point (22.5°N, 113.8°E) was 0.6 m / s (flow direction 80°), while the coastal radar data showed 0.7 m / s (flow direction 75°). The differences between the two data types were calculated: velocity difference 0.1 m / s, flow direction difference 5°. A threshold of 0.15 m / s for velocity difference and 10° for flow direction difference were set, and the data at this point did not exceed these limits. Initial weights were assigned based on the reliability of the data sources: satellite-buoy fusion data weighted at 0.6, and radar data weighted at 0.4. The Bayesian update formula was: posterior weight = prior weight × likelihood probability. It was assumed that the historical error standard deviation of the vertically corrected data was 0.05 m / s, and that of the radar data was 0.08 m / s, with the likelihood probability calculated according to a normal distribution. In the example, the likelihood value of the vertical correction data is exp(-0.1). 2 / (2×0.05 2 ))≈0.135, Radar data likelihood value = exp(-0.1 2 / (2×0.08 2 The updated weights are: vertical data weight = 0.6 × 0.135 ≈ 0.081, radar data weight = 0.4 × 0.457 ≈ 0.183, and normalized weights are 0.31 and 0.69 respectively. The weighted fused velocity is calculated as: 0.6 × 0.31 + 0.7 × 0.69 ≈ 0.66 m / s. Data points exceeding the threshold (e.g., velocity difference > 0.2 m / s) are directly discarded and interpolated. After spatiotemporal window sliding processing, a continuous three-dimensional velocity field is output.

[0021] See Figure 3 The generation of the tidal constituent ellipse element matrix based on S2 includes the following steps: S201: Based on the three-dimensional tidal reference dataset, the fundamental frequency components are extracted using the fast Fourier transform preprocessing algorithm to generate the tidal frequency characteristic spectrum. The 3D tidal baseline dataset contains hourly flow velocity data for 30 days from stations in the Bohai Bay, with a time series length of 720 hours. Fast Fourier Transform (FFT) preprocessing was used, with a frequency resolution of 0.001 cph (per week / hour). The amplitude of each frequency component was calculated. In the example, the amplitude corresponding to a frequency of 0.08 cph is 0.25 m / s, and the amplitude of a frequency of 0.04 cph is 0.18 m / s. The dominant frequency peaks were identified: the theoretical frequency of the M2 tidal constituent is 0.0805 cph, and the actual detected peak is 0.0803 cph with an amplitude of 0.48 m / s; the theoretical frequency of the K1 tidal constituent is 0.0418 cph, and the actual peak is 0.0417 cph with an amplitude of 0.22 m / s. The amplitude values ​​of characteristic tidal constituent frequencies such as O1 (0.0387 cph) and S2 (0.0833 cph) were extracted, and a spectral feature table was generated to record each tidal constituent frequency and its corresponding amplitude.

[0022] S202: Based on the tidal frequency characteristic spectrum, the constrained least squares fitting method is used to solve for the initial tidal values ​​of O1 / K1 / M2 / S2 / M4 / MS6, and the initial solution of the tidal parameters is generated. The tidal frequency characteristic spectrum provides the initial amplitude and phase of the O1 / K1 / M2 / S2 / M4 / MS6 tidal constituents. Constraints are set as follows: the amplitude of the M2 tidal constituent is 0.4-0.6 m / s, and the phase lag is 30°-60°. The least squares fitting objective function is minΣ[V o bs(t)-Σ(A i cosω i t+B i sinω i t)] 2 In the example, the observed flow velocity at t=0 is 0.52 m / s. The initial solution of the model values ​​is: M2 component A. M2 =0.45, B M2 =0.15, calculated value 0.45×cos0+0.15×sin0=0.45, residual 0.07, adjusted parameter step size 0.01, after 5 iterations A M1 =0.48, B M2 =0.18, residuals decreased to 0.02, synchronous fitting K1 tide A K1 =0.12, B K1 =0.08, finally obtaining the initial set of harmonic constants for each tidal constituent.

[0023] S203: Based on the initial solution of the tidal parameters, the major semi-axis, minor semi-axis, inclination angle and rotation direction are calculated by the complex ellipse vector iteration algorithm to generate the tidal ellipse element matrix.

[0024] The initial solution of the tidal parameters provides the harmonic constant A. i B i Calculation of the vector of a complex ellipse: semi-major axis = √(A) 2 +B 2), minor semi-axis = √(C 2 +D 2 ), where C and D are derived from the elliptic rotation matrix, in the example M2 tidal division A = 0.48, B = 0.18, Inclination angle θ = 0.5arctan(2AB / (A)) 2 -B 2 Substituting this into the equation, we get θ = 0.5arctan(2 × 0.48 × 0.18 / (0.48)). 2 -0.18 2 ))≈14°, rotation direction determination: B / A>0 is counterclockwise, in this example 0.18 / 0.48>0, marked CCW, iterative optimization: starting from the initial dip angle, adjust the ±2° range to search for the minimum velocity residual. In the example, the residual is reduced by 5% at 12°. Finally, the elliptical element matrix is ​​output to record the semi-major axis, semi-minor axis, dip angle and rotation direction parameters of each tide.

[0025] See Figure 4 The generation of shallow water corrected tidal parameter sets based on S3 includes the following steps: S301: Based on the M2 and M4 tidal constituents in the tidal constituent ellipse element matrix, the shallow water intensity factor is calculated using a nonlinear ratio model: γ=(M4 amp / M2 amp ) 2 , generating a scalar field of shallow water intensity factor; The tidal constituent ellipse element matrix provides the M2 and M4 tidal constituent amplitude parameters. Taking the grid point (31.5°N, 122°E) in the Yangtze River Estuary area as an example, the M2 tidal constituent amplitude is extracted to be 0.8 m / s, and the M4 tidal constituent amplitude is 0.25 m / s. The shallow water intensity factor γ = (M4) is calculated using a nonlinear ratio model. a mp / M2 a mp) 2 =(0.25 / 0.8) 2 ≈0.098, traversing all grid points in the study area, such as a point in the Bohai Strait where M2 amplitude is 1.2 m / s and M4 amplitude is 0.18 m / s, γ=(0.18 / 1.2) 2 =0.0225, generating a spatial distribution field of shallow water intensity factor covering the entire area.

[0026] S302: Based on the shallow water intensity factor scalar field and seabed topographic gradient data, the Legendre multinomial regression algorithm is used to fit the topographic dependence parameters α and β to generate the topographic coupling correction coefficient. The shallow water intensity factor scalar field and seafloor topographic slope data were processed collaboratively. The topographic slope data was calculated using water depth grid difference. In the example, the slope outside the Pearl River Estuary was 0.8‰, the shallow water intensity factor was 0.15, and the Legendre polynomial regression adopted the second-order form P(x)=α+βx+γx. 2x is the slope value (in‰), the polynomial order is set to 2, and the least squares objective function is minΣ(γ) o bs-(α+β·slope+γ·slope) 2 )) 2 Input 100 sets of sample data, with a slope range of 0.1‰-2.0‰ and a shallow water factor range of 0.01-0.2. Solve for the coefficients: α=0.12, β=-0.05, γ=0.008. Validation example: When the slope is 1.2‰, the predicted γ=0.12-0.05×1.2+0.008×1.44≈0.065, the measured value is 0.068, and the residual is 0.003. Output the set of terrain coupling correction coefficients.

[0027] S303: Based on the terrain coupling correction coefficient and the tidal ellipse element matrix, through the amplitude correction equation: A corr =A ori ·(1+αγ+βγ 2 Update the tidal parameters and generate a shallow water corrected tidal parameter set.

[0028] The topographic coupling correction coefficient is linked to the tidal constituent ellipse element matrix. The amplitude correction equation is Acorr=Aori·(1+αγ+βγ2). In the example, the original amplitude of the M2 tidal constituent is 0.8m / s, the shallow water intensity factor γ=0.098, the topographic parameters α=0.12, β=-0.05, and the corrected amplitude is 0.8×(1+0.12×0.098-0.05×0.098). 2 The amplitude is approximately 0.8 × 1.011 ≈ 0.809 m / s. After correction, the amplitude of the K1 tidal constituent (original amplitude 0.3 m / s, γ = 0.03) is approximately 0.3 × (1 + 0.12 × 0.03 - 0.05 × 0.0009) ≈ 0.307 m / s. The parameters are updated by traversing all tidal constituents and grid points.

[0029] See Figure 5 The optimization combination of extreme event generation based on S4 includes the following steps: S401: Based on the historical typhoon event database and the shallow water corrected tidal parameter set, the Pearson correlation matrix algorithm is used to calculate the correlation between tidal constituents and extreme current velocities, and generate the tidal constituent event correlation matrix. The historical typhoon event database contains measured maximum flow velocity data during typhoons in the Pearl River Estuary over the past decade. The shallow-water corrected tidal parameter set provides the amplitude and phase of tidal constituents such as M2 and S2. Using the Typhoon Mangkhut event as an example, the Pearson correlation coefficient r = cov(X,Y) / (σ) is calculated between the M2 tidal constituent amplitude sequence and the measured maximum flow velocity sequence. X σ YIn the example, the M2 amplitude sequence is [0.82, 0.85, 0.91] m / s (24 hours before the typhoon), the measured velocity sequence is [1.2, 1.3, 1.5] m / s, the covariance is cov = 0.028, and the standard deviation is σ. X =0.045, σ Y =0.15, r = 0.028 / (0.045×0.15)≈0.72, set a significant correlation threshold |r|>0.5, determine that the M2 tidal constituent is strongly correlated with extreme flow velocity, simultaneously calculate the correlation coefficients of K1, S2 and other tidal constituents, and generate a tidal constituent-event correlation matrix to record the correlation strength between each tidal constituent and historical typhoons.

[0030] S402: Based on the correlation matrix of tidal events, a random forest feature importance assessment is used to screen high-contribution tidal periods and generate a ranking of key tidal contributions; The correlation matrix of tidal events is input into a random forest model. 100 decision trees are set, and 70% of the samples are randomly selected from each tree. The Gini coefficient descent method is used to evaluate feature importance. In the example, when the M2 tidal event splits at a single tree node, the Gini coefficient decreases from 0.22 to 0.18, a decrease of 0.04. The K1 tidal event decreases by 0.02. The average importance score for the M2 tidal event is 0.15, for K1 it is 0.08, and for S2 it is 0.12. A high contribution threshold >0.1 is set, and the three tidal events M2, S2, and M4 are selected and sorted in descending order of score to generate a contribution ranking table.

[0031] S403: Based on the ranking of key tidal contributions, generate extreme scenario tidal combinations through Monte Carlo vector synthesis simulation, and generate optimized combinations of extreme event tidal combinations.

[0032] The key tidal constituents M2, S2, and M4 were ranked by their contribution and identified as the core tidal constituents. A Monte Carlo simulation was conducted with 1000 random trials. In each trial: the amplitude of M2 was uniformly sampled within ±20% of the baseline value of 0.8 m / s, with the initial sample at 0.96 m / s. The phase was randomly offset by +8° within the range of [-15°, +15°]. The amplitude of S2 was sampled at 0.65 m / s (baseline 0.6 m / s) with a phase offset of -5°. The amplitude of M4 was sampled at 0.28 m / s (baseline 0.25 m / s). The composite velocity V was calculated using vector synthesis. syn =√[(M2cosθ) M 2+S2cosθ S 2+M4cosθ M 4) 2 +(M2sinθ M 2+S2sinθ S 2+M4sinθ M 4) 2], instance value √[(0.96cos8°+0.65cos(-5°)+0.28cos0°) 2 +(0.96sin8°+0.65sin(-5°)+0.28sin0°) 2 The value is approximately 1.72 m / s. After statistical analysis of 1000 experimental results, the parameter combination corresponding to the 95th percentile of 1.85 m / s is taken as the extreme optimal combination.

[0033] See Figure 6 The generation of coupled shear stress field based on S5 includes the following steps: S501: Based on the M2 / S2 tidal vector in the extreme event tidal optimization combination, a spectral boundary condition generation algorithm is used to construct the wave model input field and generate a set of wave spectral boundary conditions. The extreme event tidal combination optimization provides M2 / S2 tidal vector parameters: M2 tidal semi-axis 0.85 m / s, dip angle 15°; S2 tidal semi-axis 0.6 m / s, dip angle -10°. Spectral boundary conditions are generated by inputting the tidal ellipse rotational velocity component and the east-west component U of the M2 tidal element. M2 =0.85×cos(15°)≈0.82m / s, North-South component V M2 =0.85×sin(15°)≈0.22m / s, S2 tidal constituent U S2 =0.6×cos(-10°)≈0.59m / s, V S2 =0.6×sin(-10°)≈-0.10m / s, the combined surface velocity field U=U M2 +U S2 V = V M2 +V S2 For the example grid point (22°N, 114°E), U = 0.82 + 0.59 = 1.41 m / s, V = 0.22 - 0.10 = 0.12 m / s. According to the velocity modulus |V| = √(1.41 m / s),... 2 +0.12 2 The wind speed at a height of 10 meters is approximately 1.42 m / s. Using empirical formulas, the wind speed W = 3.2 × |V| ≈ 4.54 m / s. The directional spectrum is generated using the JONSWAP spectral formula, with the main wave direction taken as the velocity vector direction θ = arctan(V / U) ≈ 4.9°, and the significant wave height H. s =0.2W 2 / g(g=9.8m / s 2 )≈0.42m, period The output wave spectrum boundary condition set includes H for each grid point s T p Main wave direction parameters.

[0034] S502: Based on the wave spectrum boundary condition set, run the third-generation wave numerical model to solve the radiation stress tensor and generate the wave radiation stress field. A third-generation wave model was driven by a set of wave spectrum boundary conditions. The grid resolution for the Pearl River Estuary region was set to 500 meters, with a time step of 10 seconds. Input boundary conditions were: significant wave height 0.42 m, peak period 5.1 s, and main wave direction 4.9°. The model solved the wave energy balance equation and calculated the radiation stress component S. xx =ρg∫[n(cos 2 θ+1)-0.5]E(f,θ)dfdθ, where ρ=1025kg / m 3 g = 9.8 m / s 2 Where n is the phase velocity and the group velocity ratio is taken as 0.8, E is the wave energy density, and θ is the wave direction angle. In this example, the single-frequency component f = 0.2 Hz has an energy density of 0.5 m³ / s. 2 / Hz, θ=5°, S xx Component = 1025 × 9.8 × [(0.8 × (cos 2 [5°+1)-0.5)×0.5]≈3500N / m, S calculated simultaneously yy S xy The components output the radiation stress tensor field across the entire region.

[0035] S503: Solving the coupled shear stress equation based on the wave radiation stress field and the shallow water modified tidal parameter set: Generate a coupled shear stress field.

[0036] The wave radiation stress field provides the tensor components, and the shallow water corrected tidal parameter set provides the tidal current vector. ρ is taken as 1025 kg / m 3 Drag coefficient C d =0.0025, V tide Let V be the velocity vector of the M2 tributary current (0.82 m / s, 0.22 m / s). wave The velocity of the wave particles is taken as 10% of the surface flow velocity, resulting in (0.141 m / s, 0.012 m / s). The first term = 1025 × 0.0025 × 0.98 2 ≈2.46N / m 2 Radiation stress divergence Instance point S xx =3500N / m, adjacent point Δx = 500 meters, Synthetic τ 耦合 ≈2.46 + 0.4 = 2.86 N / m 2 The stress field is generated by traversing the entire mesh.

[0037] See Figure 7 The generation of multi-scenario particle transport distance vectors based on S6 includes the following steps: S601: Based on the coupled shear stress field, the velocity correction model is used to update the flow velocity vector of the split current. κ is an empirical coefficient used to generate a stress-corrected velocity field; Coupled shear stress field provides τ 耦合 The value, the empirical coefficient κ, is set according to water depth stratification: 0.05 for water depth <20 meters, 0.03 for 20-50 meters, and 0.01 for >50 meters. For example, at the Pearl River Estuary grid point (22.1°N, 113.9°E), with a water depth of 18 meters, κ = 0.05, τ 耦合 =2.8N / m 2 Original M2 splitting current velocity vector (0.82 m / s, 0.22 m / s), corrected velocity: U corr =0.82+0.05×2.8≈0.96m / s, V corr =0.22+0.05×2.8≈0.36m / s, generate stress-corrected velocity field by traversing the global mesh.

[0038] S602: Calculate particle displacement using a complex harmonic integral algorithm based on a stress-corrected velocity field. Generate the original set of particle trajectories; The stress-corrected velocity field includes the harmonic constants corrected for each tidal constituent. The particle displacement is calculated using complex harmonic integration: The time integration period is taken as T = 12.42 hours (M2 tidal period), and the time step Δt = 10 minutes. For a particle at its initial position (22°N, 114°E), the velocity at t = 0 is V0 = (0.96cos0° + 0.36sin0°) ≈ 0.96 m / s, and the displacement increment ΔD... x =V x ·Δt=0.96×600≈576 meters, t=1 hour, phase increases by 29°, V1=0.96cos29°+0.36sin29°≈0.84m / s, ΔD x =0.84 × 600 = 504 meters, cumulative displacement D over 12.42 hours x =ΣΔD x ≈12.3 km, simultaneously calculate the north and south components, and output the particle trajectory coordinate sequence.

[0039] S603: Based on the original set of particle trajectories, a spring tide / mid tide / neap tide scenario separation algorithm is used to classify and statistically analyze the transport distance, generating multi-scenario particle transport distance vectors.

[0040] The original particle trajectory set is classified according to astronomical tide type: spring tide is taken on the second day after the new moon and full moon, mid-tide is taken on the first day after the first and last quarter moons, and neap tide is taken on the quadrant points of the lunar phase. An example is the trajectory displacement of a particle in the Pearl River Estuary during spring tide. The daily displacement of a mid-tide is 9.8 km, and the daily displacement of a neap tide is 7.1 km. Statistically, for 100 mass points, the average displacement is 12.5 ± 0.8 km for spring tide, 9.6 ± 0.6 km for mid-tide, and 7.0 ± 0.5 km for neap tide. Output the vector sets of mass point movement distances under the three scenarios.

[0041] See Figure 8 The generation of engineering-grade power flow characteristic parameters based on S7 includes the following steps: S701: Based on the multi-scenario particle transport distance vector, the extreme value statistical analysis algorithm is used to calculate the quantile transport distance and generate the maximum possible transport distance of the particle. The multi-scenario particle transport distance vector includes displacement data for three types of tides: spring tide, mid-tide, and neap tide. Extreme value statistical analysis uses a generalized Pareto distribution for fitting, with a return period of 50 years. The location parameter is the minimum value of the displacement sequence, the scale parameter σ = 0.2 × (75th quantile - 25th quantile), and the shape parameter ξ = 0.1. For example, the spring tide displacement sequence of 100 particles in the Pearl River Estuary [12.1, 12.3, ..., 13.0] km has a 75th quantile of 12.7 km and a 25th quantile of 12.2 km. σ = 0.2 × (12.7 - 12.2) = 0.1, ξ = 0.1. The 99th quantile is calculated as: x p = Position parameter + σ / ξ × [(n(1-p) / N) ( -ξ)-1],p=0.99,n sample size 100,N threshold exceedance 20,x p =12.0 + 0.1 / 0.1 × [(100 × 0.01 / 20)] ( -0.1)-1]≈13.5 km, outputting the maximum possible migration distance for each grid point.

[0042] S702: Based on the shallow water corrected tidal parameter set, the maximum possible tidal velocity is solved by the vector superposition extreme value search algorithm to generate the engineering extreme value velocity vector. The shallow water corrected tidal parameter set provides M2 / S2 equal tidal ellipse elements. Vector superposition extreme value search: Set the angle step size to 1°, traverse the 0°-360° direction, and calculate the composite current velocity V=√[(ΣA i cos(θ-α i )) 2 +(ΣA i sin(θ-α i )) 2 A i The semi-major axis of tidal constituents, α iInclination angle, for example, grid point M2 has a major semi-axis of 0.85 m / s and an inclination angle of 15°, and grid point S2 has a major semi-axis of 0.6 m / s and an inclination angle of -10°. When θ = 30°, V = √[(0.85cos(30-15)+0.6cos(30+10)) 2 +(0.85sin(30-15)+0.6sin(30+10)) 2 ]≈√[(0.82+0.46) 2 +(0.22+0.39) 2 =1.28m / s, the maximum value of 1.42m / s (θ = 42°) was obtained from the full-angle search, and it was marked as the engineering extreme velocity vector.

[0043] S703: Based on the engineering extreme velocity vector, the three-dimensional power flow ellipse is reconstructed using the elliptic parameterization algorithm to generate a three-dimensional elliptic dynamic topology map; An elliptic parameterization algorithm was used to input the extreme velocity vector in engineering studies to reconstruct a three-dimensional tidal current ellipse: The surface M2 tidal constituent has a major semi-axis of 0.85 m / s, a minor semi-axis of 0.35 m / s, and a dip angle of 15°; the middle layer (20 meters) has a major semi-axis of 0.70 m / s, a minor semi-axis of 0.25 m / s, and a dip angle of 18°; the bottom layer (50 meters) has a major semi-axis of 0.45 m / s, a minor semi-axis of 0.15 m / s, and a dip angle of 22°. Ten depth layers were set, and cubic splines were used for interlayer interpolation. The surface ellipse equation at the instance point (22°N, 114°E) is: (x'cos15°+y'sin15°) / 0.85 2 +(-x'sin15°+y'cos15°) / 0.35 2 =1, outputs a three-dimensional elliptical space topological mesh.

[0044] S704: Integrates the maximum possible transport distance of mass points, engineering extreme velocity vectors, and three-dimensional elliptical dynamic topology diagrams, and generates a standardized output package through the CUDA parallel compression encapsulation algorithm to generate an engineering-grade power flow characteristic parameter package.

[0045] The maximum possible transport distance of the integrated mass point is 13.5 km, the engineering extreme velocity vector is 1.42 m / s, the three-dimensional elliptical topology mesh data is used, and CUDA parallel compression is performed: 1024 thread blocks are allocated, each thread processes 100 mesh points, differential coding compression algorithm is used, the velocity difference between adjacent meshes in the example is 0.02 m / s, the difference is stored instead of the absolute value, the elliptical parameters use 16-bit floating-point precision, the compression rate is improved by 60%, and it is encapsulated into binary data packets.

[0046] The above are all preferred embodiments of this application and are not intended to limit the scope of protection of this application. Therefore, all equivalent changes made in accordance with the structure, shape and principle of this application should be covered within the scope of protection of this application.

Claims

1. A data processing method based on ocean current harmonic analysis, characterized in that... Includes the following steps: S1: Based on surface flow velocity data from satellite altimeters, horizontal flow field data from coastal radar, and vertical profile data from buoys, a spatiotemporal registration algorithm is used to align the spatiotemporal reference, and Bayesian weighted Kalman filtering is used to dynamically fuse the three types of data sources to eliminate instrument system errors and generate a three-dimensional tidal current reference dataset. S2: Based on the three-dimensional tidal reference dataset, the least squares quasi-harmonic analysis algorithm is used to solve the elliptical elements of the six tidal constituents O1, K1, M2, S2, M4, and MS4. A complex vector superposition model is introduced to calculate the amplitude and lag angle of each tidal constituent and generate the elliptical element matrix of the tidal constituents. S3: Based on the M2, M4, and MS4 tidal parameters in the tidal ellipse element matrix, the shallow water effect coefficients are constructed using the Legendre polynomial fitting algorithm. δ=α·(M4 / M2) 2 +β·(MS4 / M2) 2 α and β are generated by inversion through the topographic gradient dependence model, which corrects the tidal amplitude and generates a shallow water corrected tidal parameter set; S4: Based on the shallow water corrected tidal parameter set and the historical typhoon event database, the random forest contribution evaluation algorithm is used to dynamically identify highly correlated tidal constituents. Monte Carlo vector synthesis is used to simulate and generate tidal constituent combinations under extreme scenarios, and to generate optimized combinations of tidal constituents for extreme events. S5: Based on the M2 / S2 semi-major axis vector of the extreme event tidal optimization combination, the boundary conditions of the third-generation wave model are input, and the radiation stress-tidal shear coupling equation is adopted: Solve for the wave disturbance intensity on the tidal field to generate a coupled shear stress field; S6: Based on the coupled shear stress field and the shallow water corrected tidal parameter set, the complex vector integral algorithm is used to calculate the trajectory of tidal water particles: Output the migration distance vectors under three scenarios: spring tide, mid tide, and neap tide, and generate multi-scenario particle migration distance vectors; S7: Based on multi-scenario particle transport distance vectors, a GPU-accelerated vector synthesis algorithm is used to perform parallel calculations on millions of grids for: the maximum possible flow velocity of the tidal current, the maximum possible transport distance of particles, and the 3D tidal current elliptical visualization map, generating an engineering-grade tidal current feature parameter package.

2. The data processing method based on ocean current harmonic analysis according to claim 1, characterized in that: The generation of a 3D power flow benchmark dataset based on S1 includes the following steps: S101: Based on the sea surface height anomaly data from the satellite altimeter, a radial basis function interpolation algorithm is used for spatial resampling to generate a gridded sea surface velocity field and generate a satellite surface velocity grid. S102: Based on satellite surface velocity grid and buoy ADCP vertical profile data, a dynamic optimal interpolation method is used for vertical correction to generate a vertically corrected velocity field. S103: Based on the vertically corrected velocity field and coastal radar horizontal flow data, a three-dimensional tidal current benchmark dataset is generated by using Bayesian weighted Kalman filtering to fuse spatiotemporal outliers and remove them.

3. The data processing method based on ocean current harmonic analysis according to claim 1, characterized in that: The generation of the tidal ellipse element matrix based on S2 includes the following steps: S201: Based on the three-dimensional tidal reference dataset, the fundamental frequency components are extracted using the fast Fourier transform preprocessing algorithm to generate the tidal frequency characteristic spectrum. S202: Based on the tidal frequency characteristic spectrum, the constrained least squares fitting method is used to solve for the initial tidal values ​​of O1 / K1 / M2 / S2 / M4 / MS6, and the initial solution of the tidal parameters is generated. S203: Based on the initial solution of the tidal parameters, the major semi-axis, minor semi-axis, inclination angle and rotation direction are calculated by the complex ellipse vector iteration algorithm to generate the tidal ellipse element matrix.

4. The data processing method based on ocean current harmonic analysis according to claim 1, characterized in that: The generation of shallow water corrected tidal parameter sets based on S3 includes the following steps: S301: Based on the M2 and M4 tidal constituents in the tidal constituent ellipse element matrix, the shallow water intensity factor is calculated using a nonlinear ratio model: γ=(M4 amp / M2 amp ) 2 , generating a scalar field of shallow water intensity factor; S302: Based on the shallow water intensity factor scalar field and seabed topographic gradient data, the Legendre multinomial regression algorithm is used to fit the topographic dependence parameters α and β to generate the topographic coupling correction coefficient. S303: Based on the terrain coupling correction coefficient and the tidal ellipse element matrix, through the amplitude correction equation: A corr =A ori ·(1+αγ+βγ 2 Update the tidal parameters and generate a shallow water corrected tidal parameter set.

5. The data processing method based on ocean current harmonic analysis according to claim 1, characterized in that: The S4-based extreme event tidal optimization combination includes the following steps: S401: Based on the historical typhoon event database and the shallow water corrected tidal parameter set, the Pearson correlation matrix algorithm is used to calculate the correlation between tidal constituents and extreme current velocities, and generate the tidal constituent event correlation matrix. S402: Based on the correlation matrix of tidal events, a random forest feature importance assessment is used to screen high-contribution tidal periods and generate a ranking of key tidal contributions; S403: Based on the ranking of key tidal contributions, generate extreme scenario tidal combinations through Monte Carlo vector synthesis simulation, and generate optimized combinations of extreme event tidal combinations.

6. The data processing method based on ocean current harmonic analysis according to claim 1, characterized in that: The generation of coupled shear stress field based on S5 includes the following steps: S501: Based on the M2 / S2 tidal vector in the extreme event tidal optimization combination, a spectral boundary condition generation algorithm is used to construct the wave model input field and generate a set of wave spectral boundary conditions. S502: Based on the wave spectrum boundary condition set, run the third-generation wave numerical model to solve the radiation stress tensor and generate the wave radiation stress field. S503: Solving the coupled shear stress equation based on the wave radiation stress field and the shallow water modified tidal parameter set: Generate a coupled shear stress field.

7. The data processing method based on ocean current harmonic analysis according to claim 1, characterized in that: The generation of multi-scenario particle transport distance vectors based on S6 includes the following steps: S601: Based on the coupled shear stress field, the velocity correction model is used to update the flow velocity vector of the split current. κ is an empirical coefficient used to generate a stress-corrected velocity field; S602: Calculate particle displacement using a complex harmonic integral algorithm based on a stress-corrected velocity field. Generate the original set of particle trajectories; S603: Based on the original set of particle trajectories, a spring tide / mid tide / neap tide scenario separation algorithm is used to classify and statistically analyze the transport distance, generating multi-scenario particle transport distance vectors.

8. The data processing method based on ocean current harmonic analysis according to claim 1, characterized in that: The generation of engineering-level power flow characteristic parameters based on S7 includes the following steps: S701: Based on the multi-scenario particle transport distance vector, the extreme value statistical analysis algorithm is used to calculate the quantile transport distance and generate the maximum possible transport distance of the particle. S702: Based on the shallow water corrected tidal parameter set, the maximum possible tidal velocity is solved by the vector superposition extreme value search algorithm to generate the engineering extreme value velocity vector. S703: Based on the engineering extreme velocity vector, the three-dimensional power flow ellipse is reconstructed using the elliptic parameterization algorithm to generate a three-dimensional elliptic dynamic topology map; S704: Integrates the maximum possible transport distance of mass points, engineering extreme velocity vectors, and three-dimensional elliptical dynamic topology diagrams, and generates a standardized output package through the CUDA parallel compression encapsulation algorithm to generate an engineering-grade power flow characteristic parameter package.

Images