阅读说明：本技术 一种非均匀布阵舷侧阵声纳波束形成快速方法 (Rapid method for forming non-uniform array broadside array sonar wave beam ) 是由 楼万翔 沈文苗 傅仁琦 侯觉 于 2020-12-14 设计创作，主要内容包括：本发明公开了一种非均匀布阵舷侧阵声纳波束形成快速方法,主要过程为首先将频域波束形成每个频点的补相位求和的计算过程转化为两个DFT卷积,然后进行筛选做一个近似计算。本发明进行工程可接受的近似计算,计算量小；对阵元位置没有严格要求,克服了现有快速算法通常要均匀布阵的限制；形成的波束个数与阵元数无关,波束个数及波束指向可任意设置。(The invention discloses a method for quickly forming a beam of a non-uniform array side array sonar. The invention carries out approximate calculation acceptable for engineering, and has small calculation amount; the position of the array element is not strictly required, and the limit that the conventional fast algorithm needs to be uniformly arranged is overcome; the number of formed beams is independent of the number of array elements, and the number of the beams and the beam direction can be set arbitrarily.)

1. A method for quickly forming a sonar wave beam by a non-uniform array broadside array is characterized by comprising the following steps: firstly, converting the calculation process of the complementary phase summation of each frequency point formed by frequency domain wave beams into two DFT convolutions, namely calculating the complementary phase summation process of each frequency point formed by frequency domain wave beams by using two DFT convolutions, and then screening to perform an approximate calculation; the method specifically comprises the following steps:

1) the broadside array is set as a plane array, and the number of horizontal dimensional array elements is n₁The number of elements of the vertical dimension array is n₂The total number of array elements N ═ N₁×n₂The horizontal dimension and the vertical dimension are both non-uniformly arranged, and the coordinate of each array element is (x)_n,y_n,z_n) N is 0,1, …, N-1, wherein N is 0,1, …, N₁-1 is the first row of elements, n ═ n₁,n₁+1,,…,2n₁-1 is the second row of array elements and so on;

2) let L beams be formed, the first beam direction be theta_lAngle of pitch isThe path difference to be compensated for when the nth array element forms the first beam is delta_n,lI.e. by

In which theta is defined_lThe range of the angle is-90 degrees to 90 degrees, the heading of the carrier is-90 degrees, the heading of the stern is 90 degrees,is in the range of-90 to 90 degrees, is positive upwards and negative downwards along the horizontal plane;

3) the output of the p frequency point of the l wave beam is as follows:

in the formula w_nFor the beamforming weighting coefficients, f is the signal frequency, c is the speed of sound, the wavelength is c/f,the above formula is written as follows:

where k is the number of screening points, in the form of a fourier transform, further written as:

in the formula (I), the compound is shown in the specification,in order to perform the convolution operation,representing a discrete fourier transform, let:

obtaining:

4) to pairAnalyzing, except a few, the number of the sequences is larger, the other values are close to 0, screening the sequences, only reserving the value with larger value in the sequences, and setting reserved k to k₁,k₁+1,…,k₂Time of flightThe rest are set to 0, and have:

Technical Field

The invention relates to the field of frequency domain beam forming of a non-uniform array side array in sonar array signal processing, in particular to a rapid method for forming a non-uniform array side array sonar beam.

Background

The conventional beam forming in the frequency domain is one of the processing modules which consume the most computing resources as an essential link of the signal processing flow. In the engineering implementation process, a plurality of DSP boards or a multi-core CPU parallel processing technology is usually adopted to ensure the real-time performance of the system, so that the economic cost and the programming complexity are greatly increased. The side array sonar is one of the most common sonars for submarines, the area of the side array sonar is larger and larger in order to improve the remote warning and passive positioning capacity of the side array sonar, the number of array elements is more and more, and meanwhile, in engineering application, the array is difficult to be arranged strictly according to equal intervals. This presents a significant challenge to the implementation of conventional frequency domain beamforming.

The fast algorithms currently available are mainly the following:

one is to guide the vector (take uniform linear array as an example)Change to(where f denotes frequency, d denotes array element spacing, c denotes sound velocity, m denotes array element number, θ denotes beam angle, k denotes frequency point number, and N denotes beam number) (here, end-fire direction is 0 degree), and the calculation process is converted into FFT operation. The result of each frequency point needs to be corrected according to the frequency and the array element spacing, and the realization of broadband beam forming is relatively limited, and the method is generally only suitable for narrow bands.

And the other is based on the idea of FRFT, the beam forming result is converted into a fractional Fourier transform (FRFT) form by introducing a Brnstein formula, and then FFT is adopted for realization. The implementation process of the algorithm is not described in detail, and in the process of calculating the FRFT by using the FFT, if the FRFT is directly realized by adopting two FFTs and one FFT-1 with the same length point as the data length, the obtained result has certain difference with the direct calculation because the difference between the linear convolution and the circular convolution is not considered, and the adverse effect is generated on the subsequent processing of beam listening, target tracking identification and the like.

