Geometric random underwater acoustic channel modeling method considering non-equal sound velocity

Abstract

The invention discloses a geometric random underwater acoustic channel modeling method considering non-equal sound velocity. The method comprises the following steps: 1) constructing a shallow sea layered sound velocity profile; 2) multi-path propagation including a direct path, a sea surface reflection path and a seabed reflection path is considered; 3) calculating a propagation angle, a Doppler frequency shift, a propagation distance and a time delay parameter of a path between the receiving end and the transmitting end based on the layered sound velocity profile in combination with a geometric acoustic principle and a Snell law; and 4) establishing a channel impulse response and transmission function, calculating a time self-correlation function and a frequency correlation function at the same time, and performing analysis. According to the method, the non-equal sound velocity profile is introduced into the geometric random underwater acoustic channel model for the first time, the time-frequency fading characteristics of the channel under the equal sound velocity condition and the non-equal sound velocity condition are compared and analyzed, the physical authenticity and calculation precision of the model are improved, and reliable technical support is provided for design and performance evaluation of an underwater communication system.

Description

Technical Field

[0001] This invention relates to the field of channel modeling technology, and in particular to a geometrically random underwater acoustic channel modeling method that considers non-uniform sound speeds. Background Technology

[0002] With the rapid development of sixth-generation (6G) integrated air-space-ground-sea communication networks, underwater communication, as a crucial component for applications such as ocean sensing, marine resource exploration and development, environmental monitoring, and unmanned system collaboration, is experiencing increasing demand for research and engineering deployment. Because electromagnetic waves attenuate rapidly in seawater and have limited transmission distances, and optical links are susceptible to absorption and scattering, sound waves, with their relatively low propagation loss and long transmission distance, have become the primary means of achieving long-distance underwater communication.

[0003] However, the propagation environment of underwater acoustic channels is highly complex and significantly affected by the variability of the seawater medium, exhibiting strong time-frequency fading characteristics. Therefore, establishing a channel model that reflects the actual propagation mechanism is of great significance for the design, performance evaluation, and optimization of underwater communication systems. Especially in shallow sea environments, the speed of sound is influenced by multiple factors such as temperature, salinity, and pressure, often exhibiting a significant non-uniform speed distribution along the depth direction. This causes the propagation ray to bend, thereby altering key characteristics such as the propagation angle, distance, and time delay of the sound wave. Most existing underwater acoustic channel modeling methods are based on the assumption of isotropic profiles, failing to accurately represent the actual situation where the speed of sound varies significantly with depth in shallow seas. This leads to a large deviation between existing models and reality when describing the characteristics of underwater acoustic channels, thus limiting the accurate evaluation of communication system models in complex shallow sea environments. Summary of the Invention

[0004] The purpose of this invention is to provide a geometrically stochastic underwater acoustic channel modeling method that considers non-uniform sound speeds. This method can more realistically describe the propagation characteristics of underwater acoustic channels in shallow sea environments, overcome the shortcomings of existing technologies that ignore the variation of sound speed with depth, and thus improve the accuracy of underwater acoustic communication system model evaluation and design.

[0005] Technical Solution: To achieve the above objectives, this invention proposes a geometrically stochastic underwater acoustic channel modeling method considering non-uniform sound speeds. This method includes the following steps:

[0006] Step S1: Initialize the transceiver parameters and environmental parameters of the underwater acoustic channel model in the shallow sea scenario;

[0007] Step S2: Considering the characteristics of sound speed variation in shallow seas, construct a layered sound speed profile and calculate the propagation distance and time delay;

[0008] Step S3: Based on the constructed layered sound velocity profile and combined with the principles of geometric acoustics and Snell's law, calculate the path parameters between the transmitter and receiver, including propagation angle, Doppler frequency shift, propagation distance, and time delay.

[0009] Step S4: Calculate the large-scale fading parameters of the channel;

[0010] Step S5: Calculate the impulse response and transmission function of the channel between the transmitter and receiver based on the path parameters and large-scale fading parameters;

[0011] Step S6: Calculate the statistical characteristic parameters of the channel based on the channel transfer function.

[0012] Furthermore, step S1 specifically includes:

[0013] Step S101: Set the initial depths of the transmitter and receiver to z respectively. T and z R The horizontal distance is D;

[0014] Step S102: Set the motion parameters of the transmitter and receiver, with motion speeds v respectively. T and v R The azimuth angles of the motion are respectively and

