A real-time correction observation method for acoustic tomography of small-scale horizontal flow fields

Abstract

The present invention discloses a small-scale horizontal flow field acoustic tomography real-time correction observation method, including: arranging an acoustic transceiver system in the observation water area, using the acoustic transceiver system to collect original hydroacoustic data in the observation water area; cross-correlating the collected original hydroacoustic data to extract the direct path transmission time; performing transmission time preprocessing to eliminate abnormal values; combining the prior hydrological information and the mean of the two-way propagation time measured by the temperature and depth instrument to perform real-time station correction; performing grid division based on the station layout, calculating the length of the sound line in each grid, and constructing a coefficient matrix; setting an error threshold, continuously iterating and solving the horizontal profile flow field through the cone least squares method, and performing visual display of the horizontal flow field. The present invention can improve the accuracy of flow field measurement in small-scale water areas.

Description

Technical Field

[0001] The invention belongs to the technical field of hydrological environment monitoring, and in particular relates to a small-scale horizontal flow field acoustic tomography real-time correction observation method. Background Art

[0002] Flow field observation in small-scale waters, such as marine ranches, offshore currents, and estuaries, is closely related to the marine environment. It has very important scientific significance for the study of marine physics, chemistry, and ecology in small-scale waters, and has engineering guidance significance for marine fisheries. Therefore, the research on inventing a method for long-term, effective, and high-precision flow field observation has attracted widespread attention from scholars at home and abroad.

[0003] For example, the Chinese patent document with publication number CN109584314A discloses a method, device and electronic device for measuring the flow field on the surface of water. The position of the measuring camera is adjusted according to the position of the auxiliary calibration camera, and the measuring camera is used to collect the video of the local ripple details on the water surface; the coordinates of the large-scene video of the water surface and the video of the local ripple details on the water surface are mapped to obtain the relative coordinates of each frame of the video image of the local ripple details on the water surface in the large-scene video image of the water surface; the video of the local ripple details on the water surface is subjected to wave tracking processing to obtain the motion vector information of the wave, and then calibrated according to the three-dimensional calibration data to obtain the flow velocity vector data. However, this method can only measure the surface flow field, and cannot accurately measure the underwater flow field.

[0004] The Chinese patent document with the publication number CN 113466872 A discloses a small-scale layered horizontal two-dimensional flow field observation method, which cross-correlates the obtained raw data to obtain the sound line propagation time; high-precision sound line simulation, obtains the sound line mode, reference propagation time and time window information; multipath resolution and extraction, calculates the arrival peak propagation time and sound line length of different paths; propagation time preprocessing; vertical stratification based on sound line distribution; calculates the sound line length and propagation time of each layer, constructs a coefficient matrix, constructs a vertical layered flow field, obtains the path average flow velocity and inversion error of each layer; sets a threshold, performs iterative calculation until the inversion error meets the requirements; uses the stream function method, takes the path average flow velocity of each layer as input data, and inverts the multi-layer horizontal two-dimensional flow field; and finally performs visualization processing. Using the present invention, the inversion of multi-layer horizontal flow fields can be performed, but the error caused by station drift cannot be eliminated.

[0005] At present, there is a lack of long-term regional flow field observations in small-scale sea areas such as marine ranches, offshore currents, and estuaries. The existing flow field observation methods combined with underwater acoustic systems are relatively lacking in research on the treatment of station drift. Summary of the invention

[0006] The present invention provides a small-scale horizontal flow field acoustic tomography real-time correction observation method, which can perform real-time station correction and improve the accuracy and density of small-scale water area hydrological information acoustic tomography measurement.

[0007] A small-scale horizontal flow field acoustic tomography real-time correction observation method, comprising:

[0008] (1) Arranging an acoustic transceiver system in the observation waters, wherein the acoustic transceiver system uses four acoustic transducers to perform mutual transmission of acoustic signals;

[0009] The acoustic transceiver system is used to collect the original hydroacoustic data in the observed water area, and the temperature and depth meter is used to collect the temperature information of the layer where the acoustic transducer is located;

[0010] (2) Performing acoustic station cross-correlation on the collected raw underwater acoustic data, and extracting the transmission time corresponding to the bidirectional arrival peak by setting the time window and signal-to-noise ratio threshold;

[0011] (3) Preprocess the extracted arrival peak transmission time, remove abnormal data, define the maximum difference of the return propagation time, and remove the misidentified peaks;

[0012] (4) Combine the extracted two-way transmission time mean and the temperature and depth information of the layer where the acoustic transducer is located to calculate the spacing of the acoustic transceiver system and perform station correction;

