Sound velocity profile inversion method based on array multi-path time delay structure

By using a sound velocity profile inversion method based on array multipath delay structure, a synthetic dataset is generated and a sound velocity perturbation model is inverted using a genetic algorithm and the multipath delay structure of a deep-sea array. This solves the problem of insufficient sound velocity profile prediction accuracy in complex marine environments and achieves high-precision sound velocity profile prediction.

CN117521335BActive Publication Date: 2026-01-06NORTHWESTERN POLYTECHNICAL UNIV
View PDF 2 Cites 0 Cited by

Patent Information

Application Number
CN202311339925.4
Authority / Receiving Office
CN · China
Patent Type
Patents(China)
Current Assignee / Owner
Filing Date
2023-10-16
Publication Date
2026-01-06
Estimated Expiration
2043-10-16

AI Technical Summary

Technical Problem

Existing technologies are insufficient in predicting and inverting sound velocity profiles when facing complex marine environments such as vortices, internal waves, and strong fronts. Satellite remote sensing technology and acoustic observation methods have significant errors.

Method used

A sound velocity profile inversion method based on array multipath delay structure is adopted. By generating a synthetic dataset and using a genetic algorithm to invert the sound velocity perturbation model, and combining the multipath delay structure of a large-depth array with an empirical orthogonal decomposition algorithm, high-precision prediction of the sound velocity profile is achieved.

Benefits of technology

It improves the accuracy of sound velocity profile prediction and inversion, simplifies sound velocity profile prediction in complex marine environments, directly retrieves sound velocity profiles from acoustic observations, and reduces dependence on marine dynamic processes.

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN117521335B_ABST
    Figure CN117521335B_ABST
Patent Text Reader

Abstract

The present application relates to a kind of acoustic velocity profile inversion method based on array multi-path time delay structure, establishes acoustic velocity disturbance model, utilizes vertical displacement disturbance and sea surface acoustic velocity disturbance to characterize acoustic velocity disturbance information, realizes the nonlinear mapping of high degree of freedom disturbance structure of acoustic velocity profile to vertical displacement disturbance and sea surface acoustic velocity disturbance.Synthetic data set is established, based on empirical orthogonal decomposition algorithm and adding constraint condition, generate new synthetic data set, enhance the variability of data set, make it cover the real acoustic velocity profile of study area to the greatest extent.Based on the multi-path time delay structure between the direct signal of large depth array and the first sea surface reflection signal, construct fitness function, realize the prediction inversion of acoustic velocity profile by combining genetic algorithm, utilize the multi-path time delay arrival structure of each channel received signal of large depth array and an explosive bomb signal, can realize the prediction inversion of acoustic velocity profile, realize simple, applicable to the correction of satellite remote sensing technology assimilation prediction, improve prediction accuracy.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] This invention belongs to the fields of marine physics, marine engineering, and underwater acoustics engineering. It relates to a sound velocity profile inversion method based on an array multipath delay structure. Specifically, it involves using sound velocity samples from a public database to establish a sound velocity perturbation model and a synthetic dataset model to generate a sample set. The mean and modal basis vectors of the sample set are solved using an empirical orthogonal decomposition algorithm. Then, the sound velocity profile inversion is achieved by using a genetic algorithm and the multipath delay structure of a deep array. Background Technology

[0002] Influenced by ocean dynamic processes such as eddies, internal waves, and fronts, the high degree of freedom variation of the sound velocity profile over time significantly impacts sound propagation. Therefore, the prediction and inversion of underwater sound velocity profiles has always been a hot research topic in oceanography. Currently, the prediction and inversion of underwater sound velocity profiles mainly employ data assimilation methods based on satellite remote sensing technology and sound velocity correction methods based on acoustic observations. However, these methods still face significant technical challenges in complex ocean environments with severe sound velocity disturbances. The main reason for this is that these methods have certain limitations when dealing with complex ocean environments characterized by frequent eddies, internal waves, and strong fronts, as analyzed in detail below:

[0003] (1) Data assimilation method of satellite remote sensing technology. This method utilizes the correlation between sea surface remote sensing parameters and data characteristics of sound velocity profiles at various depths, fits the regression relationship between the two, and realizes the prediction of sound velocity profiles. There are two main methods for establishing the correlation between sea surface remote sensing parameters and sound velocity profile data characteristics: one is to establish a benchmark sound velocity profile and fit the regression relationship between satellite remote sensing parameters and sound velocity disturbance values; the other is to perform empirical orthogonal decomposition on the sound velocity profile set that changes over time, extract its mode vector and empirical orthogonal function coefficients, and fit the regression relationship between satellite remote sensing parameters and empirical orthogonal function coefficients. The prediction and inversion accuracy of sound velocity profiles is related to ocean dynamic processes, spatiotemporal resolution, and remote sensing observation accuracy. Especially when facing complex ocean environments with vortices, frequent internal waves, and strong fronts, the correlation between sea surface remote sensing parameters and extracted sound velocity profile data characteristics becomes lower, which will bring significant errors to the prediction and inversion of sound velocity profiles.