[0015] Step S103: Set shallow sea environment parameters, including seawater depth h and initial sound velocity c0 at the sea surface.

[0016] Furthermore, step S2 specifically includes:

[0017] Step S201: Construct sound velocity profiles with different gradients, and set the sound velocity gradient in each layer to be constant, and the depth of each layer to be the same. The expression for the sound velocity profile is:

[0018]

[0019] Where K represents the total number of water layers, then k = 1, 2, ..., K, c k Let g represent the speed of sound at the lower boundary of the k-th layer. k Δz represents the sound velocity gradient of the k-th layer, and Δz represents the depth of each layer.

[0020] Step S202: Calculate the propagation distance S of the sound wave through the entire water layer according to the definition of arc length:

[0021]

[0022] Among them, R k This represents the radius of curvature of the k-th sound layer. Let represent the incident angle at the lower boundary of the k-th layer and the incident angle at the lower boundary of the (k+1)-th layer. Using Snell's Law calculate;

[0023] Step S203: Calculate the propagation time T of the sound wave through the entire water layer according to the principles of geometric acoustics and Snell's law:

[0024]

[0025] Furthermore, step S3 specifically includes:

[0026] Step S301: Calculate the propagation angle of the direct path:

[0027]

[0028] Where, α L β represents the departure angle of the direct path. L Indicates the angle of arrival of the direct route. This indicates the angle of incidence of the direct path at the launching end. The angle of incidence of the direct path at the receiving end is expressed using Snell's law. Calculate, where c T and c R k represents the speed of sound at the locations of the transmitter and receiver, respectively. T and k R These represent the water layer numbers where the transmitter and receiver are located, respectively, using the formula... and calculate, This represents the floor function;

[0029] Step S302: Calculate the propagation angles corresponding to the DA and UA paths:

[0030]

[0031] in, and These represent the departure angle and arrival angle of the DA path, respectively. and These represent the departure angle and arrival angle of the UA path, respectively; s and These represent the number of surface reflections and bottom reflections experienced by the DA path, respectively. b and b represent the number of surface reflections and bottom reflections experienced by the UA path, respectively; and These represent the incident angles of the DA path at the transmitting and receiving ends, respectively. and These represent the incident angles of the UA path at the transmitting and receiving ends, respectively.

[0032] In the geometrically random underwater acoustic channel model, for the boundary reflection paths, the DA path is the path component that the last reflection occurs at the sea surface and propagates downward to the receiver, while the UA path is the path component that the last reflection occurs at the seabed and propagates upward to the receiver.

[0033] Step S303: Calculate the Doppler frequency shift f of the direct path based on the motion parameters of the transceiver. L :

[0034]

[0035] in, and Representing the maximum Doppler frequency shift at the transmitting and receiving ends, respectively, calculated as follows: and Among them, f c It is the carrier frequency;

[0036] Step S304: Based on the motion parameters and reflection path characteristics of the transceiver, calculate the Doppler frequency shift corresponding to the DA path and the UA path. and

[0037]

[0038] Step S305: Combining the sound velocity profile and geometric relationships, determine the propagation distance d of the direct path. L With delay τ L (τ):

[0039]

[0040] in, and They represent the k locations of the transmitting end, respectively. T The layer and the receiver are located in k R The radius of curvature of the layered sound ray. and They represent the k locations of the transmitting end, respectively. T The layer and the receiver are located in k R The sound velocity gradient of the layer, and These represent the direct path at layer k and the k-th layer where the transmitter is located, respectively. T The layer and the receiver are located in k R The lower boundary angle of incidence of the layer;

[0041] Step S306: Combining the sound velocity profile and geometric relationships, calculate the propagation distances of paths DA and UA, respectively. and Timely delay and

[0042]

[0043]

[0044] in, and These represent the k-th layer and the k-th layer where the DA path is located, respectively. T The layer and the receiver are located in k R The lower boundary incident angle of the layer, and These represent the k-th layer and the k-th layer where the transmitter is located, respectively. T The layer and the receiver are located in k R The incident angle at the lower boundary of the layer.

[0045] Furthermore, step S4 specifically includes:

[0046] Step S401: For the direct path, considering absorption loss and geometric spread loss, the large-scale fading parameters of the direct path are: [The large-scale fading loss L is missing from the original text]. L (d L The expression for ) is:

[0047] L L (d L ) = L a (d L )L s (d L )

