Method for calculating transmission time of multipath underwater acoustic signal based on markov model

Abstract

This invention relates to the field of underwater acoustic signal processing technology, and discloses a method for calculating the transmission time of multipath underwater acoustic signals based on a Markov model. The method includes: first, performing bandgap conditioning and correlation processing on the received underwater signal to extract the energy envelope, and mapping it to a likelihood distribution in the time delay domain to obtain multipath structure observation parameters; then, calculating the tidal period phase factor based on the signal acquisition time, and establishing a multipath state evolution model modulated by this phase factor; next, using the multipath structure observation parameters to drive the model for time-series recursion, calculating the confidence sequence of the first-arrival multipath cluster's belonging state; finally, using the confidence sequence as weights, performing a fusion estimation of the preset propagation delay corresponding to each state, and outputting a continuous transmission time. This invention can effectively suppress jump errors caused by the time-varying multipath structure due to internal tides, and significantly improve the stability and accuracy of transmission time calculation in complex shallow sea environments.

Description

Technical Field

[0001] This invention relates to the field of underwater acoustic signal processing technology, and more specifically, to a method for calculating the transmission time of multipath underwater acoustic signals based on Markov models. Background Technology

[0002] In underwater positioning, ranging, and time synchronization applications, accurately obtaining the transmission time (i.e., propagation delay) of a signal from the transmitter to the receiver is crucial. Shallow sea waveguide environments exhibit significant multipath propagation effects, and the signal received by the receiver is often a complex waveform composed of multiple superimposed paths (such as direct sound, sound reflected from the sea surface, and sound reflected from the seabed). Traditional methods for calculating transmission time typically employ matched filtering or correlation processing techniques, determining the transmission time by detecting correlation peaks or identifying the first arriving signal exceeding a preset threshold.

[0003] However, in practical shallow-sea engineering applications, the marine environment is not static. In particular, the dynamic processes such as internal tides and their derived internal solitary waves cause significant fluctuations in the seawater temperature, salinity, and depth profiles, leading to time-varying sound velocity fields. This environmental mechanism not only causes slight jitters in the transmission time of single paths but also results in drastic energy transfers in multipath arrival structures. For example, under certain tidal phases, the originally strongest direct path may attenuate to undetectable levels, while the originally weaker reflected path may amplify to become the first detectable signal.

[0004] Existing transmission time calculation schemes, whether based on peak detection with a fixed threshold or path tracking algorithms based on conventional Kalman filtering or homogeneous hidden Markov models, typically treat channel variations as irregular random interference or Gaussian white noise. These existing methods do not consider the periodic gating effect of the internal humid phase on the selection of sound propagation paths, i.e., they ignore the dominant influence of environmental conditions on which path arrives first.

[0005] Therefore, when internal tides cause the first-arrival signal to switch between different path clusters, existing algorithms often incorrectly identify this as measurement noise or abnormal jumps, resulting in phase-locking biases that are difficult to eliminate. This bias leads to a systematic error in the calculated transmission time series, which fluctuates with the tidal cycle, severely affecting the stability and accuracy of high-precision underwater navigation and long-duration communication. Summary of the Invention

[0006] This invention provides a method for calculating the transmission time of multipath underwater acoustic signals based on Markov models, which solves the technical problems mentioned in the background art.

[0007] This invention provides a method for calculating the transmission time of multipath underwater acoustic signals based on a Markov model, including:

[0008] Bandwidth conditioning and correlation processing are performed on the underwater received signal to extract the energy envelope reflecting the multipath channel impulse response.

[0009] The energy envelope is mapped to a likelihood distribution in the time delay domain, and multipath structure observation parameters reflecting the energy concentration of the first-arrival signal are extracted based on the likelihood distribution.

[0010] The tidal periodic phase factor is calculated based on the signal acquisition time, and a multipath state evolution model of the vocal tract modulated by the tidal periodic phase factor is established.

[0011] The multipath structure observation parameters are used to drive the vocal tract multipath state evolution model to perform time-series recursion and solve the confidence sequence of the first-arrival multipath cluster affiliation state;

[0012] Using the confidence sequence as weights, the preset propagation delay corresponding to each state is fused and estimated to output the continuous underwater acoustic signal transmission time.

[0013] The beneficial effects of this invention include: it can effectively solve the problem of transmission time calculation deviation caused by the time-varying multipath structure due to tidal dynamics in complex waveguide environments in shallow seas; by introducing the tidal periodic phase factor as a gating variable for state evolution and combining it with the soft observation parameters based on likelihood distribution for time-series recursion, this invention establishes the correlation between channel state switching and environmental period, thereby avoiding the jump error generated by the traditional hard decision method when path energy fluctuates, eliminating the systematic bias of tidal phase locking, and significantly improving the continuity, stability and calculation accuracy of underwater long-term communication and navigation positioning in dynamic environments. Attached Figure Description

[0014] Figure 1 This is a flowchart of the multipath underwater acoustic signal transmission time calculation method based on Markov model of the present invention;

[0015] Figure 2 This is a schematic diagram illustrating a specific implementation of the present invention. Detailed Implementation

