CN101051083B

CN101051083B - Secondary wave arrival direction estimation sonar signal processing method

Info

Publication number: CN101051083B
Application number: CN2006100727187A
Authority: CN
Inventors: 朱维庆; 刘晓东; 张方生; 陈强; 刘光军
Original assignee: Institute of Acoustics CAS
Current assignee: Institute of Acoustics CAS
Priority date: 2006-04-07
Filing date: 2006-04-07
Publication date: 2010-11-03
Anticipated expiration: 2026-04-07
Also published as: CN101051083A

Abstract

The invention discloses a two-dimensional direction of arrival estimation sonar signal processing method, comprising: 1) correcting the phase of the output signal of the primitive; 2) calculating the estimated value of the correlation function of the sonar array; 3) estimating the correlation function of the sonar array 4) Finding the first sub-array operation; 5) Finding the estimated value of the correlation function of all sub-array pairs; 6) Decomposing the eigenvalue of the estimated value of the correlation function of all sub-array pairs; 7) Finding all sub-array pairs The eigenvector estimated value of the array pair; 8) seek the estimated value of the opening angle of the target for the two directions of the sonar array, and select two estimated values of the opening angle corresponding to the minimum variance; 9) seek the estimated number of targets; 10) the second Find the sub-array operation; 11) Find the estimated value of the eigenvectors corresponding to all two-dimensional sonar arrays; 12) Find the estimated values of the opening angles of each target for the two directions of the sonar array, and select the two opening angles corresponding to the minimum variance estimated value. The invention can correctly determine the number of targets and obtain a higher-resolution target three-dimensional sound image.

Description

A Two-dimensional Direction of Arrival Estimation Sonar Signal Processing Method

技术领域technical field

本发明涉及一种二维高分辨率声纳信号处理方法，特别涉及一种利用二维波达方向(Direction of Arrival，以下简称DOA)估计获得水下目标方位信息的声纳信号处理方法。The invention relates to a two-dimensional high-resolution sonar signal processing method, in particular to a sonar signal processing method for obtaining underwater target azimuth information by using two-dimensional direction of arrival (Direction of Arrival, hereinafter referred to as DOA) estimation.

背景技术Background technique

在现有技术中，一般采用常规的波束形成技术在二维平面上形成波束，获得目标的三维声像。常规波束形成技术以付氐变换为基础，它的波束角宽为Δθ＝λ/d，其中λ是波长，d是声纳阵的最大尺寸。如果要使波束角宽Δθ减小，获得高的分辨率，通常有两种方法，一种方法是减小波长λ，频率增大，但是缩短了作用距离；另一种方法是增大声纳阵的尺寸d，但是这样会使得声纳阵造价上升，安装不便。In the prior art, conventional beamforming techniques are generally used to form beams on a two-dimensional plane to obtain a three-dimensional sound image of a target. Conventional beamforming techniques are based on the Fudi transform, and its beam angular width is Δθ=λ/d, where λ is the wavelength and d is the maximum size of the sonar array. If you want to reduce the beam angle width Δθ and obtain high resolution, there are usually two methods. One method is to reduce the wavelength λ and increase the frequency, but shorten the operating distance; the other method is to increase the sonar The size d of the array, but this will increase the cost of the sonar array and make installation inconvenient.

近年来，现代阵信号处理发展很快，利用阵的时空相关函数构成矩阵，由此矩阵提取声波的信息，包括声波入射角和振幅等，这种方法被称为波达方向(DOA)估计方法。声纳阵接收的信号中除了目标信号外，还有噪声信号，因此阵时空相关函数矩阵在泛函空间里分解成信号子空间和噪声子空间，两个子空间相互垂直。波达方向估计方法一般分两大类：一类是谱基方法，它包括噪声子空间法，又称零空间法，但在小样本、低信噪比和高信号相干性时，此类方法的性能明显下降；另一类是参量法，它包括信号子空间法，参量法的性能明显优于谱基方法。In recent years, modern array signal processing has developed rapidly. The matrix is formed by using the space-time correlation function of the array, and the information of the sound wave is extracted from the matrix, including the incident angle and amplitude of the sound wave. This method is called the direction of arrival (DOA) estimation method. . In addition to the target signal, the signal received by the sonar array also has noise signals, so the array space-time correlation function matrix is decomposed into a signal subspace and a noise subspace in the functional space, and the two subspaces are perpendicular to each other. Direction of arrival estimation methods are generally divided into two categories: one is the spectrum-based method, which includes the noise subspace method, also known as the null space method, but in the case of small samples, low signal-to-noise ratio and high signal coherence The other is the parametric method, which includes the signal subspace method, and the performance of the parametric method is obviously better than that of the spectrum-based method.

P.Kraeutner等人的专利，专利号US6130641，专利名称“Imaging methods andapparatus using model-based array signal processing”中采用现代阵信号处理中的零空间法，即噪声子空间，对换能器阵时空相关函数矩阵进行处理，获得了一维波达方向(DOA)估计的声纳阵信号处理方法，该方法比常规波束形成技术的分辨率高。The patent of P.Kraeutner et al., patent No. US6130641, the patent name "Imaging methods and apparatus using model-based array signal processing" adopts the null space method in the modern array signal processing, that is, the noise subspace, which is related to the time and space of the transducer array. The function matrix is processed to obtain a sonar array signal processing method for one-dimensional direction of arrival (DOA) estimation, which has a higher resolution than conventional beamforming techniques.

朱维庆等人的中国专利ZL01142136.3，专利名称“用于测量海底微地貌的高分辨率测深侧扫声纳系统和测量方法”和朱维庆等人的美国专利US6873570B2，专利名称“High resolution bathymetric sonar system and measuring method for measuring thephysiognomy of the seabed”采用现代阵信号处理中的信号子空间法，对换能器时空相关函数矩阵进行处理，获得了一维波达方向(DOA)估计的声纳阵信号处理方法，该方法比常规波束形成技术的分辨率高。The Chinese patent ZL01142136.3 of Zhu Weiqing et al., the patent name is "High-resolution bathymetric side-scan sonar system and measurement method for measuring seabed microtopography" and the US patent US6873570B2 of Zhu Weiqing et al., the patent name is "High resolution bathymetric sonar The system and measuring method for measuring the physiognomy of the seabed" uses the signal subspace method in the modern array signal processing to process the transducer space-time correlation function matrix to obtain the one-dimensional direction of arrival (DOA) estimated sonar array signal processing method that provides higher resolution than conventional beamforming techniques.

