CN104181499A - Ranging passive location method under azimuth angle prior condition based on linear sparse arrays - Google Patents

Abstract

The present invention provides a passive positioning method for distance measurement based on a sparse line array under the azimuth prior condition, combining the characteristics of the passive positioning multi-station system and the single-station system, by using a sparse line array in a small space On the basis of the arrangement of the observation station, the observation station is used to receive the data phase information to realize the ranging and positioning of a near-field radiation source. On the basis of careful analysis of the phase information structure of the observation data, the phase compensation is performed on the phase ambiguity part of the phase information of the observation data by using the prior azimuth information, and the array is established by using the part of the phase information of the observation data without phase ambiguity The signal model is based on the MUSIC algorithm to complete the estimation of the distance of the near-field radiation source.

Description

A range-finding passive positioning method based on sparse wiring array under the condition of azimuth angle prior

technical field

The invention relates to electronic information technology, in particular to a near-field radiation source passive positioning technology.

Background technique

Passive positioning refers to passively receiving electromagnetic wave signals from near-field radiation sources under the condition that the observation equipment itself does not radiate electromagnetic waves, and using the space in the sense of direction and distance between the observation equipment and the radiation source in a certain coordinate system With respect to the geometric relationship and its changes, by measuring the angle of arrival, time difference of arrival and/or changes in other corresponding physical quantities carried in the electromagnetic wave signal of the radiation source, the target position, velocity and acceleration of the radiation source are calculated, and the target positioning is finally realized covertly .

The passive positioning system is divided into multi-station passive positioning and single-station passive positioning.

The multi-station passive positioning method is mainly to use multiple observation stations set up in a large space to intercept the signals emitted by the target radiation source and calculate the corresponding observations, such as the arrival time difference between any two stations, etc., thus estimating The position of the target radiation source and other information. Under the multi-station passive positioning system, the time synchronization requirements between observation stations are relatively strict, and a large amount of data transmission and fusion are usually required, which makes the multi-station positioning system generally very complicated; in addition, there may be multiple observations in the multi-station positioning system Stations cannot successfully detect signals emitted by radiation sources that use low probability of intercept technologies such as adaptive beamforming, which directly leads to the failure of the multi-station passive positioning system.

The single-station passive positioning method, such as the single-station positioning method based on observations such as the direction of arrival of radiation sources, does not require complex time synchronization and data fusion between multiple observation stations. However, because the observations in the single-station method often do not fully contain the position information of the radiation source, this makes the single-station passive positioning method generally have observability problems. In order to solve the observability problem of the single-station positioning system, the single-station passive positioning method usually requires multiple measurements and long observation time intervals. In certain occasions, in order to solve the problem of observability, the single-station passive positioning system even needs to perform a specific form of motion. Therefore, the observability problem often leads to limitations in the application of various single-site passive positioning methods.

Contents of the invention

The technical problem to be solved by the present invention is to provide a more stable and effective passive positioning method for near-field radiation sources based on a batch of measurement data.

The technical solution adopted by the present invention for solving the above-mentioned technical problems is that the ranging passive positioning method based on sparse wiring array under the azimuth priori condition comprises the following steps:

Step 1. Set M observation stations on a straight line within a small range, select one observation station as No. 0 observation station, and the distance d _m between other observation points and No. 0 observation station satisfies Among them, M is the total number of observation stations in the system, M≥2, m is the number of observation stations m=0,...,M-1, λ is the working wavelength of the radiation source, θ is the azimuth arrival angle of the radiation source, R ₀ represents the radiation The distance from the source to the observation station No. 0;

Step 2. Receive the observation data at M observation stations, calculate the estimated value θ of the radiation source azimuth arrival angle, and generate the observation data matrix x(t), x(t)=a(θ,R ₀ )s(t)+ n(t), t=1,2,...,N; among them, t is the sampling time, t=1,2,...,N, N is the total sampling length, s(t) the observation station receives signal vector, n(t) is Gaussian white noise vector, a(θ,R ₀ ) is radiation source array manifold vector, $a(\mathrm{\θ},R_0)={[1,e^{-j2\mathrm{\π}{d}_1\sin \mathrm{\θ}/\mathrm{\λ}+\mathrm{j\π}{\left(d_1\cos \mathrm{\θ}\right)}^2/\left(\mathrm{\λ}{R}_0\right)},...,e^{{-j2\mathrm{\πd}}_{m-1}\sin \mathrm{\θ}/\mathrm{\λ}+\mathrm{j\π}{\left(d_{m-1}\cos \mathrm{\θ}\right)}^2/\left(\mathrm{\λ}{R}_0\right)}]}^T,$ (·) ^T is matrix transposition operation;

Step 3. Determine Φ using the estimated value θ of the angle-of-arrival for the calculated radiation source azimuth, diag represents a diagonal matrix, and the phase compensation matrix Φ is left-multiplied to the observation data matrix x(t) to complete the phase compensation;

Step 4. Find the autocovariance matrix of the data matrix obtained by phase compensation The autocovariance matrix Perform matrix eigenvalue decomposition, delete the eigenvector corresponding to the largest eigenvalue in the eigenvector matrix corresponding to the decomposed eigenvalue, and then obtain the noise subspace matrix G, and use the noise subspace matrix G to search for the MUSIC spectrum peak, thereby obtaining The distance R ₀ from the radiation source corresponding to the peak value of the criterion function to the observation station No. 0 is the estimated value of the distance of the near-field radiation source; the criterion function P _MUSIC (R ₀ ) used for MUSIC spectrum peak search is: P _MUSIC (R ₀ )=1 /||G ^H Φa(θ,R ₀ )|| ² , where ||·|| ² is the square of the 2-norm.

The invention combines the respective characteristics of the passive positioning multi-station system and the single-station system, and provides a near-field radiation source passive positioning method between the single-station and the multi-station. On the basis of arranging observation stations in the form of sparse linear arrays in a small space, the observation stations are used to receive data phase information to realize the ranging and positioning of a near-field radiation source. The invention uses the phase information of the observation data to measure the distance, which is completely different from the observations such as the time difference of arrival adopted by the traditional passive positioning method. Since the observation station array is sparsely arranged in the airspace, the array manifold vector described by the two parameters of the azimuth arrival angle and the distance of the near-field radiation source usually has a phase ambiguity problem. On the basis of careful analysis of the phase information structure of the observation data, the phase compensation is performed on the phase ambiguity part of the phase information of the observation data by using the prior azimuth information, and the array is established by using the part of the phase information of the observation data without phase ambiguity The signal model is based on the MUSIC algorithm to complete the estimation of the distance of the near-field radiation source.

Further, after the initial estimation of the distance of the near-field radiation source is completed based on the MUSIC algorithm, a high-precision maximum likelihood method is used to iteratively obtain a high-precision estimation of the distance of the near-field radiation source. That is, step 5 is also included after step 4:

Step 5: Using the estimated value of the radiation source azimuth angle of arrival θ and the estimated value of the near-field radiation source distance as initial values, according to the maximum likelihood method, use quasi-Newton iteration to obtain the radiation source azimuth arrival angle and the near-field radiation source distance. dimension parameter vector High-precision estimate of ;

The criterion function of the maximum likelihood method is:

$\stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}<span\; class="notranslate">\; =</span>{<span\; class="notranslate">\; [</span>\stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; \&theta;</span>}}<span\; class="notranslate">\; ,</span>{\stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; R</span>}}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; ]</span>}^{\mathrm{<span\; class="notranslate">\; T</span>}}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; arg</span>}\underset{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}<span\; class="notranslate">\; =</span>{<span\; class="notranslate">\; [</span>\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; ]</span>}^{\mathrm{<span\; class="notranslate">\; T</span>}}}{\mathrm{<span\; class="notranslate">\; min</span>}}\mathrm{<span\; class="notranslate">\; f</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&alpha;</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; f</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&alpha;</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; log</span>}<span\; class="notranslate">\; |</span>{\stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}\mathrm{<span\; class="notranslate">\; a</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; a</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; +</span>{\stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}{\mathrm{<span\; class="notranslate">\; I</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}<span\; class="notranslate">\; |</span><span\; class="notranslate">\; ;</span>$

is the high-precision estimate of the angle of arrival of the azimuth of the radiation source, is a high-precision estimate of the distance to the near-field radiation source, Indicates the value corresponding to α when the objective function f(α) takes the minimum value; |·| is the determinant symbol, is the estimated value of the radiation source power, is the estimated value of noise power, and I _M is the identity matrix of M×M; is the pseudo-inverse of the radiation source array manifold vector a(θ,R ₀ ), (·) ^H is the matrix conjugate transpose operation, is the sample autocovariance matrix, $\hat{R}=\underset{t=1}{\overset{N}{\mathrm{\&Sigma;}}}x\left(t\right)x^h\left(t\right)/N;$

Tr{·} is to take the matrix trace, is the vertical projection operator matrix of the orthogonal complement subspace of a(θ,R ₀ ), ${\mathrm{\&Pi;}}_a^{\&perp;}=I_m-a(\mathrm{\&theta;},R_0){\left[a^h(\mathrm{\&theta;},R_0)a(\mathrm{\&theta;},R_0)\right]}^{-1}{a}^h(\mathrm{\&theta;},R_0);$

The quasi-Newton iterative process is:

α ^k+1 = α ^k -μ _k H ^-1 (α ^k )f'(α ^k );

Among them, α ^k is the estimated value of the radiation source array manifold vector α obtained by the k-th iteration, μ _k is the iteration step size of the k-th step, H is the Hessian Hessian matrix of the objective function f(α) or its positive definite approximation Matrix, f′(α ^k ) is the gradient of the objective function; the iteration termination criterion is ||H ^-1 (α ^k )f′(α ^k )|| is less than the specified value ε and the matrix H is positive definite, or the number of iterations reaches the preset maximum number of times;

Step 6: Use the estimated value θ of the radiation source azimuth angle of arrival and the high-precision estimated value R of the near-field radiation source distance Complete radiation source location.

The beneficial effect of the present invention is that it utilizes the same space dimension as the multi-station system, has the characteristics of rapid deployment of a single-station passive positioning system, and can realize distance measurement and positioning of near-field radiation sources according to a batch of measurement data.

Description of drawings

Fig. 1 is a geometric illustration of the passive positioning problem in the present invention in a right-handed Cartesian coordinate system. Considering that M observation stations are sparsely arranged in a small area on a straight line in the space plane (the farthest distance between observation stations does not exceed a certain threshold), and a near-field radiation source T is received from the azimuth angle θ.

Fig. 2 is the statistical result of the root mean square error of ranging measured under different azimuth prior information conditions by adopting the specific example mode of the present invention, the ordinate is the root mean square value of the ranging error, and the unit is meter, and the abscissa is the azimuth θ ＝-70°,-50°,-30°,-10°,0°,20°,30°,40°,60°,70°, the unit is degree.

Fig. 3 is the statistical result of the relative ranging root mean square error measured under different azimuth prior information conditions for adopting the specific example mode of the present invention, and the ordinate is the relative ranging error root mean square value (proportional value, so there is no unit), The abscissa is the azimuth angle θ=-70°, -50°, -30°, -10°, 0°, 20°, 30°, 40°, 60°, 70°, and the unit is degree.

Detailed ways

As shown in Figure 1, select an observation station as the origin of the coordinates, record it as the No. 0 observation station, and set the line where the observation station array is located as the X-axis, and establish XY in the space plane formed by the near-field radiation source and the X-axis Right-handed Cartesian coordinate system. The angle swept by the ray OT between the near-field radiation source and the origin of the coordinates in the counterclockwise direction to the coordinate axis Y is defined as the positive azimuth arrival angle.

The used key equation of the inventive method is:

(R _m ) ² ＝(R ₀ ) ² +(d _m ) ² -2R ₀ d _m sinθ

This equation is obtained according to the law of cosines in the triangle formed by "radiation source-coordinate origin-observation station m". Among them, R _m represents the distance from the radiation source to the m-th observation station, the value of the subscript m is m=0,1,...,M-1, and M represents the number of observation stations. R ₀ represents the distance from the radiation source to the observation station No. 0. d _m represents the distance from the mth observation station to the 0th observation station, and the value of the subscript m is m=0,1,...,M-1. θ(-90°<θ≤90°) is the azimuth arrival angle of the radiation source, referred to as azimuth angle.

The used radiation source of the inventive method is to the distance of the m observation station:

${\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}<span\; class="notranslate">\; =</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; +</span>{<span\; class="notranslate">\; (</span>\frac{{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}}{{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}\frac{{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}}{{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}}\mathrm{<span\; class="notranslate">\; sin</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; /</span>\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; \&ap;</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; [</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>\frac{{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}}{{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}}\mathrm{<span\; class="notranslate">\; sin</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; +</span>\frac{\mathrm{<span\; class="notranslate">\; 1</span>}}{\mathrm{<span\; class="notranslate">\; 2</span>}}{<span\; class="notranslate">\; (</span>\frac{{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}\mathrm{<span\; class="notranslate">\; cos</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}}{{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; ]</span>$

This equation is based on the principle that x=(d _m /R ₀ ) ² -2(d _m /R ₀ ) sinθ is much smaller than 1, and the R _m obtained by the law of cosines is expanded and approximated by Taylor series. Among them, R _m represents the distance from the radiation source to the m-th observation station. R ₀ represents the distance from the radiation source to the observation station No. 0. d _m represents the distance from the mth observation station to the 0th observation station. θ represents the azimuth angle of the radiation source.

In the method of the present invention, the time delay between the spherical wave front of the radiation source from the m observation station to the coordinate origin is:

${\mathrm{<span\; class="notranslate">\; \&tau;</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}<span\; class="notranslate">\; =</span>\frac{\sqrt{{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; +</span>{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}\mathrm{<span\; class="notranslate">\; sin</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}}{\mathrm{<span\; class="notranslate">\; c</span>}}<span\; class="notranslate">\; \&ap;</span>\frac{<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}\mathrm{<span\; class="notranslate">\; sin</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}}{\mathrm{<span\; class="notranslate">\; c</span>}}<span\; class="notranslate">\; +</span>\frac{{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}\mathrm{<span\; class="notranslate">\; cos</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}{\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; c</span>}{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}}$

Among them, τ _m represents the time delay between the radiation source signal wavefront from the mth (m=1,...,M-1) observation station to the coordinate origin. Represents the square root operation. R ₀ represents the distance from the radiation source to the observation station No. 0. c stands for the speed of light. d _m represents the distance from the mth observation station to the 0th observation station. θ represents the azimuth angle of the radiation source.

The used observation station data of the inventive method is:

$\begin{array}{c}{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; e</span>}}^{{\mathrm{<span\; class="notranslate">\; j</span>}\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; \&pi;f</span>}}_{\mathrm{<span\; class="notranslate">\; c</span>}}{\mathrm{<span\; class="notranslate">\; \&tau;</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; no</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>\\ <span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; e</span>}}^{{<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; j</span>}\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; \&pi;d</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}\mathrm{<span\; class="notranslate">\; sin</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; /</span>\mathrm{<span\; class="notranslate">\; \&lambda;</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; j\&pi;</span>}{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}\mathrm{<span\; class="notranslate">\; cos</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; /</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&lambda;</span>}{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; no</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 1,2</span>}<span\; class="notranslate">\; ,</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; N</span>}\end{array}$

Among them, x _m (t) represents the output data of the mth observation station. s(t) represents the signal data received by the observation station. n _m (t) represents the space-time stationary Gaussian white noise on the mth observation station, and its distribution obeys N(0,σ ² ). λ represents the operating wavelength of the radiation source. t represents the sampling time, and the value of t is t=1,2,...,N, where N is the total length of data. The phase information in this equation consists of two terms with different conditions for no phase ambiguity. For the first item that is only related to the azimuth angle of the radiation source, there is no phase ambiguity when the value of d _m does not exceed half a wavelength; for the second item, when the value of d _m satisfies There is no phase ambiguity. In the actual positioning problem, due to the sparse arrangement of observation stations, the distance d _m from the observation station to the coordinate origin and the distance between adjacent observation stations are usually much larger than half a wavelength, so the first term usually has phase ambiguity. Under near-field conditions, only when the array spacing d _m of the observation stations satisfies There is phase ambiguity in the second term. Under near-field conditions, combining the working wavelength of the passive positioning system and reasonably setting the distance between observation stations can make the second term phase satisfy the above-mentioned unambiguous inequality condition in the sense of the entire array manifold vector in a wide near-field range, namely There is always no ambiguity for the entire array manifold manifold vector.

For those skilled in the art, it is a common technique in this field to calculate the azimuth angle according to the received signal, and there are many disclosed methods, such as phase interferometer method and spatial spectrum method, etc., for the present invention, so the existing disclosure All methods for obtaining the azimuth angle are applicable to the azimuth angle in the present invention. The focus of the present invention is to calculate the distance of the near-field radiation source, and the azimuth angle is used as a priori for calculating the distance of the near-field radiation source. The invention uses the azimuth priori knowledge to compensate the first item with phase ambiguity, then uses the part without phase ambiguity to establish a ranging model, and uses the spatial spectrum estimation MUSIC method to obtain the initial estimate of the distance of the near-field radiation source. On this basis, the high-precision maximum likelihood method is applied to perform iterative optimization around the initial estimation result based on the complete uncompensated data model to obtain a high-precision estimation of the distance of the near-field radiation source, thereby achieving the purpose of the invention. Thereby the inventive method comprises:

Step 1. Matrix processing observation data: write the observation data x _m (t) (m=0,1,...,M-1) in complex form received by M observation stations into a matrix form x(t);

Step 2. Phase compensation: use the azimuth prior information to determine the phase compensation matrix Φ, multiply Φ to the vector x(t) to the left, and perform phase compensation on the first item with phase ambiguity in each element of x(t);

Step 3. MUSIC method ranging: using the phase-compensated array data in step 2, the autocovariance matrix of the array data obtained after compensation R performs matrix eigenvalue decomposition, and deletes the eigenvector corresponding to the largest eigenvalue from the eigenvector matrix to obtain the noise subspace matrix G. Using this noise subspace matrix and according to the MUSIC method criterion, search for spectral peaks in the distance dimension , to obtain the distance estimate of the near-field radiation source, go to step 4;

Step 4. Ranging with the maximum likelihood method: use the original observation data in step 1, based on its autocovariance matrix R, and use the prior azimuth information and the estimated distance of the near-field radiation source in step 3 as the initial value, According to the principle of maximum likelihood method, the quasi-Newton iterative technique is used to obtain the high-precision estimation of the two-dimensional parameter vectors of the azimuth and distance of the near-field radiation source. Since the azimuth angle is known a priori and does not need to be considered, the maximum likelihood distance estimate is taken as the final ranging result of the near-field radiation source, that is, the inventive task is completed.

In step 1, the observation data x(t) in complex form is expressed as:

x(t)=a(θ,R ₀ )s(t)+n(t),t=1,2,...,N

Among them, x(t) is the observation data matrix of the observation station, t is the sampling time, and the value is t=1,2,...,N, a(θ,R ₀ ) is the manifold vector of the radiation source array, expressed as :

$\mathrm{<span\; class="notranslate">\; a</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>{<span\; class="notranslate">\; [</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; e</span>}}^{<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; j</span>}\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; \&pi;</span>}{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}\mathrm{<span\; class="notranslate">\; sin</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; /</span>\mathrm{<span\; class="notranslate">\; \&lambda;</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; j\&pi;</span>}{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}\mathrm{<span\; class="notranslate">\; cos</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; /</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&lambda;</span>}{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; ,</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; e</span>}}^{{<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; j</span>}\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; \&pi;d</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}}\mathrm{<span\; class="notranslate">\; sin</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; /</span>\mathrm{<span\; class="notranslate">\; \&lambda;</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; j\&pi;</span>}{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}}\mathrm{<span\; class="notranslate">\; cos</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; /</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&lambda;</span>}{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span>}<span\; class="notranslate">\; ]</span>}^{\mathrm{<span\; class="notranslate">\; T</span>}}$

Where (·) ^T is the matrix transpose operation, is the signal vector of the near-field radiation source, n(t)～N(0,σ ² I _M ) is the Gaussian white noise vector, and I _M is the identity matrix of M×M.

In step 2, the phase compensation matrix Φ is a diagonal matrix whose expression is:

$\mathrm{<span\; class="notranslate">\; \&Phi;</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; diag</span>}<span\; class="notranslate">\; \{</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; e</span>}}^{\mathrm{<span\; class="notranslate">\; j</span>}\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; \&pi;</span>}{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}\mathrm{<span\; class="notranslate">\; sin</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; /</span>\mathrm{<span\; class="notranslate">\; \&lambda;</span>}}<span\; class="notranslate">\; ,</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; e</span>}}^{\mathrm{<span\; class="notranslate">\; j</span>}\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; \&pi;</span>}{\mathrm{<span\; class="notranslate">\; d</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}}\mathrm{<span\; class="notranslate">\; sin</span>}\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; /</span>\mathrm{<span\; class="notranslate">\; \&lambda;</span>}}<span\; class="notranslate">\; \}</span>$

In step 3, the autocovariance matrix of the resulting array data after compensation Satisfy the following model:

$\stackrel{<span\; class="notranslate">\; ~</span>}{\mathrm{<span\; class="notranslate">\; R</span>}}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; E.</span>}<span\; class="notranslate">\; \{</span>\mathrm{<span\; class="notranslate">\; \&Phi;x</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; x</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; \&Phi;</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}<span\; class="notranslate">\; \}</span><span\; class="notranslate">\; =</span>{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}\mathrm{<span\; class="notranslate">\; \&Phi;a</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; a</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; \&Phi;</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}{\mathrm{<span\; class="notranslate">\; I</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}$

Where (·) ^H is the matrix conjugate transpose operation. Will Eigendecomposition, its eigenvalues are sorted as λ ₁ >λ ₂ =...=λ _M , the eigenvectors corresponding to these eigenvalues are u ₁ ,...u _k ,u _k+1 ,...u _M , noise The subspace matrix is G=[u ₂ u ₃ . . . u _M ]. The criterion function used for MUSIC peak search and distance measurement is:

P _MUSIC (R ₀ )＝1/||G ^H Φa(θ,R ₀ )|| ²

where ||·|| ² is the square of the 2-norm of the vector. When calculating, By the sample autocovariance matrix replace.

In step 4, the autocovariance matrix R of the original multiple observation data satisfies the following model:

$\mathrm{<span\; class="notranslate">\; R</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; E.</span>}<span\; class="notranslate">\; \{</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; x</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; \}</span><span\; class="notranslate">\; =</span>{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}\mathrm{<span\; class="notranslate">\; a</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; a</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}{\mathrm{<span\; class="notranslate">\; I</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}$

Among them, E{·} is the expected operation. The criterion function of the maximum likelihood method is:

$\stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}<span\; class="notranslate">\; =</span>{<span\; class="notranslate">\; [</span>\stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; \&theta;</span>}}<span\; class="notranslate">\; ,</span>{\stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; R</span>}}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; ]</span>}^{\mathrm{<span\; class="notranslate">\; T</span>}}<span\; class="notranslate">\; =</span>\underset{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}<span\; class="notranslate">\; =</span>{<span\; class="notranslate">\; [</span>\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; ]</span>}^{\mathrm{<span\; class="notranslate">\; T</span>}}}{\mathrm{<span\; class="notranslate">\; min</span>}}\mathrm{<span\; class="notranslate">\; f</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&alpha;</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\underset{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}<span\; class="notranslate">\; =</span>{<span\; class="notranslate">\; [</span>\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; ]</span>}^{\mathrm{<span\; class="notranslate">\; T</span>}}}{\mathrm{<span\; class="notranslate">\; min</span>}}\mathrm{<span\; class="notranslate">\; log</span>}<span\; class="notranslate">\; |</span>{\stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}\mathrm{<span\; class="notranslate">\; a</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; a</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; \&theta;</span>}<span\; class="notranslate">\; ,</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; +</span>{\stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}{\mathrm{<span\; class="notranslate">\; I</span>}}_{\mathrm{<span\; class="notranslate">\; m</span>}}<span\; class="notranslate">\; |</span>$

where |·| is the symbol of the determinant, the power estimation of the radiation source Noise Power Estimation ${\widehat{\mathrm{\&sigma;}}}^2=\mathrm{Tr}\left\{{\mathrm{\&Pi;}}_a^{\&perp;}\hat{R}\right\}/(m-1),$ Tr{·} is to take the matrix trace. ${\mathrm{\&Pi;}}_a^{\&perp;}=I_m-a(\mathrm{\&theta;},R_0){\left[a^h(\mathrm{\&theta;},R_0)a(\mathrm{\&theta;},R_0)\right]}^{-1}{a}^h(\mathrm{\&theta;},R_0)$ is the vertical projection operator matrix to a(θ,R ₀ ) orthogonal complement subspace. is the pseudo-inverse of a(θ,R ₀ ), the sample autocovariance matrix vector elements in is the estimated distance to the near-field radiation source.

The iterative process of the quasi-Newton iterative technique is:

α ^k+1 ＝α ^k -μ _k H ^-1 (α ^k )f′(α ^k )

where α ^k is the estimated value of α obtained by the k-th iteration. μ _k is the iteration step size of the kth step. H is the Hessian matrix of the cost function f(α) or its positive definite approximate matrix, and f′(α ^k ) is the gradient of the cost function. The iteration termination criterion is that ||H ^-1 (α ^k )f′(α ^k )|| is less than a specified value ε and H is positive definite, or the number of iterations reaches the preset maximum number.

