Reverberation Suppression Method Based on Angle-Delay-Doppler Three-Dimensional Deconvolution

Abstract

The present invention proposes a reverberation suppression method based on angle-delay-Doppler three-dimensional deconvolution, which transmits narrowband pulses to the target detection area, and then uses a set of matched filters carrying different Doppler frequency shifts to process the array element echoes to obtain a delay-Doppler two-dimensional output; conventional beamforming is performed on the delay-Doppler two-dimensional output of the array element to obtain an angle-delay-Doppler three-dimensional output; an angle-delay-Doppler three-dimensional beam map is obtained, and its absolute value is taken as the point spread function of the R-L algorithm; the absolute value and point spread function of the angle-delay-Doppler three-dimensional output are used as the input of the three-dimensional R-L algorithm, the number of iterations is selected, and the reverberation suppression result is calculated according to the deconvolution formula. The present invention uses the inherent energy distribution difference between the target echo and the strong reverberation in the three dimensions of angle, delay, and Doppler to suppress reverberation, improve the multi-dimensional resolution capability, suppress the multi-dimensional sidelobe energy leakage, and suppress the bottom and surface reverberation at the same time.

Description

Technical Field

[0001] The invention belongs to the field of array signal processing, and in particular relates to a reverberation suppression method based on angle-delay-Doppler three-dimensional deconvolution. Background Art

[0002] When detecting small moving targets in a reverberation-limited environment, high-intensity reverberation can mask the target echo. Therefore, the problem of strong reverberation suppression has attracted widespread attention.

[0003] In order to suppress the underwater reverberation by using Doppler difference, a delay-Doppler two-dimensional reverberation suppression method was proposed (Yang TC, Schindall J, Huang CF, et al, "Clutter reduction using Doppler sonar in a harbor environment," Journal of the Acoustical Society of America, 2012; 132 (5): 3053-3067.). However, this method uses a matched filter group (MFs: matched filters) to process the echo of a single array element, which has the problems of insufficient two-dimensional resolution, insufficient sidelobe suppression and low signal-to-mixing ratio. Although the delay-Doppler two-dimensional reverberation suppression method based on the Richardson-Lucy (RL) deconvolution algorithm can improve the two-dimensional resolution and sidelobe suppression capabilities to a certain extent. However, under low signal-to-mixing ratio conditions, the performance of this method will be significantly degraded (CN104780127A-A delay-Doppler RL deconvolution multipath channel estimation method-disclosed).

[0004] In order to improve the signal-to-mixing ratio by joint array processing, a method of angle-delay-Doppler three-dimensional reverberation suppression combined with conventional beamforming (CBF) was proposed (Gong Fusheng, "Joint Space-Time Three-Dimensional Estimation of Active Sonar Moving Targets", Yangtze River Information Communication, 2022; 35(04): 12-13+19.). However, this method only uses MFs and CBF to process echoes. This method still has the problems of insufficient multi-dimensional resolution and insufficient multi-dimensional sidelobe suppression. When a slow small target is close to a strong scatterer, the strong reverberation sidelobe energy leakage may mask the weak target.

[0005] Therefore, in order to simultaneously improve the multi-dimensional resolution and multi-dimensional sidelobe suppression capabilities and achieve effective suppression of strong reverberation, new solutions must be sought. Summary of the invention

[0006] In order to solve the problems of insufficient multi-dimensional resolution and multi-dimensional sidelobe suppression capabilities of traditional methods, a reverberation suppression method based on angle-delay-Doppler three-dimensional deconvolution is proposed by utilizing the inherent energy distribution differences between slow small target echoes and strong reverberation in three dimensions of angle, distance and Doppler.

[0007] The technical solution of the present invention is: a reverberation suppression method based on angle-delay-Doppler three-dimensional deconvolution, characterized in that it includes the following steps:

[0008] Step 1: Transmit pulses to the detection area, and the horizontal receiving array collects echoes;

[0009] Step 2: Calculate the delay-Doppler output of the array element echo;

[0010] Step 3: Perform conventional beamforming on the delay-Doppler output of each array element;

[0011] Step 4: Deconvolve the conventional angle-delay-Doppler three-dimensional beam output to obtain the reverberation suppression result, including the following sub-steps:

[0012] Step 4.1: Calculate the angle-delay-Doppler three-dimensional conventional beam output of the virtual target echo with no delay, no Doppler and zero angle, take its absolute value, and obtain the point spread function of the deconvolution algorithm;

[0013] Step 4.2: The absolute value of the three-dimensional beam output obtained by steps 1, 2, and 3 is used as the observation data of the iterative deconvolution algorithm;

[0014] Step 4.3: Iterate and calculate the reverberation suppression result according to the deconvolution formula.

[0015] Furthermore, the step 1 includes the following sub-steps:

[0016] Step 1.1: Initialize the detection area: define a large number of scatterers randomly distributed on the bottom and the surface, where the number of bottom scatterers is Q and their positions are different, and the number of surface scatterers is D and their positions are different; the bottom scatterers remain stationary, and the surface scatterers move randomly at a low speed; the number of moving targets is P; the detection sonar is a transmitting source and an N-element horizontal uniform linear receiving array; the underwater sound source transmits sound energy to the detection area, and the receiving array collects echoes; the echo includes the target echo and reverberation; the reverberation includes the bottom reverberation and the surface reverberation, where the bottom reverberation is the scattered wave of the bottom scatterer, and the surface reverberation is the scattered wave of the surface scatterer; the reverberation scatterer refers to the bottom scatterer and the surface scatterer; one target corresponds to one target scatterer; considering the multipath effect, the number of scattered waves from the pth target to the nth array element is M pn The number of scattered waves from the qth underwater scatterer to the nth array element is M qnThe number of scattered waves from the dth water surface scatterer to the nth array element is M dn It is assumed that one scattered wave corresponds to one virtual scatterer; virtual scatterers include target scatterers and reverberation scatterers; the number of virtual scatterers is K, which is the sum of all scattered waves including target and reverberation scatterers;

[0017] Step 1.2: Ignoring the noise, the received echo of the nth (n=1, 2, ..., N) array element is expressed as:

[0018]

[0019] Where K represents the number of virtual scatterers, that is, the total number of scattered waves; σ k represents the complex amplitude of the scattered wave of the kth virtual scatterer; s(t) represents the transmitted signal; t represents time; υ k represents the time delay of the sound wave from the sound source to the kth virtual scatterer and then to the reference receiving element; ξ k =2v k f c / c represents the Doppler frequency shift of the kth virtual scatterer; v k represents the radial velocity of the kth virtual scatterer relative to the reference receiver; f c is the center frequency of the transmitted signal; c represents the underwater sound speed; τ kn =(n-1)dsinθ k / c represents the time delay of the kth virtual scatterer echo received by the nth array element relative to the reference array element; d = c / f c / 2 represents the array element spacing; θ k It represents the azimuth angle of the kth virtual scatterer relative to the horizontal normal of the uniform linear receiving array.

[0020] Furthermore, in step 2, the calculation formula is:

[0021]

[0022] in,

[0023]

[0024] is the self-ambiguity function of the transmitted signal; cn (υ,ξ) represents the delay-Doppler output corresponding to the nth receiving array element, also known as the mutual ambiguity function; υ represents the scanning delay; ξ represents the scanning Doppler frequency shift; [] c Represents the conjugate operation.

[0025] Furthermore, in step 3, the three-dimensional angle-delay-Doppler output of the conventional beamformer is expressed as:

[0026]

[0027] Among them, sinc(θ)=sinθ / θ, u=sinθ, u k = sinθ k ,

[0028] Furthermore, in step 4.1, the angle-delay-Doppler three-dimensional beam diagram of the virtual target echo is expressed as:

[0029]

[0030] Take the absolute value of η(υ,ξ,u): PSF(υ,ξ,u)=|η(υ,ξ,u)|, PSF(υ,ξ,u) can be used as the point spread function of the deconvolution algorithm.

[0031] Furthermore, in step 4.2, the three-dimensional angle-delay-Doppler output absolute value |B(υ,ξ,u)| obtained by steps 1, 2, and 3 is used as observation data of the iterative deconvolution algorithm; and the RL iterative deconvolution algorithm is selected for processing.

[0032] Furthermore, the number of iterations is 5-50 times.

[0033] Furthermore, after obtaining the observation data and the point spread function, the RL algorithm performs the following process:

[0034]

[0035] in, is the reverberation suppression result outputted at the i-th iteration, * indicates the correlation operation, Represents a three-dimensional convolution operation. Among them,

[0036] Effects of the Invention

[0037] The technical effect of the present invention is that: in response to the problems of insufficient multi-dimensional resolution and insufficient sidelobe suppression capabilities faced by traditional methods, a reverberation suppression method based on an iterative deconvolution algorithm that jointly utilizes the three-dimensional resolution capabilities of angle, delay, and Doppler is proposed to enhance the detection capability of active sonar of slow small targets under the interference of strong adjacent reverberation.

[0038] This method transmits narrowband pulses to the target detection area, and then uses a set of matched filters with different Doppler frequency shifts to process the array element echoes to obtain a two-dimensional delay-Doppler output. Conventional beamforming is then performed on the two-dimensional delay-Doppler output of the array element to obtain a three-dimensional angle-delay-Doppler output. Secondly, the above steps are used to obtain the angle-delay-Doppler three-dimensional beam map, and its absolute value is taken as the point spread function of the RL algorithm. Finally, the absolute value and point spread function of the angle-delay-Doppler three-dimensional output are used as the input of the three-dimensional RL algorithm, the number of iterations is selected, and the reverberation suppression result is calculated according to the deconvolution formula.