但是，上述已有技术中存在的不足在于：已有的声纳信号处理方法中，只用了一维波达方向估计进行声信号处理，另一维是采用常规波束形成技术，从而使得声纳阵探测水下目标的分辨率很低；事实上，在利用声纳阵探测水下目标时，需要二维波达方向(DOA)估计的声信号处理，加上距离轴，才能获得水下目标的高分辨率三维声像。However, the disadvantages of the above-mentioned existing technologies are: in the existing sonar signal processing methods, only one-dimensional direction of arrival estimation is used for acoustic signal processing, and the other dimension is to use conventional beamforming technology, so that the sonar The resolution of detecting underwater targets with a sonar array is very low; in fact, when using a sonar array to detect underwater targets, acoustic signal processing for two-dimensional direction of arrival (DOA) estimation, plus the distance axis, is required to obtain the underwater target high-resolution three-dimensional sound image.

针对现有技术的不足，人们希望有一种二维波达方向估计的声纳信号处理方法。Aiming at the deficiencies of the existing technology, people hope to have a sonar signal processing method for two-dimensional direction of arrival estimation.

发明内容Contents of the invention

本发明的目的在于克服已有技术的不足，提供一种更高分辨率的、实际应用中具有良好性能的二维波达方向估计的声信号处理方法。The purpose of the present invention is to overcome the deficiencies of the prior art, and provide a higher-resolution acoustic signal processing method for two-dimensional direction-of-arrival estimation with good performance in practical applications.

为了达到上述目的，本发明采取如下技术方案：In order to achieve the above object, the present invention takes the following technical solutions:

一种二维波达方向估计声纳信号处理方法，包括如下步骤：A two-dimensional DOA estimation sonar signal processing method, comprising the steps of:

1)对声纳阵基元输出信号的相位进行修正，所述声纳阵为中心对称二维声纳阵；1) Correcting the phase of the output signal of the sonar array primitive, the sonar array being a centrosymmetric two-dimensional sonar array;

2)计算声纳阵相关函数估计值；2) Calculating the estimated value of the correlation function of the sonar array;

3)对上述声纳阵相关函数估计值进行特征值分解；3) Carry out eigenvalue decomposition to above-mentioned sonar array correlation function estimated value;

4)第一求子阵运算：将中心对称的二维声纳阵分成分别沿相互垂直的x和y方向的一维声纳阵，每个一维声纳阵分成多个子阵对，每个子阵对的子阵包含相同数量的基元，执行求子阵运算；4) The first sub-array operation: divide the centrosymmetric two-dimensional sonar array into one-dimensional sonar arrays along the mutually perpendicular x and y directions, each one-dimensional sonar array is divided into multiple sub-array pairs, and each sub-array The sub-arrays of the array pair contain the same number of primitives, and the sub-array operation is performed;

5)计算所有子阵对的相关函数估计值；5) Calculate the estimated value of the correlation function of all sub-array pairs;

6)对所有子阵对的相关函数估计值进行特征值分解；6) Carry out eigenvalue decomposition to the estimated value of the correlation function of all sub-array pairs;

7)计算所有子阵对的特征向量估计值；7) Calculate the estimated value of the eigenvectors of all subarray pairs;

8)计算目标对于x方向和y方向声纳阵的两个张角的估计值，挑选与最小方差对应的两个张角估计值；8) Calculate the estimated values of the target for the two opening angles of the sonar array in the x direction and the y direction, and select two estimated values of the opening angles corresponding to the minimum variance;

9)计算目标数估计，由步骤6)分别求对应x方向和y方向的目标数，如果两方向的目标数不等，则停止运算，否则继续运算；9) Calculating the target number estimate, by step 6) seeking the target numbers corresponding to the x direction and the y direction respectively, if the target numbers in the two directions are not equal, then stop the calculation, otherwise continue the calculation;

10)第二求子阵运算，将二维声纳阵展成一维，按一维声纳阵的方式选择矩阵；10) The second sub-array operation is to expand the two-dimensional sonar array into one dimension, and select the matrix in the mode of one-dimensional sonar array;

11)求所有二维声纳阵对应的特征向量估计值；11) Find the eigenvector estimates corresponding to all two-dimensional sonar arrays;

12)计算各目标对于x方向和y方向声纳阵的两个张角的估计值，挑选与最小方差对应的两个张角估计值。12) Calculate the estimated values of the two opening angles of each target for the sonar array in the x direction and the y direction, and select the two estimated values of the opening angles corresponding to the minimum variance.

在上述技术方案中，所述步骤1)具体包括：In the above technical solution, said step 1) specifically includes:

(1)由多个不相关的信号源，得到声纳阵的输出；(1) Obtain the output of the sonar array from multiple unrelated signal sources;

(2)计算等效基元间距误差；(2) Calculating the equivalent primitive spacing error;

(3)由等效基元间距误差对声纳阵基元输出信号的相位进行修正。(3) The phase of the output signal of the sonar array primitive is corrected by the error of the equivalent primitive spacing.

在上述技术方案中，步骤1)中所述中心对称二维声纳阵的形式是矩形、方框形、十字形或圆形等形状，如图1(a)-(d)所示，声纳阵的基元数在二维对称的情形下可以增减。In the above technical solution, the center-symmetrical two-dimensional sonar array in step 1) is in the form of a rectangle, a box, a cross or a circle, as shown in Figure 1(a)-(d), the acoustic The number of primitives in the nanoarray can increase or decrease in the case of two-dimensional symmetry.

与现有技术相比，本发明的优点在于：Compared with the prior art, the present invention has the advantages of:

1)可以正确判定目标数；1) The number of targets can be correctly determined;

2)利用等效基元间距误差对声纳阵基元输出信号的相位进行修正；2) Correct the phase of the output signal of the sonar array primitives by using the equivalent primitive spacing error;

3)本发明采用多子阵的处理方法，提高目标方位估计精度；3) The present invention adopts the processing method of multiple sub-arrays to improve the target orientation estimation accuracy;

4)获得水下目标的三维声像分辨率更高。4) The resolution of the three-dimensional audio image of the underwater target is higher.

附图说明Description of drawings

图1(a)表示本发明的矩形声纳阵示意图；Fig. 1 (a) represents the rectangular sonar array schematic diagram of the present invention;

图1(b)表示本发明的方框形声纳阵示意图；Fig. 1 (b) represents the block shape sonar array schematic diagram of the present invention;

图1(c)表示本发明的十字形声纳阵示意图；Fig. 1 (c) shows the cross-shaped sonar array schematic diagram of the present invention;

图1(d)表示本发明的圆形声纳阵示意图；Fig. 1 (d) represents the circular sonar array schematic diagram of the present invention;

图2表示本发明实施例示意图；Fig. 2 shows the schematic diagram of the embodiment of the present invention;

图3表示本发明实施例示意图；Fig. 3 shows the schematic diagram of the embodiment of the present invention;

图4本发明一实施例将十字形声纳阵展开后的子阵选择示意图。Fig. 4 is a schematic diagram of sub-array selection after the cross-shaped sonar array is deployed according to an embodiment of the present invention.

具体实施方式Detailed ways

下面结合附图和具体实施方式对本发明作进一步详细描述：Below in conjunction with accompanying drawing and specific embodiment the present invention is described in further detail:

能够再现本发明的中心对称的二维声纳阵形式有很多种，而且声纳阵元的数目不应该受到限制，以二维中心对称的方式可以增减声纳阵元数目。There are many forms of centrosymmetric two-dimensional sonar arrays that can reproduce the present invention, and the number of sonar array elements should not be limited, and the number of sonar array elements can be increased or decreased in a two-dimensional centrosymmetric manner.

本实施例以图2中的十字阵为例，来详细说明本发明的二维波达方向估计声纳信号处理方法。本发明对所有二维中心对称的二维声纳阵均适用，本领域技术人员根据实施例的描述，将本方法应用到其他形式的二维中心对称的声纳阵，是可以胜任的。In this embodiment, the cross array in FIG. 2 is taken as an example to describe the two-dimensional DOA estimation sonar signal processing method of the present invention in detail. The present invention is applicable to all two-dimensional centrosymmetric two-dimensional sonar arrays, and those skilled in the art can apply the method to other forms of two-dimensional centrosymmetric sonar arrays according to the description of the embodiments.

如图1(a)-(d)所示的四种声纳阵：矩形、方框形、十字形和圆形声纳阵，四种声纳阵在x方向和y方向的基元数最多为5个，可分辨4个相干目标。并且，在两个方向最多基元数相等的情形下，十字形状的二维声纳阵的基元数是最少的。如图2所示，本实施例的十字形二维声纳阵在xy平面内，共有9个基元，5号基元为基准点，3、4、5、6、7号基元顺序排在x轴上，构成x方向的一维声纳阵；1、2、5、8、9号基元顺序排在y轴上，构成y方向的一维声纳阵。图2中声纳阵前方有目标A，它的波矢量为K＝k(sinφcosθ sinφsinθ cosφ)^T，φ称为俯仰角，θ称为方位角，k为波数，同时也得到：Four sonar arrays as shown in Figure 1(a)-(d): rectangular, box-shaped, cross-shaped and circular sonar arrays, the four sonar arrays have the largest number of primitives in the x and y directions There are 5, and 4 coherent targets can be distinguished. Moreover, under the condition that the maximum number of primitives in the two directions is equal, the number of primitives of the cross-shaped two-dimensional sonar array is the least. As shown in Figure 2, the cross-shaped two-dimensional sonar array of this embodiment has 9 primitives in the xy plane, the 5th primitive is the reference point, and the 3rd, 4th, 5th, 6th, and 7th primitives are arranged in sequence On the x-axis, a one-dimensional sonar array in the x-direction is formed; primitives 1, 2, 5, 8, and 9 are sequentially arranged on the y-axis, forming a one-dimensional sonar array in the y-direction. In Figure 2, there is a target A in front of the sonar array, and its wave vector is K=k(sinφcosθ sinφsinθ cosφ) ^T , φ is called the pitch angle, θ is called the azimuth angle, k is the wave number, and at the same time:

sinφcosθ＝sinα (1a)sinφcosθ＝sinα (1a)

sinφsinθ＝sinβ (1b)sinφsinθ=sinβ (1b)

称α为x轴声纳阵的张角，以顺时针为正方向；称β为y轴声纳阵的张角，以逆时针为正方向，则K＝k(sinαcosθ sinβ cosθ)^T。Call α the opening angle of the x-axis sonar array, with clockwise as the positive direction; call β the opening angle of the y-axis sonar array, and take the counterclockwise as the positive direction, then K=k(sinαcosθ sinβ cosθ) ^T .

本发明波达方向估计方法解决的技术问题之一是目标数的正确判定。如图3所示，当存在多目标，并且它们之间的连线与x轴或者y轴平行时，会存在目标数的不正确的判定。在图3中，十字形声纳阵布设与图2相同，存在目标1和目标2两个目标，它们的连线与y轴平行，在xy平面的投影关于x轴对称，角α₁和α₂分别是两个目标对x轴的张角，角β₁和β₂分别是两个目标对y轴的张角。由图3可见，α₁和α₂相等，β₁和β₂大小相等，方向相反，由此得出的结果是对x轴的声纳阵只存在一个目标，对y轴的声纳阵存在两个目标。对于上述情况，如果用直接二维波达方向估计方法，则对x轴和y轴声纳阵的张角均为不等的二个值，也即对x轴的声纳阵多出了一个不存在的α₃。这显然是错误的。本发明的解决方法是，首先，将x轴和y轴的均匀线阵分别用一维波达方向估计方法估计张角α和β，如果发现α和β的个数不相等，则此计算即为最终结果。否则，采用二维波达方向估计方法，获得α和β精度较高的估计。One of the technical problems solved by the direction of arrival estimation method of the present invention is the correct determination of the number of targets. As shown in FIG. 3 , when there are multiple targets and the lines between them are parallel to the x-axis or the y-axis, there will be an incorrect determination of the number of targets. In Figure 3, the layout of the cross-shaped sonar array is the same as that in Figure 2. There are two targets, Target 1 and Target 2, and their connection lines are parallel to the y-axis. The projection on the xy plane is symmetrical about the x-axis, and the angles α ₁ and α ₂ are the opening angles of the two targets to the x-axis, and angles β ₁ and β ₂ are the opening angles of the two targets to the y-axis. It can be seen from Figure 3 that α ₁ and α ₂ are equal, β ₁ and β ₂ are equal in size and opposite in direction, and the result obtained from this is that there is only one target for the sonar array on the x-axis, and there is a target for the sonar array on the y-axis. two goals. For the above situation, if the direct two-dimensional direction of arrival estimation method is used, the opening angles of the x-axis and y-axis sonar arrays are both unequal two values, that is, there is one more for the x-axis sonar array α ₃ does not exist. This is clearly wrong. The solution of the present invention is, at first, use the one-dimensional DOA estimation method to estimate the opening angles α and β for the uniform line arrays of the x-axis and y-axis respectively, if it is found that the numbers of α and β are not equal, then this calculation is for the final result. Otherwise, the two-dimensional direction of arrival estimation method is used to obtain estimates of α and β with higher accuracy.

波达方向估计是非线性估计，对硬件有较高的要求，主要是对接收通道的一致性有较高的要求。影响通道一致性的主要因素有几个方面：(1)声纳阵基元相位和幅度的误差，其中主要是相位误差。(2)基元间距的不一致性，它产生基元接收信号的相位误差。(3)入射波阵面的相位起伏。上述三种因素均可形式上等效为基元间距d_l产生了误差Δd_l，在本实施例中，d_l是第l个阵元与第5个阵元即中心阵元的间距。本发明波达方向估计方法解决的技术问题之二是如何由等效基元间距误差对声纳阵基元输出信号的相位进行修正。在本实施例中，构造目标函数：Direction of Arrival estimation is a nonlinear estimation, which has high requirements on hardware, mainly on the consistency of the receiving channel. There are several main factors affecting channel consistency: (1) Phase and amplitude errors of sonar array elements, mainly phase errors. (2) The inconsistency of the cell spacing, which produces the phase error of the cell received signal. (3) The phase fluctuation of the incident wavefront. The above three factors can be formally equivalent to an error Δd _l generated by the distance d _l of the elements. In this embodiment, d _l is the distance between the lth array element and the fifth array element, that is, the central array element. The second technical problem solved by the DOA estimation method of the present invention is how to correct the phase of the output signal of the sonar array element by the error of the equivalent element spacing. In this example, construct the objective function:

$F f ((Δ Δ {d d}_{l l},, Δθ Δθ)) = = min min {| | | | θ θ - - {θ θ}_{A A} | | | |}_{F f}^{22} - - - - - - ((22))$

此处，对于图3中的十字阵，l＝1……8，θ_A＝[θ_A1……θ_Ah]，它是理论入射角，h是测量次数，Δθ是θ的误差角，||·||F是Frobenius范数。找到满足上式的Δd_l和Δθ，则基元的相位误差为：Here, for the cross matrix in Figure 3, l=1...8, θ _A =[θ _A1 ...θ _Ah ], it is the theoretical incident angle, h is the number of measurements, Δθ is the error angle of θ, || ·||F is the Frobenius norm. Find Δd _l and Δθ that satisfy the above formula, then the phase error of the primitive is:

Δφ_l＝kΔd_lsinθ m＝1，……8 (3)此Δφ_l加在相应通道信号的相位上，进行补偿。Δφ _l = kΔd _l sinθ m = 1, ... 8 (3) This Δφ _l is added to the phase of the corresponding channel signal for compensation.

下面的步骤详细说明了如图3所示的十字形声纳阵的二维波达方向估计声纳信号处理方法。当然，本领域技术人员根据下述步骤对其他形式二维中心对阵的声纳阵采用二维波达方向估计声纳信号处理，是可以胜任的。步骤包括：The following steps detail the sonar signal processing method for the two-dimensional DOA estimation of the cross-shaped sonar array shown in FIG. 3 . Of course, those skilled in the art are competent to use the two-dimensional direction of arrival estimation sonar signal processing for other forms of two-dimensional centrally arrayed sonar arrays according to the following steps. Steps include:

1).设有M个不相关的信号源s，得声纳阵输出x，由下式表示1). There are M uncorrelated signal sources s, and the output x of the sonar array is obtained, which is expressed by the following formula

x＝A(α，β)s+n (4)其中x＝A(α, β)s+n (4) where

x＝[x₁ x₂ ……x_L]^T (5a)x=[x ₁ x ₂ ... x _L ] ^T (5a)

A＝(α，β)＝[a(α₁，β₁)……a(α_M，β_M)] (5b)A=(α,β)=[a(α ₁ ,β ₁ )...a(α _M ,β _M )] (5b)

α＝[α₁……α_M]^T (5c)α=[α ₁ ... α _M ] ^T (5c)

β＝[β₁……β_M]^T (5d)β=[β ₁ ... β _M ] ^T (5d)

$a a (({α α}_{i i},, {β β}_{i i})) = = g g (({α α}_{i i},, {β β}_{i i})) [[{e e}^{- - j j (({L L}_{22} - - 11)) kd kd sin sin {β β}_{i i} / / 22} . . . . . . 11,, {e e}^{- - jkd jkd sin sin {β β}_{i i}},, {e e}^{j j (({L L}_{11} - - 11)) kd kd sin sin {α α}_{i i} / / 22},, . . . . . . {e e}^{jkd jkd sin sin {α α}_{i i}},,$

$11,, {e e}^{- - jkd jkd sin sin {α α}_{i i}},, . . . . . .,, {e e}^{- - j j (({L L}_{11} - - 11)) kd kd sin sin {α α}_{i i} / / 22},, {e e}^{jkd jkd sin sin {β β}_{i i}},, . . . . . .,, {e e}^{j j (({L L}_{22} - - 11)) kd kd sin sin {β β}_{i i} / / 22} {]]}^{T T} - - - - - - ((55 e e))$

s＝[s₁……s_M]^T (5f)s=[s ₁ ...s _M ] ^T (5f)

n＝[n₁……n_L]^T (5g)n=[n ₁ ... n _L ] ^T (5g)

其中n是零均值的空间白噪声，g(α_i，β_i)是基元的指向性系数，L是基元数，L₁和L₂分别是x方向和y方向的基元数，k是波数，d为基元间距。where n is the spatial white noise with zero mean value, g(α _i , β _i ) is the directivity coefficient of the primitive, L is the number of primitives, L ₁ and L ₂ are the number of primitives in the x direction and y direction respectively, k is the wave number, and d is the element spacing.

2).求等效基元间距误差Δd_m。目标函数为：2). Calculate the equivalent element spacing error Δd _m . The objective function is:

$F f (({Δd Δd}_{m m},, Δθ Δθ)) = = min min {| | | | θ θ - - {θ θ}_{A A} | | | |}_{F f}^{22},, m m = = 11 . . . . . . L L - - 11 - - - - - - ((66))$

其中θ_A＝[θ_A1……θ_Ah]，它是理论声波入射角，h是测量次数，θ是测得的声波入射角，Δθ是θ的误差角，L是声纳阵基元数，||·||_F是Frobenius范数。求得满足上式的Δd_m和Δθ。Where θ _A = [θ _A1 ... θ _Ah ], it is the theoretical sound wave incidence angle, h is the number of measurements, θ is the measured sound wave incidence angle, Δθ is the error angle of θ, L is the number of sonar array primitives, ||·|| _F is the Frobenius norm. Find Δd _m and Δθ that satisfy the above formula.

3).求声纳阵基元输出信号x_l的相位修正Δφ_l；3). Find the phase correction Δφ _l of the output signal x _l of the sonar array element;

Δφ_l＝kΔd_lsinθ l＝1……L-1 (7)Δφ _l = kΔd _l sinθ l = 1...L-1 (7)

将Δφ_l加在x_l的相位上。Add Δφ _l to the phase of x _l .

4).求声纳阵相关函数估计值即4). Find the estimated value of the correlation function of the sonar array Right now

$\overset{^^}{R R} = = \frac{11}{N N} {Σ Σ}_{n no = = 11}^{N N} {x x}_{n no} {x x}_{n no}^{H h} - - - - - - ((88))$

5).对声纳阵相关函数估计值进行特征值分解：5). Estimated value of sonar array correlation function Perform an eigenvalue decomposition:

$\overset{^^}{R R} = = A A {\overset{^^}{R R}}_{s the s} {A A}^{H h} + + {σ σ}^{22} I I = = {\overset{^^}{U u}}_{s the s} {\overset{^^}{Σ Σ}}_{s the s} {\overset{^^}{U u}}_{s the s}^{H h} + + {\overset{^^}{U u}}_{n no} {\overset{^^}{Σ Σ}}_{n no} {\overset{^^}{U u}}_{n no}^{H h} - - - - - - ((99))$