Converting the calculation process of the complementary summation of each frequency point formed by linear array frequency domain wave beams into linear convolution operation, calculating the result after M (array element number) points by using two times of FFT and one time of IFFT (inverse FFT), constructing a new sequence by adopting data translation, calculating the result of the previous M points by using one time of FFT and one time of IFFT, splicing the two results to obtain the result of all directions, and finally summing the squares of the result modules of each frequency point to obtain the wave beam forming result. However, this method requires the array to be arranged at exactly equal intervals and the number of beams must be greater than two array elements.

Disclosure of Invention

The invention aims to overcome the defects in the prior art and provide a method for quickly forming a non-uniformly arrayed broadside array sonar wave beam.

The purpose of the invention is achieved by the following technical scheme: firstly, converting a calculation process of complementary phase summation of each frequency point formed by frequency domain wave beams into two DFT convolutions, namely calculating the complementary phase summation process of each frequency point formed by the frequency domain wave beams by using two DFT convolutions, and then screening to perform approximate calculation; the method comprises the following steps:

1) the broadside array is set as a plane array, and the number of horizontal dimensional array elements is n₁The number of elements of the vertical dimension array is n₂The total number of array elements N ═ N₁×n₂The horizontal dimension and the vertical dimension are both non-uniformly arranged, and the coordinate of each array element is (x)_n,y_n,z_n) N is 0,1, …, N-1, wherein N is 0,1, …, N₁-1 is the first row of elements, n ═ n₁,n₁+1,,…,2n₁-1 is the second row of array elements and so on;

2) let L beams be formed, the first beam direction be theta_lAngle of pitch isThe path difference to be compensated for when the nth array element forms the first beam is delta_n,lI.e. by

In the formula [ theta ]_lThe range of the angle is-90 degrees to 90 degrees, the heading of the carrier is-90 degrees, the heading of the stern is 90 degrees,is in the range of-90 to 90 degrees, is positive upwards and negative downwards along the horizontal plane;

3) the output of the p frequency point of the l wave beam is as follows:

in the formula w_nFor the beamforming weighting coefficients, f is the signal frequency, c is the speed of sound, the wavelength is c/f,the above formula is written as follows:

where k represents the number of screening points, in the form of a fourier transform, further written as:

in the formula (I), the compound is shown in the specification,for convolution operations, let:

obtaining:

4) to pairAnalyzing, except a few, the number of the sequences is larger, the other values are close to 0, screening the sequences, only reserving the value with larger value in the sequences, and setting reserved k to k₁,k₁+1,…,k₂Time of flightThe rest are set to 0, and have:

the invention has the beneficial effects that: the invention firstly converts the calculation process of the complementary phase summation of each frequency point formed by frequency domain wave beams into two DFT convolutions, then carries out screening to carry out approximate calculation, carries out approximate calculation acceptable by engineering, and has small calculation amount; the position of the array element is not strictly required, and the limit that the conventional fast algorithm needs to be uniformly arranged is overcome; the number of formed beams is independent of the number of array elements, and the number of the beams and the beam direction can be set arbitrarily.

Drawings

FIG. 1 is a schematic view of a broadside array according to the present invention.

Fig. 2 is a diagram comparing the beamforming results of the present invention and the conventional algorithm at azimuth 0 ° and frequency 1.5 kHz.

FIG. 3 is a diagram illustrating the beamforming results of the conventional algorithm of the present invention at an azimuth angle of-45 ° and a frequency of 1.5 kHz.

Fig. 4 is a diagram comparing the beamforming results of the present invention and the conventional algorithm at azimuth 0 ° and frequency 0.9 kHz.

FIG. 5 is a diagram illustrating the beamforming results of the conventional algorithm of the present invention at an azimuth angle of-45 ° and a frequency of 0.9 kHz.

Fig. 6 is a diagram comparing the beamforming results of the present invention and the conventional algorithm at azimuth 0 ° and frequency 0.3 kHz.

FIG. 7 is a diagram illustrating the beamforming results of the conventional algorithm of the present invention at an azimuth angle of-45 ° and a frequency of 0.3 kHz.

Detailed Description

The invention will be described in detail below with reference to the following drawings:

as shown in fig. 1, the method for quickly forming the beam of the non-uniformly arrayed broadside array sonar mainly comprises the following steps:

1) the broadside array is set as a plane array, and the number of horizontal dimensional array elements is n₁The number of elements of the vertical dimension array is n₂The total number of array elements N ═ N₁×n₂The horizontal dimension and the vertical dimension are both non-uniformly arranged, and the coordinate of each array element is (x)_n,y_n,z_n) N is 0,1, …, N-1, wherein N is 0,1, …, N₁-1 is the first row of elements, n ═ n₁,n₁+1,,…,2n₁-1 is the second row of array elements and so on;