[0039] The method of the present invention can jointly utilize the inherent energy distribution differences between the target echo and the strong reverberation in three dimensions of angle, time delay, and Doppler to suppress reverberation, improve multidimensional resolution, suppress three-dimensional sidelobe energy leakage, and simultaneously suppress bottom reverberation and surface reverberation, thereby effectively improving the detection performance of the sonar system.

[0040] The basic principle and implementation scheme of the present invention have been verified by computer numerical simulation, and the results show that: compared with the traditional method, the method proposed in the present invention jointly utilizes the inherent energy distribution differences between the moving target echo and the strong reverberation in three dimensions of angle, delay and Doppler to suppress reverberation, effectively improves the multi-dimensional resolution capability of the detection sonar for slow small target echoes and strong reverberation, and at the same time suppresses the multi-dimensional sidelobe energy leakage of the strong reverberation, thereby improving the signal-to-mixture ratio of the target unit. BRIEF DESCRIPTION OF THE DRAWINGS

[0041] Figure 1 This is a schematic diagram of the detection sonar working scene;

[0042] Figure 2 is a flowchart of the entire process of the present invention;

[0043] Figure 3 The absolute value result of the delay-Doppler two-dimensional output of the 96th array element echo in the implementation example;

[0044] Figure 4 is the result obtained by deconvolution processing the absolute value of the delay-Doppler two-dimensional output of the 96th array element echo in the implementation example;

[0045] Figure 5 In the implementation example, the absolute value of the traditional angle-delay-Doppler three-dimensional output is a range-velocity slice at the target angle;

[0046] Figure 6 An angle-range slice of the absolute value of the conventional angle-delay-Doppler three-dimensional output at the target speed in the implementation example;

[0047] Figure 7Angle-velocity slices of the absolute value of the traditional angle-delay-Doppler three-dimensional output at the target distance in the implementation example;

[0048] Figure 8 For the implementation example Figure 5 The result obtained by deconvolution of the distance-velocity slice;

[0049] Fig. 9 The output of the three-dimensional deconvolution process in the implementation example is a range-velocity slice at a target angle;

[0050] Fig.10 An angle-distance slice at the target speed for the output of the three-dimensional deconvolution process in the implementation example;

[0051] Fig.11 It is the angle-velocity slice at the target distance of the output of the three-dimensional deconvolution processing in the implementation example. DETAILED DESCRIPTION

[0052] See also Figure 1-Figure 11 , the main contents of the present invention are:

[0053] (1) Use array element domain data to calculate the traditional angle-delay-Doppler three-dimensional output. The transmitting source emits acoustic energy to the target detection scene, and then the horizontal uniform linear receiving array collects the echo. The echo contains bottom reverberation, surface reverberation and target echo. According to the center frequency of the transmitted signal and the target velocity prior, a set of matched filters with different Doppler frequency shifts are designed. The matched filter group is used to process the echo of each array element. Each array element echo corresponds to a frame of delay-Doppler output. Then the sector to be scanned is discretized, and the array manifold matrix is ​​designed according to the array parameters and the transmission signal parameters. The scanning vectors are selected one by one, and conventional beamforming is performed on the delay-Doppler output of each array element, and finally the traditional three-dimensional beam output is obtained.

[0054] The details of this step include:

[0055] 1) Modeling of array element echo

[0056] The detection scene in this embodiment is as follows Figure 1As shown: a fixed-base transmitting and receiving sonar detects P slow-moving small targets in a reverberation-limited environment; a large number of scatterers are randomly distributed on the bottom and the surface of the water, the bottom scatterers remain stationary, and the surface scatterers move randomly at a low speed. In this embodiment, the detection sonar is simplified to a transmitting source and an N-element horizontal uniform linear array (ULA: uniform linear array). The sound source transmits sound energy to the detection scene, and the receiving array collects echoes. The echo includes the target echo and reverberation. The reverberation includes the bottom reverberation and the surface reverberation. In the array element echo model, the bottom reverberation is regarded as a scattered wave from Q bottom scatterers with different positions, and the surface reverberation is regarded as a scattered wave from D surface scatterers with different positions. The bottom reverberation is considered to have no Doppler effect. The surface reverberation has a Doppler effect. The moving target echo is considered to have a certain Doppler frequency shift. Therefore, there is a large difference in energy distribution between the bottom reverberation and the target echo in the angle dimension, delay dimension, and Doppler dimension. The difference between the water surface reverberation and the target echo in the Doppler dimension is small, but there is a large difference in energy distribution in the delay dimension and the angle dimension. For far-field targets and scatterers, the Doppler differences on each array element are ignored. Therefore, in a reverberation-limited environment, ignoring noise and considering multipath effects, the received echo of the nth (n=1,2,...,N) array element can be expressed as