[0016] The subject matter described herein will now be discussed with reference to exemplary embodiments. It should be understood that these embodiments are discussed only to enable those skilled in the art to better understand and implement the subject matter described herein, and changes may be made to the function and arrangement of the elements discussed without departing from the scope of this specification. Various processes or components may be omitted, substituted, or added as needed in the examples. Furthermore, features described in some examples may be combined in other examples.

[0017] like Figure 1 As shown, the method for calculating the transmission time of multipath underwater acoustic signals based on the Markov model includes:

[0018] Bandwidth conditioning and correlation processing are performed on the underwater received signal to extract the energy envelope reflecting the multipath channel impulse response.

[0019] The energy envelope is mapped to a likelihood distribution in the time delay domain, and multipath structure observation parameters reflecting the energy concentration of the first-arrival signal are extracted based on the likelihood distribution.

[0020] The tidal periodic phase factor is calculated based on the signal acquisition time, and a multipath state evolution model of the vocal tract modulated by the tidal periodic phase factor is established.

[0021] The multipath structure observation parameters are used to drive the vocal tract multipath state evolution model to perform time-series recursion and solve the confidence sequence of the first-arrival multipath cluster affiliation state;

[0022] Using the confidence sequence as weights, the preset propagation delay corresponding to each state is fused and estimated to output the continuous underwater acoustic signal transmission time.

[0023] Preferably, the underwater received signal undergoes band conditioning and correlation processing to extract the energy envelope reflecting the multipath channel impulse response, including:

[0024] Set system sampling rate With the upper bound of maximum delay Constructing a discrete time-delay grid ,in , , This represents the number of time-delay grid points;

[0025] Using bandpass filters For the first Underwater received signal per cycle Perform filtering:

[0026] ;

[0027] in, This represents the convolution operation;

[0028] For the filtered signal Amplitude standardization is performed to obtain a standardized signal. :

[0029] ;

[0030] in For frame length, To prevent zero constant;

[0031] Using a pre-set pilot sequence The energy envelope of the multipath channel impulse response is obtained by performing cross-correlation on the standardized signal and taking the modulus value. :

[0032] ;

[0033] ;

[0034] in This is a time-delayed response.

[0035] The system sampling rate is the sampling frequency when collecting underwater received signals, preferably between 10 kHz and 100 kHz, to match the typical bandwidth of underwater acoustic signals and ensure that key frequency components are not lost after signal sampling.

[0036] The maximum time delay upper bound is a preset maximum range of propagation time, preferably 1 to 10 seconds, calculated from the maximum geometric distance between the transmitter and receiver and the minimum sound speed of seawater, ensuring coverage of all possible multipath propagation times.

[0037] A bandpass filter is a device that performs frequency domain filtering on received signals. It is preferably a finite-length unit impulse response filter or an infinite-length unit impulse response filter with a passband frequency range of 1 kHz to 20 kHz and a filter order of 16 to 64 to match the operating frequency of underwater acoustic detection signals and filter out external low-frequency noise and high-frequency interference.

[0038] The underwater received signal in the nth cycle is the superimposed underwater acoustic signal collected by the receiver in the nth transmission cycle, which can be obtained by an underwater acoustic transducer (hydrophone).

[0039] The observation window frame length is the length of the time window for energy calculation of the filtered signal, preferably 1024 to 4096 points, to match the preset pilot sequence length and ensure complete coverage of the pilot signal duration.

[0040] The zero-prevention constant is a small positive number that prevents calculation errors when the signal energy is zero. It is preferably between 10 to the power of -6 and 10 to the power of -4, so as not to affect the amplitude range after the signal is standardized, and to effectively prevent calculation errors with a denominator of zero.

[0041] The preset pilot sequence is a known reference signal that is pre-set and transmitted by the transmitter. It is preferably a pseudo-random sequence or an orthogonal sequence with a length of 1024 to 4096 points to have good autocorrelation characteristics, so that the receiver can identify the signal arrival time through cross-correlation calculation.

[0042] First, the received signal is frequency-domain filtered using a bandpass filter to retain effective components that match the pilot signal frequency, while filtering out low-frequency current noise and high-frequency electromagnetic interference from the marine environment. Next, the root-mean-square amplitude of the filtered signal within the observation window is calculated, and this is used as a benchmark to normalize the signal amplitude, eliminating the impact of amplitude fluctuations across different transmission periods. Finally, a cross-correlation operation is performed between the normalized signal and a preset pilot sequence. The magnitude of the cross-correlation result highlights the energy peak of the multipath signal, generating the energy envelope of the multipath channel impulse response. For example, when ship navigation noise is mixed into the received signal, the bandpass filter can remove the noise frequency components, the normalization process ensures that the signal amplitudes of different periods are on the same order of magnitude, and the cross-correlation operation can accurately capture the energy peak corresponding to the arrival time of each path signal.

[0043] The introduction of a zero-prevention constant during amplitude normalization addresses the extreme situation where signal energy may attenuate to near zero in shallow sea environments, preventing anomalies where the denominator in the normalization calculation would be zero. For example, when an internal solitary wave causes extremely low received signal energy in a certain period, the zero-prevention constant ensures that the normalization operation is performed normally, preventing subsequent energy envelope extraction failures due to calculation errors.

