High-resolution internal wave sound velocity field inversion method based on propagation time delay

Abstract

The invention discloses a high-resolution internal wave sound velocity field inversion method based on propagation time delay, relates to the technical field of ocean observation, and aims to solve the problem that fine changes of a sound velocity field under internal wave driving cannot be effectively captured in the prior art. And cumulative observation results at different times are equivalent to multiple groups of virtual nodes and virtual observation networks at the same moment, so that high-resolution inversion of an internal wave sound velocity field is realized, and fine changes of the sound velocity field driven by internal waves are effectively captured. The defect that satellite observation lacks deep information is overcome, and full-water-depth sound velocity field observation can be achieved.

Description

TECHNICAL FIELD

[0001] The application relates to the technical field of ocean observation, in particular to a high-resolution internal wave sound speed field inversion method based on propagation time delay. BACKGROUND

[0002] As a key submesoscale dynamic process in the ocean, the internal wave causes significant interference to underwater high-precision positioning, navigation and detection due to the spatial and temporal variation of the sound speed field. Accurate acquisition of the distribution characteristics of the internal wave is an important basis for realizing ocean environment monitoring and disaster warning. However, the existing ocean observation technology has obvious limitations in acquiring the fine structure of the sound speed field driven by the internal wave: the spatial resolution of the buoy observation is insufficient and the position is fixed, which is difficult to depict the spatial details of the internal wave; the anchor observation can obtain continuous data at a fixed point, but the cost is high, the layout is complex, and the coverage is limited; satellite remote sensing can only reflect the sea surface information and cannot obtain the sound speed variation in the water body. Although the ocean acoustic tomography technology as a non-intrusive observation method has gradually attracted attention, the existing sound speed inversion research is mostly focused on the equivalent sound speed profile or idealized layered model, and the high-resolution inversion of the internal wave sound speed field with spatial variation still faces challenges. Traditional inversion methods usually rely on dense observation networks or sufficient prior sound speed data (such as CTD profiles), and in the actual scene of the open sea and deep sea, the observation nodes are sparse and the prior information is scarce, which makes it difficult for these methods to effectively reconstruct the fine structure of the complex sound speed field caused by the internal wave.

[0003] Although the acoustic observation using the sea surface and seabed nodes has the advantages of flexible layout and relatively low cost, and the direct sound propagation time delay is sensitive to sound speed variation, under the condition of limited observation data and insufficient prior information, how to realize the high-resolution and high-precision inversion of the internal wave sound speed field is still a key technical problem that has not been overcome in the current field of ocean acoustics. The existing technology cannot effectively capture the fine variation of the sound speed field driven by the internal wave, which restricts the development of ocean environment perception and warning capability. SUMMARY

[0004] The purpose of the present application is to provide a high-resolution internal wave sound speed field inversion method based on propagation time delay to solve the problem that the existing technology cannot effectively capture the fine variation of the sound speed field driven by the internal wave.

[0005] The technical scheme adopted by the present application to solve the above technical problems is:

[0006] A high-resolution internal wave sound speed field inversion method based on propagation time delay, comprising the following steps:

[0007] Step 1: acoustic sensors are arranged on the sea surface and the sea bottom of the sea area to be observed, each acoustic sensor is taken as a node, i.e. a sea surface node and a sea bottom node, the sea surface node and the sea bottom node are located in the same plane, and the sea bottom node is taken as a transmitting end and the sea surface node is taken as a receiving end to form an observation network, and the observation network is provided with a CTD chain;

[0008] Step 2: each sea bottom node and each sea surface node is a sound ray, all sound rays form an observation network, and in internal wave motion, the observation network is observed as a whole according to a set interval time and obtains sound propagation time delay, and the cumulative observation times are M;

[0009] Step 3: all nodes are taken as observation nodes, M groups of equivalent virtual nodes and virtual observation sound rays are obtained through M times of measurement according to the relative position change between the internal wave motion and the observation nodes, and the equivalent virtual nodes and the virtual observation sound rays are taken as nodes and sound rays;

[0010] Step 4: a spatial resolution is set, and the sea area to be observed is discretized to obtain a two-dimensional grid;