2) let L beams be formed, the first beam direction be theta_lAngle of pitch isThe path difference to be compensated for when the nth array element forms the first beam is delta_n,lI.e. by

In the formula [ theta ]_lThe range of the angle is-90 degrees to 90 degrees, the heading of the carrier is-90 degrees, the heading of the stern is 90 degrees,is in the range of-90 to 90 degrees, is positive upwards and negative downwards along the horizontal plane;

3) the output of the p frequency point of the l wave beam is as follows:

in the formula w_nFor beamforming weighting coefficients, f is the signal frequency, c is the speed of sound, and the wavelength isλ＝c/f，For the output of the p frequency point of the n array element, the above formula is rewritten as follows:

where k represents the number of screening points, in the form of a fourier transform, further written as:

in the formula (I), the compound is shown in the specification,for convolution operations, let:

obtaining:

4) to pairAnalyzing, except that a few numbers are large, the other values are close to 0, screening the values, and setting the reserved k as k₁,k₁+1,…,k₂Time of flightThe other is set to 0, and the approximate pairErrors caused by the beamformed output are negligible, so there are:

in the formula k₂-k₁+1<N, therefore, the calculation amount can be greatly saved by adopting the formula. Screening points k₁,k₂The specific values of (a) are related to the beam number, the operating frequency, the number of array elements, etc., but the values can be calculated in advance.

The calculation process of the invention is as follows:

1) calculating the beam forming compensation coefficient of the No. l beam of the p frequency point

2) Careful analysisOnly the larger value in the sequence is reserved to form a new compensation coefficient sequencek′∈[k₁，k₂]；

3) Calculating compensation coefficients of all preformed beams of all frequency points, and storing the compensation coefficients in a database;

4) to the received array element signal r_n(t) performing FFT to obtain a frequency domain signal

5) The frequency domain signal is subjected to an FFT,

6) computing the beam domain frequency domain output as

The performance comparison analysis of the conventional algorithm and the rapid algorithm of the invention is as follows:

the conventional algorithm is to complement the phase of the array element data of each frequency point and then sum the phase compensated results of the array elements. The multiplication of complex numbers in the process of beam forming calculation is the key influencing the speed of the algorithm, and the influence of the complex number multiplication is mainly considered here. Assuming that the number of array elements is N, the number of beams is L, and the number of frequency points is p, the frequency domain conventional beamforming complex multiplication is p × N × L. The fast algorithm adds n-point FFT's one time, so the number of complex multiplications is p × n/2 × log₂(n)+p*L*(k₂-k₁+1) it can be seen that the number of multiplications of the fast algorithm is independent of the number of array elements. Comparing the conventional algorithm with the fast algorithm, we can get:

from the above formula, it can be seen that the performance of the algorithm is improved regardless of the frequency point number, in practice, complex addition and some other operations are involved, and the actual operation speed may be slightly different. The above formula is rewritten as:

as can be seen from the above formula, the fast algorithm is suitable for matrixes with a large number of array elements or beams.

Example (b): the space broadside array is composed of four linear arrays, each linear array is composed of 140 linear arrays which are arranged in a non-uniform mode, the total number of the linear arrays is 560, and the shape of the linear arrays is shown in figure 1.

181 wave beams are preformed, the screening points in the rapid algorithm are 61, the wave beam forming results under different angles and different frequencies of the conventional algorithm and the rapid algorithm are shown in figure 2, the frequencies are 1.5kHz, 0.9kHz and 0.3kHz respectively, and the angles are 0 degree to 45 degrees respectively. As can be seen from FIG. 2, there is a slight error in the beam results of the fast algorithm and the conventional algorithm at different angles, which is caused by the approximate calculation of the screening pointsAnd (4) obtaining the finished product. The errors of the two methods are very small, the target directions of the two methods are completely consistent, the width of a main lobe and the height of a side lobe are completely consistent, and the engineering use is not influenced. Comparing the operation time of different array elements and different beam numbers of the two algorithms, wherein the operation times are calculated according to the number of timesAnd (4) calculating the formula. And carrying out simulation calculation by utilizing matlab 7.11.0 to obtain the operation time.

TABLE 1 comparison of operation times and time of different array elements, pre-forming 181 wave beams, frequency point number 1

TABLE 2 comparison of the number of different beams, array element number 280(1/2 array), frequency point number 1, and the number of times of operation and time of the two methods

As can be seen from tables 1 and 2, the operation times and operation time of the fast algorithm are less than those of the conventional algorithm. When the number of array elements is multiplied, the operation speed of the fast algorithm is obviously superior to that of the conventional algorithm.

It should be understood that equivalent substitutions and changes to the technical solution and the inventive concept of the present invention should be made by those skilled in the art to the protection scope of the appended claims.