[0004] (2) Correction method based on acoustic observations. This method uses traditional correlation models or neural network models to fit the regression relationship between the sound velocity profile and acoustic observations such as sound propagation and multipath structure. It then corrects the prediction based on satellite remote sensing assimilation, reducing the error accumulated over time. This method is closely related to the degree of correlation between the sound velocity profile and acoustic observations, as well as ocean dynamic processes. Especially in complex ocean environments with vortices, frequent internal waves, and strong fronts, the high degree of freedom perturbation of the sound velocity profile can cause high degree of freedom changes in acoustic observations, making it difficult to find a one-to-one mapping relationship between the sound velocity profile and acoustic quantities, which can introduce significant errors into the prediction and inversion of the sound velocity profile.

[0005] In short, data assimilation methods based on satellite remote sensing and correction methods based on acoustic observations generally lead to significant deviations in the results of sound velocity profile prediction and inversion in sea areas with severe sound velocity disturbances. Therefore, it is necessary to seek new principles and technical approaches. Summary of the Invention

[0006] Technical problems to be solved

[0007] To avoid the shortcomings of existing technologies, this invention proposes a sound velocity profile inversion method based on an array multipath delay structure, which can quantify sound velocity disturbances from the perspective of the model and is particularly suitable for sound velocity profile inversion under deep-sea sound velocity disturbance conditions.

[0008] Technical solution

[0009] A method for inverting sound velocity profiles based on an array multipath delay structure, characterized by the following steps:

[0010] Step 1: Generate several sound velocity profiles (SSPs) within the study area using publicly available historical Argo data and publicly available WoA18 data. s Using the base sample as the mean, the average value of the resulting sound speed profile is used as the reference sound speed profile ssp. ref ;

[0011] For depths above 1000m, historical Argo data was used;

[0012] For depths below 1000m, use woa18 data;

[0013] Using the reference sound speed profile ssp ref The data is used as input to the sound speed perturbation model to obtain the perturbation dataset ssp. dis ;

[0014] Step 2: Perturb the dataset ssp dis The data is used as input to the empirical orthogonal decomposition algorithm to obtain a perturbation dataset expressed in terms of the mean sound velocity profile, mode vector and corresponding coefficients, resulting in Formula 5 as follows;

[0015]

[0016] Among them: ssp mean It is the mean of the perturbation dataset, Φ k It is the k-th mode vector of the empirical orthogonal decomposition, α k For the coefficients corresponding to the k-th mode vector in a selected set of empirical orthogonal decomposition coefficients α, ssp dis This is the perturbation dataset corresponding to the selected set of empirical orthogonal decomposition coefficients;

[0017] Calculate the corresponding coefficients α of the modal basis vectors k mean μ k and standard deviation σ k From a uniform distribution U(μ) k +2σ k <= α k <=μ k +2σ k Randomly select a group α from ) k The coefficients are used to generate the corresponding sound velocity profiles using Formula 5.

[0018] For the new set of sound speed profiles, enforce physical constraints:

[0019] -γΔc min,j ≤c syn,j -c mean,j ≤γΔc max,j

[0020] γ≥1,j=1,2,...,m

[0021] c mean,j c represents the average speed of sound in the perturbed dataset at depth j. syn,j Let Δc be the speed of sound in the synthetic dataset at depth j. max,j and Δc min,j represents the maximum and minimum sound speed perturbations of the perturbation dataset compared to its average sound speed at depth j, respectively, where γ is the tolerance factor;

[0022] The synthesized datasets obtained after constraints and filtering are combined to form the final synthesized dataset ssp. syn ;

[0023] The synthetic dataset is then subjected to EOF empirical orthogonal decomposition again using the following formula 7:

[0024]

[0025] Where: ssp' mean It is the mean of the synthetic dataset, Φ k ' is the k-th modality vector obtained by performing empirical orthogonal decomposition on the synthetic dataset, αk ' represents the coefficient corresponding to the k-th mode vector in a selected set of empirical orthogonal decomposition coefficients α', ssp syn The selected set of empirical orthogonal decomposition coefficients corresponds to the synthetic dataset;

[0026] Step 3: Use only the coefficient α k The first three coefficients α1', α2', and α3' are used as the inversion parameters for the genetic algorithm. The range of variation of coefficients α1', α2', and α3' is calculated, and a set of coefficients α1', α2', and α3' is arbitrarily selected within the range of empirical orthogonal decomposition coefficients. The sound velocity profile that maps one-to-one with it is found in the synthetic dataset using Formula 7, and then the multipath delay arrival structure dly(a) corresponding to the sound velocity profile is calculated using the Bellhop acoustic ray model. k ',i,j);

[0027] The constructed fitness function is the objective function f(α) k Formula 8 is represented by ').