[0011] Step 5: nodes beyond the sea area to be observed are deleted, and the sound rays corresponding to the nodes are also deleted, then the intercepts of the intersection points of each sound ray and the two-dimensional grid are obtained based on the sound velocity profile measured by the CTD chain to form a sound ray path matrix A, wherein the row vectors of the sound ray path matrix A are the intercept vectors of the grids passed by the sound rays, and the intercepts in the grids not passed by the sound rays are 0;

[0012] Step 6: an observation equation set of time delay-internal wave sound velocity field is established according to the sound ray path matrix A and the corresponding sound propagation time delay;

[0013] Step 7: the observation equation set is solved by using a least square method under the constraint of the sound velocity profile to obtain a slowness field of the observation sea area, and the reciprocal of the slowness field is taken to obtain an internal wave sound velocity field of the observation area.

[0014] Further, in step 1, the acoustic sensors arranged on the sea surface and the acoustic sensors arranged on the sea bottom are arranged in pairs, and the range of the acoustic sensors is less than 8km.

[0015] Further, in step 5, the intercepts of the intersection points of each sound ray and the two-dimensional grid are obtained by a sound ray tracing algorithm.

[0016] Further, the sound velocity profile is expressed as:

[0017]

[0018] wherein, temperature is represented by T, salinity is represented by S, and depth is represented by D.

[0019] Furthermore, the specific steps of step 6 are as follows:

[0020] First, by mapping the sound ray to the propagation time delay, the following equation is formed:

[0021] ,

[0022] in, Indicates the first The propagation delay of the observed sound ray, Represents the ray path matrix A of the ray path matrix. The first voice in Intercept within each grid cell For the first The slowness value within each grid, i.e., the reciprocal of the speed of sound. Indicates the number of grid cells in a two-dimensional grid;

[0023] Using the propagation time delay and ray path matrix A corresponding to all sound rays, the observation equations for the time delay-internal wave sound velocity field are obtained, expressed as:

[0024] ,

[0025] in, This represents the observed direct sound delay vector. Let represent the column vector of the slow-degree field to be determined.

[0026] Furthermore, the objective function for solving the observation equations using the least squares method in step 7 is expressed as:

[0027] ,

[0028] in, This is the regularization parameter.

[0029] Furthermore, the slow-degree field of the observed sea area is represented as follows:

[0030] ,

[0031] Where I represents the identity matrix, This represents the slowness value within each grid cell.

[0032] Furthermore, the number of grid cells in the two-dimensional grid The following conditions must be met:

[0033] ,

[0034] in, and These represent the number of sound rays measured in a single instance and the cumulative number of observations, respectively.

[0035] In a two-dimensional mesh, the mesh cells are uniformly distributed horizontally but non-uniformly distributed vertically.

[0036] Furthermore, the horizontal position of the equivalent virtual node is represented as follows:

[0037] ,

[0038] in, Indicates in The equivalent virtual node horizontal position at time. Indicates in Horizontal position of the time node This represents the horizontal velocity of the internal wave at the observation section, with the direction being the positive direction of the horizontal coordinate axis.

[0039] Furthermore, the specific steps in step 5 for obtaining the intercept of each sound ray intersecting the two-dimensional grid based on the sound velocity profile measured by the CTD chain are as follows:

[0040] Using the sound velocity profile measured by the CTD chain as the estimated profile, the sound ray under the sound velocity profile is calculated by the fourth-order Runge-Kutta method and the bisection method. The sound ray intercept in the two-dimensional grid is obtained by calculating the intersection of the sound ray and the two-dimensional grid.

[0041] The beneficial effects of this invention are:

[0042] This application introduces the physical characteristics of the conformal nature of the sound velocity field under internal wave motion, equating the accumulated observation results at different times to multiple sets of virtual nodes and a virtual observation network at the same moment, thereby achieving high-resolution inversion of the internal wave sound velocity field and effectively capturing the fine changes in the sound velocity field driven by internal waves. This overcomes the deficiency of satellite observations in lacking deep information and enables full-depth sound velocity field observation.

[0043] This application has the ability to invert internal wave velocity fields in a sparse node networking mode on the sea surface and seabed. It does not require dense equipment deployment, which reduces the observation cost and operation difficulty. It is suitable for areas with limited observation conditions such as the open sea and deep sea.