其中

是信号相关函数的估计值，σ²是噪声方差的估计值，

和

分别是信号的特征向量和特征值的估计值，和

分别是噪声的特征向量和特征值的估计值；上标H表示共轭转置运算；in

is the estimated value of the signal correlation function, ^σ2 is the estimated value of the noise variance,

and

are the estimated eigenvectors and eigenvalues of the signal, respectively, and

are the estimated values of the eigenvectors and eigenvalues of the noise, respectively; the superscript H represents the conjugate transpose operation;

6).第一求子阵运算。将二维十字声纳阵分成x和y方向的一维声纳阵，相应的基元数分别为L₁和L₂，相应的输出为w和v。执行下列求子称运算，得到子阵对：x方向 y方向

输出w由下式表示6). The first sub-array operation. Divide the two-dimensional cross sonar array into one-dimensional sonar arrays in x and y directions, the corresponding numbers of primitives are L ₁ and L ₂ , and the corresponding outputs are w and v. Execute the following sub-scale operations to obtain the sub-array pair: x direction y direction

The output w is represented by

w＝A(α)s+n (12)w＝A(α)s+n (12)

其中in

$w w = = {[[{w w}_{11} {w w}_{22} . . . . . . . . . . . . {w w}_{{L L}_{11}}]]}^{T T} - - - - - - ((1313 a a))$

A(α)＝[a(α₁)……a(α_M)] (13b)A(α)=[a(α ₁ )...a(α _M )] (13b)

α＝[α₁……α_M]^T (13c)α=[α ₁ ... α _M ] ^T (13c)

$a a (({α α}_{i i})) = = g g (({α α}_{i i})) [[{e e}^{- - (({L L}_{11} - - 11)) kd kd sin sin {α α}_{11} / / 22},, . . . . . .,, {e e}^{jkd jkd sin sin {α α}_{i i}},, 11,, {e e}^{- - jk jk sin sin {α α}_{i i}},, . . . . . .,, {e e}^{- - j j (({L L}_{11} - - 11)) kd kd sin sin {α α}_{i i} / / 22} {]]}^{T T} - - - - - - ((1313 d d))$

s＝[s₁……s_M]^T (13e)s=[s ₁ ...s _M ] ^T (13e)

$n no = = {[[{n no}_{11} . . . . . . . . . . . . {n no}_{{L L}_{11}}]]}^{T T} - - - - - - ((1313 f f))$

输出v由下式表示The output v is represented by

v＝A(β)s+n (14)v＝A(β)s+n (14)

其中in

$v v = = {[[{v v}_{11} {v v}_{22} . . . . . . . . . . . . {v v}_{{L L}_{22}}]]}^{T T} - - - - - - ((1515 a a))$

A(β)＝[a(β₁)……a(β_M)] (15b)A(β)=[a(β ₁ )...a(β _M )] (15b)

β＝[β₁……β_M]^T (15c)β=[β ₁ ... β _M ] ^T (15c)

$a a (({β β}_{i i})) = = g g (({β β}_{i i})) [[{e e}^{- - j j (({L L}_{22} - - 11)) kd kd sin sin {β β}_{i i} / / 22},, . . . . . .,, {e e}^{- - jkd jkd sin sin {β β}_{i i}},, 11,, {e e}^{jkd jkd sin sin {β β}_{i i}},, . . . . . .,, {e e}^{- - j j (({L L}_{22} - - 11)) kd kd sin sin {β β}_{i i} / / 22} {]]}^{T T} - - - - - - ((1515 d d))$

s＝[s₁……s_M]^T (15e)s=[s ₁ ... s _M ] ^T (15e)

$n no = = {[[{n no}_{11} . . . . . . . . . . . . {n no}_{{L L}_{22}}]]}^{T T} - - - - - - ((1515 f f))$

由图2可见L₁和L₂可分别为4、3、2，得到多个子阵对。以x方向为例，子阵对包含的基元的形式可为：L＝4，2子阵，由基元3^#到6^#和4^#到7^#构成；L＝3，3子阵，由基元3^#到5^#，4^#到6^#和5^#到7^#构成；L＝2，4子阵，由基元3^#4^#、4^#5^#、5^#6^#和6^#7^#构成。更为详细的技术内容已记载在本申请人的在先申请中，申请号为200510085511.9，名称为“一种高分辨率测深侧扫声纳系统和信号处理方法”。对于y方向采取类似计算。It can be seen from Figure 2 that L ₁ and L ₂ can be 4, 3, and 2 respectively, and multiple sub-array pairs can be obtained. Taking the x direction as an example, the form of the primitives contained in the sub-array pair can be: L=4, 2 sub-arrays, composed of primitives 3 ^# to 6 ^# and 4 ^# to 7 ^# ; L=3, 3 sub-arrays, Consists of primitives 3 ^# to 5 ^# , 4 ^# to 6 ^# and 5 ^# to 7 ^# ; L=2, 4 subarrays, consisting of primitives 3 ^# 4 #, 4 ^# 5 ^# ^, 5 ^# 6 ^# and 6 ^# 7 ^# composition. More detailed technical content has been recorded in the applicant's previous application, the application number is 200510085511.9, and the title is "a high-resolution depth-sounding side-scan sonar system and signal processing method". Similar calculations are taken for the y direction.

7).求所有子阵对相关函数估计值，即7). Find the estimated value of the correlation function of all subarray pairs, namely

${\overset{^^}{R R}}_{W W} = = \frac{11}{N N} {Σ Σ}_{n no = = 11}^{N N} {ww ww}^{H h} - - - - - - ((1616))$

${\overset{^^}{R R}}_{V V} = = \frac{11}{N N} {Σ Σ}_{n no = = 11}^{N N} {v v}_{n no} {v v}_{n no}^{H h} - - - - - - ((1717))$

8).对所有子阵对的

Figure DEST_PATH_GA20191100200610072718701D00013

和

Figure DEST_PATH_GA20191100200610072718701D00014

进行特征值分解8). For all sub-array pairs

and

Do eigenvalue decomposition

${\overset{^^}{R R}}_{11} = = A A ((α α)) {\overset{^^}{R R}}_{s the s 11} A A {((α α))}^{H h} + + {σ σ}^{22} I I = = {\overset{^^}{U u}}_{s the s 11} {\overset{^^}{Σ Σ}}_{s the s 11} {U u}_{s the s 11}^{H h} + + {\overset{^^}{U u}}_{n no} {\overset{^^}{Σ Σ}}_{n no} {\overset{^^}{U u}}_{n no} - - - - - - ((1818))$

${\overset{^^}{R R}}_{22} = = A A ((β β)) {\overset{^^}{R R}}_{s the s 22} A A {((β β))}^{H h} + + {σ σ}^{22} I I = = {\overset{^^}{U u}}_{s the s 22} {\overset{^^}{Σ Σ}}_{s the s 22} {\overset{^^}{U u}}_{s the s 22} + + {\overset{^^}{U u}}_{n no} {\overset{^^}{Σ Σ}}_{n no} {\overset{^^}{U u}}_{n no} - - - - - - ((1919))$