[0028]

[0029] Where n represents the number of channels in the array; dlyrl(i,j) is the multipath delay structure corresponding to the acoustic observation data;

[0030] A genetic algorithm is used to estimate the parameters α1', α2', and α3' to be inverted, and then the mean ssp of the synthetic dataset is combined. mean 'And the first three modal vectors Φ1', Φ2' and Φ3' of the empirical orthogonal decomposition of the synthetic dataset, substitute the parameters into Equation 7 to realize the inversion of the sound speed profile, that is, to find the sound speed profile with the lowest fitness in the synthetic dataset.

[0031] Several sound speed profiles (ssp) were generated using empirical formulas for sound speed. s As a base sample.

[0032] The sound speed disturbance models are: Model A, in which the vertical movement of seawater lifts the isotherms to the ocean surface, thus thinning the mixed layer; and Model B, in which there is no upward movement of isotherms when there is a strong temperature and salinity anomaly in the upper ocean, indicating a warm upper ocean characteristic.

[0033] The sonic perturbation A model: When perturbation A occurs, the depth structure model of the vertical displacement of the feature center is as follows:

[0034]

[0035] Where z is the depth of the feature center, ζ0 is the vertical displacement perturbation, and D z To set the vertical scale of the perturbation. It is a normalization factor, and the occurrence of this perturbation produces upper ocean enhancement features;

[0036] Calculate sound speed disturbances using nonlinear equations:

[0037]

[0038] Where, γ a =0.016s -1 This is the adiabatic sound velocity gradient. The second term in square brackets is determined to be the sound velocity at depth z of the stationary block after it has been adiabatically raised to a new depth.

[0039] The sound speed perturbation B model: When perturbation B occurs, the depth structure model of the sound speed at the feature center is as follows:

[0040]

[0041] Where z is the depth of the feature center z0 = 2.5 * h, δc0 is the sound speed perturbation value, and h is the thickness of the mixing layer at the perturbation center.

[0042] The reference sound speed profile ssp ref The data is used as input to the sound speed perturbation model to obtain the perturbation dataset ssp. dis At that time, the surface sound velocity disturbance δc0 was set to a range of (-1,1) m / s with an interval of 0.1 m / s, and the vertical displacement disturbance ζ0 was set to a range of (-50,90) m with an interval of 1 m.

[0043] The specific process of using a genetic algorithm to estimate the parameters α1', α2', and α3' to be inverted in step 3 is as follows:

[0044] (a) Set the fitness function threshold to η and the maximum number of iterations W0: Pre-encode the parameters α1', α2', and α3' to be inverted using binary codes. Each seed of the parameter corresponds to a unique binary code. Based on the set number of initial seeds N, generate N sets of parameters to be inverted as initial seeds. The pre-encoding results of α1', α2', and α3' are respectively represented by (A1, A2, ..., A...). N ), (B1,B2,...,B N (C1,C2,...,C) and (C1,C2,...,C) N This means that each initial seed corresponds to a unique set of parameters α. k ';

[0045] (b) Next, the constructed fitness function formula 8 is used to calculate each initial seed, i.e., each selected set of parameters α. k The corresponding fitness value is used as the initial fitness value, and W = 1 is set.

[0046] (c) Then, based on the selected parameter α... k The corresponding fitness value is used to select, crossover, and mutate the first generation of seeds to produce the next generation of seeds;

[0047] (d) Calculate the fitness of the next generation seed according to fitness function formula 8. If the fitness function does not reach the set threshold f(a) k If f(a) > η or the number of iterations has not reached the algorithm-defined number of iterations W < W0, let W = W + 1, and return to step (c) until f(a) > η. k Until ') <= η or W = W0.

[0048] The range of variation of the coefficients α1', α2', and α3': the maximum and minimum values ​​of the inversion parameters α1', α2', and α3' (α 1min a 1max ), (α ni2m α 2max ), (α ni3m α 3max ) represents the range of variation of the coefficients α1', α2', and α3' corresponding to the first three basis vectors of the calculated mode.

[0049] The value of the tolerance factor γ variable is 1.

[0050] An electronic device, characterized in that it includes a processor and a memory, wherein the processor is used to execute a computer program stored in the memory to implement the data migration method of the sound velocity profile inversion method.

[0051] Beneficial effects

[0052] This invention proposes a sound velocity profile inversion method based on array multipath delay structure. A sound velocity perturbation model is established, and vertical displacement perturbation and sea surface sound velocity perturbation are used to characterize sound velocity perturbation information. This achieves a nonlinear mapping from the high-degree-of-freedom perturbation structure of the sound velocity profile to vertical displacement perturbation and sea surface sound velocity perturbation. A synthetic dataset is proposed. Based on an empirical orthogonal decomposition algorithm and with added constraints, a new synthetic dataset is generated, enhancing the dataset's variability and maximizing its coverage of the true sound velocity profile of the study area. Based on the multipath delay structure between the direct signal and the primary sea surface reflection signal from a deep-depth array, a fitness function is constructed, and a genetic algorithm is used to achieve the prediction and inversion of the sound velocity profile. This method utilizes the multipath delay arrival structure of the signals received from each channel of the deep-depth array and a single explosive bomb signal to achieve the prediction and inversion of the sound velocity profile. It is simple to implement and suitable for correcting satellite remote sensing assimilation predictions, thereby improving prediction accuracy.

[0053] Compared to existing technologies, the beneficial effects are reflected in:

[0054] (1) This invention uses a sound speed perturbation model and a method to build a synthetic dataset to construct a dataset, which covers the real sound speed profile of the study area to a greater extent compared with other methods, and improves the prediction and inversion accuracy of the sound speed profile.

[0055] (2) This invention uses a genetic algorithm to estimate the sound velocity profile based on the multipath time delay arrival structure of the direct wave and the first sea surface reflected wave of the deep receiving array. This multipath structure is sensitive to changes in the sound velocity profile and has a unique mapping relationship with the sound velocity profile, thus having high estimation accuracy.

[0056] (3) This invention does not require modeling of complex ocean dynamic processes and ocean meteorological processes. It directly retrieves the sound velocity profile from acoustic observations, making it simple to implement. Attached Figure Description

[0057] Figure 1 (a) Schematic diagram of the experimental area (the rectangular area is the research sea area, the array deployment position is (18.00°N, 129.50°E, the sound source position is 22.46 km south of the array deployment position)); (b) publicly available Argo data sound velocity profile of the experimental area, the red line is the reference sound velocity profile obtained by averaging; (c) publicly available WOA18 data sound velocity profile of the acoustic observation experimental location.

[0058] Figure 2 (a) Sound speed profile samples of the perturbation dataset established using the sound speed perturbation model; the red line represents the reference sound speed profile; (b) EOF in the figure. i (c) Mode vectors of the i-th order empirical orthogonal function after empirical orthogonal decomposition of the perturbation dataset; (d) Synthetic dataset sound velocity profile samples, with the red line representing the reference sound velocity profile; i It is the mode vector of the i-th order empirical orthogonal function of the synthetic dataset after empirical orthogonal decomposition.

[0059] Figure 3 Genetic Algorithm Structure Flowchart

[0060] Figure 4 : Schematic diagram of the genetic algorithm precoding, selection, crossover and mutation process; (a) Schematic diagram of precoding N groups of selected coefficients using empirical orthogonal function coefficient α1' as an example; (b) Schematic diagram of the genetic algorithm selecting the next generation of parents based on the current generation fitness value; (c) Schematic diagram of the genetic algorithm performing crossover on the two parents selected in the Wth generation to obtain the W+1th generation; (d) Schematic diagram of the genetic algorithm performing mutation on the W+1th generation sample obtained from the crossover process.

[0061] Figure 5 : Convergence image of sound velocity profile prediction and inversion using genetic algorithm.

[0062] Figure 6 A comparison chart of predicted inversion results, mean values ​​of WOA18 and Argo data, and measured acoustic experimental data.

[0063] Figure 7 Comparison of the multipath time delay structure corresponding to the predicted sound velocity profile and the multipath time delay structure of the actual acoustic experimental data. Detailed Implementation

[0064] The present invention will now be further described in conjunction with the embodiments and accompanying drawings:

[0065] Based on the sound velocity profile inversion method using array multipath delay structure, several sound velocity profiles (SSPs) are generated within the study sea area using publicly available Argo historical data and WOA18 data, employing empirical sound velocity formulas. s Using these as the base samples, the mean of the resulting sound speed profile is then taken as the reference sound speed profile ssp. ref Based on the reference sound velocity profile, a sound velocity perturbation model is used to establish a perturbation dataset. This model divides the sound velocity perturbation into two parts: vertical displacement perturbation and ocean surface sound velocity perturbation. Next, an empirical orthogonal decomposition algorithm and defined constraints are used to establish a synthetic dataset. The synthetic dataset is then subjected to another empirical orthogonal decomposition to obtain the range of values ​​for the parameters to be estimated. The Bellhop acoustic ray model is used to calculate the array multipath delay structure corresponding to each set of parameters to be estimated. Finally, a fitness function is constructed based on the array multipath delay structure observed in acoustic experiments, and a genetic algorithm is used to estimate the parameters, thereby further realizing the inversion of the sound velocity profile. The process consists of the following three steps.

[0066] Step 1: Select the experimental sea area (e.g., Figure 1 (a) shows that the publicly available historical Argo data generates several sound speed profiles (SSP) using empirical formulas for sound speed. s As the base sample, the sound velocity profile above 1000m was calculated using publicly available Argo data from the experimental sea area, while the profile below 1000m was extended using WOA18 data from the acoustic experiment observation location. The average of these data was then calculated, and the resulting average sound velocity profile was used as the reference sound velocity profile (SSP). ref .

[0067] Based on the reference sound speed profile ssp ref A sound speed perturbation model was established, and a perturbation dataset was generated. The sound speed perturbation in the ocean is divided into two parts. Perturbation A is the vertical movement of seawater that lifts the isotherms to the ocean surface, thereby thinning the mixed layer. Perturbation B simulates the warm characteristics of the upper ocean without the isotherms moving upward. This situation may occur when there are strong temperature and salinity anomalies in the upper ocean.

[0068] When disturbance A occurs, the depth structure model of the vertical displacement of the feature center is as follows:

[0069]

[0070] Where z is the depth of the feature center, ζ0 is the vertical displacement perturbation, and D z To set the vertical scale of the perturbation. It is a normalization factor, and the occurrence of this perturbation produces upper ocean enhancement features.

[0071] Calculate sound speed disturbances using nonlinear equations:

[0072]

[0073] Where, γ a =0.016s -1 This is the adiabatic sound velocity gradient. The second term in square brackets is determined to be the sound velocity at depth z of the stationary block after it has been adiabatically raised to a new depth.

[0074] When disturbance B occurs, the depth structure model of the sound velocity at the feature center is as follows:

[0075]

[0076] Where z is the depth of the feature center, δc0 is the sound speed perturbation value, and h is the thickness of the mixing layer at the perturbation center.

[0077] z0 = 2.5 * h (4)

[0078] The surface acoustic velocity perturbation δc0 is set to range from (-1, 1) m / s with an interval of 0.1 m / s, and the vertical displacement perturbation ζ0 is set to range from (-50, 90) m with an interval of 1 m. A perturbation dataset ssp is generated using the acoustic velocity perturbation model. dis

[0079] Step 2: Establish a synthetic dataset based on the empirical orthogonal decomposition algorithm to enhance the variability of the dataset and make it cover the real sound speed profile of the region to a large extent. First, the perturbation dataset is decomposed into the mean of the sound speed profile, the mode vector and the corresponding coefficients using the empirical orthogonal decomposition algorithm, as shown in Equation (5).

[0080]

[0081] ssp mean It is the mean of the perturbation dataset, Φ k It is the k-th mode vector of the empirical orthogonal decomposition, α k For the coefficients corresponding to the k-th mode vector in a selected set of empirical orthogonal decomposition coefficients α, ssp disThis is the perturbation dataset corresponding to the selected set of empirical orthogonal decomposition coefficients. The mean and standard deviation (μ1, σ1), (μ2, σ2), and (μ3, σ3) of the coefficients α1, α2, and α3 corresponding to the first three basis vectors of the modes are calculated. To generate the synthetic dataset, it is necessary to obtain the perturbation dataset from the uniform distribution U(μ1, σ2, σ3). k +2σ k <= α k <=μ k +2σ k Randomly select a group α from ) k The coefficient is then used to generate the corresponding sound speed profile using equation (5), by using this arbitrary α k The coefficient set, while preserving the EOF basis vector Φ k Without changing α, a new set of sound speed profiles can be generated. k It is randomly sampled, therefore the synthesized sound speed profile dataset ssp syn This may be physically impractical. Therefore, physical constraints are enforced on these generated synthetic datasets as shown in Equation (6).

[0082]

[0083] c mean,j c represents the average speed of sound in the perturbed dataset at depth j. syn,j Let Δc be the speed of sound in the synthetic dataset at depth j. max,j and Δc min,j These represent the maximum and minimum sound speed perturbations of the perturbed dataset compared to its average sound speed at depth j, respectively. γ is a tolerance factor, which also provides corresponding constraints on the gradient of the synthetic dataset, ensuring its smoothness. In this patent, the value of the γ variable is set to 1. After constraints and filtering, the synthetic dataset is obtained, and then the perturbed datasets are merged to form the final synthetic dataset ssp. syn .

[0084] Then, the synthetic dataset was subjected to EOF empirical orthogonal decomposition again using equation (7). According to the results of the second EOF empirical orthogonal decomposition, it can be reflected that when the first 3 modal basis vectors are selected, the modal variance contribution rate accounts for 99.98%. The first 3 modal basis vectors can represent the data feature quantity of the entire synthetic dataset. The coefficients α corresponding to the first 3 modal basis vectors k 'These are the coefficients to be inverted.'

[0085]

[0086] At this time, ssp' mean It is the mean of the synthetic dataset, Φ k ' is the k-th modality vector obtained by performing empirical orthogonal decomposition on the synthetic dataset, α k' represents the coefficient corresponding to the k-th mode vector in a selected set of empirical orthogonal decomposition coefficients α', ssp syn The selected set of empirical orthogonal decomposition coefficients corresponds to a synthetic dataset. The maximum and minimum values ​​(α1', α2', and α3') of the coefficients α1', α2', and α3' corresponding to the first three basis vectors of the modes are calculated. 1min a 1max ), (α 2min α 2max ), (α 3min α 3max This allows us to obtain the maximum and minimum ranges of sound speed perturbation variation in the synthetic dataset, within which a set of EOF coefficients α is selected. k This allows it to cover any sound speed profile within the entire synthetic dataset.

[0087] Step 3: Select the relative time delay difference between the direct wave and the first sea surface reflected wave as the independent variable to retrieve the sound speed profile. Within the range of the empirical orthogonal decomposition coefficients obtained in Step 2, arbitrarily select a set of coefficients α. k ', the sound velocity profile that maps one-to-one with it can be found in the synthetic dataset using equation (7), and then the multipath arrival structure dly(a) corresponding to the sound velocity profile can be calculated using the Bellhop acoustic ray model. k ',i,j). When j=1, dly(a k ',i,j) represents the selected set of parameters α k The relative time delay of the direct wave in the i-th channel of the array is calculated using the Bellhop model. When j=2, dly(a k ',i,j) represents the selected set of parameters α k The relative time delay of the first sea surface reflected wave in the i-th channel of the array is calculated using the Bellhop model. Let dlyrl(i,j) be the multipath delay structure corresponding to the acoustic observation data. Similarly, when j=1, dlyrl(i,j) represents the relative time delay of the direct wave in the i-th channel of the acoustic observation data; when j=2, dlyrl(i,j) represents the relative time delay of the first sea surface reflected wave in the i-th channel of the acoustic observation data. Therefore, the constructed fitness function (objective function) f(α) k As shown in equation (8):

[0088]

[0089] Where n represents the number of channels in the array.

[0090] After constructing the fitness function, a genetic algorithm is used to estimate the parameters α1', α2', and α3' to be inverted. Finally, the mean ssp of the synthetic dataset is combined. mean'And the first three modal vectors Φ1', Φ2' and Φ3' of the empirical orthogonal decomposition of the synthetic dataset, substitute the parameters into equation (7) to realize the inversion of the sound speed profile, that is, to find the sound speed profile with the lowest fitness in the synthetic dataset. The specific process of parameter estimation by the genetic algorithm is as follows:

[0091] (a) Set the fitness function threshold to η and the maximum number of iterations W0. Pre-encode the parameters α1', α2', and α3' to be inverted using binary codes. Each seed of a parameter corresponds to a unique binary code. Based on the set initial seed number N, generate N sets of parameters to be inverted as initial seeds. The pre-encoding results of α1', α2', and α3' are respectively represented by (A1, A2, ..., A...). N ), (B1,B2,...,B N (C1,C2,...,C) and (C1,C2,...,C) N This means that each initial seed corresponds to a unique set of parameters α. k '.

[0092] (b) Then, the constructed fitness function (8) is used to calculate each initial seed, i.e., each set of selected parameters α. k The corresponding fitness value is used as the initial fitness value, and W = 1 is set.

[0093] (c) Then, based on the selected parameter α... k The corresponding fitness value is used to select, crossover, and mutate the first generation of seeds to produce the next generation of seeds.

[0094] (d) Calculate the fitness of the next generation seed according to the fitness function (8). If the fitness function does not reach the set threshold f(a) k If f(a) > η or the number of iterations has not reached the algorithm-defined number of iterations W < W0, let W = W + 1 and return to step (c) until f(a) > η. k Until ') <= η or W = W0.

[0095] Further explanation with reference to the attached diagram:

[0096] Figure 1 (a) The area of ​​the acoustic observation experiment, the deployment location of the array, and the location of the sound source are given; (b) Publicly available historical argo data and their mean (reference sound velocity profile) are given; (c) Publicly available woa18 data sound velocity profile of the acoustic experiment observation location is given.

[0097] Figure 2 (a) gives the following Figure 1 (a) Using the reference sound speed profile in the model as a benchmark, a perturbation dataset of sound speed profile samples is established using the sound speed perturbation model; (b) provides the sound speed profile samples of the perturbation dataset. Figure 2(a) The first three mode vectors obtained by empirical orthogonal decomposition of the perturbation dataset account for 99.96% of the variance; (c) A sample of the sound speed profile of the synthetic dataset is given based on the coefficient statistical characteristics obtained by empirical orthogonal decomposition of the perturbation dataset; (d) A sample of the sound speed profile of the synthetic dataset is given. Figure 2 (c) The first three modal vectors obtained by empirical orthogonal decomposition of the synthetic dataset account for 99.98% of the variance contribution rate of the first three modal vectors. The synthetic dataset has 7126 sound velocity profile samples, while the perturbation dataset has 2961 sound velocity profile samples. This method increases the number of sound velocity profile samples by more than 100%.

[0098] Figure 3 A flowchart for using a genetic algorithm to predict and invert sound velocity profiles is presented. First, the initial number of seeds is set to 250, and the maximum number of iterations is 15. The parameters to be inverted are pre-encoded. Then, the fitness function is used to calculate the fitness of each initial seed. Based on the fitness of each seed, selection, crossover and mutation processes are performed to generate the next generation of seeds. The iteration stops when the fitness value reaches a set threshold or the number of iterations reaches the set maximum number.

[0099] Figure 4 The following diagrams illustrate the specific processes of selection, crossover, and mutation in a genetic algorithm: (a) A diagram illustrating the pre-encoding of N groups of selected coefficients, using the empirical orthogonal function coefficients α1' as an example, where each group of empirical orthogonal function coefficients corresponds to a unique binary code. (b) A diagram illustrating the pre-encoding of N groups of selected coefficients according to the initial seed (A1, A2, ..., A...). N A pie chart was plotted using the corresponding fitness values ​​as proportions. Two parental samples A were selected from the disk through two random samplings. i and A j This is the selection process. (c) shows the crossover process of the genetic algorithm, which... Figure 4 (b) The binary codes of the two parent samples selected are partially crossed over according to a pre-set probability, which is the crossover process, thus generating two new parameters α1'. (d) shows the mutation process of the genetic algorithm, which... Figure 4 (c) The process by which the generated new seed causes mutations in the binary code portion according to a pre-set mutation probability is called the mutation process.

[0100] Figure 5 The fitness convergence graph of the sound velocity profile prediction and inversion using the genetic algorithm is given as the number of iterations increases. The fitness of the genetic algorithm converges when the number of iterations reaches the 10th round.

[0101] Figure 6The paper presents a comparison chart of the predicted inversion results, the mean values ​​of the WOA18 and Argo data, and the acoustic experimental measured data. The root mean square error (RMSE) between the sound velocity profile prediction and inversion results using the method proposed in this patent and the actual sound velocity profile observed in the acoustic experiment is 1.3 m / s. The RMSE between the mean value of the disclosed Argo data sound velocity profile (reference sound velocity profile) and the actual sound velocity profile observed in the acoustic experiment is 3.2 m / s. The RMSE between the disclosed WOA18 sound velocity profile data and the actual sound velocity profile observed in the acoustic experiment is 6.48 m / s.

[0102] Figure 7 A comparison diagram is provided between the multipath time delay structure corresponding to the predicted sound velocity profile and the multipath time delay structure of the actual acoustic experimental data.

[0103] Table 1 presents the statistical properties of the empirical orthogonal function coefficients α1, α2, and α3 obtained by empirical orthogonal decomposition of the perturbation dataset.

[0104] Table 1 Statistical characteristics of EOF decomposition coefficients in the perturbation dataset

[0105]

[0106] Table 2 shows the maximum and minimum values ​​of the empirical orthogonal function coefficients α1', α2', and α3' obtained by performing empirical orthogonal decomposition on the synthetic dataset, which represent the range of values ​​for the parameters to be inverted by the genetic algorithm.

[0107] Table 2 EOF decomposition coefficients of the synthetic dataset

[0108]

[0109] The method proposed in this invention has achieved significant implementation results in typical examples. The sound velocity profile inversion method based on array multipath delay structure has superior performance. It does not require modeling of complex ocean dynamic and ocean meteorological processes, and directly predicts and inverts the sound velocity profile from acoustic observations, making it simple to implement.

Claims

1. A sound velocity profile inversion method based on an array multi-path time delay structure, characterized in that The steps are as follows: Step 1: Generate several sound speed profiles in the study area using the published argo historical data and published woa18 data As the base sample, the mean value is obtained, and the average sound speed profile is obtained as the reference sound speed profile ; Above 1000m depth, argo historical data is used; Below 1000m, woa18 data is used; with reference sound speed profile data as input to a sound speed perturbation model, resulting in a perturbation dataset ; Step 2: Perturb the dataset The data is used as input to the empirical orthogonal decomposition algorithm to obtain a perturbation dataset expressed in terms of the mean sound velocity profile, mode vector and corresponding coefficients, resulting in Formula 5 as follows; wherein: is the mean of the perturbation data set, is the i-th modal vector of the empirical orthogonal decomposition, k is the i-th coefficient of the selected set of empirical orthogonal decomposition coefficients, k is the i-th coefficient corresponding to the i-th modal vector, k is the perturbation data set corresponding to the selected set of empirical orthogonal decomposition coefficients; calculating the mean and the standard deviation of the i-th coefficient , randomly drawing a set of coefficients from a uniform distribution U( ), and generating the corresponding sound speed profile using equation 5;​​​ For a new set of sound speed profiles, physical constraints are enforced: c denotes the average sound velocity of the perturbed data set at a depth of c denotes the sound velocity of the synthetic data set at a depth of c denotes the sound velocity of the synthetic data set at a depth of c denotes the sound velocity of the synthetic data set at a depth of c denotes the sound velocity of the synthetic data set at a depth of c denotes the sound velocity of the synthetic data set at a depth of c denotes the sound velocity of the synthetic data set at a depth of c denotes the sound velocity of the synthetic data set at a depth of The synthetic datasets that have been constrained and screened are brought together to form a final synthetic dataset ; The EOF empirical orthogonal decomposition is performed again on the synthetic dataset using the following equation 7: in: It is the mean of the synthetic dataset. The first step is to perform empirical orthogonal decomposition on the synthetic dataset. k One modal vector, For a selected set of empirical orthogonal decomposition coefficients The Middle k The coefficients corresponding to each modal vector The selected set of empirical orthogonal decomposition coefficients corresponds to the synthetic dataset; Step 3: only use the coefficients of the first 3 orders , and as the parameters to be inverted by the genetic algorithm, calculate the variation ranges of the coefficients , and , and randomly select a set of coefficients , and from the range of the empirical orthogonal decomposition coefficients; find the sound velocity profile that is one-to-one mapped with the set of coefficients in the synthetic data set by formula 7, and then calculate the multi-path time delay arrival structure corresponding to the sound velocity profile using the bellhop acoustic ray model ; The fitness function of the construction, i.e. the objective function For equation 8: wherein, represents the number of channels of the array; is a multi-path time delay structure corresponding to the acoustic observation data; The parameters to be inverted are estimated using a genetic algorithm , and and then substituted into equation 7 using the mean of the synthetic dataset and the first three mode vectors of the empirical orthogonal decomposition of the synthetic dataset , and to achieve the inversion of the sound speed profile, i.e. finding the sound speed profile in the synthetic dataset that has the lowest fitness.

2. The method of claim 1, wherein: Empirical sound speed formulae are used to generate a number of sound speed profiles as the base sample.

3. The method of claim 1, wherein: The sound speed perturbation model is: the vertical movement of seawater lifts the isotherm to the ocean surface, thereby thinning the mixed layer thickness sound speed perturbation A model and the sound speed perturbation B model without isotherm lifting when strong temperature and salinity anomalies appear in the upper layer of the ocean.

4. The method of claim 3, wherein: The sound speed perturbation A model: when perturbation A occurs, the depth structure of the characteristic center vertical displacement is modeled as: wherein, is the depth of the characteristic center, is the vertical displacement perturbation, is the vertical scale of the set perturbation, is a normalization factor, the occurrence of which generates a superimposed feature of the upper ocean; The sound speed perturbation is calculated using a nonlinear equation: where, is the adiabatic sound speed gradient, the second term in brackets is determined as the sound speed at depth after the adiabatic uplift to the new depth of the stationary block.

5. The method of claim 3, wherein: The sound speed perturbation B model: when perturbation B occurs, the depth structure of the characteristic center sound speed is modeled as: wherein depth of the feature center , is a surface speed perturbation value, is a mixed layer thickness of the perturbation center.

6. The method of claim 1, wherein: the reference sound speed profile data as input to the sound speed perturbation model, resulting in a perturbed data set when the surface sound speed perturbation is set range (-1, 1) m / s, interval 0.1 m / s, vertical displacement perturbation range (-50, 90) m, interval 1 m.

7. The method of claim 1, wherein: The step 3 uses genetic algorithm to inverse the parameters to be inverted , and The specific process of inverse estimation is as follows: (a) Set the fitness function threshold to and maximum number of iterations Parameters for inversion , and Each parameter is pre-encoded separately, using binary code. Each seed parameter corresponds to a unique binary code, depending on the set initial number of seeds. ,generate The parameters to be inverted are used as the initial seed. , and The precoding results are respectively used , and This means that each initial seed corresponds to a unique set of parameters. ; (b) Next, the fitness function of the formulation 8 is used to calculate each initial seed, i.e. each set of parameters selected The corresponding fitness value is taken as the initial fitness value, and let ; (c) thereafter according to each selected set of parameters corresponding fitness values, the initial seeds are selected, crossed and mutated to produce the next generation of seeds; (d) Calculate the fitness of the next generation of seeds according to fitness function formula 8, if the fitness function does not reach the set threshold or the number of iterations does not reach the number of iterations defined by the algorithm , let return to step (c) until or .

8. The method of claim 1, wherein: the coefficients , and the range of variation of the coefficients , and the maximum and minimum values of the coefficients , , , , , are the coefficients calculated for the first three modal basis vectors , and the range of variation of the coefficients.

9. The method of claim 1, wherein: The tolerance factor The value of the variable is 1.

10. An electronic device, comprising: A data migration method including a processor and a memory, the processor being configured to implement the steps of the sound speed profile inversion method of any one of claims 1 to 9 when executing a computer program stored in the memory.

Citation Information

Patent Citations

  • Sound velocity profile inversion method based on empirical orthogonal function method

    CN113218493A

  • Sound velocity profile inversion method based on inverted multi-beam echometer

    US20210231800A1