[0013] (5) In the two-dimensional horizontal profile, combined with the station layout, horizontal grid division is performed, the coordinates of the grid center point are calculated, the stream function is Taylor expanded, and the corresponding known coefficient matrix and the coefficient matrix to be solved are established;

[0014] (6) Perform two-dimensional horizontal grid flow field inversion by using the conical least squares method combined with the L-curve method, and constrain the error threshold; if the error exceeds the threshold, return to the flow field inversion until the requirements are met;

[0015] (7) Visualize the two-dimensional horizontal flow field obtained by inversion.

[0016] Furthermore, in step (1), the acoustic transceiver system includes a float, an acoustic transducer, a temperature-salinity-depth profiler and a data acquisition system;

[0017] The data acquisition system is located on the shore and connected to the acoustic transducers; four acoustic transducers are arranged at the same depth underwater, one end of which is connected to a buoy moored under the water surface and the other end is anchored at the bottom of the water.

[0018] In the acoustic transceiver system, each acoustic station uses the same-order M sequence with the same frequency but different modes, and adopts the same-transmit and same-receive mode for acoustic signals, ensuring that each pair of acoustic transducers receives a two-way acoustic signal after each signal is transmitted.

[0019] The specific process of step (2) is as follows:

[0020] The collected original underwater acoustic data are cross-correlated between two sound stations. By setting the time window and signal-to-noise ratio threshold, the transmission time of the first arrival peak, that is, the transmission time of the correlation peak corresponding to the direct path, is distinguished and extracted.

[0021] The specific process of step (4) is as follows:

[0022] (4-1) Based on the extracted two-way propagation time, calculate the mean two-way propagation time of each path during each signal transmission Among them, t + is the forward transmission time, t - is the back propagation time;

[0023] (4-2) Calculate the reference sound velocity C based on the temperature of the layer where the transducer is located measured by the temperature depth meter 0 , calculate the corrected distance between the two stations

[0024] (4-3) The station position is corrected by solving the intersection point of the station position line, and the formula is as follows:

[0025]

[0026]

[0027]

[0028]

[0029]

[0030]

[0031] Among them, (x n ,y n ) is the initial coordinate point of the station, (Δx n ,Δy n ) is the position drift of station n, (x n +Δx n ,y n +Δy n ) is the corrected station coordinate; since |Δx n -Δx m |<<|x n -x m |,|Δy n -Δy m |<<|y n -y m|, perform Taylor expansion on the above equation, and set station 1 close to the shore as a fixed point, let Δx 1 =Δy 1 =0 After simplifying, we get:

[0032]

[0033]

[0034]

[0035]

[0036]

[0037]

[0038] By solving the above equations, the station correction can be completed.

[0039] In step (5), when performing horizontal grid division, each acoustic station is projected to the horizontal section where the acoustic transducer depth is located, and r grids are divided. The sound line propagates in the grid. The specific parameters and formulas are as follows:

[0040] For each pair of sound stations, the direct path is:

[0041]

[0042] After Taylor expansion, we get

[0043]

[0044] Among them, l ij is the length of the i-th sound line passing through the j-th grid, C 0 is the sound velocity of the layer where the acoustic transducer is located, V j represents the projection of the flow velocity of the jth grid along the propagation path, represents the forward and reverse reference transmission time along the ith sound ray, Δt i is the two-way propagation time difference.

[0045] In step (5), the specific process of establishing the corresponding known coefficient matrix and the coefficient matrix to be solved is:

[0046] Introducing stream functions Among them, the unknown coefficient matrix Known coefficient matrix y) is the coordinate of the point to be measured; L x , L y is the length and width of the inversion area; N x , Ny are the number of grids in the inversion area along the x and y directions, respectively.

[0047] Calculate the north component of velocity based on the stream function and the east component of velocity

[0048] In step (6), the specific process of performing two-dimensional horizontal grid flow field inversion is as follows:

[0049] Substitute the north component u and the east component v of the velocity into get Written as a matrix equation y = Ex + n; where x = D m is the quantity to be solved, n is the solution error, y=Δt i is the measured time difference of the sound transmission between the two stations. is the transformation matrix;

[0050] The cone least squares method is used to solve the above matrix equation. By setting the penalty function J = n T n+a 2 x T x=(y-Ex) T (y-Ex)+α 2 x T x, minimize the penalty function to get the expected optimal solution of x Then the stream function expression is obtained, and the flow field in the observation area can be obtained by substituting the coordinates of the point to be solved.

[0051] Compared with the prior art, the present invention has the following beneficial effects:

[0052] 1. The present invention only requires a few observation stations to perform underwater acoustic tomography high-resolution imaging inversion, without the need to set up multiple fixed-point flow field monitoring equipment.

[0053] 2. The present invention adopts a mutual return transmission mode of acoustic signals, with a design of simultaneous transmission and reception of signals, to distinguish and extract multiple sound lines passing through the monitoring area in a small scale range, and obtain accurate acoustic signal transmission time.

[0054] 3. The present invention uses a real-time position correction algorithm to perform real-time correction of the position of the acoustic transceiver system, thereby improving the accuracy and reliability of flow field inversion.

[0055] 4. The present invention obtains two-dimensional flow field information through the stream function method, and can effectively improve the accuracy of small-scale flow field observation through the cyclic iteration method. BRIEF DESCRIPTION OF THE DRAWINGS

[0056] Figure 1 This is a flow chart of a small-scale horizontal flow field acoustic tomography real-time correction observation method of the present invention;

[0057] Figure 2 This is a schematic diagram of the arrangement of observation stations in a preferred embodiment of the present invention;

[0058] Figure 3 A schematic diagram of a station correction algorithm in a preferred embodiment of the present invention;

[0059] Figure 4 A comparison diagram of the standing distance before and after correction in a preferred embodiment of the present invention;

[0060] Figure 5 It is a result diagram of horizontal cross-section flow field between four stations in a preferred embodiment of the present invention;

[0061] Figure 6 This is a comparison diagram of the flow rate before and after correction and the ADCP flow rate in a preferred embodiment of the present invention. DETAILED DESCRIPTION

[0062] The present invention is further described in detail below in conjunction with the accompanying drawings and embodiments. It should be pointed out that the embodiments described below are intended to facilitate the understanding of the present invention and do not have any limiting effect on the present invention.

[0063] like Figure 1 As shown, a small-scale horizontal flow field acoustic tomography real-time correction observation method includes the following steps:

[0064] Step 1: Arrange the acoustic transceiver system in the observation waters and collect the original data of the observation waters.

[0065] In this embodiment, the acoustic transceiver system includes a buoy, an acoustic transducer, a temperature-salinity-depth profiler, a data acquisition system, etc. Four acoustic transducers S1, S2, S3, and S4 are fixed underwater in the water area to be observed to perform mutual transmission of acoustic signals. Among them, the acoustic signal adopts a simultaneous transmission and reception mode, and adopts the same-order M sequence with the same frequency but different modes to ensure that each pair of transducers can receive a two-way acoustic signal after each signal is transmitted.

[0066] Acoustic transducers S1, S2, S3, and S4 can receive or transmit acoustic signals with a center frequency of 50kHz. Figure 2 The bottom surface mooring method shown is fixed underwater. Specifically, the lower end of the acoustic transducer is connected to a weight through a rope and anchored to the bottom of the water, and the upper end is connected to a buoy through a rope. The transducer is connected to the data acquisition system through a cable. The data acquisition system, battery, etc. are placed in the buoy, and the cable is always kept in a relaxed state. The temperature of the layer where the transducer is located can be collected by CTD (temperature-salinity-depth instrument) or TD (temperature-depth instrument).

[0067] Step 2: Perform two-station cross-correlation on the collected original observation data, and distinguish and extract the transmission time of the first arriving signal, that is, the transmission time of the correlation peak corresponding to the direct path, by setting an appropriate time window and signal-to-noise ratio threshold.

[0068] Step 3: Preprocessing of propagation time to check system errors. First, correct the system errors and remove abnormal data. Define the maximum difference of the return propagation time, remove abnormal points, and obtain high-quality observation data.

[0069] Step 4: When performing position correction, calculate the reference sound velocity through the real-time temperature information obtained by TD and the mean of the two-way propagation time extracted in step 3 to perform position correction.

[0070] C(T,S,D)=1448.96+4.591T0.05304T 2 +2.734×10 4 T 3 +1.340(S35)+1.630×10 2 D+1.675×10 7 D 2 1.025×10 2 T(S35)7.139×10 -13 TD 3

[0071] C, T, D, and S are sound velocity, temperature, depth, and salinity, respectively.

[0072] The corrected station spacing

[0073] like Figure 3 As shown, station 1 bound to the shore is set as the reference point, then the distance between the two stations is:

[0074]

[0075]

[0076]

[0077]

[0078]

[0079]