[0057]

[0058] Where p, q, d, and l represent the pth target, the qth underwater scatterer, the dth surface scatterer, and the lth multipath path, respectively. M represents the number of multipath scattered waves returned from the target or scatterer received by the array element. σ represents the complex amplitude of the multipath scattered wave. α = (c + v) / (cv) represents the Doppler coefficient. Here, c is the speed of sound in water, and v is the radial velocity of the moving target or surface scatterer relative to the receiving element along a certain path. If it is assumed that the target or surface scatterer moves horizontally at a radial velocity v0, the radial velocity corresponding to each multipath path is expressed as is the angle between the multipath and the horizontal plane. μ represents the time delay of the sound wave from the sound source to the target or scatterer along a certain path and then to the receiving element. s(t) represents the narrowband transmission signal.

[0059] Since the transmitted pulse is a narrowband signal and the speed of the target and the water surface scatterer is small, the Doppler effect can be simply regarded as a Doppler frequency shift. At the same time, under far-field conditions, the delay difference of a certain path arrival wave on each array element can be simplified to a phase shift form. Assuming that each multipath arrival wave corresponds to only one virtual scatterer, and the virtual scatterers include the target, the bottom scatterer and the water surface scatterer, then equation (1) is simplified to

[0060]

[0061] Where K represents the number of virtual scatterers, that is, the total number of multipath arrival waves. k represents the kth virtual scatterer. υ represents the delay of the sound wave from the sound source to the virtual scatterer and then to the reference receiving element. ξ=2vf c / c represents Doppler frequency shift. c is the center frequency of the transmitted signal. kn =(n-1)dsinθ k / c represents the delay of the kth virtual scatterer echo received by the nth array element relative to the reference receiving element. d = c / f c / 2 represents the array element spacing. θ represents the azimuth angle of the virtual scatterer relative to the horizontal normal of the ULA.

[0062] It should be pointed out that in a shallow waveguide environment and soft bottom conditions, the Doppler difference of the multipath paths corresponding to the same far-field target is small, and the radial velocity is concentrated near v0, that is, the Doppler difference between the target multipath echo and the bottom reverberation can still be distinguished.

[0063] 2) Calculation of the element echo delay-Doppler two-dimensional output:

[0064] This step follows the traditional delay-Doppler two-dimensional underwater reverberation suppression method. This method mainly uses the Doppler difference between the moving target echo and the underwater reverberation, and directly calculates the mutual ambiguity function of the echo of the nth array element:

[0065]

[0066] in,[·] c represents the conjugate operation. The fuzzy plane is obtained by a set of matched filters with different Doppler frequency shifts. Substituting equation (2) into equation (3), we get

[0067]

[0068] in,

[0069]

[0070] is the self-ambiguity function of the transmitted signal.

[0071] At this point, the delay-Doppler output calculation of the array element echo is completed.

[0072] 3) Regarding the calculation of the traditional angle-delay-Doppler three-dimensional output, the specific contents involved are as follows:

[0073] In this implementation example, the moving target echo and the strong reverberation do not overlap in the angle-delay-Doppler joint domain, that is, the moving target echo and the strong reverberation have inherent resolution potential in the three-dimensional feature domain.

[0074] By observing equation (4), we can find that the phase difference 2πf in the delay-Doppler two-dimensional output of each element echo is c τ kn is not used. This phase difference depends on the azimuth angle. Therefore, the delay-Doppler output of each element echo is processed by CBF to redistribute the energy to the angle dimension. Here is the output of a conventional beamformer

[0075]

[0076] in, represents the complex weight value on the nth receiving element. Further rearranging formula (6) yields

[0077]

[0078] Among them, sinc(θ)=sinθ / θ, u=sinθ, At this point, the traditional angle-delay-Doppler three-dimensional output calculation is completed. After this processing, the angle-dimensional resolution of the target echo and reverberation is utilized, and the signal-to-mixing ratio of the target echo unit is improved. However, when the slow small target echo is adjacent to a strong reverberation scatterer, even if there is an inherent energy distribution difference between the two in the three-dimensional joint domain, the sidelobe energy leakage of the strong scatterer masks the target echo, reflecting the inadequacy of the multidimensional resolution and multidimensional sidelobe suppression capabilities faced by the traditional method when limited by strong reverberation.

[0079] Transform equation (7) to obtain the following form

[0080]

[0081] in,

[0082]

[0083] is the angle-delay-Doppler three-dimensional beam diagram. δ(υ,ξ,u) represents the three-dimensional Dirac function. Represents a 3D convolution operation.