[0048] Among them, L a (·) represents absorption loss, L s (·) represents the geometric spread loss, and the absorption loss L a The expression for (d) is:

[0049]

[0050] Where d represents the total propagation distance, α(f) represents the absorption coefficient, and is calculated using the formula: Where f represents the carrier frequency;

[0051] Geometric spread loss L s The expression for (d) is:

[0052] L s (d)=d -β / 2

[0053] Where β represents the expansion coefficient and is usually set to 2;

[0054] Step S402: For the DA and UA paths, considering absorption loss, geometric spread loss, and bottom reflection loss as large-scale fading parameters, the large-scale fading losses of the DA and UA paths are then calculated. and The expressions are as follows:

[0055]

[0056] Among them, L b (·) represents the bottom reflection loss, and its expression is:

[0057]

[0058] in, ρ represents the angle of incidence at the lower boundary of the k-th water layer. k and ρ b Let c represent the densities of seawater and seabed at the lower boundary of the k-th water layer, respectively. b This represents the speed of sound on the seabed at the lower boundary of the k-th water layer.

[0059] Furthermore, step S5 specifically includes:

[0060] Step S501: Divide the channel impulse response from the transmitter to the receiver into direct path and non-direct path. The general expression for the channel impulse response h(t,τ) is:

[0061] h(t,τ)=h L (t,τ)+h N (t,τ)

[0062] Among them, h L (t,τ) represents the channel impulse response of the direct path, h N (t,τ) represents the channel impulse response of a non-direct path, where h L The expression for (t,τ) is:

[0063]

[0064] Among them, a L The channel gain of the direct path is expressed as follows: K R Represents the Rice factor; θ L The phase shift representing the direct path is defined as a random variable that follows a uniform distribution in (0, 2π].

[0065] Step S502, Impulse response h of non-direct path N (t,τ) is divided into two types based on the type of arrival path: DA path and UA path. The two parts are summed according to a certain ratio to obtain the complete expression of the non-direct path impulse response:

[0066] h N (t,τ)=h DA (t,τ)+h UA (t,τ)

[0067] Among them, the impulse responses h of the DA and UA paths DA (t,τ) and h UA The expression for (t,τ) is as follows:

[0068]

[0069] The channel gain expression for the DA path is: The expression for the channel gain of the UA path is: And η DA and η UA Used to balance the contributions of paths DA and UA to the total power, satisfying η DA +η UA =1; N S N represents the maximum number of times the DA path contacts the sea surface. B This indicates the maximum number of times the UA pathway contacts the seabed. Each DA pathway and UA pathway is represented by... and Composed of micro-scattered rays, and Let represent the phase shifts of the DA and UA paths, respectively, and define them as random variables that follow a uniform distribution in (0, 2π].

[0070] Step S503: Perform a Fourier transform on the constructed channel impulse response to obtain the expression for the corresponding channel transfer function:

[0071] H(t,f)=H L (t,f)+H N (t,f)

[0072]

[0073] H N (t,f)=H DA (t,f)+H UA (t,f)

[0074]

[0075] Where H(t,f) represents the overall channel transfer function, H L (t,f) and H N (t,f) represent the channel transfer functions for the direct path and the indirect path, respectively; H DA (t,f) and H UA (t,f) represent the channel transmission functions of the DA path and the UA path, respectively.

[0076] Furthermore, step S6 specifically includes:

[0077] Step S601: Based on the established channel transmission function, calculate the time autocorrelation function of the channel, the expression of which is:

[0078]

[0079]

[0080] Where Δt represents the time interval, I0(κ) is the zeroth-order modified Bessel function of the first kind, and κ represents the concentration level. and express and The average value;

[0081] Step S602: Based on the established channel transmission function, calculate the channel frequency correlation function, the expression of which is:

[0082]

[0083] Where Δf represents the frequency interval.

[0084] Beneficial effects: Compared with the prior art, the technical solution of the present invention has the following beneficial technical effects:

[0085] This invention incorporates the characteristics of sound speed variation with water depth in shallow sea environments into geometric stochastic underwater acoustic channel modeling. This allows for a more accurate reflection of the influence of the vertical sound speed gradient on ray bending, propagation path length, and time delay characteristics, effectively overcoming the problem of significant deviations between existing models based on the assumption of constant sound speed and actual environments. By constructing layered sound speed profiles and combining geometric acoustics principles and Snell's law for multipath analysis, this invention can more accurately determine the propagation angle, Doppler frequency shift, propagation distance, and time delay parameters for each path. Furthermore, it establishes the channel impulse response and transmission function, while simultaneously calculating the time autocorrelation function and frequency correlation function. The proposed method can flexibly adjust model parameters according to different sound speed profiles in shallow sea environments, exhibiting strong adaptability and generalization ability, significantly improving the physical realism and computational accuracy of the model. Attached Figure Description

[0086] Figure 1 This is a flowchart illustrating a geometric stochastic underwater acoustic channel modeling method considering non-uniform sound speeds provided in this embodiment;

[0087] Figure 2 This is a schematic diagram of a geometrically random underwater acoustic channel model considering non-uniform sound speeds provided in this embodiment.

[0088] Figure 3This is a schematic diagram of the channel time autocorrelation function at different times and different sound speed profiles provided in this embodiment;

[0089] Figure 4 This is a schematic diagram of the channel frequency correlation function for different frequencies and different sound speed profiles provided in this embodiment. Detailed Implementation

[0090] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0091] Example 1

[0092] See Figures 1-4 This embodiment provides a geometrically stochastic underwater acoustic channel modeling method considering non-uniform sound speeds, specifically including the following steps:

[0093] Step S1: Initialize the transceiver parameters and environmental parameters of the underwater acoustic channel model in the shallow sea scenario;

[0094] Specifically, in this embodiment, step S1 includes:

[0095] Step S101: Set the initial depths of the transmitter and receiver to z respectively. T and z R The horizontal distance is D;

[0096] Step S102: Set the motion parameters of the transmitter and receiver, with motion speeds v respectively. T and v R The azimuth angles of the motion are respectively and

[0097] Step S103: Set shallow sea environment parameters, including seawater depth h and initial sound velocity c0 at the sea surface.

[0098] Step S2: Considering the characteristics of sound speed variation in shallow seas, construct a layered sound speed profile and calculate the propagation distance and time delay;

[0099] Specifically, in this embodiment, step S2 includes:

[0100] Step S201: Construct sound velocity profiles with different gradients, and set the sound velocity gradient in each layer to be constant, and the depth of each layer to be the same. The expression for the sound velocity profile is:

[0101]

[0102] Where K represents the total number of water layers, then k = 1, 2, ..., K, c k Let g represent the speed of sound at the lower boundary of the k-th layer. k Δz represents the sound velocity gradient of the k-th layer, and Δz represents the depth of each layer.

[0103] Step S202: Calculate the propagation distance S of the sound wave through the entire water layer according to the definition of arc length:

[0104]

[0105] Among them, R k This represents the radius of curvature of the k-th sound layer. Let represent the incident angle at the lower boundary of the k-th layer and the incident angle at the lower boundary of the (k+1)-th layer. Using Snell's Law calculate;

[0106] Step S203: Calculate the propagation time T of the sound wave through the entire water layer according to the principles of geometric acoustics and Snell's law:

[0107]

[0108] Step S3: Based on the constructed layered sound velocity profile and combined with the principles of geometric acoustics and Snell's law, calculate the path parameters between the transmitter and receiver, including propagation angle, Doppler frequency shift, propagation distance, and time delay.

[0109] Specifically, in this embodiment, step S3 includes:

[0110] Step S301: Calculate the propagation angle of the direct path:

[0111]

[0112] Where, α L β represents the departure angle of the direct path. L Indicates the angle of arrival of the direct route. This indicates the angle of incidence of the direct path at the launching end. The angle of incidence of the direct path at the receiving end is expressed using Snell's law. Calculate, where c T and c R k represents the speed of sound at the locations of the transmitter and receiver, respectively. T and k R These represent the water layer numbers where the transmitter and receiver are located, respectively, using the formula... and calculate, This represents the floor function;

[0113] Step S302: Calculate the propagation angles corresponding to the DA and UA paths:

[0114]

[0115] in, and These represent the departure angle and arrival angle of the DA path, respectively. and These represent the departure angle and arrival angle of the UA path, respectively; s and These represent the number of surface reflections and bottom reflections experienced by the DA path, respectively. b and b represent the number of surface reflections and bottom reflections experienced by the UA path, respectively; and These represent the incident angles of the DA path at the transmitting and receiving ends, respectively. and These represent the incident angles of the UA path at the transmitting and receiving ends, respectively.