其中各量同公式(8)，只是下标s1和s2分别表示x和y方向的子阵。The quantities are the same as formula (8), except that the subscripts s1 and s2 represent the sub-arrays in the x and y directions respectively.

9).计算所有子阵对的特征向量估计值：9). Calculate the estimated value of the eigenvectors of all subarray pairs:

${\overset{^^}{U u}}_{s the s 1111} = = {J J}_{11} {\overset{^^}{U u}}_{s the s 11}$ ${\overset{^^}{U u}}_{s the s 1212} = = {J J}_{22} {\overset{^^}{U u}}_{s the s 11} - - - - - - ((2020))$

${\overset{^^}{U u}}_{s the s 21 twenty one} = = {J J}_{33} {\overset{^^}{U u}}_{s the s 22}$ ${\overset{^^}{U u}}_{s the s 22 twenty two} = = {J J}_{44} {\overset{^^}{U u}}_{s the s 22} - - - - - - ((21 twenty one))$

10).求所有张角的估计值

和

10). Find the estimated value of all opening angles

and

在约束条件 $({\hat{U}}_{s 12} + Δ {\hat{U}}_{s 12}) &Element; range ({\hat{U}}_{s 11} + Δ {\hat{U}}_{s 11})$ 下：in constraints $({\hat{u}}_{the s 12} + Δ {\hat{u}}_{the s 12}) &Element; range ({\hat{u}}_{the s 11} + Δ {\hat{u}}_{the s 11})$ Down:

Figure DEST_PATH_GA20191100200610072718701D000114

在约束条件 $({\hat{U}}_{s 21} + Δ {\hat{U}}_{s 21}) &Element; range ({\hat{U}}_{s 22} + Δ {\hat{U}}_{s 22})$ 下：in constraints $({\hat{u}}_{the s twenty one} + Δ {\hat{u}}_{the s twenty one}) &Element; range ({\hat{u}}_{the s twenty two} + Δ {\hat{u}}_{the s twenty two})$ Down:

Figure DEST_PATH_GA20191100200610072718701D000116

其中和

Figure DEST_PATH_GA20191100200610072718701D000118

分别是

Figure DEST_PATH_GA20191100200610072718701D000119

和

Figure DEST_PATH_GA20191100200610072718701D000120

的误差。in and

respectively

and

error.

11)比较所有子阵对的

和

Figure DEST_PATH_GA20191100200610072718701D000122

挑选与最小方差对应的

和

Figure DEST_PATH_GA20191100200610072718701D000124

11) Compare all subarray pairs

and

Pick the one that corresponds to the smallest variance

and

12)目标数估计，由步骤8)分别求目标的个数

Figure DEST_PATH_GA20191100200610072718701D000125

i＝1…I，

Figure DEST_PATH_GA20191100200610072718701D000126

i＝1…J，如果I≠J，则停止运算，否则继续运算；12) estimate the number of targets, by step 8) ask the number of targets respectively

i=1...I,

i=1...J, if I≠J, stop the operation, otherwise continue the operation;

13)第二求子阵运算，将图2中的十字声纳阵展开成一维，按一维声纳阵类似的方式得选择矩阵，获得子阵对；13) the second seeks the sub-array operation, expands the cross sonar array in Fig. 2 into one-dimensional, obtains the sub-array pair by selecting the matrix in a similar manner to the one-dimensional sonar array;

本实施例中，x方向选择矩阵：In this embodiment, the x-direction selection matrix:

${J J}_{DX DX 11} = = [\begin{matrix} 00 & 00 & 11 & 00 & 00 & 00 & 00 & 00 & 00 \\ 00 & 00 & 00 & 11 & 00 & 00 & 00 & 00 & 00 \\ 00 & 00 & 00 & 00 & 11 & 00 & 00 & 00 & 00 \\ 00 & 00 & 00 & 00 & 00 & 11 & 00 & 00 & 00 \end{matrix}] - - - - - - ((24 twenty four))$

${J J}_{DX DX 22} = = [\begin{matrix} 00 & 00 & 00 & 11 & 00 & 00 & 00 & 00 & 00 \\ 00 & 00 & 00 & 00 & 11 & 00 & 00 & 00 & 00 \\ 00 & 00 & 00 & 00 & 00 & 11 & 00 & 00 & 00 \\ 00 & 00 & 00 & 00 & 00 & 00 & 11 & 00 & 00 \end{matrix}] - - - - - - ((2525))$

本实施例中，y方向选择矩阵：In this embodiment, the y direction selection matrix:

${J J}_{DY Dy 11} = = [\begin{matrix} 11 & 00 & 00 & 00 & 00 & 00 & 00 & 00 & 00 \\ 00 & 11 & 00 & 00 & 00 & 00 & 00 & 00 & 00 \\ 00 & 00 & 00 & 00 & 11 & 00 & 00 & 00 & 00 \\ 00 & 00 & 00 & 00 & 00 & 00 & 00 & 11 & 00 \end{matrix}] - - - - - - ((2626))$

${J J}_{DY Dy 22} = = [\begin{matrix} 00 & 11 & 00 & 00 & 00 & 00 & 00 & 00 & 00 \\ 00 & 00 & 00 & 00 & 11 & 00 & 00 & 00 & 00 \\ 00 & 00 & 00 & 00 & 00 & 00 & 00 & 11 & 00 \\ 00 & 00 & 00 & 00 & 00 & 00 & 00 & 00 & 11 \end{matrix}] - - - - - - ((2727))$

与公式(24)、(25)、(26)和(27)对应的子阵选择示意图见图4。See Figure 4 for a schematic diagram of sub-array selection corresponding to formulas (24), (25), (26) and (27).

本领域技术人员清楚，由图4可知，与一维情况一样有多种子阵对选择方式，例如x方向，2子阵，由基元3^#到6^#和4^#到7^#构成；3子阵，由基元3^#到5^#，4^#到6^#和5^#到7^#构成；4子阵，由基元3^#4^#、4^#5^#、5^#6^#和6^#7^#构成；y方向类似计算；得出与公式(25)、(26)、(27)和(28)对应的选择矩阵。It is clear to those skilled in the art that, as can be seen from Fig. 4, there are multiple sub-array pair selection modes as in the one-dimensional case, such as the x direction, 2 sub-arrays, consisting of primitives 3 ^# to 6 ^# and 4 ^# to 7 ^# ; 3 sub-arrays Array, composed of primitives 3 ^# to 5 ^# , 4 ^# to 6 ^# and 5 ^# to 7 ^# ; 4 sub-arrays, composed of primitives 3 ^# 4 ^# , 4 ^# 5 ^# , 5 ^# 6 ^# and 6 ^# 7 ^# Composition; similar calculation in the y direction; obtain the selection matrix corresponding to the formulas (25), (26), (27) and (28).