[0044] This application, through multi-time-point cumulative observation processing combined with regularization constraints of a single measured profile, has strong anti-noise interference capabilities and still possesses good inversion capabilities even with time delay errors of less than 100μs and node position errors of less than 25cm. Attached Figure Description

[0045] Figure 1 This is the overall flowchart of this application;

[0046] Figure 2 A schematic diagram showing the distribution of seabed acoustic observation nodes and grid division.

[0047] Figure 3 Schematic diagram of the equivalent virtual observation node and observation network under internal wave motion Figure 1 ;

[0048] Figure 4 Schematic diagram of the equivalent virtual observation node and observation network under internal wave motion Figure 2 ;

[0049] Figure 5 A schematic diagram illustrating the principle of the observation equation for observing the sound ray path and corresponding time delay;

[0050] Figure 6 A schematic diagram of the true sound velocity field of high-resolution internal waves;

[0051] Figure 7 A schematic diagram for calculating the sound velocity field for high-resolution internal waves;

[0052] Figure 8 A schematic diagram of the error in the high-resolution internal wave sound velocity field;

[0053] Figure 9 A schematic diagram of the true sound velocity field of high-resolution internal waves under Gaussian noise.

[0054] Figure 10 This is a schematic diagram for calculating the sound velocity field of high-resolution internal waves in the presence of Gaussian noise. Detailed Implementation

[0055] It should be noted that, where there is no conflict, the various embodiments disclosed in this application can be combined with each other.

[0056] Specific Implementation Method 1: This implementation method describes a high-resolution internal wave sound velocity field inversion method based on propagation time delay, including:

[0057] Step 1: Deploy multiple pairs of acoustic sensor nodes and one CTD chain on the sea surface and seabed of the area to be observed. The surface and seabed nodes can communicate with each other. Multiple pairs of surface and seabed acoustic nodes are located on the same plane, with a range of less than 8 km. The seabed nodes act as transmitters, and the surface nodes as receivers, forming an observation network. The sound speed is calculated from temperature, salinity, and depth data measured via a single CTD chain. The specific calculation formula is as follows:

[0058]

[0059] in the formula Temperature (degrees Celsius). Salinity (parts per thousand). Depth (meters);

[0060] Step 2: Based on node locations, multiple direct acoustic measurement paths from the sea surface to the seabed are formed. The number of acoustic ray paths in a single measurement network is N, where... , and This represents the number of receiving nodes and transmitting nodes. During the internal wave motion, the observation network performs an overall observation at regular intervals and obtains the sound propagation time delay. The cumulative number of observations is M.

[0061] Since the number of observed sound rays N generated in a single observation is usually much smaller than the number of internal wave sound velocity field parameters under high spatial resolution conditions, the resulting set of observation equations is rank deficient and cannot correctly invert the sound velocity field parameters. This is an ill-posed problem in the inversion of the internal wave sound velocity field and a key issue affecting the error of the inversion results.

[0062] To address the ill-posed problem of internal wave sound velocity field inversion mentioned above, this invention utilizes the shape-preserving characteristic of the sound velocity field in internal wave motion, namely... By using a cumulative processing method based on multiple observations, it is possible to generate [data] within M observations. By observing multiple acoustic lines and increasing the observation density, a high-resolution internal wave velocity field can be derived without requiring extensive prior hydrological information as a constraint.

[0063] Step 3: Based on the characteristics of internal wave motion and the changing relative positions of the internal wave and the observation node, multiple sets of equivalent virtual nodes and virtual observation sound lines are formed through repeated measurements. It is assumed that the horizontal velocity of the internal wave at the observation cross-section is... The direction is the positive direction of the horizontal coordinate axis. The horizontal position of the time node is , then Equivalent virtual node horizontal position at time The virtual observation sound ray translation method is the same;