[0116] In the geometrically random underwater acoustic channel model, for the boundary reflection paths, the DA path is the path component that the last reflection occurs at the sea surface and propagates downward to the receiver, while the UA path is the path component that the last reflection occurs at the seabed and propagates upward to the receiver.

[0117] Step S303: Calculate the Doppler frequency shift f of the direct path based on the motion parameters of the transceiver. L :

[0118]

[0119] in, and Representing the maximum Doppler frequency shift at the transmitting and receiving ends, respectively, calculated as follows: and Among them, f c It is the carrier frequency;

[0120] Step S304: Based on the motion parameters and reflection path characteristics of the transceiver, calculate the Doppler frequency shift corresponding to the DA path and the UA path. and

[0121]

[0122] Step S305: Combining the sound velocity profile and geometric relationships, determine the propagation distance d of the direct path. L With delay τ L (τ):

[0123]

[0124] in, and They represent the k locations of the transmitting end, respectively. T The layer and the receiver are located in k R The radius of curvature of the layered sound ray. and They represent the k locations of the transmitting end, respectively. T The layer and the receiver are located in k R The sound velocity gradient of the layer, and These represent the direct path at layer k and the k-th layer where the transmitter is located, respectively. T The layer and the receiver are located in k R The lower boundary angle of incidence of the layer;

[0125] Step S306: Combining the sound velocity profile and geometric relationships, calculate the propagation distances of paths DA and UA, respectively. and Timely delay and

[0126]

[0127]

[0128] in, and These represent the k-th layer and the k-th layer where the DA path is located, respectively. T The layer and the receiver are located in k R The lower boundary incident angle of the layer, and These represent the k-th layer and the k-th layer where the transmitter is located, respectively. t The layer and the receiver are located in k R The incident angle at the lower boundary of the layer.

[0129] Step S4: Calculate the large-scale fading parameters of the channel;

[0130] Specifically, in this embodiment, step S4 includes:

[0131] Step S401: For the direct path, considering absorption loss and geometric spread loss, the large-scale fading parameters of the direct path are: [The large-scale fading loss L is missing from the original text]. L (d L The expression for ) is:

[0132] L L (d L ) = L a (d L )L a (d L )

[0133] Among them, L a (·) represents absorption loss, L s (·) represents the geometric spread loss, and the absorption loss L a The expression for (d) is:

[0134]

[0135] Where d represents the total propagation distance, α(f) represents the absorption coefficient, and is calculated using the formula: Where f represents the carrier frequency;

[0136] Geometric spread loss L s The expression for (d) is:

[0137] L s (d)=d -β / 2

[0138] Where β represents the expansion coefficient and is usually set to 2;

[0139] Step S402: For the DA and UA paths, considering absorption loss, geometric spread loss, and bottom reflection loss as large-scale fading parameters, the large-scale fading losses of the DA and UA paths are then calculated. and The expressions are as follows:

[0140]

[0141] Among them, L b (·) represents the bottom reflection loss, and its expression is:

[0142]

[0143] in, ρ represents the angle of incidence at the lower boundary of the k-th water layer. k and ρ b Let c represent the densities of seawater and seabed at the lower boundary of the k-th water layer, respectively. b This represents the speed of sound on the seabed at the lower boundary of the k-th water layer.

[0144] Step S5: Calculate the impulse response and transmission function of the channel between the transmitter and receiver based on the path parameters and large-scale fading parameters;

[0145] Specifically, in this embodiment, step S5 includes:

[0146] Step S501: Divide the channel impulse response from the transmitter to the receiver into direct path and non-direct path. The general expression for the channel impulse response h(t,τ) is:

[0147] h(t,τ)=hL (t,τ)+h N (t,τ)

[0148] Among them, h L (t,τ) represents the channel impulse response of the direct path, h N (t,τ) represents the channel impulse response of a non-direct path, where h L The expression for (t,τ) is:

[0149]

[0150] Among them, a L The channel gain of the direct path is expressed as follows: K R Represents the Rice factor; θ L The phase shift representing the direct path is defined as a random variable that follows a uniform distribution in (0, 2π].

[0151] Step S502, Impulse response h of non-direct path N (t,τ) is divided into two types based on the type of arrival path: DA path and UA path. The two parts are summed according to a certain ratio to obtain the complete expression of the non-direct path impulse response:

[0152] h N (t,τ)=h DA (t,τ)+h UA (t,τ)