[0080] Where (x n ,y n ) is the initial coordinate point of the station, (Δx n ,Δy n) is the position drift of station n, (x n +Δx n ,y n +Δy n ) is the corrected station coordinate. Since |Δx n -Δx m |<<|x n -x m |,|Δy n -Δy m |<<|y n -y m |, Taylor expand the above equation. Since station 1 is set as a fixed point, let Δx 1 =Δy 1 =0 After simplifying, we get:

[0081]

[0082]

[0083]

[0084]

[0085]

[0086]

[0087] Solving the above formula, we can get the corrected station position.

[0088] Step 5: Use the stream function method to perform horizontal grid division at each moment, and construct the known coefficient matrix and the coefficient matrix to be solved. The specific process is as follows:

[0089] Introducing stream functions Among them, the unknown coefficient matrix Known coefficient matrix

[0090] The north component of velocity and the east component of velocity Substitution Since Vdl=udx+vdy, we get

[0091]

[0092] The above equation can be written as a matrix equation y = Ex + n. Where x = {D m} is the quantity to be solved, n is the solution error, y={Δt i} is the measured time difference of the sound transmission between the two stations, is the transformation matrix.

[0093] Step 6: Apply the tapered least squares method to solve the above matrix equation, by setting the penalty function J = n T n+a 2 x T x=(y-Ex) T (y-Ex)+α 2 x T The values ​​of x,α are determined by bounding the expected error within a set threshold and are updated in real time during the experiment to track the dynamic environment.

[0094] Minimize the penalty function to get the expected optimal solution of x Then we can get the expressions of u and v, and substitute them into the coordinates of the points to be solved to find the flow field in the observation area.

[0095] Step 7: Finally, visualize the two-dimensional flow field.

[0096] In order to verify the effect of the present invention, an observation experiment was conducted on a certain water area of ​​the Huangcai Reservoir in Changsha using the present invention to obtain two-dimensional horizontal flow field information, specifically:

[0097] It includes the following devices: four sets of 50kHZ frequency acoustic transducers, four sets of acoustic transceiver systems, floats, buoy weights, TD, ADCP, etc.

[0098] Step 2 is used to distinguish and extract the transmission time of the first arriving signal between each pair of stations, that is, the transmission time of the correlation peak corresponding to the direct path.

[0099] The transmission time preprocessing in step 3 is used to remove outliers in the transmission time and obtain more accurate observation data.

[0100] Using the real-time station correction method in step 4, the reference sound velocity at the depth to be measured is calculated based on the hydrological information measured by the temperature and depth instrument. Figure 4 The actual station spacing shown is corrected.

[0101] After the grid is divided, the length of the sound line passing through each grid and the two-way propagation time difference of the sound line are calculated according to the reference sound speed. The corresponding program of this method is written by MATLAB, and the horizontal temperature information of the two-dimensional horizontal profile is inverted according to steps 5 and 6. Figure 5 Shown are the two-dimensional flow field results of the horizontal profile at 9 times in the observed water area. Figure 6 This is the comparison between the inverted velocity and the ADCP velocity at a certain point before and after correction.

[0102] The small-scale horizontal profile flow field with real-time station position correction more accurately reflects the state and changing trend of the underwater horizontal flow field during the observation period. Compared with single-point observation and surface observation methods, the observation area is wider and the accuracy is improved compared with that before correction, which shows the effectiveness and accuracy of this method in observing horizontal profile flow fields in small-scale waters.

[0103] The embodiments described above provide a detailed description of the technical solutions and beneficial effects of the present invention. It should be understood that the above are only specific embodiments of the present invention and are not intended to limit the present invention. Any modifications, supplements and equivalent substitutions made within the scope of the principles of the present invention should be included in the protection scope of the present invention.

Claims

1. A real-time correction observation method for acoustic tomography of small-scale horizontal flow fields. It is characterized in that include: (1) Arranging an acoustic transceiver system in the observation waters, wherein the acoustic transceiver system uses four acoustic transducers to perform mutual transmission of acoustic signals; The acoustic transceiver system is used to collect the original hydroacoustic data in the observed water area, and the temperature and depth meter is used to collect the temperature information of the layer where the acoustic transducer is located; (2) Performing station cross-correlation on the collected raw underwater acoustic data, and extracting the transmission time corresponding to the bidirectional arrival peak between the acoustic stations by setting the time window and signal-to-noise ratio threshold; (3) Preprocess the extracted arrival peak transmission time, remove abnormal data, define the maximum difference of the return propagation time, and remove the misidentified peaks; (4) Combine the extracted two-way transmission time mean and the temperature and depth information of the layer where the acoustic transducer is located to calculate the spacing of the acoustic transceiver system and perform station correction; (5) In the two-dimensional horizontal profile, combined with the station layout, horizontal grid division is performed, the coordinates of the grid center point are calculated, the stream function is Taylor expanded, and the corresponding known coefficient matrix and the coefficient matrix to be solved are established; When performing horizontal grid division, each acoustic station is projected to the horizontal section where the acoustic transducer depth is located, and r grids are divided. The sound line propagates within the grid. The specific parameters and formulas are as follows: For each pair of sound stations, the direct path is: After Taylor expansion, we get Among them, l ij is the length of the i-th sound line passing through the j-th grid, C 0 is the sound velocity of the layer where the acoustic transducer is located, V j represents the projection of the flow velocity of the jth grid along the propagation path, represents the forward and reverse reference transmission time along the ith sound ray, Δt i is the two-way propagation time difference; (6) Perform two-dimensional horizontal grid flow field inversion by using the conical least squares method combined with the L-curve method, and constrain the error threshold; if the error exceeds the threshold, return to the flow field inversion until the requirements are met; (7) Visualize the two-dimensional horizontal flow field obtained by inversion.