[0064] Step 4: To determine the required spatial resolution for the inversion results, the observed sea area is discretized into a two-dimensional grid consisting of several rows and columns in both the horizontal and vertical directions. Each two-dimensional grid comprises multiple cells. The grid is horizontally uniformly divided, but vertically non-uniformly divided, with a denser grid on the surface and a sparser grid on the bottom. The grid cells are then numbered. The number of grid cells in each grid is... ,in , and The number of rows and columns represents the number of grid cells. To ensure the accuracy of subsequent algorithm calculations, the following conditions must be met. ,in and This refers to the number of sound rays measured in a single instance and the cumulative number of observations in step two;

[0065] Step 5: Based on the sound velocity profile measured by the CTD chain, use the ray tracing algorithm to calculate the intercept of each sound ray with the grid in multiple observation paths, forming a sound ray path matrix A. The row vectors of the path matrix are the intercept vectors of the grids through which the sound ray passes, and the grid number position of the row vector is the sound ray intercept in the corresponding grid. The intercept is 0 in the grids through which the sound ray has not passed.

[0066] Specifically, step five includes:

[0067] For all nodes, including equivalent virtual nodes, nodes outside the observation area and their corresponding ray propagation delays are removed to obtain all observed ray equations and propagation delays involved in the calculation. The measured sound velocity profile is used as the estimated profile, and the fourth-order Runge-Kutta method and the bisection method are used to calculate the observed ray paths under the sound velocity profile. The ray intercept within the grid is calculated by calculating the intersection of the ray and the grid.

[0068] Step Six: Based on the path matrix and the corresponding propagation delay, establish the observation equations for the time delay-internal wave velocity field; after mapping the sound ray path and propagation delay, the following equations are formed:

[0069] ;

[0070] in Indicates the first The propagation delay of the observed sound ray, Represents the first in matrix A The first voice in Intercept within each grid cell For the first The slowness value within each grid, i.e., the reciprocal of the speed of sound. Combining all observed sound ray paths and propagation time delay equations to form a set of observation equations, the time delay-internal wave speed of sound field observation equation set is as follows:

[0071] ;

[0072] in To observe the direct sound time delay vector, It is the column vector of the slow-degree field to be determined.

[0073] Step 7: Under the constraint of the measured sound velocity profile, solve for the slowness field of the observation area using the least squares method, and obtain the internal wave sound velocity field of the observation area by taking the reciprocal. Applying regularization constraints under the known profile, the objective function can be expressed as:

[0074] ;

[0075] in This is the sound velocity profile obtained from CTD measurements. Here is the regularization parameter. The slowness value within each grid is solved using the least squares method, and the solution matrix expression is:

[0076] ;

[0077] Where I is the identity matrix, the internal wave velocity field is obtained by taking the reciprocal of the slowness. This invention is based on the direct sound propagation delay information between sparse seabed nodes and utilizes the characteristics of ocean internal wave propagation. Through multiple measurement accumulation and processing, the ill-posedness of the number of parameters in high-resolution velocity field inversion, which is much greater than the number of observation information, is greatly reduced or even resolved. This significantly reduces the dependence on prior hydrological information, while also ensuring the robustness of the inversion method under the regularization constraint of a single measurement profile.

[0078] Example

[0079] The simulation environment, i.e. the true data, was based on the results of the LLC4320 South China Sea high-resolution numerical model. The simulation survey area was selected as 20.01°N, 118.24°E to 118.30°E, with a horizontal distance of 8km and a depth of 2700m. The total time length was 4 weeks, totaling 672 hours, with a simulation data time interval of 1 hour. The horizontal motion velocity of the internal waves was approximately 3m / s.

[0080] Step 1) Deploy 5 pairs of acoustic sensor nodes and 1 CTD chain on the sea surface and seabed of the simulation area. The distance between the nodes is 1.9km and the starting and ending range is 8km. The seabed nodes are the transmitting end and the sea surface nodes are the receiving end, forming an observation network. Select a sound velocity profile in the observation area as the real profile obtained by CTD measurement.

[0081] Step 2) Based on the node location, the number of sound ray paths in a single measurement network is 25. The observation network performs an overall observation every certain period of time. In the true value environment, the sound ray paths between the sea surface and seabed nodes are traced using the fourth-order Runge-Kutta method and the bisection method. The sound propagation delay of the corresponding sound ray is calculated as the propagation delay observation value in the simulation. The cumulative number of observations for each inversion in the simulation is 15.