[0153] Among them, the impulse responses h of the DA and UA paths DA (t,τ) and h UA The expression for (t,τ) is as follows:

[0154]

[0155] The channel gain expression for the DA path is: The expression for the channel gain of the UA path is: And η DA and η UA Used to balance the contributions of paths DA and UA to the total power, satisfying η DA +η UA =1; N S N represents the maximum number of times the DA path contacts the sea surface. B This indicates the maximum number of times the UA pathway contacts the seabed. Each DA pathway and UA pathway is represented by... and Composed of micro-scattered rays, and Let represent the phase shifts of the DA and UA paths, respectively, and define them as random variables that follow a uniform distribution in (0, 2π].

[0156] Step S503: Perform a Fourier transform on the constructed channel impulse response to obtain the expression for the corresponding channel transfer function:

[0157] H(t,f)=H L (t,f)+H N (t,f)

[0158]

[0159] H N (t,f)=H Ds (t,f)+H UA (t,f)

[0160]

[0161] Where H(t,f) represents the overall channel transfer function, H L (t,f) and H N (t,f) represent the channel transfer functions for the direct path and the indirect path, respectively; H Da (t,f) and H UA (t,f) represent the channel transmission functions of the DA path and the UA path, respectively.

[0162] Step S6: Calculate the statistical characteristic parameters of the channel based on the channel transfer function;

[0163] Specifically, in this embodiment, step S6 includes:

[0164] Step S601: Based on the established channel transmission function, calculate the time autocorrelation function of the channel, the expression of which is:

[0165]

[0166] Where Δt represents the time interval, I0(κ) is the zeroth-order modified Bessel function of the first kind, and κ represents the concentration level. and express and The average value;

[0167] Step S602: Based on the established channel transmission function, calculate the channel frequency correlation function, the expression of which is:

[0168]

[0169] Where Δf represents the frequency interval.

[0170] To verify the practical value of the geometric stochastic underwater acoustic channel modeling method considering non-uniform sound speeds provided in this embodiment, two simulation experiments were conducted in this embodiment:

[0171] Experiment 1

[0172] The model is designed with the transmitter and receiver moving at a speed of 3 m / s and an initial distance of 1 km. The azimuth angles of the transmitter and receiver are set accordingly. and The time autocorrelation function of the channel was simulated and analyzed at 180° and 0°, respectively, under different time points and different sound velocity profiles. The results are as follows: Figure 3 As shown, compared with the isotropic profile, the time autocorrelation function decays faster under the non-isotropic profile condition, and the channel's time autocorrelation gradually decreases with the increase of the sound speed gradient, verifying the significant impact of the sound speed gradient on the channel's time correlation characteristics. Furthermore, the value of the time autocorrelation function varies significantly at different times, verifying the effectiveness and rationality of the established model in characterizing time-domain nonstationarity.

[0173] Experiment 2

[0174] The model is designed with the transmitter and receiver moving at a speed of 3 m / s and an initial distance of 1 km. The azimuth angles of the transmitter and receiver are set accordingly. and Simulations were performed on the frequency correlation function of the channel at 180° and 0°, under different sound velocity profiles at t=0, and the results are as follows. Figure 4 As shown, compared with the isotropic profile, the frequency correlation function decays faster under the non-isotropic profile condition, and the frequency correlation of the channel gradually decreases with the increase of the sound speed profile gradient. This verifies that the sound speed gradient has a significant impact on the channel frequency domain correlation, thus further proving the rationality and accuracy of the established model in describing the frequency domain characteristics.

[0175] It should be understood that the sequence number of each step in the above embodiments does not imply the order of execution. The execution order of each process should be determined by its function and internal logic, and should not constitute any limitation on the implementation process of the embodiments of this application.

[0176] The above-described embodiments are only used to illustrate the technical solutions of this application, and are not intended to limit them. Although this application has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features. Such modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of this application, and should all be included within the protection scope of this application.

Claims

1. A method for modeling geometrically stochastic underwater acoustic channels considering non-uniform sound speeds, characterized in that, The method includes the following steps: Step S1: Initialize the transceiver parameters and environmental parameters of the underwater acoustic channel model in the shallow sea scenario; Step S2: Considering the characteristics of sound speed variation in shallow seas, construct a layered sound speed profile and calculate the propagation distance and time delay; Step S3: Based on the constructed layered sound velocity profile and combined with the principles of geometric acoustics and Snell's law, calculate the path parameters between the transmitter and receiver, including propagation angle, Doppler frequency shift, propagation distance, and time delay. Step S4: Calculate the large-scale fading parameters of the channel; Step S5: Calculate the impulse response and transmission function of the channel between the transmitter and receiver based on the path parameters and large-scale fading parameters; Step S6: Calculate the statistical characteristic parameters of the channel based on the channel transfer function.