2. The small-scale horizontal flow field acoustic tomography real-time correction observation method according to claim 1, It is characterized in that In step (1), the acoustic transceiver system includes a float, an acoustic transducer, a temperature-salinity-depth profiler and a data acquisition system; The data acquisition system is located on the shore and connected to the acoustic transducers; four acoustic transducers are arranged at the same depth underwater, one end of which is connected to a buoy moored under the water surface and the other end is anchored at the bottom of the water.

3. The small-scale horizontal flow field acoustic tomography real-time correction observation method according to claim 1, It is characterized in that In step (1), in the acoustic transceiver system, each acoustic station uses the same-order M sequence with the same frequency but different modes, and adopts the same-transmit and same-receive mode for acoustic signals, ensuring that each pair of acoustic transducers receives a two-way acoustic signal after each signal is transmitted.

4. The small-scale horizontal flow field acoustic tomography real-time correction observation method according to claim 1, It is characterized in that The specific process of step (2) is as follows: The collected original underwater acoustic data are cross-correlated between two sound stations. By setting the time window and signal-to-noise ratio threshold, the transmission time of the first arrival peak, that is, the transmission time of the correlation peak corresponding to the direct path, is distinguished and extracted.

5. The small-scale horizontal flow field acoustic tomography real-time correction observation method according to claim 1, It is characterized in that The specific process of step (4) is as follows: (4-1) Based on the extracted two-way propagation time, calculate the mean two-way propagation time of each path during each signal transmission Among them, t + is the forward transmission time, t - is the back propagation time; (4-2) Calculate the reference sound velocity C based on the temperature of the layer where the transducer is located measured by the temperature depth meter 0 , calculate the corrected distance between the two stations (4-3) The station position is corrected by solving the intersection point of the station position line, and the formula is as follows: Among them, (x n ,y n ) is the initial coordinate point of the station, (Δx n ,Δy n ) is the position drift of station n, (x n +Δx n ,y n +Δy n ) is the corrected station coordinate; since |Δx n -Δx m |<<|x n -x m |,|Δy n -Δy m |<<|y n -y m |, perform Taylor expansion on the above equation, and set station 1 close to the shore as a fixed point, let Δx 1 =Δy 1 =0 After simplifying, we get: By solving the above equations, the station correction can be completed.

6. The small-scale horizontal flow field acoustic tomography real-time correction observation method according to claim 1, It is characterized in that In step (5), the specific process of establishing the corresponding known coefficient matrix and the coefficient matrix to be solved is: Introducing stream functions Among them, the unknown coefficient matrix Known coefficient matrix (x, y) is the coordinate of the point to be measured; L x , L y is the length and width of the inversion area; N x , N y are the number of grids in the inversion area along the x and y directions, respectively; Calculate the north component of velocity based on the stream function and the east component of velocity 7. The small-scale horizontal flow field acoustic tomography real-time correction observation method according to claim 6, It is characterized in that In step (6), the specific process of performing two-dimensional horizontal grid flow field inversion is as follows: Substitute the north component u and the east component v of the velocity into get Written as a matrix equation y = Ex + n; where x = D m is the quantity to be solved, n is the solution error, y=Δt i is the measured time difference of the sound transmission between the two stations. is the transformation matrix; The cone least squares method is used to solve the above matrix equation. By setting the penalty function J = n T n+α 2 x T x=(y-Ex) T (y-Ex)+α 2 x T x, minimize the penalty function to get the expected optimal solution of x Then the stream function expression is obtained, and the coordinates of the point to be solved are substituted to obtain the flow field in the observation area.

Images