[0084] By observing equation (8), we can find that it is the wide main lobe and high side lobe of η(υ,ξ,u) that blur the inherent three-dimensional joint domain differences between virtual scatterers. Therefore, in order to restore the inherent distinguishability of each dimension, the convolution effect caused by η(υ,ξ,u) must be eliminated or suppressed. The beam pattern represented by equation (9) is easy to obtain. Under the premise of having an angle-delay-Doppler three-dimensional beam pattern, the deconvolution algorithm can be used to suppress the convolution effect in equation (8), improve the three-dimensional resolution performance, restore the inherent characteristic domain differences between the target echo and the reverberation, and effectively suppress the three-dimensional sidelobe energy leakage.

[0085] (2) Perform iterative deconvolution on the absolute value of the traditional angle-delay-Doppler three-dimensional output to obtain the reverberation suppression result. First, calculate the angle-delay-Doppler three-dimensional output of the virtual target echo with no delay, no frequency shift, and zero angle, that is, the angle-delay-Doppler three-dimensional beam diagram, take its absolute value, and obtain the point spread function of the deconvolution algorithm. Then use the absolute value of the traditional angle-delay-Doppler three-dimensional output as the observation data of the RL algorithm. Finally, iteratively calculate according to the deconvolution formula to obtain the reverberation suppression result.

[0086] The details of this step include:

[0087] There are many deconvolution algorithms used in engineering practice. However, considering the robustness of the algorithm, the RL iterative deconvolution algorithm is selected in this implementation example. Because under actual non-ideal conditions, environmental uncertainty and system errors always exist and are difficult to ignore. This places high demands on the robustness of the algorithm. Studies have shown that the robustness of the RL algorithm is comparable to that of CBF and MF, which are recognized as robust algorithms.

[0088] Since the RL algorithm has non-negative constraints, the observation data and the point spread function are non-negative. Here, the absolute value of the traditional three-dimensional output |B(υ,ξ,u)| is used as the observation data of the deconvolution algorithm. At this time, there is approximately the relationship shown in formula (10).

[0089]

[0090] in,

[0091] PSF(υ,ξ,u)=|η(υ,ξ,u)| (11)

[0092] It is the absolute value processing of the angle-delay-Doppler three-dimensional beam pattern, which serves as the point spread function (PSF) of the deconvolution algorithm.

[0093]

[0094] It represents the original distribution of the virtual scatterer in the angle-delay-Doppler three-dimensional feature domain, that is, the result of removing the convolution effect. It can be seen from formula (12) that the original distribution of the target and the underwater scatterer and the surface scatterer in the three-dimensional feature domain can be regarded as the weighted sum of the 3-D Dirac functions carrying different position information. Its shape is similar to the high-resolution low-sidelobe nail-board function. If the original distribution of the target and the underwater scatterer and the surface scatterer in the three-dimensional feature domain does not overlap, the two will show obvious distinguishability. Therefore, the original distribution represented by formula (12) is the ideal reverberation suppression result that the present invention wants to obtain.

[0095] However, it should be pointed out that the approximate relationship presented in formula (10) shows that |B(υ,ξ,u)| cannot be perfectly expressed as the three-dimensional convolution of PSF and the original distribution, and there is a certain approximate error. Therefore, the deconvolution process cannot completely eliminate the convolution effect. A small amount of side lobes remain in the processing result. Therefore, in order to avoid error amplification in the iterative deconvolution process, too high an iteration number cannot be selected in practical applications. In this embodiment, 10-20 times are selected. Experimental verification shows that this deconvolution process can significantly suppress the convolution effect, while improving the multi-dimensional resolution and multi-dimensional sidelobe suppression capabilities.

[0096] After obtaining the observed data and the point spread function, the RL algorithm will perform the process described by equation (13): iteratively estimate the original distribution ρ(υ,ξ,u).

[0097]

[0098] in, is the reverberation suppression result of the i-th iteration, and * indicates the correlation operation. In fact, the iterative estimation process described by equation (13) can be accelerated by fast Fourier transform. At this point, the complete angle-delay-Doppler three-dimensional deconvolution reverberation suppression processing flow is established.

[0099] (3) The angle-delay-Doppler three-dimensional deconvolution processing results under the reverberation background are given by computer numerical simulation, which proves that the method proposed in the present invention has better reverberation suppression performance than the traditional method. This method combines the multi-domain energy distribution difference between the slow small target echo and the strong reverberation to effectively improve the multi-dimensional resolution ability and effectively suppress the three-dimensional sidelobe energy leakage of the strong reverberation, which is more suitable for the detection of slow small targets under the restriction of strong reverberation.

[0100] 1) Signal transmission and echo collection:

[0101] The sonar system adopts a one-element transmission and a 96-element horizontal ULA reception. The ULA adopts a half-wavelength array. The transmission signal takes a CW pulse as an example, with a center frequency of 100kHz and a pulse width of 20ms. The sampling frequency is 250kHz. The speed of sound in water is 1500m / s. Considering the far-field conditions of the shallow waveguide, the bottom is soft. The horizontal detection range is 500m-600m and the width is 100m. Without loss of generality, a single ideal point target is set at a relative distance of 50m. The target approaches the sonar horizontally at a uniform speed of 1m / s. The target azimuth is set to 4°. At the same time, in order to facilitate the demonstration of the reverberation suppression effect of this method, strong bottom scatterers are set at relative distances of 20m, 45m and 80m, marked as 1, 2, 3, and the azimuths are set to -15°, 0° and 25° respectively; a strong surface scatterer is set at a relative distance of 60m, the azimuth is set to 8°, and the radial speed is set to 0.75m / s. A large number of static scatterers are randomly and evenly distributed on the bottom of the water; a large number of dynamic scatterers are randomly and evenly distributed on the surface of the water, and the radial velocity of the scatterers obeys the zero-mean Gaussian distribution. The amplitude of the simulated reverberation sequence obeys the zero-mean Gaussian distribution, and the phase obeys the [0,2π] uniform distribution. The signal-to-mixing ratio in the moving target echo band is set to -16.09dB.

[0102] 2) Delay-Doppler two-dimensional output of array element echo:

[0103] According to the formula υ=2r / c, the horizontal distance r is estimated. The distance scanning interval is 0.036m. According to the formula ξ=2vf c / c, and estimate the radial velocity v. The velocity scanning range is ±10m / s, and the interval is 0.05m / s. It should be noted that in the given reverberation suppression results, the time delay υ is converted to distance, the Doppler frequency shift ξ is converted to velocity, and the sine value u is converted to angle. Therefore, the time delay-Doppler output is converted to distance-velocity output, the u-time delay output is converted to angle-distance output, and the u-Doppler output is converted to angle-velocity output. A set of matched filters carrying different velocity information is designed according to the center frequency of the transmitted signal and the velocity scanning point. The matched filter group is used to scan the distance and velocity of each array element echo in turn to obtain a two-dimensional output of distance-velocity. Figure 3 The absolute value of the range-velocity output of the 96th element echo is given. The position of the target echo is marked by a circle, the position of the strong bottom scatterer generating reverberation is marked by a box, the position of the strong surface scatterer generating reverberation is marked by an ellipse, and the zero velocity line is marked by a dotted line. As can be seen from the figure, the target echo cannot be observed at all within the circle. The sidelobe energy leakage of the strong bottom reverberation and the strong surface reverberation has a serious masking effect on the target echo.

[0104] Since the traditional two-dimensional output of range-velocity superimposes the energy in different directions, the echoes of four strong scatterers at different directions all appear on the same range-velocity plane. In this case, the signal-to-mixing ratio of the unit where the target echo is located is very low, and the detection performance of the sonar system is severely limited. At this time, if the absolute value of the self-blur function is used as the point spread function, and the deconvolution algorithm is used to directly process the absolute value of the range-velocity output of the 96th array element echo, the reverberation suppression result is as follows: Figure 4 As shown in Figure 2, the number of iterations is 15. Figure 4 , we can see that the target cannot be observed at the circle mark. This shows that the range-velocity reverberation suppression method based on deconvolution fails under low signal-to-mixture ratio conditions.

[0105] 3) Angle-Delay-Doppler 3D output:

[0106] In order to improve the signal-to-confusion ratio by jointly utilizing the inherent resolution potential of the angle dimension, the range-velocity output should be processed by CBF. Set the horizontal sector scanning range to u∈[0.5,0.5], and divide it evenly at intervals of 0.005 to obtain non-uniform discrete scanning angles. Design the array manifold according to the array parameters and signal parameters. Select the scanning vectors one by one, perform weighted summation on the range-velocity output of each array element, and obtain the traditional three-dimensional angle-range-velocity output B(υ,ξ,u). For ease of display, the three-dimensional output absolute value |B(υ,ξ,u)| is sliced ​​according to the position parameters of the target. Figure 5 The distance-velocity slice of |B(υ,ξ,u)| at the target angle is given. As can be seen from the figure, after CBF processing, relative to Figure 3 , the energy output at the circle mark is improved. At the same time, the first and third scatterers do not appear in this slice. This shows that the CBF processing partially filters out the reverberation energy in other directions in the angle domain, so the angle-dimensional resolution potential of the target echo and the strong scatterer is utilized. However, it is worth noting that since the azimuth angles of the second scatterer and the strong scatterer on the water surface differ from the target by only 4°, there is still high sidelobe energy leakage in the angle dimension generated by the second scatterer and the strong scatterer on the water surface in this slice. This result shows that the strong reverberation sidelobe energy leakage of CBF needs to be further suppressed. Observe again Figure 5 ,It can be seen that the slow small target echo is seriously confused with the strong sidelobe energy leakage of the bottom scatterer and the surface scatterer, and even merges with the reverberation main lobe, resulting in a large number of false alarms, making it difficult to accurately estimate the target's position parameters.