2. The method for modeling a geometrically stochastic underwater acoustic channel considering non-uniform sound speeds according to claim 1, characterized in that, Step S1 specifically includes: Step S101: Set the initial depths of the transmitter and receiver to z respectively. T and z R The horizontal distance is D; Step S102: Set the motion parameters of the transmitter and receiver, with motion speeds v respectively. T and v R The azimuth angles of the motion are respectively and Step S103: Set shallow sea environment parameters, including seawater depth h and initial sound velocity c0 at the sea surface.

3. The method for modeling a geometrically stochastic underwater acoustic channel considering non-uniform sound speeds according to claim 2, characterized in that, Step S2 specifically includes: Step S201: Construct sound velocity profiles with different gradients, and set the sound velocity gradient in each layer to be constant, and the depth of each layer to be the same. The expression for the sound velocity profile is: Where K represents the total number of water layers, then k = 1, 2, ..., K, c k Let g represent the speed of sound at the lower boundary of the k-th layer. k Δz represents the sound velocity gradient of the k-th layer, and Δz represents the depth of each layer. Step S202: Calculate the propagation distance S of the sound wave through the entire water layer according to the definition of arc length: Among them, R k This represents the radius of curvature of the k-th sound layer. Let represent the incident angle at the lower boundary of the k-th layer and the incident angle at the lower boundary of the (k+1)-th layer. Using Snell's Law calculate; Step S203: Calculate the propagation time T of the sound wave through the entire water layer according to the principles of geometric acoustics and Snell's law:

4. The method for modeling a geometrically stochastic underwater acoustic channel considering non-uniform sound speeds according to claim 3, characterized in that, Step S3 specifically includes: Step S301: Calculate the propagation angle of the direct path: Where, α L β represents the departure angle of the direct path. L Indicates the angle of arrival of the direct route. This indicates the angle of incidence of the direct path at the launching end. The angle of incidence of the direct path at the receiving end is expressed using Snell's law. Calculate, where c T and c R k represents the speed of sound at the locations of the transmitter and receiver, respectively. T and k R These represent the water layer numbers where the transmitter and receiver are located, respectively, using the formula... and calculate, This represents the floor function; Step S302: Calculate the propagation angles corresponding to the DA and UA paths: in, and These represent the departure angle and arrival angle of the DA path, respectively. and These represent the departure angle and arrival angle of the UA path, respectively; s and These represent the number of surface reflections and bottom reflections experienced by the DA path, respectively. b and b represent the number of surface reflections and bottom reflections experienced by the UA path, respectively; and These represent the incident angles of the DA path at the transmitting and receiving ends, respectively. and These represent the incident angles of the UA path at the transmitting and receiving ends, respectively. In the geometrically random underwater acoustic channel model, for the boundary reflection paths, the DA path is the path component that the last reflection occurs at the sea surface and propagates downward to the receiver, while the UA path is the path component that the last reflection occurs at the seabed and propagates upward to the receiver. Step S303: Calculate the Doppler frequency shift f of the direct path based on the motion parameters of the transceiver. L : in, and Representing the maximum Doppler frequency shift at the transmitting and receiving ends, respectively, calculated as follows: and Among them, f c It is the carrier frequency; Step S304: Based on the motion parameters and reflection path characteristics of the transceiver, calculate the Doppler frequency shift corresponding to the DA path and the UA path. and Step S305: Combining the sound velocity profile and geometric relationships, determine the propagation distance d of the direct path. L With delay τ L (τ): in, and They represent the k locations of the transmitting end, respectively. T The layer and the receiver are located in k R The radius of curvature of the layered sound ray. and They represent the k locations of the transmitting end, respectively. T The layer and the receiver are located in k R The sound velocity gradient of the layer, and These represent the direct path at layer k and the k-th layer where the transmitter is located, respectively. T The layer and the receiver are located in k R The lower boundary angle of incidence of the layer; Step S306: Combining the sound velocity profile and geometric relationships, calculate the propagation distances of paths DA and UA, respectively. and Timely delay and in, and These represent the k-th layer and the k-th layer where the DA path is located, respectively. T The layer and the receiver are located in k R The lower boundary incident angle of the layer, and These represent the k-th layer and the k-th layer where the transmitter is located, respectively. T The layer and the receiver are located in k R The incident angle at the lower boundary of the layer.