[0082] Step 3) Based on the characteristics of internal wave motion, and according to the changes in the relative positions of the internal wave and the observation node, multiple sets of equivalent virtual nodes and virtual observation sound lines are formed, as shown in the schematic diagram. Figure 3 ;

[0083] Step 4) Discretize the sea area to be observed. Divide it into 4 grids on average in the horizontal direction with a grid spacing of 2km. Divide it into 74 uneven layers in the vertical direction. The surface sound velocity changes more frequently, so the grid is denser in the surface layer and sparser in the bottom layer. The minimum grid spacing is 1.5m. Therefore, the number of grid units in the simulation area is 74×4, which means there are 296 unknown parameters in the sound velocity field inversion. These parameters are numbered from 1 to 296.

[0084] Step 5) For all nodes, including equivalent virtual nodes, remove nodes that are outside the observation area and their corresponding ray propagation delays; Based on one known sound speed profile from Step 1), calculate the intercepts between the ray and the grid in multiple observations to form a ray path matrix A, where the row vectors of the path matrix are the intercept vectors of the grids through which the ray passes, and the grid number position of the row vector is the ray intercept within the corresponding grid. The intercept is 0 in grids where the ray has not passed.

[0085] Step 6) Based on the path matrix and the corresponding propagation delay, establish the observation equations for the time delay-internal wave velocity field; after mapping the sound ray path and propagation delay, the following equations are formed:

[0086] ;

[0087] in Indicates the first The propagation delay of the observed sound ray, Represents the first in matrix A The first voice in Intercept within each grid cell For the first The slowness value within each grid, i.e., the reciprocal of the speed of sound. Combining all observed sound ray paths and propagation time delay equations to form a set of observation equations, the time delay-internal wave speed of sound field observation equation set is as follows:

[0088] ;

[0089] in To observe the direct sound time delay vector, It is the column vector of the slow-degree field to be determined.

[0090] Step 7) Under the constraint of the measured sound velocity profile, solve for the slowness field of the observation area using the least squares method, and obtain the internal wave sound velocity field of the observation area by taking the reciprocal. Under the known profile, regularization constraints are applied, and the objective function can be expressed as:

[0091] ;

[0092] in This is the sound velocity profile obtained from CTD measurements. The regularization parameter is set to 5 in the simulation. The slowness value within each grid is solved using the least squares method, and the solution matrix expression is as follows:

[0093] ;

[0094] Where I is the identity matrix, the internal wave velocity field is obtained by taking the reciprocal of the slowness. Under noise-free conditions, the internal wave velocity field inversion error is less than 3.6 m / s. By increasing the number of cumulative observations, higher spatial resolution inversion results can be obtained. This application utilizes the propagation characteristics of ocean internal waves and employs a method of cumulative processing through multiple measurements to significantly reduce or even eliminate the ill-posedness of the number of parameters in high-resolution velocity field inversion, which is much greater than the number of observed information. This greatly reduces the dependence on prior hydrological information, while also ensuring the robustness of the inversion method under the regularization constraint of a single measurement profile.

[0095] It should be noted that the specific embodiments are merely explanations and illustrations of the technical solution of the present invention and should not be used to limit the scope of protection. Any modifications made in accordance with the claims and specification of the present invention that are only partial should still fall within the protection scope of the present invention.

Claims