14)计算所有二维十字声纳阵对应的特征向量估计值：14) Calculate the eigenvector estimates corresponding to all two-dimensional cross sonar arrays:

${\overset{^^}{U u}}_{DX DX 11} = = {J J}_{DX DX 11} {\overset{^^}{U u}}_{s the s}$ ${\overset{^^}{U u}}_{DX DX 22} = = {J J}_{DX DX 22} {\overset{^^}{U u}}_{s the s} - - - - - - ((2828))$

${\overset{^^}{U u}}_{DY Dy 11} = = {J J}_{DY Dy 11} {\overset{^^}{U u}}_{s the s}$ ${\overset{^^}{U u}}_{DY Dy 22} = = {J J}_{DY Dy 22} {\overset{^^}{U u}}_{s the s} - - - - - - ((2929))$

15)求所有张角的估计值

和 15) Find the estimated value of all opening angles

and

在约束条件 $({\hat{U}}_{DX 2} + Δ {\hat{U}}_{DX 2}) &Element; range ({\hat{U}}_{DX 1} + Δ {\hat{U}}_{DX 1})$ 下in constraints $({\hat{u}}_{DX 2} + Δ {\hat{u}}_{DX 2}) &Element; range ({\hat{u}}_{DX 1} + Δ {\hat{u}}_{DX 1})$ Down

在约束条件 $({\hat{U}}_{DY 2} + Δ {\hat{U}}_{DY 2}) &Element; range ({\hat{U}}_{DY 1} + Δ {\hat{U}}_{DY 1})$ 下in constraints $({\hat{u}}_{Dy 2} + Δ {\hat{u}}_{Dy 2}) &Element; range ({\hat{u}}_{Dy 1} + Δ {\hat{u}}_{Dy 1})$ Down

16)比较所有子阵对的

和

选择与最小方差对应的

和

16) Compare all subarray pairs

and

Choose the one that corresponds to the minimum variance

and

本领域人员熟知，对于不同形式和不同数目的声纳阵采用本发明的实际应用来说，具体公式会有不同的形式，但实质是相同的。It is well known to those skilled in the art that for the practical application of the present invention for different forms and different numbers of sonar arrays, the specific formulas may have different forms, but the essence is the same.

最后所应说明的是，以上实施例仅用以说明本发明的技术方案而非限制。尽管参照实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，对本发明的技术方案进行修改或者等同替换，都不脱离本发明技术方案的精神和范围，其均应涵盖在本发明的权利要求范围当中。Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention rather than limit them. Although the present invention has been described in detail with reference to the embodiments, those skilled in the art should understand that modifications or equivalent replacements to the technical solutions of the present invention do not depart from the spirit and scope of the technical solutions of the present invention, and all of them should be included in the scope of the present invention. within the scope of the claims.

Claims

1. a secondary wave arrival direction estimution sonar signal processing method comprises the steps:

1) phase place of sonar battle array primitive output signal is revised, described sonar battle array is a center Symmetrical Two-dimentional sonar battle array;

2) calculate sonar battle array related function estimated value;

3) above-mentioned sonar battle array related function estimated value is carried out characteristic value decomposition;

4) first ask the submatrix computing;

5) calculate the right related function estimated value of all submatrixs;

6) the right related function estimated value of all submatrixs is carried out characteristic value decomposition;

7) calculate the right proper vector estimated value of all submatrixs;

8) calculate the estimated value of target, select two the subtended angle estimated values corresponding with minimum variance for two subtended angles of x direction and y direction sonar battle array;

9) calculate number of targets and estimate, ask the number of targets of corresponding orthogonal two coordinate axis x directions and y direction respectively, if the number of targets of two directions does not wait, then stop computing, otherwise continue computing by step 6);

10) second ask the submatrix computing;

11) ask all two-dimentional sonar battle array characteristic of correspondence vector estimated values;

12) calculate the estimated value of each target, select two the subtended angle estimated values corresponding with minimum variance for two subtended angles of x direction and y direction sonar battle array.

2. according to the described secondary wave arrival direction estimution sonar signal processing method of claim 1, it is characterized in that described step 1) is revised the phase place of sonar battle array primitive output signal, specifically comprises the steps:

(1) by a plurality of incoherent signal sources, obtains the output of sonar battle array;

(2) calculate equivalent primitive interval error;

(3) by equivalent primitive interval error the phase place of sonar battle array primitive output signal is revised.

3. according to claim 1 or 2 described secondary wave arrival direction estimution sonar signal processing methods, it is characterized in that, in described step 4), centrosymmetric two-dimentional sonar battle array is divided into respectively one dimension sonar battle array along orthogonal x and y direction, it is right that each one dimension sonar battle array is divided into a plurality of submatrixs, the submatrix that each submatrix is right comprises the primitive of equal number, carries out and asks the submatrix computing.

4. according to the described secondary wave arrival direction estimution sonar signal processing method of claim 3, it is characterized in that, in described step 10),, press the mode selection matrix of one dimension sonar battle array two-dimentional sonar battle array generate one dimension.

5. according to the described secondary wave arrival direction estimution sonar signal processing method of claim 1, it is characterized in that, the form of the Symmetrical Two-dimentional sonar battle array of center described in the step 1) is rectangle, square frame shape, cruciform or circle, and the primitive number of sonar battle array can increase or reduce under the situation of two-dimensional symmetric.

6. secondary wave arrival direction estimution sonar signal processing method, for the centrosymmetric cruciform sonar of two dimension battle array, step comprises:

1). for M incoherent signal source s arranged, get sonar battle array output x, be expressed from the next:

x＝A(α，β)s+n

X=[x wherein ₁x ₂X _L] ^T

A＝(α，β)＝[a(α ₁，β ₁)……a(α _M，β _M)]

α＝[α ₁……α _M] ^T

β＝[β ₁……β _M] ^T

a (α_{i}, β_{i}) = g (α_{i}, β_{i}) [e^{- j (L_{2} - 1) kd \sin β_{i} / 2} \cdot \cdot \cdot 1, e^{- jkd \sin β_{i}}, e^{j (L_{1} - 1) kd \sin α_{i} / 2}, \cdot \cdot \cdot e^{jkd \sin α_{i}},

1, e^{- j kd \sin α_{i}}, \cdot \cdot \cdot, e^{- j (L_{1} - 1) kd \sin α_{i} / 2}, e^{j kd \sin β_{i}}, \cdot \cdot \cdot, e^{j (L_{2} - 1) kd \sin β_{i} / 2}]^{T}

s＝[s ₁……s _M] ^T

n＝[n ₁……n _L] ^T

Wherein n is the space white noise of zero-mean, g (α _i, β _i) be the directive property coefficient of primitive, L is the primitive number, L ₁And L ₂Be respectively the primitive number of x direction and y direction, k is a wave number, and d is the primitive spacing;