[0107] Figure 6 This gives the angle-distance slice of |B(υ,ξ,u)| at the target velocity. Figure 6, it can be found that due to the slow speed of the target, the strong side lobes of the velocity dimension of the bottom reverberation have leaked into this velocity channel, resulting in a high reverberation background level. The angle dimension side lobes of the strong bottom scatterers even mask the small target echo. At the same time, there is a large amount of water surface reverberation in this velocity channel. The angle dimension side lobes of the strong water surface scatterers have a masking effect on the target echo.

[0108] Figure 7 The angle-velocity slices of |B(υ,ξ,u)| at the target distance are given. It can be clearly seen from the figure that there is a strong reverberation energy output in the 0° direction and the 8° direction. Its essence is the sidelobe energy leakage of the second scatterer and the strong scatterer on the water surface in the distance dimension. At this time, the target echo at the circle mark is submerged in the strong sidelobe energy leakage. Therefore, in general, the strong sidelobe energy leakage of the reverberation generated by the strong scatterer in each dimension is very likely to mask the slow small target echo, resulting in the inability of the sonar system to detect small targets.

[0109] 4) Three-dimensional deconvolution processing:

[0110] To illustrate the advantages of three-dimensional joint deconvolution reverberation suppression, we first use the absolute value of the self-blur function as a point spread function to directly slice the distance-velocity at the target angle (such as Figure 5 ) is deconvoluted and iterated 15 times. The reverberation suppression result is as follows Figure 8 As shown in the figure, even though the second scatterer and the strong scatterer on the water surface are not in the same azimuth as the small target, there is a high energy output at the box mark and the ellipse mark. The reason is that the angular sidelobe energy leakage of the scatterers in other azimuths is not effectively suppressed, resulting in a decrease in the signal-to-mixing ratio in the angular slice.

[0111] Next, the absolute value of the three-dimensional output of angle-distance-velocity is jointly deconvolved. First, the angle-delay-Doppler three-dimensional beam diagram of the virtual target echo with no delay, no frequency shift and zero angle is calculated, and its absolute value is taken to obtain the point spread function of the deconvolution algorithm. Then, the three-dimensional output of the echo |B(υ,ξ,u)| is used as the observation data of the deconvolution algorithm. The number of iterations i=15 is selected. |B(υ,ξ)u| is used as the initial estimate of the original distribution ρ(υ,ξ,u) of the virtual scatterer, and then the deconvolution formula is used for iterative calculation to finally obtain the reverberation suppression result.

[0112] Fig. 9 Give the target angle after 3D deconvolution It can be clearly seen from this figure that, except for the target echo, there is no other energy output outside the zero-velocity reverberation zone. Figure 8, the angular sidelobe energy output of the second scatterer at the box mark is effectively suppressed, and the intensity of the target echo relative to the scatterer is improved. There is no energy output at the ellipse mark, that is, the angular sidelobe energy leakage of the strong scatterer on the water surface is effectively suppressed. This result shows that after the joint deconvolution processing, the signal-to-mixing ratio of the target echo unit is effectively improved. Figure 5 By comparison, it can be seen that after using the three-dimensional joint deconvolution processing, the sidelobe energy leakage of the bottom reverberation and the surface reverberation is significantly suppressed, the reverberation background level of the entire slice is significantly reduced, and the target echo is no longer confused with the reverberation. Figure 4 , the target echo is well preserved.

[0113] Fig.10 The target velocity after three-dimensional deconvolution is given. Angle-distance slice. From this result, we can see that relative to Figure 6 , the reverberation background level of the slice is significantly reduced. The velocity-dimensional reverberation sidelobe energy of the strong bottom scatterer leaking into the target velocity channel is significantly suppressed. At the same time, the sidelobe energy leakage of the angle-dimensional and distance-dimensional strong surface scatterers is also effectively suppressed, and the target echo is no longer masked by the strong reverberation sidelobe.

[0114] Fig.11 The given value is the target distance after 3D deconvolution. Angle-velocity slice. It can be clearly seen from the figure that compared with Figure 7 , the overall reverberation background level of the slice is significantly reduced. At the same time, the deconvolution processing of the distance dimension effectively suppresses the sidelobe energy leakage of the second scatterer and the strong scatterer on the water surface in the target distance channel. Figure 7 , the energy output ratio between the target echo and the second scatterer is effectively improved. The strong scatterer on the water surface at the ellipse mark has no sidelobe energy leakage in the range dimension.

[0115] In summary, the proposed method effectively utilizes the multi-domain energy distribution difference between the moving target echo and the strong reverberation, and can suppress the bottom reverberation and the surface reverberation at the same time. Even the array processing used in it has a certain inhibitory effect on white noise.