5. The method for modeling a geometrically stochastic underwater acoustic channel considering non-uniform sound speeds according to claim 4, characterized in that, Step S4 specifically includes: Step S401: For the direct path, considering absorption loss and geometric spread loss, the large-scale fading parameters of the direct path are: [The large-scale fading loss L is missing from the original text]. L (d L The expression for ) is: L L (d L )=L a (d L )L s (d L ) Among them, L a (·) represents absorption loss, L s (·) represents the geometric spread loss, and the absorption loss L a The expression for (d) is: Where d represents the total propagation distance, α(f) represents the absorption coefficient, and is calculated using the formula: Where f represents the carrier frequency; Geometric spread loss L s The expression for (d) is: L s (d)=d -β / 2 Where β represents the expansion coefficient and is usually set to 2; Step S402: For the DA and UA paths, considering absorption loss, geometric spread loss, and bottom reflection loss as large-scale fading parameters, the large-scale fading losses of the DA and UA paths are then calculated. and The expressions are as follows: Among them, L b (·) represents the bottom reflection loss, and its expression is: in, ρ represents the angle of incidence at the lower boundary of the k-th water layer. k and ρ b Let c represent the densities of seawater and seabed at the lower boundary of the k-th water layer, respectively. b This represents the speed of sound on the seabed at the lower boundary of the k-th water layer.

6. The method for modeling a geometrically stochastic underwater acoustic channel considering non-uniform sound speeds according to claim 5, characterized in that, Step S5 specifically includes: Step S501: Divide the channel impulse response from the transmitter to the receiver into direct path and non-direct path. The general expression for the channel impulse response h(t,τ) is: h(t,τ)=h L (t,τ)+h N (t,τ) Among them, h L (t,τ) represents the channel impulse response of the direct path, h N (t,τ) represents the channel impulse response of a non-direct path, where h L The expression for (t,τ) is: Among them, a L The channel gain of the direct path is expressed as follows: K R Represents the Rice factor; θ L The phase shift representing the direct path is defined as a random variable that follows a uniform distribution in (0, 2π]. Step S502, Impulse response h of non-direct path N (t,τ) is divided into two types based on the type of arrival path: DA path and UA path. The two parts are summed according to a certain ratio to obtain the complete expression of the non-direct path impulse response: h N (t,τ)=h DA (t,τ)+h UA (t,τ) Among them, the impulse responses h of the DA and UA paths DA (t,τ) and h UA The expression for (t,τ) is as follows: The channel gain expression for the DA path is: The expression for the channel gain of the UA path is: And η DA and η UA Used to balance the contributions of paths DA and UA to the total power, satisfying η DA +η UA =1; N S N represents the maximum number of times the DA path contacts the sea surface. B This indicates the maximum number of times the UA pathway contacts the seabed. Each DA pathway and UA pathway is represented by... and Composed of micro-scattered rays, and Let represent the phase shifts of the DA and UA paths, respectively, and define them as random variables that follow a uniform distribution in (0, 2π]. Step S503: Perform a Fourier transform on the constructed channel impulse response to obtain the expression for the corresponding channel transfer function: H(t,f)=H L (t,f)+H N (t,f) H N (t,f)=H DA (t,f)+H UA (t,f) Where H(t,f) represents the overall channel transfer function, H L (t,f) and H N (t,f) represent the channel transfer functions for the direct path and the indirect path, respectively; H DA (t,f) and H UA (t,f) represent the channel transmission functions of the DA path and the UA path, respectively.

7. The method for modeling a geometrically stochastic underwater acoustic channel considering non-uniform sound speeds according to claim 6, characterized in that, Step S6 specifically includes: Step S601: Based on the established channel transmission function, calculate the time autocorrelation function of the channel, the expression of which is: Where Δt represents the time interval, I0(κ) is the zeroth-order modified Bessel function of the first kind, and κ represents the concentration level. and express and The average value; Step S602: Based on the established channel transmission function, calculate the channel frequency correlation function, the expression of which is: Where Δf represents the frequency interval.

Images