[0044] The specific parameters of the bandpass filter are as follows: If a finite-length unit impulse response filter is selected, the passband frequency range is set to 1 kHz to 20 kHz, the stopband attenuation is not less than 40 dB, and the filter order is 32; if an infinite-length unit impulse response filter is selected, the passband ripple is not more than 0.5 dB, the stopband attenuation is not less than 60 dB, and the order is 16, to ensure effective filtering of non-target frequency components.

[0045] The value of the observation window frame length must be consistent with the preset pilot sequence length. For example, when the pilot sequence length is 2048 points, the observation window frame length should be set to 2048 points to ensure that the pilot signal can be completely covered when calculating the root mean square amplitude, and to avoid energy calculation deviation due to mismatch in window length.

[0046] The specific value of the zero-prevention constant is 10 to the power of negative 5. This value is much smaller than the typical amplitude after signal standardization (usually between 0.1 and 1), and will not have a substantial impact on the standardization result. At the same time, it can reliably avoid the calculation anomaly of zero denominator.

[0047] The preset pilot sequence is a pseudo-random sequence (such as an m-sequence) with a length of 2048 points. The autocorrelation function of this sequence exhibits a sharp peak when the time delay is zero and approaches zero at other time delays. This makes it easier for the receiver to accurately identify the arrival time of each multipath signal through cross-correlation calculations, thereby reducing mutual interference between multipaths.

[0048] Preferably, mapping the energy envelope to a likelihood distribution in the time-delay domain includes:

[0049] Set temperature coefficient ;

[0050] Based on the energy envelope of the multipath channel impulse response Calculate the first The cycle is in the first Likelihood distribution of the time delay domain at each time delay grid point :

[0051] ;

[0052] in The total number of time-delay grid points. It is an exponential function.

[0053] The temperature coefficient is a key parameter for adjusting the concentration of the likelihood distribution in the time delay domain. It is preferably between 0.1 and 10 to balance the concentration of the distribution with the ability to resist interference. This ensures that the contribution of strong energy multipaths is highlighted, while avoiding the loss of effective information of weak energy multipaths due to excessive concentration of the distribution.

[0054] The time-delay domain likelihood distribution transforms the energy envelope of the multipath channel impulse response into a probabilistic distribution, reflecting the reliability of the multipath signal corresponding to each time-delay grid point.

[0055] By introducing a temperature coefficient and applying exponential mapping and normalization, the energy envelope is transformed into a time-delay domain likelihood distribution. This soft processing method differs from traditional hard-decision methods. The specific process is as follows: First, the energy envelope is multiplied by the temperature coefficient to amplify or reduce energy differences. Then, the exponential function value is calculated for the product at each time-delay grid point to strengthen the weight of high-energy points. Subsequently, the exponential function values ​​of all grid points are summed to obtain a normalization factor. Finally, the exponential function values ​​of each point are divided by this factor to obtain a likelihood distribution with a probability sum of 1. For example, when the energy difference of multipath signals is large, a temperature coefficient of 5 can make the likelihood values ​​corresponding to high-energy multipaths more prominent; when the multipath energy distribution is uniform, a temperature coefficient of 1 can evenly retain the probability contribution of each multipath, avoiding the omission of weak-energy multipath information.

[0056] The specific adjustment principle of the temperature coefficient is as follows: when the signal-to-noise ratio (SNR) of the received signal is higher than 20 dB, a value of 5 to 10 is preferred to enhance the distribution concentration and highlight the main path signal; when the SNR is between 10 and 20 dB, a value of 1 to 5 is preferred to balance the contributions of the main and secondary paths; when the SNR is lower than 10 dB or when multipath is dense, a value of 0.1 to 1 is preferred to expand the distribution coverage and avoid losing weak multipath information. For example, in a calm shallow sea environment with a high SNR, the temperature coefficient is set to 8; when internal waves are active and multipath is dense, the temperature coefficient is set to 0.5 to ensure that the likelihood distribution can fully reflect the multipath structure.

[0057] Preferably, multipath structure observation parameters reflecting the energy concentration of the first-arrival signal are extracted based on the likelihood distribution, including:

[0058] Calculate the centroid of the time delay distribution using the following formula. :

[0059] ;

[0060] Calculate the time delay distribution diffusion using the following formula. :

[0061] ;

[0062] Set the early arrival weighted decay factor The weighted energy moment of the early arrival signal is calculated according to the following formula. :

[0063] ;

[0064] The multipath structure observation parameters are obtained by combining them. :

[0065] ;

[0066] in, For the first The delay value of each delay grid point. For the first The likelihood distribution value of the time delay domain for each period, This represents the total number of time-delay grid points.

[0067] The early arrival weighted attenuation factor is an adjustment parameter that highlights the energy contribution of the first arrival signal and suppresses the interference of the late arrival signal. It is preferably between 1 and 10 to adapt to the time interval characteristics of multipath signals, ensuring that the weight of the first arrival signal is dominant, while not completely ignoring the information of subsequent effective multipath signals.

[0068] The centroid of the time delay distribution is a statistical measure that reflects the location of the energy concentration of a multipath signal. It is obtained by weighted summation of the time delay values ​​and the corresponding likelihood distribution values ​​at each time delay grid point.

[0069] The time delay distribution diffusion is a statistical measure that reflects the degree of energy dispersion of multipath signals. It is obtained by weighted summation of the squared differences between each time delay grid point and the centroid through the likelihood distribution values.

[0070] The early arrival signal weighted energy moment is a statistical measure that highlights the proportion of energy of the first arrival signal. It is obtained by weighting the time delay by the early arrival attenuation factor and then summing the likelihood distribution values.

[0071] The multipath structure observation parameters are a combination of parameters that comprehensively characterize the energy accumulation characteristics of the first-arrival signal. They consist of the centroid of the time delay distribution, the spread of the time delay distribution, and the weighted energy moment of the early-arrival signal.

[0072] By combining statistical measures across three dimensions, a multipath structure observation parameter is constructed to comprehensively capture the energy distribution characteristics of multipath signals, distinguishing it from traditional single-dimensional observation indicators. Specifically, the process is as follows: the centroid of the delay distribution locates the core delay position where energy is concentrated, reflecting the approximate arrival interval of the first-arrival signal; the diffusion of the delay distribution reflects the dispersion of multipath signals, indicating the density of multipath events in the channel; and the weighted energy moment of the early-arrival signal strengthens the weight of the first-arrival signal through an attenuation factor, highlighting the need for first-arrival detection. For example, when internal humidity causes the first-arrival cluster to switch from a direct path to a reflection path, the centroid of the delay distribution will shift accordingly, and the diffusion may change due to the redistribution of multipath energy. The weighted energy moment of the early-arrival signal can quickly respond to energy changes in the first-arrival path. The combination of these three parameters can fully characterize the switching process, providing a comprehensive basis for subsequent state identification.

[0073] The specific adjustment principle of the early arrival weighted attenuation factor is as follows: when the time interval between the first and second arrival paths of a multipath signal is greater than 10 milliseconds, a value of 5 to 10 is preferred to strengthen the weight of the first arrival signal; when the time interval is between 3 and 10 milliseconds, a value of 3 to 5 is used to balance the contributions of the first and second arrival signals; when the time interval is less than 3 milliseconds or when multipaths overlap densely, a value of 1 to 3 is used to avoid excessive suppression of adjacent effective multipaths. For example, when the multipath interval is large during long-distance transmission in shallow waters, the attenuation factor is set to 8; when the multipaths are dense at close range in shallow waters, the attenuation factor is set to 2 to ensure accurate capture of the energy concentration characteristics of the first arrival signal.

[0074] Preferably, the tidal period phase factor is calculated based on the signal acquisition time, including:

[0075] Set semi-diurnal tidal cycle Calculate the tidal angular frequency :

[0076] ;

[0077] Based on the current signal acquisition time The tidal period phase factor is calculated using the following formula. :

[0078] ;

[0079] in, The preset reference start time, The preset reference initial phase, This represents a modulo operation that constrains a numerical value to the range of zero to twice the value of pi.

[0080] The semi-diurnal tidal period is the standard period of the semi-diurnal tidal period in ocean tides. It is 12.42 hours in length and is an internationally recognized ocean tidal constant determined based on global ocean tidal observation data.

[0081] Tidal angular frequency is a parameter that describes the rate of change of the tidal cycle and is calculated from the semi-diurnal tidal cycle.

[0082] The signal acquisition time is the specific time point at which the receiver acquires the underwater acoustic signal of the current cycle. It can be obtained through a high-precision clock module (such as a GPS timing module) at the receiver.

[0083] The reference start time is a preset time reference point for calculating the tidal phase, preferably the system startup time or the local tidal observation start time, so as to unify the time calculation reference and ensure the continuity of the phase factor.

[0084] The reference initial phase is a preset tidal phase starting offset, preferably 0 radians, to simplify the calculation process and ensure the regularity of tidal phase changes over time.

[0085] The tidal period phase factor is a key variable for quantifying the periodicity of tides. It reflects the current position of the tide in the period and is obtained by modulo operation from angular frequency, time difference, and initial phase.

[0086] The tidal periodicity is transformed into a quantifiable phase factor, which serves as a gating variable for the subsequent state evolution model. The process of linking environmental periodicity with the signal processing model is as follows: First, the internationally recognized 12.42-hour semi-diurnal tidal period is used to calculate the tidal angular frequency. Then, the time difference between the acquisition time and the reference start time is obtained, and the angular frequency is multiplied by this difference to obtain the phase component that varies linearly with time. Subsequently, the reference initial phase is superimposed to compensate for the time reference deviation. Finally, a modulo operation between 0 and 6.28 is performed on the superposition result to constrain the phase within one period, ensuring the periodicity of the phase factor. For example, if the reference start time is set to the 0 point of system startup and the reference initial phase is 0 radians, when the acquisition time is 6.21 hours (half of the semi-diurnal tidal period), the angular frequency is 2π divided by 12.42, and the calculated phase component is π radians. After modulo operation, the phase factor is π radians, accurately corresponding to the phase state of the tidal half-cycle.

[0087] The reference start time setting needs to be synchronized with the receiver clock. The system startup time or the start time of the tide forecast should be selected first. For example, if the system starts at 0:00 on a certain day, this should be used as the reference start time to ensure that the time difference calculation benchmark of all collection times is consistent.

[0088] The preferred value for the initial reference phase is 0 radians. If it is necessary to adapt to the tidal phase offset of a specific sea area in actual application, it can be adjusted to a specific value between 0 and 6.28 based on the tidal observation data of that sea area. For example, if the tidal initial phase offset of a certain sea area is 0.314 radians, then the initial reference phase is set to 0.314 radians.

[0089] The specific implementation of the modulo operation is as follows: divide the superposition result by 6.28 and take the remainder as the final phase factor. For example, if the superposition result is 7.536, then the remainder of 7.536 divided by 6.28 is 1.256, and the corresponding phase factor is 1.256 radians, ensuring that the result is between 0 and 6.28.

[0090] Preferably, establishing a multipath state evolution model of the vocal tract modulated by the tidal periodic phase factor includes:

[0091] Define a discrete hidden state representing the category to which the first-reach multipath cluster belongs. ,in The total number of states;

[0092] The phase factor during the tidal cycle is calculated using the following formula. Below, from the previous moment in hidden state Transition to the hidden state at the current moment Nonnormalized transfer intensity value :

[0093] ;

[0094] in, For baseline transfer intensity, and These are the phase gating coefficients;

[0095] The hidden state transition probability matrix that dynamically varies with the tidal period phase factor is calculated using the following formula. :

[0096] ;

[0097] in, It is an exponential function. For summation index.

[0098] The total number of discrete hidden states is the number of categories to which the first-reach multipath cluster belongs, preferably 2 to 4. The typical cluster classification of shallow sea multipath propagation usually includes direct path clusters, surface reflection path clusters, and bottom reflection path clusters. Too many clusters will increase computational complexity, while too few clusters will not be able to cover key multipath clusters.

[0099] The baseline transition strength is the strength parameter of the basic transition between hidden states, preferably between 0 and 1, to ensure the rationality of the initial state transition and avoid model imbalance caused by excessively high or low baseline transition probabilities.

[0100] The phase gating coefficient is a parameter that adjusts the degree of influence of the tidal period phase on the state transition. It is preferably between -1 and 1 to balance the intensity of phase modulation, so as to reflect the gating effect of the tidal phase without causing excessive fluctuations in the transition probability.

[0101] The unnormalized transfer intensity value is the state transfer intensity that has not undergone normalization processing, and is obtained by superimposing the baseline transfer intensity and the phase-weighted component.

[0102] The hidden state transition probability matrix is ​​a matrix that describes the transition probabilities between different hidden states. Its elements change dynamically with the tidal period phase factor to ensure that the sum of the probabilities in each row is 1.

[0103] The hidden states are directly defined as the category to which the first-reach multipath cluster belongs, rather than the generalized channel state. A non-homogeneous state evolution model is constructed by modulating the transfer intensity using the cosine and sine components of the tidal period phase. The specific process is as follows: First, discrete hidden states are defined based on the characteristics of shallow sea multipaths. For example, when K=3, state 1 represents the direct path cluster, state 2 represents the surface reflection path cluster, and state 3 represents the bottom reflection path cluster. Then, for any transition between two states, the baseline transfer intensity is superimposed with the cosine-weighted and sine-weighted components of the tidal phase to form a phase-variable non-normalized transfer intensity. Finally, an exponential normalization operation is performed on all non-normalized transfer intensities corresponding to each initial state to ensure that the sum of each row of the transfer probability matrix is ​​1. For example, when the tidal phase is 0 radians, the cosine term is 1, the sine term is 0, and the non-normalized transfer intensity is... The corresponding state transition probability increases accordingly; when the phase is π radians, the cosine term is -1, and the transition strength is... The transition probability is adjusted accordingly to make the state transition conform to the tidal cycle.

[0104] The total number of discrete hidden states needs to be determined based on the multipath cluster observation results of the specific application sea area. For example, in shallow seas with short-range transmission, there are few multipath clusters, so a value of 2 is sufficient; in long-range transmission, there are abundant multipath clusters, so a value of 3 or 4 is appropriate. For example, in a certain shallow sea area, where direct clusters and surface reflection clusters mainly exist, the total number of discrete hidden states is set to 2. The baseline transfer intensity and phase gating coefficient are initialized as follows: the baseline transfer intensity is initialized to 0.5, and the phase gating coefficient is initialized to 0 to ensure the stability of the initial state of the model; subsequently, iterative updates are made using the expectation-maximization algorithm on offline training data until the model fitting effect meets the requirements.

[0105] The specific implementation of exponential normalization is as follows: calculate the exponential function value for the non-normalized transition strength value of all target states j corresponding to each initial state i, and then divide the exponential function value of each target state j by the sum of the exponential function values ​​of all target states to obtain the transition probability of each target state under the initial state, ensuring that the sum of the probabilities of each row is 1.

[0106] Preferably, the multipath structure observation parameters are used to drive the vocal tract multipath state evolution model to perform time-series recursion, and the confidence sequence of the first-arrival multipath cluster affiliation state is calculated, including:

[0107] Initialize the confidence sequence at time 1 ;

[0108] For the At each time point, the confidence sequence at the current time is recursively calculated using the following formula. :

[0109] ;

[0110] Perform normalization:

[0111] ;

[0112] in, The multipath structure observation parameters at the current moment are... In the hidden state and the tidal period phase factor The observed emission probability density is as follows. The state of the previous moment Confidence level, Let be the hidden state transition probability matrix.

[0113] The initial confidence level is the confidence level of the state to which each first-arrival multipath cluster belongs at the start of the time-series recursion. It is preferably uniformly distributed. For example, when the total number of hidden states is 3, the initial confidence level of each state is 0.333, so as to ensure that the initial weights of each state are fair when there is no initial observation information and to avoid model bias.

[0114] The observation emission probability density is the probability density of the occurrence of multipath structure observation parameters given the hidden state and tidal period phase factor. It is preferably a Gaussian distribution, so that the statistical characteristics of multipath structure observation parameters conform to the normal distribution law of continuous random variables and are suitable for the quantitative characteristics of observation data.

[0115] The state prediction probability is the predicted state value at the current time moment, which is obtained based on the confidence level of the previous time moment and the hidden state transition probability matrix, and is generated by a weighted summation operation.

[0116] The current confidence sequence is the confidence sequence of the state of each first-arrival multipath cluster at the current time. It is obtained by multiplying the state prediction probability and the observation emission probability density and then normalizing it, and reflects the current confidence level of each state.

[0117] By combining multipath structure observation parameters with a latent state transition probability matrix modulated by tidal phase, a time-series recursive logic of initial confidence initialization, time-by-time prediction, observation update, and normalization is used to calculate the confidence sequence of the first-arrival multipath cluster's belonging state. Specifically, the process is as follows: First, a uniformly distributed initial confidence is calculated based on the observation emission probability at the initial time. Then, at each time step, the confidence sequence from the previous time step is weighted and summed with the dynamic transition probability matrix to obtain a state prediction probability independent of the current observation. Next, the observation emission probability density of each latent state is calculated using the current multipath structure observation parameters to measure the degree of matching between the current observation and each state. Finally, the state prediction probability is multiplied by the observation emission probability density, and normalization is applied to ensure that the sum of the confidences of each state is 1, yielding the confidence sequence for the current time step. For example, if the total number of hidden states is 2, and the confidence of state 1 is 0.8 and the confidence of state 2 is 0.2 at the previous time step, and the transition matrix corresponding to the tidal phase has a probability of 0.7 for state 1 to state 1 and 0.3 for state 2, and a probability of 0.4 for state 2 to state 1 and 0.6 for state 2, then the weighted sum gives the predicted probability of state 1 as 0.8 × 0.7 + 0.2 × 0.4 = 0.64 and the predicted probability of state 2 as 0.8 × 0.3 + 0.2 × 0.6 = 0.36. If the current observation parameter is multiplied by the emission probability density of state 1 (0.9) and the emission probability density of state 2 (0.1), we get 0.64 × 0.9 = 0.576 and 0.36 × 0.1 = 0.036. After normalization, the confidence of state 1 is 0.576 ÷ (0.576 + 0.036) = 0.941 and the confidence of state 2 is 0.059, reflecting the state affiliation supported by the current observation.

[0118] The observed emission probability density is specifically set as a Gaussian distribution, with its mean being the benchmark value of the multipath structure observation parameters in the corresponding hidden state superimposed with the tidal phase fine-tuning. The covariance matrix is ​​a 3×3 positive definite matrix, with the diagonal elements corresponding to the variances of the centroid of the time delay distribution, the spread of the time delay distribution, and the weighted energy moment of the early arrival signal, respectively. The off-diagonal elements are the covariances between the parameters, which can be calculated by the weighted covariance of the offline training data.

[0119] The initial confidence level is calculated as follows: when the total number of hidden states is K, the initial confidence level of each state is set to 1÷K. For example, when K=3, the initial confidence level is 0.333, ensuring that the initial state is unbiased. If there are initial observation parameters, the initial observation parameters can be substituted into the emission probability density of each state and then normalized to obtain the initial confidence level.

[0120] The specific implementation of normalization is as follows: calculate the sum of the state prediction probability × observation emission probability density corresponding to all hidden states, and then divide the product of each hidden state by the sum to obtain the confidence of each state at the current time, ensuring that the sum of the confidence of all states is 1, thus guaranteeing the probability consistency of the recursive process.

[0121] Preferably, using the confidence sequence as weights, the preset propagation delay corresponding to each state is fused and estimated to output the continuous underwater acoustic signal transmission time, including:

[0122] For each hidden state The local mapping propagation delay corresponding to the hidden state is calculated using the following formula. :

[0123] ;

[0124] in, For a predetermined mapping coefficient vector, The observation parameters of the multipath structure are as follows. Based on the basic deviation value, Indicates the transpose operation;

[0125] The transmission time of the continuous underwater acoustic signal is calculated using the following formula. :

[0126] ;

[0127] in, For the corresponding hidden state in the confidence sequence The value, This represents the total number of hidden states.

[0128] The mapping coefficient vector is a linear projection coefficient that transforms the multipath structure observation parameters into local propagation delays in each hidden state. It is preferably between -10 and 10 to adapt to the numerical range of the multipath structure observation parameters and ensure that the delay result after linear projection is within the actual multipath propagation time interval.

[0129] The base deviation value is the calibration offset of the local mapping propagation delay in each hidden state, preferably 0 to 1 second, to compensate for the system deviation generated during linear projection and to fit the actual propagation delay characteristics of different multipath clusters.

[0130] The local mapping propagation delay is an estimate of the propagation delay corresponding to a single hidden state, obtained by superimposing the multipath structure observation parameters with the basic deviation value through linear projection.

[0131] The continuous underwater acoustic signal transmission time is the final output transmission time result, which is obtained by weighted summation of the local mapping propagation delays of all hidden states through a confidence sequence.

[0132] An independent local mapping propagation delay calculation method is designed for each hidden state belonging to the first-reach multipath cluster, and a fusion estimation is performed using the confidence level of each state as the weight to avoid transmission time jumps caused by hard decisions in a single state. The specific process is as follows: For each hidden state, a dedicated mapping coefficient vector and a basic deviation value are pre-determined. The multipath structure observation parameters are linearly projected using the mapping coefficient vector to transform the three-dimensional observation parameters into one-dimensional delay components. The basic deviation value is then superimposed to complete the calibration, obtaining the local mapping propagation delay of that state. Subsequently, all hidden states are traversed, and the local mapping propagation delay of each state is multiplied by the corresponding confidence level. All product results are then summed to obtain the final continuous transmission time. For example, hidden state 1 corresponds to the direct path cluster with a confidence level of 0.7 and a local mapping propagation delay of 0.5 seconds; hidden state 2 corresponds to the surface reflection path cluster with a confidence level of 0.3 and a local mapping propagation delay of 0.6 seconds; after fusion calculation, the transmission time is 0.7×0.5+0.3×0.6=0.53 seconds. Even if there is a subsequent state switch, the gradual change in confidence level will make the transmission time transition smoothly without any jumps.

[0133] The mapping coefficient vector has a dimension of 3, which corresponds one-to-one with the three components of the multipath structure observation parameters (time delay distribution centroid, time delay distribution diffusion, and early arrival signal weighted energy moment). It is initialized with small random numbers between -0.1 and 0.1. Subsequently, it is fitted and optimized by the weighted least squares algorithm of offline training data to ensure that the projection result accurately matches the multipath time delay characteristics of the corresponding hidden state.

[0134] The initial value of the basic deviation is 0. If there is equipment delay or channel offset in actual application, it can be adjusted to a specific value between 0 and 1 second according to historical transmission time data. For example, for a certain hidden state corresponding to the bottom reflection path cluster, the basic deviation value is set to 0.05 seconds after actual measurement and calibration.

[0135] The specific implementation of the weighted summation operation is as follows: first, obtain the confidence level of each hidden state at the current time and the corresponding local mapping propagation delay, calculate the product of the confidence level of each hidden state and the local mapping propagation delay, and then sum all the product results. The sum is the final continuous underwater acoustic signal transmission time.

[0136] like Figure 2 As shown, Figure 2 This exhibit showcases a multi-path underwater acoustic signal transmission scenario in a shallow sea environment, including key elements such as the transmitter, receiver (hydrophone), sea surface, and seabed. It clearly presents three signal transmission paths: the direct path (transmitter directly to receiver), the surface reflection path (reflected by the sea surface to reach receiver), and the bottom reflection path (reflected by the seabed to reach receiver). It also marks the sound speed disturbance area caused by internal tides, intuitively demonstrating the environmental impact of internal tides on underwater acoustic signal transmission.

[0137] The embodiments of this example have been described above. However, this example is not limited to the specific implementation methods described above. The specific implementation methods described above are merely illustrative and not restrictive. Those skilled in the art can make many other forms based on the guidance of this example, and all of them are within the protection scope of this example.

Claims

1. A method for calculating the transmission time of multipath underwater acoustic signals based on a Markov model, characterized in that, include: Bandwidth conditioning and correlation processing are performed on the underwater received signal to extract the energy envelope reflecting the multipath channel impulse response. The energy envelope is mapped to a likelihood distribution in the time delay domain, and multipath structure observation parameters reflecting the energy concentration of the first-arrival signal are extracted based on the likelihood distribution. The tidal periodic phase factor is calculated based on the signal acquisition time, and a multipath state evolution model of the vocal tract modulated by the tidal periodic phase factor is established. The multipath structure observation parameters are used to drive the vocal tract multipath state evolution model to perform time-series recursion and solve the confidence sequence of the first-arrival multipath cluster affiliation state; Using the confidence sequence as weights, the preset propagation delay corresponding to each state is fused and estimated to output the continuous underwater acoustic signal transmission time.