Example

This embodiment takes 7 observation stations sparsely arranged on a straight line, and the working wavelength of the near-field radiation source is 3 cm and the distance from the observation station is 30.5 kilometers as an example, that is, M=7, λ=3cm, R ₀ =30.5km. With the leftmost observation station on the straight line as the coordinate origin, and the straight line where the observation station is located as the abscissa axis, a right-handed rectangular coordinate system is established. The distances between the remaining six observation stations and the coordinate origin are d ₁ =10m, d ₂ =20m, d ₃ =35m, d ₄ =60m, d ₅ =83m, d ₆ =100m. The total length of observation data received in this example is N=1000, and the signal-to-noise ratio is 15dB. Investigate that the near-field radiation source is incident on the observatory at 10 azimuth angles θ=-70°, -50°, -30°, -10°, 0°, 20°, 30°, 40°, 60°, and 70° The distance measuring performance of the present invention is used in the array. For each azimuth, 200 Monte-Carlo experiments are performed to calculate the root mean square error of ranging and the root mean square error of relative ranging.

The embodiment process is as follows:

Step 1. In the 10 azimuths investigated, starting from the azimuth θ=-70°, using the selected current azimuth as a priori, constructing a phase compensation matrix Φ;

Step 2. Build the sample autocorrelation matrix of the array received data in complex form:

$\stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; R</span>}}<span\; class="notranslate">\; =</span>\underset{\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 1</span>}}{\overset{\mathrm{<span\; class="notranslate">\; 1000</span>}}{\mathrm{<span\; class="notranslate">\; \&Sigma;</span>}}}\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; x</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; /</span>\mathrm{<span\; class="notranslate">\; 1000</span>}$

Step 3. Calculate the sample autocovariance matrix after phase compensation:

$\stackrel{<span\; class="notranslate">\; ^</span>}{\stackrel{<span\; class="notranslate">\; ~</span>}{\mathrm{<span\; class="notranslate">\; R</span>}}}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; \&Phi;</span>}\stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; R</span>}}{\mathrm{<span\; class="notranslate">\; \&Phi;</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}$

Make a matrix The eigendecomposition of , a matrix composed of all other eigenvectors except the eigenvector corresponding to its largest eigenvalue, is used as an estimate of the noise subspace matrix, and the symbol is denoted as

Step 4. Combine the current azimuth value with the matrix Substituting into the ranging criterion function of the MUSIC method, search for spectral peaks in the distance dimension, and the search step is 1 km. Take the corresponding distance value under the spectral peak as the ranging estimated value of the near-field radiation source;

Step 5. Use the current azimuth value and the range estimation value obtained in step 4 as the initial value of the iteration, and substitute it into the criterion function of the maximum likelihood method to perform quasi-Newton iteration. In the quasi-Newton iterative process, μ _k =(1/2) ^k , ε=10 ^-6 , and the preset maximum number of iterations is 50 times. After the iterative process terminates, the estimated values of the azimuth and distance two-dimensional parameter vectors are obtained. Only care about distance parameter estimates Take it as the final estimate of the distance to the near-field radiation source;

Step 6. Repeat steps 2-5, do 200 Monte-Carlo simulation experiments, and calculate the ranging root mean square error (the smaller the value, the higher the ranging accuracy) and the relative ranging root mean square error (the closer the value is to 0, The higher the ranging accuracy), the calculation formula is:

in is the final estimated value of the near-field radiation source distance in the lth Monte-Carlo simulation experiment;

Step 7. Select the next azimuth as the current azimuth, and go to step 1 until all 10 azimuths have been inspected.

Fig. 2 is the statistical result of the root mean square error of ranging measured under different azimuth prior information conditions by adopting the specific example mode of the present invention, the ordinate is the root mean square value of the ranging error, and the unit is meter, and the abscissa is the azimuth θ ＝-70°,-50°,-30°,-10°,0°,20°,30°,40°,60°,70°, the unit is degree.

Fig. 3 is the statistical result of the relative ranging root mean square error measured under different azimuth prior information conditions for adopting the specific example mode of the present invention, and the ordinate is the relative ranging error root mean square value (proportional value, so there is no unit), The abscissa is the azimuth angle θ=-70°, -50°, -30°, -10°, 0°, 20°, 30°, 40°, 60°, 70°, and the unit is degree.

It can be seen from the figure that, for a near-field radiation source at a distance of 30.5 kilometers, the root mean square of the ranging errors under the prior conditions of each azimuth angle under consideration are all on the order of hundreds of meters, and the relative ranging errors are all on the order of zero percent. The order of magnitude of several, these two evaluation indexes have all illustrated the effectiveness of the present invention.

Claims (2)

1. the range finding passive location method based on bare cloth linear array under the priori conditions of position angle, is characterized in that, comprises the following steps:

Step 1. is arranged on M research station on interior among a small circle straight line, and selecting a research station is No. 0 research station, the distance d of other observation station and No. 0 research station _mmeet wherein, M is research station sum in system, M>=2, and m is research station numbering m=0 ..., M-1, λ is radiation source operation wavelength, θ is the radiation source orientation angle of arrival, R ₀represent the distance of No. 0 research station of radiation source to the;

M research station received observation data by step 2., calculates the estimated values theta of the radiation source orientation angle of arrival, generates observation data matrix x (t), x (t)=a (θ, R ₀) s (t)+n (t), t=1,2 ..., N; Wherein, t is sampling instant, t=1, and 2 ..., N, N is sampling total length, the signal vector that s (t) research station receives, n (t) is white Gaussian noise vector, a (θ, R ₀) be radiant array stream shape vector, $a(\mathrm{\&theta;},R_0)={[1,e^{-j2\mathrm{\&pi;}{d}_1\sin \mathrm{\&theta;}/\mathrm{\&lambda;}+\mathrm{j\&pi;}{\left(d_1\cos \mathrm{\&theta;}\right)}^2/\left(\mathrm{\&lambda;}{R}_0\right)},...,e^{{-j2\mathrm{\&pi;d}}_{M-1}\sin \mathrm{\&theta;}/\mathrm{\&lambda;}+\mathrm{j\&pi;}{\left(d_{M-1}\cos \mathrm{\&theta;}\right)}^2/\left(\mathrm{\&lambda;}{R}_0\right)}]}^T,$ () ^tfor matrix transpose computing;

The prior imformation θ that step 3. utilization is calculated the radiation source orientation angle of arrival determines Φ, diag represents diagonal matrix, and phase compensation matrix Φ premultiplication is completed to phase compensation to observation data matrix x (t);

Step 4. is asked the auto-covariance matrix of phase compensation the data obtained matrix by auto-covariance matrix carry out Eigenvalue Decomposition, after deleting the corresponding proper vector of eigenvalue of maximum in the corresponding eigenvectors matrix of eigenwert obtaining in decomposition, obtain noise subspace matrix G, utilize noise subspace matrix G to carry out MUSIC spectrum peak search, thereby obtain the distance R of No. 0 research station of the corresponding radiation source to the of criterion function peak value ₀for near-field thermal radiation spacing is from estimated value; MUSIC spectrum peak search criterion function P used _mUSIC(R ₀) be: P _mUSIC(R ₀)=1/||G ^hΦ a (θ, R ₀) || ², wherein, || || ²for 2-norm square.

2. the range finding passive location method based on bare cloth linear array under the priori conditions of position angle as claimed in claim 1, is characterized in that, also comprises step 5 after step 4:

Step 5: the priori value θ of the use radiation source orientation angle of arrival and near-field thermal radiation spacing as initial value, according to maximum likelihood method criterion, adopt the plan Newton iteration acquisition radiation source orientation angle of arrival and near-field thermal radiation spacing from bidimensional parameter vector from estimated value high precision estimated value;

The criterion function of maximum likelihood method is:

$\widehat{\mathrm{\&alpha;}}={[\widehat{\mathrm{\&theta;}},{\hat{R}}_0]}^T=\arg \underset{\mathrm{\&alpha;}={[\mathrm{\&theta;},R_0]}^T}{\min }f\left(\mathrm{\&alpha;}\right),f\left(\mathrm{\&alpha;}\right)=\log |{\widehat{\mathrm{\&sigma;}}}_s^{2}a(\mathrm{\&theta;},R_0)a^H(\mathrm{\&theta;},R_0)+{\widehat{\mathrm{\&sigma;}}}^2{I}_M|;$

for the high precision estimated value of the radiation source orientation angle of arrival, for near-field thermal radiation spacing from high precision estimated value, represent the value of corresponding α in the time that objective function f (α) gets minimum value; || be determinant symbol, for power of radiation source estimated value, for noise power estimation value, I _mfor the unit matrix of M × M; radiant array stream shape vector a (θ, R ₀) pseudoinverse, () ^hfor the computing of Matrix Conjugate transposition, for sample auto-covariance matrix, $\hat{R}=\underset{t=1}{\overset{N}{\mathrm{\&Sigma;}}}x\left(t\right)x^H\left(t\right)/N;$

tr{} is for getting trace of a matrix, for a (θ, R ₀) the vertical projection operator matrix of orthocomplement subspace, ${\mathrm{\&Pi;}}_a^{\&perp;}=I_M-a(\mathrm{\&theta;},R_0){\left[a^H(\mathrm{\&theta;},R_0)a(\mathrm{\&theta;},R_0)\right]}^{-1}{a}^H(\mathrm{\&theta;},R_0);$

Intending Newton iteration process is:

α ^k+1＝α ^k-μ _kH ^-1(α ^k)f′(α ^k)；

Wherein, α ^kthe radiant array stream shape vector α estimated value that k step iteration obtains, μ _kbe the iteration step length of k step, H is the black plug Hessian matrix of objective function f (α) or the approximate matrix of its positive definite, f ' (α ^k) be the gradient of objective function; Iteration stop criterion is || H ^-1(α ^k) f ' (α ^k) || be less than and specify numerical value ε and matrix H positive definite, or iterations reaches default maximum times;

Step 6: utilize radiation source orientation angle of arrival θ and near-field thermal radiation spacing from high precision estimated value complete radiation source location.