[0116] According to the implementation examples, it can be seen that the method proposed in the present invention can more effectively restore the inherent energy distribution difference between the moving target echo and the strong reverberation in the three-dimensional joint feature domain compared with the traditional method, can significantly improve the multi-dimensional resolution capability at the same time, can significantly suppress the multi-dimensional sidelobe energy leakage at the same time, and effectively improve the signal-to-mixture ratio of the target unit, and is more suitable for slow small target detection under the limitation of strong reverberation.

Claims

1. A reverberation suppression method based on angle-delay-Doppler three-dimensional deconvolution, characterized in that: The following steps are involved: Step 1: Transmit pulses to the detection area, and the horizontal receiving array collects echoes; It includes the following sub-steps: Step 1.1: Initialize the detection area: define a large number of scatterers randomly distributed on the bottom and the surface, where the number of bottom scatterers is Q and their positions are different, and the number of surface scatterers is D and their positions are different; the bottom scatterers remain stationary, and the surface scatterers move randomly at a low speed; the number of moving targets is P; the detection sonar is a transmitting source and an N-element horizontal uniform linear receiving array; the underwater sound source transmits sound energy to the detection area, and the receiving array collects echoes; the echo includes the target echo and reverberation; the reverberation includes the bottom reverberation and the surface reverberation, where the bottom reverberation is the scattered wave of the bottom scatterer, and the surface reverberation is the scattered wave of the surface scatterer; the reverberation scatterer refers to the bottom scatterer and the surface scatterer; one target corresponds to one target scatterer; considering the multipath effect, the number of scattered waves from the pth target to the nth array element is M pn The number of scattered waves from the qth underwater scatterer to the nth array element is M qn The number of scattered waves from the dth water surface scatterer to the nth array element is M dn It is assumed that one scattered wave corresponds to one virtual scatterer; virtual scatterers include target scatterers and reverberation scatterers; the number of virtual scatterers is K, which is the sum of all scattered waves including target and reverberation scatterers; Step 1.2: Ignoring the noise, the received echo of the nth (n=1, 2, ..., N) array element is expressed as: Where K represents the number of virtual scatterers, that is, the total number of scattered waves; σ k represents the complex amplitude of the scattered wave of the kth virtual scatterer; s(t) represents the transmitted signal; t represents time; υ k represents the time delay of the sound wave from the sound source to the kth virtual scatterer and then to the reference receiving element; ξ k =2v k f c / c represents the Doppler frequency shift of the kth virtual scatterer; v k represents the radial velocity of the kth virtual scatterer relative to the reference receiver; f c is the center frequency of the transmitted signal; c represents the underwater sound speed; τ kn =(n-1)dsinθ k / c represents the time delay of the kth virtual scatterer echo received by the nth array element relative to the reference array element; d = c / f c / 2 represents the array element spacing; θ k represents the azimuth angle of the kth virtual scatterer relative to the horizontal normal of the uniform linear receiving array; Step 2: Calculate the delay-Doppler output of the array element echo; the calculation formula is: in, is the self-ambiguity function of the transmitted signal; cn (υ,ξ) represents the delay-Doppler output corresponding to the nth receiving array element, also known as the mutual ambiguity function; υ represents the scanning delay; ξ represents the scanning Doppler frequency shift; [] c represents the conjugate operation; Step 3: Perform conventional beamforming on the delay-Doppler output of each array element; the three-dimensional angle-delay-Doppler output of the conventional beamformer is expressed as: Where, sinc(θ) = sinθ / θ, u = sinθ, u k =sinθ k , Step 4: Deconvolve the conventional angle-delay-Doppler three-dimensional beam output to obtain the reverberation suppression result, including the following sub-steps: Step 4.1: Calculate the angle-delay-Doppler three-dimensional conventional beam output of the virtual target echo with no delay, no Doppler and zero angle, take its absolute value, and obtain the point spread function of the deconvolution algorithm; the angle-delay-Doppler three-dimensional beam diagram of the virtual target echo is expressed as: Take the absolute value of η(υ,ξ,u): PSF(υ,ξ,u)=η(υ,ξ,u)|, PSF(υ,ξ,u) can be used as the point spread function of the deconvolution algorithm; Step 4.2: The absolute value of the three-dimensional beam output obtained by steps 1, 2, and 3 is used as the observation data of the iterative deconvolution algorithm; and the RL iterative deconvolution algorithm is selected for processing; the RL algorithm performs the following process: in, is the reverberation suppression result outputted at the i-th iteration, * indicates the correlation operation, represents a three-dimensional convolution operation, where Step 4.3: Iterate and calculate the reverberation suppression result according to the deconvolution formula.

2. The reverberation suppression method based on angle-delay-Doppler three-dimensional deconvolution according to claim 1, characterized in that: The number of iterations is 5-50 times.

Images