Wherein, α _iAnd β _iBe respectively the subtended angle of i target to x axle and y axle;

2). calculate equivalent primitive interval error Δ d _m, objective function is:

F (Δ d_{m}, Δθ) = \min {| | θ {- θ}_{A} | |}_{F}^{2}

m＝1……L-1

θ wherein _A=[θ _A1θ _Ah], θ _ABe theoretical sound wave incident angle, h measures number of times, and θ is the sound wave incident angle that records, and Δ θ is the error angle of θ, and L is a sonar battle array primitive number, || || _FBe the Frobenius norm, calculate the Δ d that satisfies following formula _mWith Δ θ;

3). calculate sonar battle array primitive output signal x _lPhase place correction Δ φ _l

Δφ _l＝kΔd _lsinθ l＝1……L-1

With Δ φ _lBe added in x _lPhase place on;

4). calculate sonar battle array related function estimated value

, promptly

\hat{R} = \frac{1}{N} Σ_{n = 1}^{N} x_{n} x_{n}^{H}

5). to sonar battle array related function estimated value

Carry out characteristic value decomposition:

\hat{R} = A {\hat{R}}_{s} A^{H} + σ^{2} I = {\hat{U}}_{s} {\hat{Σ}}_{s} {\hat{U}}_{s}^{H} + {\hat{U}}_{n} {\hat{Σ}}_{n} {\hat{U}}_{n}^{H}

Wherein

Be the estimated value of signal correction function, σ ²Be the estimated value of noise variance,

With

Be respectively the proper vector of signal and the estimated value of eigenwert,

With It is respectively the estimated value of characteristics of noise vector sum eigenwert; Subscript H represents the conjugate transpose computing;

6). first asks the submatrix computing; Two-dimentional cross sonar battle array is divided into the one dimension sonar battle array of x and y direction, and corresponding primitive number is respectively L ₁And L ₂, be output as w and v accordingly, carry out the following son of asking and claim computing, it is right to obtain submatrix:

The x direction

J_{1} = {[\begin{matrix} I_{L_{1} - 1} & 0 \end{matrix}]}_{(L_{1} - 1) \times L_{1}},

J_{2} = {[\begin{matrix} 0 & I_{L_{1} - 1} \end{matrix}]}_{(L_{1} - 1) \times L_{1}}

The y direction

J_{3} = {[\begin{matrix} I_{L_{2} - 1} & 0 \end{matrix}]}_{(L_{2} - 1) \times L_{2}},

J_{4} = {[\begin{matrix} 0 & I_{L_{2} - 1} \end{matrix}]}_{(L_{2} - 1) \times L_{2}}

Output w is expressed from the next:

w＝A(α)s+n

Wherein:

w = {[\begin{matrix} w_{1} & w_{2} & \cdot \cdot \cdot \cdot \cdot \cdot & w_{L_{1}} \end{matrix}]}^{T}

A(α)＝[α(α ₁)……α(α _M)]

α＝[α ₁……α _M] ^T

a (α_{i}) = g (α_{i}) {[e^{j (L_{1} - 1) kd \sin α_{i} / 2}, \cdot \cdot \cdot, e^{jkd \sin α_{i}}, 1, e^{- jk \sin α_{i}}, \cdot \cdot \cdot, e^{- j (L_{1} - 1) kd \sin α_{i} / 2}]}^{T}

s＝[s ₁......s _M] ^T

n = {{[n}_{1} \cdot \cdot \cdot \cdot \cdot \cdot n_{L_{1}}]}^{T}

Output v is expressed from the next:

v＝A(β)s+n

Wherein

v = {[\begin{matrix} v_{1} & v_{2} & \cdot \cdot \cdot \cdot \cdot \cdot & v_{L_{2}} \end{matrix}]}^{T}

A(β)＝[a(β ₁)……a(β _M)]

β＝[β ₁……β _M] ^T

a (β_{i}) = g (β_{i}) {[e^{- j (L_{2} - 1) kd \sin β_{i} / 2}, \cdot \cdot \cdot, e^{- jkd \sin β_{i}}, 1, e^{jkd \sin β_{i}}, \cdot \cdot \cdot, e^{- j (L_{2} - 1) kd \sin β_{i} / 2}]}^{T}

s＝[s ₁……s _M] ^T

n = {{[n}_{1} \cdot \cdot \cdot \cdot \cdot \cdot n_{L_{2}}]}^{T}

7). calculate all submatrix pair correlation function estimated values, formula is as follows:

{\hat{R}}_{W} = \frac{1}{N} Σ_{n = 1}^{N} w w^{H}

{\hat{R}}_{V} = \frac{1}{N} Σ_{n = 1}^{N} v_{n} v_{n}^{H}

8). right to all submatrixs

With

Carry out characteristic value decomposition

{\hat{R}}_{1} = A (α) {\hat{R}}_{s 1} A {(α)}^{H} + σ^{2} I = {\hat{U}}_{s 1} {\hat{Σ}}_{s 1} U_{s 1}^{H} + {\hat{U}}_{n} {\hat{Σ}}_{n} {\hat{U}}_{n}

{\hat{R}}_{2} = A (β) {\hat{R}}_{s 2} A {(β)}^{H} + σ^{2} I = {\hat{U}}_{s 2} {\hat{Σ}}_{s 2} {\hat{U}}_{s 2} + {\hat{U}}_{n} {\hat{Σ}}_{n} {\hat{U}}_{n}

Wherein subscript s1 and s2 represent the submatrix of x and y direction respectively;

9). calculate the right proper vector estimated value of all submatrixs:

{\hat{U}}_{s 11} = J_{1} {\hat{U}}_{s 1}

{\hat{U}}_{s 12} = J_{2} {\hat{U}}_{s 1}

{\hat{U}}_{s 21} = J_{3} {\hat{U}}_{s 2}

{\hat{U}}_{s 22} = J_{4} {\hat{U}}_{s 2}

10). calculate the estimated value of all subtended angles

With

In constraint condition

({\hat{U}}_{s 12} + Δ {\hat{U}}_{s 12}) &Element; range ({\hat{U}}_{s 11} + Δ {\hat{U}}_{s 11})

Down:

In constraint condition

({\hat{U}}_{s 21} + Δ {\hat{U}}_{s 21}) &Element; range ({\hat{U}}_{s 22} + Δ {\hat{U}}_{s 22})

Down:

Wherein

With

Be respectively

With Error;

11) relatively all submatrixs are right

With , select corresponding with minimum variance

With

12) number of targets is estimated; Ask the number of target respectively by step 8) , i=1 ... I,

, i=1 ... J if I ≠ J then stops computing, otherwise continues computing;

13) second ask the submatrix computing, cross sonar battle array is launched into one dimension, get selection matrix by the mode of one dimension sonar battle array, it is right to obtain submatrix;

14) calculate all two-dimentional cross sonar battle array characteristic of correspondence vector estimated values: