Underwater vehicle positioning method based on fitting sound velocity profile

Abstract

The embodiment of the invention provides an underwater vehicle positioning method based on a fitting sound velocity profile, and the method comprises the steps: forming a long baseline positioning system through a plurality of observation stations disposed underwater, and enabling an underwater vehicle to be located in a measurement range of the long baseline positioning system; an acoustic signal is transmitted to the observation station through the underwater vehicle, time delay measurement data is acquired, and the time delay measurement data refers to the difference value between the time when the observation station receives the acoustic signal and the time when the underwater vehicle transmits the acoustic signal; constructing a long baseline positioning system model based on time delay intersection; and performing iterative calculation on the long baseline positioning system model according to the time delay measurement data to determine the position of the underwater vehicle. According to the method, the G-N algorithm is improved, optimal position parameter estimation can be obtained through iteration, meanwhile, the proposed Gaussian kernel B spline algorithm can better fit the bending part of the sound velocity curve through parameter optimization, and a high-precision continuous sound velocity profile can be obtained. According to the technical scheme, the positioning precision of the underwater vehicle can be improved, and rapid development of the underwater Internet of Things can be promoted.

Description

Technical Field

[0001] The present invention relates to the technical field of underwater Internet of Things, and particularly relates to a positioning method for an underwater vehicle based on a fitted sound velocity profile. Background Art

[0002] The Internet of Underwater Things (IoUT), as an emerging Internet of Things, plays an important role in environmental monitoring, resource development, and information communication. The positioning of an underwater vehicle is the core of the IoUT network, and the positioning accuracy affects the network performance, coverage effect, connection security, etc. of the IoUT. With the development of ocean development activities, high-precision underwater positioning has become the main research direction, with important academic value and application value.

[0003] Electromagnetic waves attenuate extremely fast in water, and the global satellite navigation system cannot be directly used for the navigation and positioning of underwater vehicles. Sound waves can propagate over long distances underwater, and the navigation and positioning in the IoUT network rely on a positioning system based on underwater sound, such as an ultra-short baseline (USBL) positioning system, a short baseline (SBL) positioning system, and a long baseline (LBL) positioning system. The USBL positioning system realizes navigation and positioning by measuring the distance and angle from the underwater vehicle to the USBL system. When the underwater vehicle is at a long distance, it cannot meet the requirements of high-precision positioning. The SBL system and the LBL system have similar positioning principles, but since they are usually arranged at the bottom of a ship and have a short baseline, the long-distance positioning accuracy is low. The LBL system is usually arranged on the seabed. Since its baseline length is in the range of hundreds of meters to thousands of meters, it can achieve high-precision positioning within a large measurement range, and thus is widely used in the positioning of underwater vehicles.

[0004] During the positioning process, the propagation time of the sound signal from the underwater vehicle to the LBL system can be measured with high precision. The straight-line distance from the underwater vehicle to the LBL system can be obtained by multiplying the propagation time by the underwater sound speed. Further, by using a non-linear iterative method to solve the LBL positioning model, the optimal estimation of the position of the underwater vehicle can be obtained.

[0005] However, due to the inhomogeneity of the seawater medium, the underwater sound speed continuously changes, and the propagation path of the sound signal is a curve. In practical applications, the true straight-line distance from the underwater vehicle to the LBL system cannot be simply obtained by multiplying the sound signal propagation time by the underwater sound speed. Therefore, the spatio-temporal variation of the sound speed structure is the key factor affecting the positioning performance of the LBL system. Currently, according to different processing ideas, the underwater sound speed processing methods can be mainly divided into four types: constant sound speed method, effective sound speed method, equivalent sound speed profile method, and sound signal tracking method.

[0006] The constant sound speed method simplifies the underwater sound speed to a constant value. According to the sound speed profile, the average sound speed is calculated as the underwater sound speed, and the straight-line distance from the submersible to the LBL system is calculated based on the propagation time of the sound signal, and finally the position estimation of the submersible is obtained. However, this sound speed calculation method is relatively rough, and the sound speed error in the positioning model is large, resulting in low positioning accuracy. In addition, the Kalman filtering algorithm based on uncertain least squares proposed by assuming that the underwater sound speed is an unknown constant estimates the position of the submersible and the underwater sound speed at the same time. This method is still based on the assumption that the sound speed is constant, so the accuracy improvement is limited.

[0007] The effective sound speed is defined as the ratio of the straight-line distance between the sound signal emission point and the reception point to the propagation time. This method constructs an effective sound speed table before underwater positioning, finds the effective sound speed according to the time delay data during positioning, and then obtains the distance from the positioning system to the submersible, and finally calculates the position of the submersible. Although the effective sound speed method has high positioning accuracy and fast calculation speed, it requires the depth data of the submersible, so its application range in IoUT is limited.

[0008] The equivalent sound speed profile method approximates the actual complex sound speed profile as a sound speed profile with a constant gradient, which simplifies the calculation process. The introduction of the equivalent sound speed profile simplifies the underwater sound signal propagation trajectory to an arc, improving the calculation speed of the underwater submersible position. Although the equivalent sound speed profile method improves the calculation efficiency, it loses accuracy when processing the sound speed profile, and at the same time, it also requires measuring the depth of the submersible.

[0009] The sound signal tracking method is based on ray acoustics theory and combines the sound speed profile to realize the tracking of the sound signal propagation trajectory. The sound signal tracking method has a rigorous theory and can achieve high-precision underwater positioning in IoUT. Although this method does not require the known incident angle of the sound signal, its calculation efficiency is slow, time-consuming, and when the incident angle of the sound signal is large, the positioning accuracy is poor. In addition, since the actual sound speed profile is a continuous function related to depth, simplifying the sound speed within the layer to a constant gradient change does not conform to the actual sound speed change, so the positioning accuracy of IoUT will be lost.

[0010] In summary, the existing technologies cannot construct a high-precision positioning model based on the propagation time of sound signals, and it is difficult to obtain the precise positioning of the submersible. Therefore, how to accurately position the submersible is a problem that needs to be solved. Summary of the Invention

[0011] An embodiment of the present invention provides an LBL positioning method and system for underwater Internet of Things that fits the sound speed profile to construct a high-precision positioning model according to the propagation time of sound signals.

[0012] To achieve the above object, an embodiment of the present invention provides a positioning method for a submersible based on a fitted sound velocity profile, including: forming a long baseline positioning system by several measuring stations arranged underwater, and positioning the submersible within the measurement range of the long baseline positioning system; transmitting acoustic signals from the submersible to the measuring stations, and obtaining time delay measurement data, where the time delay measurement data refers to the difference between the time when the measuring station receives the acoustic signal and the time when the submersible transmits the acoustic signal; constructing a long baseline positioning system model based on time delay intersection; and performing iterative calculations on the long baseline positioning system model according to the time delay measurement data to determine the position of the submersible.

[0013] Further, the long baseline positioning system model based on time delay intersection is specifically:

[0014] t i = f(c(z), S i , X), i = 1,..., n;

[0015] where c(z) is the underwater sound velocity profile, z represents the underwater depth, S i = [x i , y i , z i T , i = 1,..., n are the positions of the measuring stations, X = [x, y, z] T is the position of the submersible, and the time delay measurement data is t i = [t1, t2,..., t n T .

[0016] Further, the iterative calculation of the long baseline positioning system model according to the time delay measurement data to determine the position of the submersible includes:

[0017] Calculating an initial value X0 of the position of the submersible according to the traditional positioning model of the prior art;

[0018] Determining an initial value k = 1 of the number of iterations k, the maximum value k of the number of iterations Xmax and the iteration termination threshold ε X ;

[0019] S1. Determining the current value of the effective sound velocity c d ;

[0020] S2. Calculating the current value of the acoustic signal propagation time;

[0021] S3. Obtaining the current value of the time delay residual Δt according to the time delay measurement data and the current value of the acoustic signal propagation time;

[0022] S4. According to the effective sound velocity c d ​​The current value of and the current value of the time delay residual Δt are used to determine the current value of the distance residual ΔR;

[0023] S5. Let k = k + 1;

[0024] S6. Update the current value X of the submersible position according to the iteration formula and the current value of the distance residual ΔR k ;

[0025] S7. If ||X k - X k-1 || > ε X and k < k Xmax , then return to step S1.

[0026] Furthermore, the effective sound speed calculation formula is:

[0027]

[0028] The following formula is used to obtain the time delay residual Δt:

[0029] Δt = [Δt1,..., Δt n T

[0030] = [t1,..., t n T - [t c1 ,..., t cn T ;

[0031] The following formula is used to obtain the distance residual ΔR:

[0032] ΔR = [c d1 Δt1,..., c dn Δt n T ;

[0033] The iteration formula is specifically:

[0034] X k = X k-1 + λ(J T J) -1 J T ΔR, λ ∈ [0, 1],

[0035] where λ is the step size.

[0036] Furthermore, the calculation method of the sound signal propagation time is:

[0037] Construct the horizontal distance expression of the sound signal propagation trajectory;

[0038] ​​​​Using the bisection method and the Romberg algorithm, solve the horizontal distance expression of the sound signal propagation trajectory to obtain the incident angle of the sound signal.

[0039] Calculate the propagation time of the sound signal according to the incident angle of the sound signal.

[0040] Among them, the horizontal distance expression of the sound signal propagation trajectory is specifically:

[0041]

[0042] In the formula, l is the horizontal distance of the sound signal propagation trajectory, z1 and z2 are the depths of the submarine vehicle emitting the sound signal and the corresponding sound signal receiving station respectively, (x1, y1) is the position coordinate of the submarine vehicle emitting the sound signal, (x2, y2) is the position coordinate of the sound signal receiving station, c0 is the sound speed value at the sound signal emission point, and θ0 is the initial incident angle.

[0043] Furthermore, in order to better fit the curved part of the sound speed curve and improve the positioning accuracy, a continuous sound speed profile can be used to replace the underwater sound speed profile c(z).

[0044] Furthermore, the continuous sound speed profile is obtained by optimizing the Gaussian kernel B-spline σ of the underwater sound speed profile, specifically including:

[0045] Determine the initial value k = 1 of the iteration number k, the maximum value k of the iteration number Bmax and the iteration termination threshold ε B ;

[0046] Divide the depth interval of the underwater sound speed profile equally through m spline nodes, and the spline nodes are T1,..., T M ;

[0047] Set the initial value σ0 of the parameter σ, and obtain the initial value Φ0 of the fitting accuracy through the spline nodes and the initial value σ0 of the parameter σ.

[0048] S11. Update the matrix Ξ;

[0049] S12. Calculate the Jacobian matrix J(σ) with respect to the parameter σ, and then obtain the Jacobian matrix J(τ) with respect to the parameter τ according to the Jacobian matrix J(σ);

[0050] S13. Let k = k + 1;

[0051] S14. Determine the current value τ of the parameter τ according to the Jacobian matrix J(τ) and the optimization process iteration calculation formula k ;

[0052] S15. Update the current value σ of the parameter σ according to the formula τ = lnσ through the current value τ of the parameter τ k ; k ;

[0053] S16. Calculate the continuous sound speed profile according to the matrix Ξ;

[0054] S17. Update the current value Φ of the fitting accuracy according to the current value σ of the parameter σ k ; k ;

[0055] S18. If |Φ k - Φ k-1 | > ε B and k < k Bmax , then return to step S11;

[0056] Determine the final and the corresponding continuous sound speed profile

[0057] After obtaining the continuous sound speed profile , it can be used to replace the underwater sound speed profile c(z) in the foregoing formula.

[0058] Furthermore, the form of the iterative calculation formula for the optimization process is:

[0059] τ k = τ k-1 + λ(J(τ) T J(τ)) -1 J(τ) T e

[0060] = τ k-1 + λ(J(τ) T J(τ)) -1 J(τ) T (I - Ξ[Ξ T Ξ] -1 Ξ T )c, λ ∈ [0, 1];

[0061] The formula for determining the fitting accuracy is:

[0062]

[0063] Among them, the method for determining the residual e is:

[0064] e = c - Ξβ = (I - Ξ[Ξ T Ξ] -1 Ξ T )c;

[0065] The method for determining the Jacobian matrix J(σ) is:

[0066]

[0067] The method for determining the Jacobian matrix J(τ) is as follows:

[0068] J(τ) = σ · J(σ);

[0069] Continuous sound speed profile The calculation method is as follows:

[0070] χ G = ΞΞ + c.

[0071] The above technical solution has the following beneficial effects:

[0072] The technical solution of this application improves the G-N algorithm, calculates the distance residual based on the effective sound speed, adjusts the distance residual calculated by the effective sound speed during the iteration process, and the optimal position parameter estimation can be obtained through iteration; at the same time, the Gaussian kernel B-spline is proposed. By optimizing the parameters, it can better fit the curved part of the sound speed curve. Under the condition of the same nodes, it can greatly improve the fitting accuracy of discrete data and effectively avoid the Ronge phenomenon at the data edge. Therefore, according to the Gaussian kernel B-spline fitting discrete sound speed profiles, high-precision continuous sound speed profiles can be obtained. The above means are all beneficial to improving the positioning accuracy of the submersible, and thus beneficial to promoting the rapid development of the underwater Internet of Things. BRIEF DESCRIPTION OF THE DRAWINGS

[0073] In order to more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the following will briefly introduce the drawings required for the description of the embodiments or the prior art. Obviously, the drawings in the following description are only some embodiments of the present invention. For those of ordinary skill in the art, other drawings can be obtained based on these drawings without creative efforts.

[0074] Figure 1 is a flowchart of a submersible positioning method based on fitting sound speed profiles according to an embodiment of the present invention;

[0075] Figure 2 is a schematic diagram of a long baseline positioning system;

[0076] Figure 3 is a schematic diagram of discrete underwater sound speed profiles in the simulation experiment of the present invention;

[0077] Figure 4 is a schematic diagram of the positions of measurement stations and the trajectory of the submersible in the simulation experiment of the present invention. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0078] Next, the technical solutions in the embodiments of the present invention will be clearly and completely described in conjunction with the accompanying drawings in the embodiments of the present invention. Obviously, the described embodiments are only a part of the embodiments of the present invention, rather than all the embodiments. All other embodiments obtained by those of ordinary skill in the art based on the embodiments of the present invention without creative efforts shall fall within the protection scope of the present invention.

[0079] Next, a specific embodiment is used to detail the derivation and implementation process of this method.

[0080] I. Time-delay intersection positioning model based on acoustic signal tracking

[0081] 1. Traditional LBL positioning model

[0082] The LBL long baseline positioning system in IoUT is fixed on the seabed as shown in Figure 2 . It includes multiple measurement stations. An acoustic signal transmitter is installed outside the underwater submersible. The long baseline positioning system is time-synchronized with the acoustic signal transmitter. During the underwater positioning process, the acoustic signal transmitter X emits an acoustic signal, and the long baseline positioning system receives the acoustic signal. According to the acoustic signal emission time, the one-way propagation time t from the beacon X to the measurement station S i , i = 1,..., n can be obtained. i .

[0083] According to formula (1), the slant range R from the submersible X to the measurement station S i can be calculated. i .

[0084] R i = ct i , (1)

[0085] where c is the underwater sound speed.

[0086] According to the sound speed profile data, the geometric mean is calculated according to formula (2) as the underwater sound speed c:

[0087]

[0088] where z s , z e are the depths of the starting point and the ending point of the acoustic signal respectively.

[0089] Assume that the positions of the measurement stations of the long baseline positioning system in IoUT are S i = [x i , y i , z i T , i = 1,..., n, and the position of the underwater submersible is X = [x, y, z] T ​, the time delay measurement value of the long baseline system is \(t = [t_1,t_2,\cdots,t n T . According to the geometric position relationship between the submersible and the LBL, we have:

[0090]

[0091] According to formulas (1) and (3), the LBL system positioning model is obtained as:

[0092]

[0093] Since the nonlinear degree of the positioning model (4) is relatively high and it is difficult to obtain the analytical solution of the position parameters, an iterative algorithm can be used to obtain the numerical solution of the position parameters. When using the G - N (Gauss - Newton) algorithm to solve, assume that the position of the submersible at the \(k\) - th iteration is \(X k , the iterative formula of model (4) is:

[0094] X k+1 =X k +(λ / 2) k (J T J) -1 J T ΔR, λ∈[0,1]. (5)

[0095] In iteration (5), \(J\) is the Jacobian matrix, which is calculated according to formula (6).

[0096]

[0097] According to the position \(X\) of the submersible k , the calculated distance value \(R c :

[0098] R c =[||S1 - X k ||,\cdots,||S n - X k ||] T . (7)

[0099] In iteration (5), the distance residual ΔR is calculated according to formula (8):

[0100] ΔR = [R1,R 2, \cdots,R n T - [R c1 ,R c2 ,\cdots,R cn T . (8)

[0101] ​​​Although the traditional LBL positioning model does not require the depth of the submersible to be known in advance, the simplification of the speed of sound to a constant value leads to large systematic errors in the distance data, resulting in limited positioning accuracy of the submersible. Therefore, it is necessary to build a time-delay rendezvous positioning model based on the acoustic signal tracking theory to avoid systematic errors caused by inaccurate distance calculation.

[0102] 2. Time-delay rendezvous positioning model

[0103] Since the underwater sound speed varies significantly with time and space, the change in sound speed causes the sound signal to refract along the propagation direction, resulting in a curved propagation path of the sound signal.

[0104] Since the underwater sound speed does not change much in the horizontal direction in a short period of time, this paper assumes that the underwater sound speed is a function related to depth:

[0105] c=c(z). (9)

[0106] The incident angle of the acoustic signal at point A is θ, which is the angle between the acoustic signal and the vertical direction. By intercepting a sufficiently small trajectory segment at point A, the horizontal propagation distance dl can be obtained:

[0107] dl=tanθdz. (10)

[0108] At depth z ι The speed of sound is c ι , the incident angle is θ ι , according to Snell's law of refraction:

[0109]

[0110] Where p is a constant.

[0111] According to formulas (10) and (11), we can get:

[0112]

[0113] Where θ0 is the initial incident angle.

[0114] Assume that the acoustic signal transmitting point (the sender) and receiving point (the receiver) are P1 = [x1, y1, z1] T and P2 = [x2, y2, z2] T By integrating formula (12), we can get the horizontal distance of the acoustic signal propagation trajectory as:

[0115]

[0116] Where z1 and z2 are the depths of P1 and P2 respectively.

[0117] Similarly, the variation in the propagation distance of the acoustic signal, \(dR\), can be obtained as follows:

[0118]

[0119] Integrating the above equation gives the propagation distance of the acoustic signal as:

[0120]

[0121] According to Equation (15), the propagation time of the acoustic signal is:

[0122]

[0123] When the acoustic signal emission point and the receiving points \(P1\) and \(P2\) are determined, according to Equation (13), we have:

[0124]

[0125] The horizontal propagation distance \(l\) of the acoustic signal increases monotonically with \(\theta_0\). Since \(z1\) and \(z2\) are known, the initial incident angle \(\theta_0\) of the acoustic signal is obtained using the bisection method according to Equation (17), and then the propagation time \(t\) of the acoustic signal can be calculated according to Equation (16). The specific algorithm is shown in Table 1 below:

[0126]

[0127] Table 1. Flow of Algorithm 1

[0128] Note that in the art, it is well known that the symbols representing the propagation time of the acoustic signal (also known as the propagation delay) include \(t\) and \(t\) c , where \(t\) is the theoretically derived value, and \(t\) c is the calculated value obtained by substituting the fitted sound speed profile into Equation (16). When calculating the delay residual \(\Delta t\) later, \(t\) c is compared with the delay measurement data.

[0129] Algorithm 1 shows that when the start and end points \(P1\) and \(P2\) of the acoustic signal are determined, the propagation time \(t\) of the acoustic signal between \(P1\) and \(P2\) can be uniquely determined. Therefore, the delay data \(t\) can be written as an implicit function of the positions of two underwater points \(P1\), \(P2\) and the sound speed profile \(c = c(z)\):

[0130] \(t = f(c(z), P1, P2). (18)\)

[0131] The implicit positioning model of the underwater LBL system based on time-delay intersection constructed in this paper is:

[0132] \(t\) i \(= f(c(z), S\) i , X), i = 1,..., n. (19)\)

[0133] The traditional average sound velocity positioning model needs to convert the time delay measurement value into distance data, and then solve the position parameters of the underwater vehicle. When calculating the distance data, the sound velocity error will be introduced, resulting in limited positioning accuracy. The time delay intersection positioning model constructed in this paper directly estimates the position parameters of the underwater vehicle according to the time delay data, without the need to calculate the distance data, thus avoiding the introduction of the sound velocity system error into the model. The numerical iterative algorithm for solving the time delay intersection model will be given below.

[0134] 3. Underwater positioning iterative algorithm

[0135] Since it is necessary to integrate formulas (13) and (16) in the iterative algorithm, the integral solution algorithm is given first. In this paper, the Romberg integration algorithm is used to solve the horizontal distance and time delay data. The Romberg integration formula is also called the successive halving acceleration method. It is a method for accelerating the calculation of integrals constructed on the basis of the trapezoidal formula, Simpson's formula, and Cotes formula. As an extrapolation algorithm, it improves the integral accuracy without increasing the amount of calculation. The Romberg integration formula has a fast convergence speed and high iterative accuracy, and is suitable for solving integrals in a finite interval

[31] .

[0136] Assume that the function to be integrated is g = g(ξ), and the integration interval is [α, β]. Divide the interval φ equally, and the trapezoidal formula can be written as

[32] :

[0137]

[0138] where h is the interval length, h = (β - α) / φ.

[0139] The extrapolation formula is

[33] :

[0140]

[0141] where η = 2, 3, 4 are the Simpson formula, Cotes formula, and Romberg formula respectively.

[0142] Furthermore, formula (21) can be written as:

[0143]

[0144] The Romberg integration sequence can be represented by Table 2:

[0145]

[0146]

[0147] Table 2 Romberg integration sequence

[0148] Formula for solving using Romberg integration algorithm The specific steps are shown in Table 3 as follows:

[0149]

[0150] Table 3. Flow of Algorithm 2

[0151] The IoUT long-baseline time-delay intersection positioning model constructed in this paper is a non-linear implicit model, which can be solved by intelligent optimization algorithms such as genetic algorithm, particle swarm algorithm, simulated annealing algorithm, differential evolution algorithm, etc. However, the intelligent optimization algorithms have a slow solving speed and require manual participation in adjusting the iteration parameters. In the traditional average sound speed positioning model, numerical iteration algorithms are often used to solve the position parameters, such as G-N iteration, L-M iteration (Levenberg-Marquardt). However, the distance residual (8) needs to be calculated in the iteration algorithm, and the residual of the positioning model in this paper is the time-delay residual, which cannot be directly used. Here, the distance residual ΔR can be calculated based on the effective sound speed, and an improved G-N iteration solution is given to solve the position of the underwater vehicle.

[0152] The effective sound speed is defined as the ratio of the straight-line distance between P1 and P2 to the propagation time, and the calculation formula is:

[0153]

[0154] Assume that the position of the underwater vehicle at the k-th iteration is X k , and according to formula (23), the effective sound speed of the underwater vehicle to the LBL system in the IoUT can be calculated as:

[0155]

[0156] According to the position X of the underwater vehicle in the iteration k , the calculated value t of the sound signal propagation time-delay can be obtained c :

[0157] t c =[f(c(z),S1,X k ),...,f(c(z),S n ,X k )] T . (25)

[0158] Furthermore, the time-delay residual Δt can be obtained:

[0159] Δt=[Δt1,...,Δt n T

[0160] =[t1,...,t n T -[t c1 ​​,...,t cn T . (26)

[0161] According to the effective sound speed c d and the time delay residual Δt, the distance residual ΔR can be obtained as follows:

[0162] ΔR = [c d1 Δt1,..., c dn Δt n T . (27)

[0163] The improved G-N iteration formula is as follows:

[0164] X k = X k-1 + λ(J T J) -1 J T ΔR, λ ∈ [0, 1]. (28)

[0165] Compared with the intelligent optimization algorithm, the improved G-N iteration algorithm takes less time while ensuring the accuracy. The specific steps are shown in Table 4:

[0166]

[0167]

[0168] Table 4. Flow of Algorithm 3

[0169] Remark 1: Although the calculation method of the above iteration correction amount is the same as the iteration form of the traditional average sound speed positioning model in form, there are great differences in the calculation method of the distance residual. The average sound speed positioning model assumes that the underwater sound speed is a constant value, and the error is large when calculating the distance residual; while the improved G-N algorithm calculates the distance residual based on the effective sound speed, and adjusts the effective sound speed to calculate the distance residual during the iteration process, and the optimal position parameter estimation can be obtained through iteration.

[0170] Remark 2: In the iteration formula of the improved G-N algorithm, the effective sound speed is used as the sound speed. Due to a small deviation between the position of the submersible and the true value during the iteration, this effective sound speed is actually an approximation of the true value of the effective sound speed. Therefore, there is a small deviation in the iteration correction amount, and the calculation of the iteration step size will lose some accuracy, resulting in a relatively large number of iterations. However, the iteration of the positioning model will finally converge to the position parameter corresponding to the minimum residual.

[0171] Remark 3: Compared with the in-layer constant gradient ray tracing, the positioning method proposed in this paper is based on the ray tracing of the continuous sound speed profile. Therefore, the positioning accuracy is higher when the ray incident angle is larger.

[0172] ​​According to the positioning process of Algorithm 3, in the time-delay intersection positioning model, the key to achieving high-precision positioning is to obtain a continuous high-precision sound speed profile. In the actual IoUT navigation and positioning experiment, the sound speed measurement device samples at equal intervals to obtain a discrete sound speed profile, which cannot be directly used in the positioning method of this paper. When the accuracy of the sound speed measurement device is fixed, improving the fitting accuracy of the discrete sound speed profile is an important way to improve the positioning accuracy. Therefore, the method of using Gaussian kernel B-spline to fit the discrete sound speed profile data can be used to obtain a high-precision continuous sound speed profile.

[0173] II. Gaussian Kernel B-Spline Fitting

[0174] To obtain a continuous sound speed profile, polynomials, B-splines, etc. can be used to fit the measured sound speed profile. When using polynomials to fit discrete sound speed profile data, it is assumed that it satisfies the polynomial law. However, in actual measurements, the underwater sound speed profile has local characteristics at different times and spaces, making it difficult to accurately fit with a single polynomial. Moreover, high-order polynomials may exhibit the Runge phenomenon, with severe jitter at the data endpoints, resulting in large fitting errors. B-spline functions can reduce the Runge phenomenon, and the fitted curve is continuous. Therefore, they are often used for discrete data fitting. The commonly used B-spline function is the cubic standard B-spline. Below, the cubic standard B-spline is briefly introduced, and the Gaussian kernel B-spline is proposed to address its deficiencies to improve the fitting accuracy of discrete sound speed data.

[0175] 1. Cubic Standard B-Spline Fitting

[0176] Assume that the depth interval of the sound speed profile is z ∈ [α, β], and the depth interval is equally divided. Assume the equal division is:

[0177] Λ: α ≤ T1 < T2 <... < T M ≤ β, (29)

[0178] where T1,..., T M are the knots.

[0179] The expression of the cubic standard B-spline basis function is as follows:

[0180]

[0181] Let ζ = (z - T j ) / ΔT, and substituting it into (30) gives:

[0182]

[0183] where ΔT is the knot spacing.

[0184] Denote the curve constructed by the cubic standard B-spline function as S C (z), and its expression is:

[0185]

[0186] where is the spline coefficient.

[0187] Let

[0188]

[0189] The curve S C (z) can be written as:

[0190] S C (z) = B C (z)β C .(34)

[0191] Assume that the discrete sound speed data is:

[0192] c ι = c(z ι ), ι = 1, 2,..., m, α = z1 < z2... < z m = β. (35)

[0193] Using a cubic standard B-spline to fit the sound speed data c(t ι ), we can obtain:

[0194]

[0195] Let

[0196] c = [c(z1), c(z2),..., c(z n )] T

[0197] e C = [e C (z1), e C (z2), …, e C (z n )] T .(37)

[0198] Write the cubic standard B-spline basis function at each depth in matrix form:

[0199]

[0200] Equation (36) can be written in vector form as:

[0201] c = S C + e C = Ξ C β C + e C . (39)

[0202] According to the least squares method, the spline coefficients can be calculated as follows:

[0203] β C = ((Ξ C )) T Ξ C ) -1 (Ξ C ) T c(40)

[0204] After 3 - time standard B - spline fitting, the obtained continuous sound - speed profile is:

[0205] χ C (z) = Ξ C (z)((Ξ C )) T Ξ C ) -1 (Ξ C ) T c.(41)

[0206] When the sound - speed data changes greatly, the fitting error of the 3 - time standard B - spline is large. If the fitting accuracy is to be improved, only by increasing the number of nodes M. Therefore, this paper proposes the Gaussian kernel B - spline. When the number of nodes is fixed, by optimizing σ, the representation error of the spline is reduced, and the fitting accuracy of the sound - speed profile is improved.

[0207] 2. Gaussian kernel B - spline fitting

[0208] The basis function of the Gaussian kernel B - spline is, in form, the probability density function of a normal distribution with an expected value of 0 and a variance of σ 2 The basis function B G (ζ,σ) is as follows:

[0209]

[0210] Its derivative is:

[0211]

[0212] The numerical table of the Gaussian kernel B - spline is as follows:

[0213]

[0214] Table 5. Numerical table of the Gaussian kernel B - spline (σ = 1)

[0215] Compared with the 3rd - order standard B - spline, the Gaussian kernel B - spline can reduce the representation error and improve the fitting accuracy by adjusting the parameter σ. For data with relatively stable changes, when σ takes a smaller value, the Gaussian kernel B - spline basis function is relatively 'flat', and it can better fit the data with stable changes; for data with relatively drastic changes, when σ takes a larger value, the Gaussian kernel B - spline basis function is relatively'steep', and it can better fit the data near the feature points. Under the requirement of the same accuracy, the Gaussian kernel B - spline requires fewer knot numbers, saves the number of parameters to be estimated, improves the modeling efficiency, and reduces the computational complexity.

[0216] Construct a curve S using the Gaussian kernel B - spline basis function G (z, σ), and its expression is:

[0217]

[0218] Let:

[0219]

[0220] The curve S G (z) can be written as:

[0221] S G (z, σ) = B G β G . (46)

[0222] Using the Gaussian kernel B - spline to fit the sound - speed data c(z ι ), we can obtain:

[0223]

[0224] Let:

[0225]

[0226] Similarly, according to the least - squares method, the coefficients of the Gaussian kernel B - spline can be obtained:

[0227] β G = ((Ξ G ) T Ξ G ) -1 (Ξ G ) T c. (49)

[0228] According to the measured discrete underwater sound - speed profile, after fitting with the Gaussian kernel B - spline, the obtained continuous sound - speed profile is:

[0229] χ G (z, σ) = Ξ G (z, σ)((Ξ G ) T ΞG ) -1 (Ξ G ) T c.(50)

[0230] Substituting the continuous sound speed data into the time-delay intersection positioning model (19), the LBL time-delay positioning model based on the Gaussian kernel B-spline fitting sound speed profile can be obtained as follows:

[0231] t = f(χ G (z, σ), S i , X), i = 1,..., n. (51)

[0232] When the number of nodes is fixed, the value of σ seriously affects the fitting accuracy of the sound speed profile, thereby restricting the positioning accuracy of the time-delay intersection model (51). Therefore, how to improve the fitting accuracy of the discrete sound speed profile by optimizing σ is another key content of this study.

[0233] 3. Gaussian kernel B-spline σ optimization method

[0234] For the convenience of representation, in the following, B, B, Ξ, S, β, V are used to replace B G , B G , Ξ G , S G , β G , χ G . The parameter estimation of the Gaussian kernel B-spline fitting model (50) can be reduced to a non-linear optimization problem, that is, to find the parameters β, σ such that:

[0235]

[0236] where Q(β, σ) = ||c - Ξβ|| 2 =(c - Ξβ) T (c - Ξβ).

[0237] The optimization problem (52) is in an M + 1-dimensional space. Since the spline coefficient β is a parameter related to σ, in order to minimize Q(β, σ), it is desired to reduce the dimension of the optimization problem (52) to a 1-dimensional space with only the parameter σ.

[0238] Denote the projection matrix of the space generated by the column vectors of Ξ as P(σ), and P(σ) is also an idempotent matrix:

[0239] P(σ) = Ξ[Ξ T Ξ] -1 Ξ T . (53)

[0240] Let:

[0241] H(σ) = I - P. (54)

[0242] Then we have:

[0243] HΞ = 0. (55)

[0244] Consider the following non - linear optimization problem to find such that:

[0245]

[0246] where:

[0247] Π(σ)= ||Hc|| 2 = c T Hc. (57)

[0248] The optimization problem (56) only involves the parameter σ, so it is an optimization problem in a one - dimensional space.

[0249] Theorems 1 and 2 prove that the optimization problems (52) and (56) are equivalent.

[0250] Theorem 1 (Equal minimum values)

[0251]

[0252] Theorem 2 (Equivalent optimal sets)

[0253] Suppose makes Q(β,σ) reach the minimum, then must make Π(σ) reach the minimum; conversely, if makes Π(σ) reach the minimum, let then must make Q(β,σ) reach the minimum.

[0254] From Theorems 1 and 2, to make Q(β,σ) reach the minimum, we only need to consider the non - linear optimization problem (56) and optimize σ.

[0255] Since (56) is a constrained optimization problem, for computational convenience, it is transformed into an unconstrained optimization problem.

[0256] Theorem 3: In the Gaussian kernel B - spline, σ satisfies σ > 0. Let τ = lnσ, then τ ∈ R, and τ and σ are in one - to - one correspondence.

[0257] According to Theorem 3, the constrained optimization problem is transformed into the following unconstrained optimization problem:

[0258]

[0259] To solve the problem (56) using non - linear optimization methods, we need to calculate the Jacobian matrix J(τ) of V with respect to τ.

[0260] Denote the generalized inverse of matrix Ξ as Ξ + , the solution of model (52) can be written as

[0261] The fitted sound speed can be obtained as:

[0262] V = ΞΞ + c.(60)

[0263] The Jacobian matrix J(σ) of V with respect to σ can be written as:

[0264]

[0265] According to the definition, the partial derivative of the Gaussian kernel B-spline basis function with respect to σ is:

[0266]

[0267] Furthermore, it can be obtained that:

[0268]

[0269] The following gives the process of

[0270] Lemma 1: Denote P Ξ = ΞΞ + , Ξ P = Ξ + Ξ, Ξ P ⊥ = I - Ξ + Ξ, there is:

[0271] D( Ξ P) = Ξ + D(Ξ) Ξ P ⊥ +(Ξ + D(Ξ) Ξ P ⊥ ) T , (64)

[0272] where D(A) represents the partial derivative of matrix A with respect to σ.

[0273] Theorem 4: The partial derivative of matrix Ξ + with respect to σ is:

[0274]

[0275] According to formulas (61) and (65), the Jacobian matrix J(σ) of V with respect to σ is:

[0276]

[0277] According to Theorems 3 and 4, the Jacobian matrix J(τ) of V with respect to τ can be obtained as follows:

[0278] J(τ) = σ · J(σ). (67)

[0279] According to formula (49), the residual e can be obtained as follows:

[0280] e = c - Ξβ = (I - Ξ[Ξ T Ξ] -1 Ξ T )c. (68)

[0281] According to the G - N iteration, the iteration formula in the optimization process can be obtained as follows:

[0282] τ k = τ k-1 + λ(J(τ) T J(τ)) -1 J(τ) T e

[0283] = τ k-1 + λ(J(τ) T J(τ)) -1 J(τ) T (I - Ξ[Ξ T Ξ] -1 Ξ T )c, λ ∈ [0, 1]. (69) The data fitting accuracy is calculated according to the root mean square error:

[0284]

[0285] By optimizing the parameter σ in the Gaussian kernel B - spline, the fitting error of the discrete sound speed data can be reduced. The specific process is shown in Table 6.

[0286]

[0287] Table 6. Flow of Algorithm 4

[0288] In summary, in the IoUT navigation and positioning experiment, first, the discrete sound speed data is fitted based on the Gaussian kernel B - spline to obtain a continuous sound speed profile (50), which is substituted into the time - delay intersection positioning algorithm (51), and the position parameters of the submersible are obtained through the improved G - N iteration (28).

[0289] III. Simulation Experiments

[0290] To verify the effectiveness of the time-delay intersection positioning algorithm based on Gaussian kernel B-spline fitting of the sound speed profile proposed in this paper, an IoUT submersible positioning scenario is constructed for numerical simulation verification. The time-delay intersection positioning algorithm proposed in this paper first needs to fit the discrete sound speed profile and then solve the position parameters of the submersible. Therefore, this simulation is divided into two parts. The first part is the fitting simulation of the discrete sound speed profile to obtain a high-precision continuous sound speed profile. The second part is to obtain high-precision submersible position parameters based on the continuous sound speed profile.

[0291] 1. Fitting simulation of discrete sound speed

[0292] Assume that the random error of sound speed measurement is σ c = 0.3 m / s, and the underwater discrete sound speed profile is as Figure 3 shown. When fitting the discrete sound speed profile, the 10th-degree polynomial (F1), the 3rd-degree standard B-spline (F2), and the Gaussian kernel B-spline (F3) proposed in this paper are used for fitting respectively. By comparing the fitting results of the three methods, the accuracy of the Gaussian kernel B-spline fitting method proposed in this paper is illustrated.

[0293] Set the number of spline nodes and select nodes at equal intervals in the depth range. Using the above three methods to fit the discrete sound speed profile, it can be seen from the results that the fitting effects of F1 and F2 are very poor at the curve bending positions, seriously deviating from the sound speed data, and there are relatively serious Runge phenomena at both ends of the data. The fitting effect of F3 is better, relatively fitting the sound speed data, and the Runge phenomenon at both ends of the data is relatively slight.

[0294] After fitting the residuals of the three methods, it can be seen that the fitting error of F2 is the largest, and the error is very large at the positions where the sound speed profile bends, that is, around depths of 35 m, 44 m, and 56 m. The fitting error of F1 is slightly better than that of F2, but the fitting accuracy is also poor at the bending positions. The fitting method F3 proposed in this paper has the smallest fitting error and the fitting residuals are relatively uniform, indicating that the fitting accuracy is higher at the bending positions.

[0295] The fitting accuracy is calculated according to the formula, and the fitting accuracies of the above three methods are shown in Table 7. In the case of 5 nodes, the fitting accuracy of the 3rd-degree standard B-spline is the worst. When using the Gaussian kernel B-spline to fit the sound speed profile, the parameters need to be optimized. Set the initial value of the Gaussian kernel B-spline parameters as, and optimize the parameters according to Algorithm 4 to obtain the optimal value. The fitting accuracy of the Gaussian kernel B-spline is improved from the initial 0.3309 to 0.3085.

[0296]

[0297] Table 7. Fitting accuracy table of discrete sound speed profile

[0298] According to the above simulation results, it can be obtained that the accuracy of traditional polynomial and cubic standard B-spline fitting of discrete sound velocity data is limited. If the fitting accuracy is to be improved, the degree of the polynomial or the number of nodes needs to be increased. The Gaussian kernel B-spline proposed in this paper can better fit the curved part of the sound velocity curve by optimizing parameters. Under the condition of the same number of nodes, it can greatly improve the fitting accuracy of discrete data and effectively avoid the Ronge phenomenon at the data edges. Therefore, according to the Gaussian kernel B-spline fitting of the discrete sound velocity profile, a high-precision continuous sound velocity profile can be obtained.

[0299] 2. Simulation of AUV positioning

[0300] To verify the effectiveness of the underwater time-delay positioning method proposed in this paper, a numerical simulation of the positioning of an offshore IoUT AUV can be constructed. As Figure 4 shown, five long baseline stations are arranged on the seabed in the IoUT, and the distance between the stations is about 500 m. The AUV moves within the measurement range of the long baseline system. The time synchronization between the AUV and the long baseline system is achieved. During the positioning experiment, the AUV emits acoustic signals, and the long baseline system receives the acoustic signals. The true time delay can be obtained according to the emission time of the acoustic signals.

[0301] Assume that the position parameters of the AUV at time are. Random errors are added to the true time delay and the station positions. Assume that the random time-delay error is, and the station site error is. The random errors are added to the true values to obtain the simulation measurement data.

[0302] The above simulation is based on the Gaussian kernel B-spline fitting of the discrete sound velocity profile to obtain a continuous sound velocity profile. Based on the continuous sound velocity profile, the time-delay intersection positioning model proposed in this paper is used to calculate the position parameters of the AUV. In addition, the traditional average sound velocity positioning model and the in-layer constant gradient ray-tracing positioning model are also used as comparisons in the simulation. By comparing the positioning accuracies of the three models, the accuracy of the positioning model proposed in this paper is illustrated. Denote the average sound velocity positioning model, the in-layer constant gradient ray-tracing positioning model, and the time-delay intersection positioning model of this paper as M1, M2, and M3 respectively, where the initial iteration values of models M2 and M3 are the calculation results of the traditional model M1.

[0303] The positioning model M2 is based on constant gradient sound velocity tracking, and the sound velocity varies linearly within each layer; the model M3 of this paper realizes ray tracing based on the continuous sound velocity profile through fitting the discrete sound velocity.

[0304] Due to the inaccuracy of the average sound velocity in the traditional model M1, the positioning errors at different positions vary greatly, and the overall error is the largest. The overall error of model M2 is relatively small, but due to the too large acoustic ray incident angle when the AUV is too far from the station position, the positioning error is large at some positions. The positioning error of model M3 is the smallest, and the overall error is relatively stable.

[0305] Let the positioning accuracies in the x, y, and z directions be Δ x , Δ y , Δ z , Δ x , Δ y , Δ z Calculated according to the mean absolute error (MAE), the position accuracy Δ X Calculated according to formula (71):

[0306]

[0307] The error statistical results are shown in Table 8. It can be seen that compared with models M1 and M2, the proposed model M3 in this paper has the highest positioning accuracy. In addition, compared with the x and y directions, the error in the depth direction z is larger, which is due to the larger horizontal range and smaller depth of the LBL measurement system, resulting in a poorer observation geometry in the depth direction.

[0308]

[0309] Table 8. Positioning accuracies of various methods

[0310] According to the above simulation results, it can be obtained that: the traditional average sound speed positioning model simplifies the underwater sound speed to a constant value, resulting in a large positioning error. Although the constant gradient ray tracing positioning model improves the positioning accuracy, it simplifies the underwater sound speed of each layer to a linear change, which does not conform to the actual change of the underwater sound speed. Therefore, the improvement of the positioning accuracy is limited, and when the ray incident angle is large, the positioning error is large. The time-delay intersection positioning model constructed in this paper is based on the continuous sound speed profile, which is more in line with the actual ray propagation situation. Therefore, the positioning accuracy is higher, and because the ray tracing is based on the continuous sound speed profile, the model in this paper can maintain a high positioning accuracy under the condition of a large incident angle.

[0311] According to the foregoing specific embodiments and simulation experiments, the conclusion is drawn that: due to the influence of the complex underwater environment, the underwater sound speed changes continuously, resulting in a curved propagation path of the sound signal, which affects the positioning accuracy of the IoUT network. To reduce the influence of the sound speed change on the underwater positioning accuracy, assuming that the sound speed only changes with depth in the measurement range, this paper constructs a long baseline time-delay intersection positioning model based on the ray tracing of the continuous sound speed profile, and proposes an improved numerical iteration algorithm based on the effective sound speed. Since this positioning model needs to obtain the continuous underwater sound speed profile, this method proposes to fit the discrete sound speed profile with Gaussian kernel B-spline. At the same time, since the fitting accuracy of this method is related to the parameters, this method gives a parameter optimization algorithm to obtain a high-precision continuous sound speed profile by optimizing the parameters.

[0312] The simulation results show that, compared with polynomial fitting and cubic standard B-spline fitting, the Gaussian kernel B-spline proposed by this method can obtain a high-precision continuous sound speed profile by optimizing parameters. Compared with the average sound speed positioning and the in-layer constant gradient ray tracing positioning, the ray tracing positioning model based on Gaussian kernel B-spline proposed by this method has higher accuracy; in addition, in the case of large incident angles, the constant gradient ray tracing method has a large positioning error, while this method can still maintain high accuracy. Therefore, the positioning model proposed by this method can provide a theoretical basis and technical support for the navigation and positioning of underwater vehicles in the IoUT network.

[0313] However, the long baseline system time-delay intersection positioning model proposed by this method is highly dependent on the accuracy of the underwater sound speed. The accuracy of the underwater sound speed is affected by the accuracy of the measurement equipment and the complex underwater environment. Since the accuracy of the measurement equipment can reach a very high level, the influence of the complex environment is the main factor affecting the sound speed accuracy. Even at the same location, the sound speed can change significantly within a day. Therefore, in order to reduce the positioning error of the underwater vehicle caused by the change of the sound speed over time, it is necessary to measure the sound speed within a short period before the positioning experiment.

[0314] In order to enable any person skilled in the art to implement or use the present invention, the above-described disclosed embodiments have been described. For those skilled in the art; various modification methods of these embodiments are obvious, and the general principles defined in this application can also be applied to other embodiments without departing from the spirit and protection scope of the present disclosure. Therefore, the present disclosure is not limited to the embodiments given in this application, but is consistent with the broadest scope of the principles and novel features disclosed in this application.

[0315] The above-described specific implementation manners further elaborate on the purpose, technical solutions, and beneficial effects of the present invention. It should be understood that the above is only the specific implementation manner of the present invention and is not used to limit the protection scope of the present invention. Any modifications, equivalent replacements, improvements, etc. made within the spirit and principle of the present invention shall be included in the protection scope of the present invention.

Claims

1. A positioning method for a submersible based on a fitted sound velocity profile, characterized in that Including: A long baseline positioning system is composed of several underwater stations, and the submersible is located within the measurement range of the long baseline positioning system; The submersible emits acoustic signals to the stations and obtains time delay measurement data, where the time delay measurement data refers to the difference between the time when the station receives the acoustic signal and the time when the submersible emits the acoustic signal; Construct a long baseline positioning system model based on time delay intersection; Perform iterative calculations on the long baseline positioning system model according to the time delay measurement data to determine the position of the submersible.

2. The method for positioning an underwater vehicle based on a fitted sound velocity profile according to claim 1, wherein The long baseline positioning system model based on time delay intersection is specifically: t i = f(c(z), S i , X), i = 1, ..., n; Among them, c(z) is the underwater sound speed profile, z represents the underwater depth, S i = [x i , y i , z i T , i = 1,..., n are the measurement station positions, X = [x, y, z] T is the position of the submarine vehicle, and the time delay measurement data is t i = [t1, t2,..., t n T .​​ 3. The method for positioning a submersible based on a fitted sound velocity profile according to claim 2, wherein The iterative calculation of the long baseline positioning system model according to the time delay measurement data to determine the position of the submersible includes: Calculate the initial value X0 of the submersible position; Determine the initial value k = 1 of the iteration number k, the maximum value k of the iteration number Xmax and the iteration termination threshold ε X ; S1. Determine the current value of the effective sound speed c d ; S2. Calculate the current value of the acoustic signal propagation time; S3. Obtain the current value of the time delay residual Δt according to the time delay measurement data and the current value of the acoustic signal propagation time; S4. Determine the current value of the range residual ΔR based on the current value of the effective sound speed c d and the current value of the time delay residual Δt. S5. Let k = k + 1; S6. Update the current value X of the position of the submersible according to the iteration formula and the current value of the distance residual ΔR k ; S7. If ||X k -X k-1 || > ε X and k < k Xmax , then return to step S1.

4. The method for positioning an underwater vehicle based on a fitted sound velocity profile according to claim 3, wherein The calculation formula of the effective sound speed is: The following formula is used to obtain the time delay residual Δt: Δt = [Δt1,..., Δt n T ​ =[t1,...,t n T -[t c1 ,...,t cn T ;​​ The following formula is used to obtain the distance residual ΔR: ΔR = [c d1 Δt1,..., c dn Δt n T ;​ The iterative formula is specifically: X k = X k-1 + λ(J T J) -1 J T ΔR, λ ∈ [0, 1], where λ is the step size.

5. The method for positioning an underwater vehicle based on a fitted sound velocity profile according to claim 4, wherein The calculation method of the acoustic signal propagation time is: Construct a horizontal distance expression of the acoustic signal propagation trajectory; Use the bisection method and the Romberg algorithm to solve the horizontal distance expression of the acoustic signal propagation trajectory to obtain the acoustic signal incident angle; Calculate the acoustic signal propagation time according to the acoustic signal incident angle; where the horizontal distance expression of the acoustic signal propagation trajectory is specifically: In the formula, l is the horizontal distance of the acoustic signal propagation trajectory, z1 and z2 are the depths of the submersible emitting the acoustic signal and the corresponding station receiving the acoustic signal respectively, (x1, y1) is the position coordinate of the submersible emitting the acoustic signal, (x2, y2) is the position coordinate of the station receiving the acoustic signal, c0 is the sound speed value at the acoustic signal emission point, and θ0 is the initial incident angle.

6. The method for positioning an underwater vehicle based on a fitted sound velocity profile according to claim 2, wherein, The underwater sound speed profile c(z) is a continuous sound speed profile 7. The method for positioning an underwater vehicle based on a fitted sound velocity profile according to claim 6, wherein The continuous sound speed profile is obtained by optimizing the underwater sound speed profile with a Gaussian kernel B-spline σ, which specifically includes: Determine the initial value k = 1 of the number of iterations k, the maximum value k of the number of iterations Bmax and the iteration termination threshold ε B ; The depth interval of the underwater sound speed profile is equally divided by m spline nodes, and the spline nodes are T1,..., T M ; Set the initial value σ0 of the parameter σ, and obtain the initial value Φ0 of the fitting accuracy through the spline nodes and the initial value σ0 of the parameter σ; S11. Update the matrix Ξ; S12. Calculate the Jacobian matrix J(σ) of the parameter σ, and then obtain the Jacobian matrix J(τ) of the parameter τ according to the Jacobian matrix J(σ); S13. Let k = k + 1; S14. Determine the current value τ of the parameter τ according to the Jacobian matrix J(τ) and the iterative calculation formula of the optimization process k ; S15. Update the current value σ of parameter σ according to the current value τ of the parameter τ k k ;​ S16. Calculate the continuous sound speed profile according to the matrix Ξ; S17. Update the current value Φ of the fitting accuracy according to the current value σ of the parameter σ k k ;​ S18. If |Φ k - Φ k-1 | > ε B and k < k Bmax , then return to step S11; Determine the final and the corresponding continuous sound velocity profile 8. The method for positioning an underwater vehicle based on a fitted sound velocity profile according to claim 7, characterized in that The form of the iterative calculation formula of the optimization process is: τ k = τ k-1 + λ(J(τ) T J(τ)) -1 J(τ) T e = τ k-1 + λ(J(τ) T J(τ)) -1 J(τ) T (I - Ξ[Ξ T Ξ] -1 Ξ T )c, λ ∈ [0, 1]; where I is the identity matrix; The determination formula of the fitting accuracy is: where the determination method of the residual e is: e = c - Ξβ = (I - Ξ[Ξ T Ξ] -1 Ξ T )c; The determination method of the Jacobian matrix J(σ) is: Among them, Ξ P = Ξ + Ξ, Ξ P ⊥ = I - Ξ + Ξ,, P Ξ = ΞΞ + , The determination method of the Jacobian matrix J(τ) is: J(τ) = σ · J(σ); Continuous sound velocity profile The calculation method is as follows: χ G = ΞΞ + c.

Images