2. The method for calculating the transmission time of multipath underwater acoustic signals based on a Markov model according to claim 1, characterized in that, The underwater received signal undergoes band conditioning and correlation processing to extract the energy envelope reflecting the multipath channel impulse response, including: Set the system sampling rate and the upper bound of the maximum delay to construct a discrete delay grid; The acquired underwater received signal was filtered in the frequency domain using a bandpass filter to obtain the filtered signal. Calculate the root mean square amplitude of the filtered signal within the observation window, and perform amplitude normalization on the signal based on the root mean square amplitude to obtain the normalized signal; The standardized signal is cross-correlated with a preset pilot sequence, and the magnitude of the cross-correlation result is calculated to generate the energy envelope of the multipath channel impulse response.

3. The method for calculating the transmission time of multipath underwater acoustic signals based on a Markov model according to claim 1, characterized in that, Mapping the energy envelope to a likelihood distribution in the time-delay domain includes: Set a temperature coefficient to adjust the concentration of the probability distribution; The exponential function value is calculated using the product of the energy envelope of the multipath channel impulse response and the temperature coefficient as the independent variable. The sum of the exponential function values ​​at all time-delay grid points is used as the normalization factor; Calculate the ratio of the exponential function value to the normalization factor at each time delay grid point to generate the likelihood distribution of the time delay domain.

4. The method for calculating the transmission time of multipath underwater acoustic signals based on a Markov model according to claim 3, characterized in that, Based on the likelihood distribution, multipath structure observation parameters reflecting the energy concentration of the first-arrival signal are extracted, including: The centroid of the time delay distribution is obtained by summing the products of the time delay value corresponding to each time delay grid point and the likelihood distribution value of the time delay domain at that time delay grid point. Calculate the square of the difference between the time delay value corresponding to each time delay grid point and the centroid of the time delay distribution, and use the likelihood distribution value of the time delay domain to perform a weighted sum of the square of the difference to obtain the time delay distribution diffusion degree; Set an early arrival weighted attenuation factor, calculate the exponential function value with the negative of the product of the early arrival weighted attenuation factor and the time delay value corresponding to each time delay grid point as the independent variable, and use the likelihood distribution value in the time delay domain to perform a weighted summation of the exponential function value to obtain the early arrival signal weighted energy moment; The centroid of the time delay distribution, the spread of the time delay distribution, and the weighted energy moment of the early arrival signal are combined to form the observation parameters of the multipath structure.

5. The method for calculating the transmission time of multipath underwater acoustic signals based on a Markov model according to claim 1, characterized in that, The tidal period phase factor is calculated based on the signal acquisition time, including: Set a semi-diurnal tidal period and calculate the tidal angular frequency based on the semi-diurnal tidal period; Obtain the current signal acquisition time and calculate the time difference between the signal acquisition time and the preset reference start time; Calculate the product of the tidal angular frequency and the time difference, superimpose a preset reference initial phase, and perform a modulo operation on the superposition result with a modulus of twice pi to obtain the tidal period phase factor.

6. The method for calculating the transmission time of multipath underwater acoustic signals based on a Markov model according to claim 5, characterized in that, Establish a multipath state evolution model of the vocal tract modulated by the tidal periodic phase factor, including: Define a set of discrete hidden states representing the category to which the first-reach multipath cluster belongs; For any two hidden states, the nonnormalized transition intensity value is calculated by superimposing the cosine-weighted component and the sine-weighted component of the tidal period phase factor on the baseline transition intensity. An exponential normalization operation is performed on the non-normalized transition intensity values ​​corresponding to all hidden states to generate a hidden state transition probability matrix that dynamically changes with the tidal period phase factor.

7. The method for calculating the transmission time of multipath underwater acoustic signals based on a Markov model according to claim 6, characterized in that, The multipath structure observation parameters are used to drive the vocal tract multipath state evolution model for time-series recursion, solving the confidence sequence of the first-arrival multipath cluster affiliation state, including: Calculate the initial confidence level based on the observed launch probability at the initial moment; At each subsequent moment, the observation emission probability density of each hidden state is calculated using the multipath structure observation parameters; The state prediction probability is obtained by weighted summing of the confidence sequence from the previous time step and the hidden state transition probability matrix. The state prediction probability is multiplied by the observed emission probability density and normalized to obtain the confidence sequence of the first-arrival multipath cluster attribution state at the current time.

8. The method for calculating the transmission time of multipath underwater acoustic signals based on a Markov model according to claim 7, characterized in that, Using the confidence sequence as weights, the preset propagation delay corresponding to each state is fused and estimated to output the continuous underwater acoustic signal transmission time, including: For each hidden state, the multipath structure observation parameters are linearly projected using a pre-determined mapping coefficient vector and the basic deviation value is superimposed to generate the local mapping propagation delay corresponding to that hidden state. The continuous underwater acoustic signal transmission time is obtained by performing a weighted summation operation on the local mapping propagation delay corresponding to all hidden states using the confidence sequence.

Images