1. A high-resolution internal wave velocity field inversion method based on propagation time delay, characterized in that... Includes the following steps: Step 1: Deploy acoustic sensors on the sea surface and seabed of the area to be observed, and treat each acoustic sensor as a node, namely a sea surface node and a seabed node. The sea surface node and the seabed node are located on the same plane, and the seabed node is the transmitter and the sea surface node is the receiver to form an observation network. The observation network is equipped with a CTD chain. Step 2: Each seabed node and each sea surface node constitute a sound ray. All sound rays form an observation network. During the internal wave motion, the observation network conducts overall observations according to the set interval time and obtains the sound propagation time delay. The cumulative number of observations is M. Step 3: Treat all nodes as observation nodes. Based on the change in the relative position between the internal wave motion and the observation nodes, obtain M sets of equivalent virtual nodes and virtual observation sound lines through M measurements, and use the equivalent virtual nodes and virtual observation sound lines as nodes and sound lines. Step 4: Set the spatial resolution and discretize the sea area to be observed to obtain a two-dimensional grid; Step 5: Delete the nodes that are outside the sea area to be observed and delete the corresponding sound ray. Then, based on the sound speed profile measured by the CTD chain, obtain the intercept of each sound ray and the two-dimensional grid, thus forming the sound ray path matrix A. The row vector of the sound ray path matrix A is the intercept vector of the grid through which the sound ray passes. The intercept in the grid through which the sound ray has not passed is 0. Step 6: Based on the sound ray path matrix A and its corresponding sound propagation time delay, establish the observation equations for the time delay-internal wave sound velocity field; Step 7: Under the constraint of the sound velocity profile, solve the observation equations using the least squares method to obtain the slow velocity field of the observation sea area, and take the reciprocal to obtain the internal wave sound velocity field of the observation area.

2. The high-resolution internal wave velocity field inversion method based on propagation time delay according to claim 1, characterized in that... In step 1, the acoustic sensors deployed on the sea surface and the acoustic sensors deployed on the seabed are arranged in pairs, and the starting and ending range of the acoustic sensors is less than 8 km.

3. The high-resolution internal wave velocity field inversion method based on propagation time delay according to claim 1, characterized in that... In step 5, the intercept of each sound ray at the intersection point with the two-dimensional grid is determined by a sound ray tracking algorithm.

4. The high-resolution internal wave velocity field inversion method based on propagation time delay according to claim 1, characterized in that... The sound velocity profile is represented as follows: in, Indicates temperature. Indicates salinity. Indicates depth.

5. The high-resolution internal wave velocity field inversion method based on propagation time delay according to claim 4, characterized in that... The specific steps of step 6 are as follows: First, by mapping the sound ray to the propagation time delay, the following equation is formed: ， in, Indicates the first The propagation delay of the observed sound ray, Represents the ray path matrix A of the ray path matrix. The first voice in Intercept within each grid cell For the first The slowness value within each grid, i.e., the reciprocal of the speed of sound. Indicates the number of grid cells in a two-dimensional grid; Using the propagation time delay and ray path matrix A corresponding to all sound rays, the observation equations for the time delay-internal wave sound velocity field are obtained, expressed as: ， in, This represents the observed direct sound delay vector. Let represent the column vector of the slow-degree field to be determined.

6. The high-resolution internal wave velocity field inversion method based on propagation time delay according to claim 5, characterized in that... The objective function for solving the observation equations using the least squares method in step 7 is expressed as: ， in, This is the regularization parameter.

7. The high-resolution internal wave velocity field inversion method based on propagation time delay according to claim 6, characterized in that... The slow field of the observed sea area is represented as follows: ， Where I represents the identity matrix, This represents the slowness value within each grid cell.

8. The high-resolution internal wave velocity field inversion method based on propagation time delay according to claim 1, characterized in that... The number of grid cells in the two-dimensional grid The following conditions must be met: ， in, and These represent the number of sound rays measured in a single instance and the cumulative number of observations, respectively. In a two-dimensional mesh, the mesh cells are uniformly distributed horizontally but non-uniformly distributed vertically.

9. The high-resolution internal wave velocity field inversion method based on propagation time delay according to claim 1, characterized in that... The horizontal position of the equivalent virtual node is represented as follows: ， in, Indicates in The equivalent virtual node horizontal position at time. Indicates in Horizontal position of the time node This represents the horizontal velocity of the internal wave at the observation section, with the direction being the positive direction of the horizontal coordinate axis.

10. The high-resolution internal wave velocity field inversion method based on propagation time delay according to claim 1, characterized in that... The specific steps for obtaining the intercept of each sound ray at the intersection point of the two-dimensional grid based on the sound velocity profile measured by the CTD chain in step 5 are as follows: Using the sound velocity profile measured by the CTD chain as the estimated profile, the sound ray under the sound velocity profile is calculated by the fourth-order Runge-Kutta method and the bisection method. The sound ray intercept in the two-dimensional grid is obtained by calculating the intersection of the sound ray and the two-dimensional grid.

Images