CN103926572B - A kind of clutter suppression method of battle array radar self adaptation subspace, airborne anon-normal side - Google Patents

Abstract

The invention belongs to Radar Technology field, disclose the clutter suppression method of battle array radar self adaptation subspace, a kind of airborne anon-normal side, it first passes through the clutter power spectrum of matching echo data estimation with theoretical clutter curve to estimate the unknown structure parameters of clutter Two dimensional Distribution curve, the process of matching adopt LTS method improve the robustness of parameter estimation, then utilize the structure parameters that estimates to calculate clutter subspace, then the orthogonal complement space that data are corresponding to clutter subspace is carried out projection carry out clutter reduction.Its simulation experiment result shows that the method for parameter estimation of the present invention has robustness and accuracy, utilizes and estimates that the weight vector that parameter calculates can obtain good clutter recognition effect.

Description

A kind of clutter suppression method of battle array radar self adaptation subspace, airborne anon-normal side

Technical field

The invention belongs to Radar Technology field, relate to the clutter suppression method that battle array radar self adaptation subspace, a kind of airborne anon-normal side is provided, specifically by estimating that the unknown structure parameters of clutter Two dimensional Distribution curve calculates the clutter subspace of distance unit to be detected, then the orthogonal complement space that data are corresponding to clutter subspace is carried out projection carry out clutter reduction.

Background technology

Airborne radar is mostly in down looks duty, inevitably to observe many land clutters, simultaneously because aircraft motion, clutter spectrum is broadening significantly, and traditional pulse Doppler technology is no longer valid.Space-time adaptive processes (spacetimeadaptiveprocessing, STAP) technology can effectively suppress land clutter, improves the detection performance of moving target.STAP technology generally chooses the adjacency unit of distance unit to be detected as training sample to calculate adaptive weight, when training sample meets independent same distribution (independentidenticallydistributed, during condition IID), along with the performance increasing self-adaptive processing of number of samples converges on optimum gradually.In reality, airborne radar is in order to realize all direction searching, generally set up multiple array antenna at fuselage simultaneously, certain azimuth scan is responsible for by each antenna, and now axial the and carrier aircraft velocity attitude of each aerial array will exist certain angle, occurs that array is that anon-normal is sidelong situation about putting.Anon-normal side battle array structure can cause that clutter angle doppler characterization is with distance change, makes training sample no longer meet independent same distribution condition, has had a strong impact on the estimated accuracy of covariance matrix, makes space-time adaptive process performance reduce.

For the non-stationary problem (namely clutter is apart from non-stationary) that clutter angle doppler characterization changes with distance, Chinese scholars has pointed out multiple solution.Zatman proposes derivative updating method (derivative-basedupdating, DBU), and the method supposes that weight vector is the linear function of distance, but when clutter apart from non-stationary very violent time, assumed above no longer set up.Melvin and Davis proposes the method (angledopplercompensation, ADW) of angle Doppler effect correction, and main-lobe clutter center is compensated by the method, and the compensation to sidelobe clutter is not then very good.Ries and Lapierre etc. propose the compensation method (registration-basedcompensation based on registration, RBC), the method is based on clutter distribution curve, the clutter power spectrum of different distance unit is carried out registration process to the non-stationary problem of the distance solving clutter, but during registration, the mapping relations of different distance unit curve are difficult to determine.Friedlander etc. propose the algorithm based on clutter subspace, the method directly utilizes the structure parameters of clutter Two dimensional Distribution curve and calculates the clutter subspace apart from unit to be detected, then the orthogonal complement space that echo data is corresponding to clutter subspace is carried out projection to curb clutter component.The method is without using training sample, it is to avoid sample is with the non-stationary problem of distance change.But the calculating of clutter subspace needs to utilize the structure parameters of clutter Two dimensional Distribution curve, and the structure parameters of clutter Two dimensional Distribution curve and radar system parameters, aircraft motion parameter etc. are relevant, and in reality, these parameters are not necessarily all known.

Summary of the invention

Clutter for airborne anon-normal side battle array radar needs to utilize this problem of structure parameters of clutter Two dimensional Distribution curve apart from non-stationary and subspace method, it is an object of the invention to provide the clutter suppression method of battle array radar self adaptation subspace, a kind of airborne anon-normal side, it is not necessary to utilize training sample and clutter recognition excellent performance.

In order to achieve the above object, the present invention is achieved by the following technical solutions.

The clutter suppression method of battle array radar self adaptation subspace, a kind of airborne anon-normal side, it is characterised in that comprise the following steps:

Step 1, sets up airborne anon-normal side battle array radar signal model, and wherein, radar operation wavelength is λ, launches M pulse in coherent processing inteval, and pulse recurrence frequency is f_r, aerial array is the uniform line-array of N number of array element composition, and array element distance is d, and carrier aircraft height is H, and carrier aircraft speed is v, and the axial angle of carrier aircraft velocity attitude and aerial array is ψ, the angle of pitch of ground clutter scattering object and azimuth respectively θ and, the oblique distance of carrier aircraft and ground clutter scattering object is R；

Its echo-signal can be expressed as

$X=x_c+n=\underset{l=1}{\overset{N_c}{\mathrm{\Σ}}}{\mathrm{\α}}_l{u}_l+n$ l=1,2…N_c

In formula, n is noise component(s), x_cFor clutter component, N_cFor clutter block number independent in ground clutter ring that oblique distance is R, α_lFor the amplitude of clutter the l clutter block of ring, and separate between different clutter block, u_lSteering vector during for clutter the l clutter block of ring empty；

The clutter Two dimensional Distribution curve of its clutter ring is

${\left(\frac{\mathrm{\λ}}{d}\right)}^2{w}_s^2+{\left(\frac{f_r\mathrm{\λ}}{2v}\right)}^2{w}_d^2-\left(\frac{{\mathrm{\λ}}^2{f}_r\cos \mathrm{\ψ}}{\mathrm{vd}}\right)w_s{w}_d={\cos }^2\mathrm{\θ}{\sin }^2\mathrm{\ψ}$

The simplified function expression formula of above formula

C=Γ(η(λ,d,f_r,v,ψ,θ))

In formula, C is clutter curve, and Γ is the operator being calculated clutter Two dimensional Distribution curve by parameter vector, and η is the parameter vector of clutter Two dimensional Distribution curve, wherein in the parameter vector η of clutter Two dimensional Distribution curve, and λ, d, f_rFor radar system known parameters, v, ψ, θ are unknown quantity；

Step 2, carries out angle Doppler's plane discretization to echo-signal, then does weighting two-dimension fourier transform spatially and temporally and obtain fourier spectra P_b, then to fourier spectra P_bTake absolute value square, obtain correspondence clutter power spectrum；

First thresholding, the grid cell of filtering clutter Two dimensional Distribution extra curvature are set, utilize the grid cell structure search submatrix on the clutter Two dimensional Distribution curve obtainedAccording to echo-signal complex magnitude γ on angle Doppler's plane correspondence grid cell, there is sparse characteristic, pass through Corresponding Sparse Algorithm $\begin{array}{c}\min {\left|\right|\mathrm{\γ}\left|\right|}_1,\\ s.t.{\left|\right|X-\widetilde{\mathrm{\Φ\γ}}\left|\right|}_2^2\≤\mathrm{\ϵ}\end{array}$ Solve the complex magnitude γ of grid cell on clutter Two dimensional Distribution curve, in formula | | | |_pRepresent L_pNorm, ε represents the limits of error that Corresponding Sparse Algorithm solves, for setting value；

According to the complex magnitude γ of grid cell on clutter Two dimensional Distribution curve, the second thresholding is set, relatively small magnitude on filtering clutter Two dimensional Distribution curve | γ_i| secondary grid cell, obtain the Power Spectrum Distribution typical grid unit in corresponding space-time two-dimensional plane distribution of reflection clutter；

Step 3, the clutter Two dimensional Distribution curve to clutter ring

${\left(\frac{\mathrm{\λ}}{d}\right)}^2{w}_s^2+{\left(\frac{f_r\mathrm{\λ}}{2v}\right)}^2{w}_d^2-\left(\frac{{\mathrm{\λ}}^2{f}_r\cos \mathrm{\ψ}}{\mathrm{vd}}\right)w_s{w}_d={\cos }^2\mathrm{\θ}{\sin }^2\mathrm{\ψ}$

Carrying out conversion to arrange, being write as vector form is y=a^Tβ, β=[1/v in formula²cosψ/vsin²ψcos²θ] T, a=[(f_rλw_d/2)²-λ₂f_rw_sw_d/d-1]^T,Wherein in β and η, unknown parameter v, ψ, θ are relevant, a, known parameters λ, d, f in y and η_rAnd the space-time two-dimensional frequency of clutter is relevant；T represents that transposition operates.

The space-time two-dimensional frequency of P2 the typical grid unit more than the second thresholding then step 2 obtained and configuration known parameter substitute into y=a^TIn a, the y of β, it is possible to obtain below equation group

$\hat{y}=\hat{A}\mathrm{\β}+\mathrm{\ξ}$

The measured value of P2 grid cell correspondence y in formula $\hat{y}={\left[\begin{array}{ccc}{\hat{y}}_1& ...& {\hat{y}}_P\end{array}\right]}^T$ For the response vector of P2 × 1 dimension, the measured value of a that P2 typical grid unit is corresponding $\hat{A}={\left[\begin{array}{ccc}{\hat{a}}_1& ...& {\hat{a}}_{N_p}\end{array}\right]}^T$ For the measurement matrix of P2 × 3 dimension, ξ is the error vector of P2 × 1 dimension, for setting value；

Step 4, adopts the equation group that minimum intercepting least square method solution procedure 3 obtainsObtain the estimated value of βAnd then solve obtain the structure parameters v of clutter Two dimensional Distribution curve, estimated value that ψ, θ are corresponding

Step 5, utilizes the structure parameters of the clutter Two dimensional Distribution curve estimatedWith known parameters λ, d, f_r, calculate steering vector matrix during the corresponding clutter sky of clutter Two dimensional Distribution curveOrthogonal Subspaces projection operatorWith filter weights vector W；

Step 6, constructs corresponding space-time filtering device according to filter weights vector W, it is suppressed that the clutter in echo-signal X, obtains the final doppler spectral after current clutter block clutter recognition: Z=W^HX。

Specifically, kth Doppler channel filtering weight vector W is utilized_kConstruct corresponding space-time filtering device, it is suppressed that the clutter in echo-signal X, obtain current clutter block, kth Doppler's passage is output as:K=1,2 ... K, K are Doppler's port number, and H represents conjugate transposition operation；Obtain the final doppler spectral after current clutter block clutter recognition: Z=[Z₁,Z₂,…,Z_K]^T, T represents that transposition operates.

The feature of technical scheme and further improvement is that:

(1) the concrete sub-step of step 2 is:

2a) by angle Doppler's plane discretization, grid cell number corresponding to Doppler frequency is N_d, grid cell number corresponding to spatial domain frequency is N_s；

Steering vector searching matrix Φ when setting up empty based on the grid cell in clutter block, concrete form $\mathrm{\Φ}=[u(w_{d,1},w_{s,1}),...,u(w_{d,i},w_{s,j}),...u(w_{d,N_d},w_{s,N_s})],$ Wherein u (w_d,i,w_s,j) in i=1 ... N_d;j=1,…N_sDenotation coordination be (i, the space-time two-dimensional frequency of grid cell j), namelyIn formulaAmass for Kronecker, wherein u_t(w_{D, i}) and u_s(w_s,j) respectively coordinate is (i, the time domain steering vector of grid cell j) and spatial domain steering vector in clutter block；

Then echo-signal X is expressed as X=Φ γ,For echo-signal complex magnitude on angle Doppler's plane correspondence grid cell, wherein T represents that transposition operates；

2b) angle Doppler's plane discretization echo-signal being done weighting two-dimension fourier transform spatially and temporally, its form isP in formula_bFor fourier spectra, t_wFor weight coefficient spatially and temporally,Amass for Hadamard, obtain fourier spectra P_b；

2c) to the fourier spectra P obtained after weighting spatially and temporally two-dimension fourier transform_bTake absolute value square, obtain the clutter power spectrum of correspondence, and it carried out the first Threshold detection, detect power P1 the grid cell more than the first threshold T H1, the space-time two-dimensional frequency of its correspondence is $\left[\begin{array}{ccc}u(w_{d,1}^1,w_{s,1}^1)& ...& u(w_{d,P1}^1,w_{s,P1}^1)\end{array}\right]$ Subscript 1 represents by the first thresholding TH1；And build a search submatrix $\widetilde{\mathrm{\Φ}}=\left[\begin{array}{ccc}u(w_{d,1}^1,w_{s,1}^1)& ...& u(w_{d,P1}^1,w_{s,P1}^1)\end{array}\right],$ It is a subset of searching matrix Φ；

2d) echo-signal complex magnitude γ on angle Doppler's plane correspondence grid cell has sparse characteristic, passes through Corresponding Sparse Algorithm $\begin{array}{c}\min {\left|\right|\mathrm{\γ}\left|\right|}_1,\\ s.t.{\left|\right|X-\widetilde{\mathrm{\Φ\γ}}\left|\right|}_2^2\≤\mathrm{\ϵ}\end{array}$ Solve the complex magnitude γ on estimation unit, in formula | | | |_pRepresent L_pNorm, ε represents the limits of error that Corresponding Sparse Algorithm solves, for setting value；

After 2e) estimating the amplitude of complex magnitude γ on grid cell, according to amplitude | γ_i| >=TH2, i=1 ... N_dN_s, in formula, TH2 is the second threshold value set, relatively small magnitude on filtering clutter Two dimensional Distribution curve | γ_i| secondary grid cell, obtain the Power Spectrum Distribution typical grid unit in corresponding space-time two-dimensional plane distribution of reflection clutter；

The space-time two-dimensional frequency that P2 typical grid unit 2f) detecting is corresponding is

$\left[\begin{array}{ccc}u(w_{d,1}^2,w_{s,1}^2)& ...& u(w_{d,P2}^2,w_{s,P2}^2)\end{array}\right],$ Subscript 2 represents by the second thresholding TH2；Typical grid unit, in the distribution of space-time two-dimensional plane, reflects the Power Spectrum Distribution of clutter.

(2) the concrete sub-step of step 4 is:

4a) adopt minimum intercepting least square method solving-optimizing problem

$\begin{array}{c}\min \underset{i=1}{\overset{P}{\mathrm{\Σ}}}{z}_i{({\hat{y}}_i-{\hat{a}}_i\mathrm{\β})}^2\\ s.t.z^{T}e=Q\\ z_i\∈\left\{\mathrm{0,1}\right\}\end{array}---\left(7\right)$

Z in formula_iRepresent flag bit, z_i=1 represents that these data are normal point, z_i=0 represents that these data are noise spot, and e represents complete 1 vector, and Q is the number of nonzero element in z, and Q < P；

4b) optimization problem above is obtained the estimated value of βAnd then calculating obtains estimated value corresponding to v, ψ, θ respectivelyIts computing formula is as follows

$\hat{v}=1/\sqrt{{\widehat{\mathrm{\&beta;}}}_{\mathrm{LTS}}\left(1\right)}$

$\widehat{\mathrm{\&psi;}}=\arccos ({\widehat{\mathrm{\&beta;}}}_{\mathrm{LTS}}\left(2\right)/\sqrt{{\widehat{\mathrm{\&beta;}}}_{\mathrm{LTS}}\left(1\right)})$

$\widehat{\mathrm{\&theta;}}=\arccos \sqrt{{\widehat{\mathrm{\&beta;}}}_{\mathrm{LTS}}\left(3\right)/(1-{({\widehat{\mathrm{\&beta;}}}_{\mathrm{LTS}}\left(2\right)/\sqrt{{\widehat{\mathrm{\&beta;}}}_{\mathrm{LTS}}\left(1\right)})}^2)}$

(3) the concrete sub-step of step 5 is:

5a) the estimation structure parameters that will obtainWith known parameters λ, d, f_r, substitute into normalization Doppler frequencyWith normalization spatial frequencyIn, obtain the time domain steering vector of the l estimation corresponding to clutter block $u_t\left({\hat{w}}_{d,l}\right)={\left[\begin{array}{cc}1& e^{j2\mathrm{\&pi;}{\hat{w}}_{d,l}}...e^{j2\mathrm{\&pi;}(M-1){\hat{w}}_{d,l}}\end{array}\right]}^T$ With spatial domain steering vector $u_s\left({\hat{w}}_{s,l}\right)={\left[\begin{array}{cc}1& e^{j2\mathrm{\&pi;}{\hat{w}}_{s,l}}...e^{j2\mathrm{\&pi;}(N-1){\hat{w}}_{s,l}}\end{array}\right]}^T,$ During the estimation that then the l clutter block now is corresponding empty, steering vector is ${\hat{V}}_{c,l}=u_t\left({\hat{w}}_{d,l}\right)\&CircleTimes;u_s\left({\hat{w}}_{s,l}\right),$ Steering vector matrix during the clutter sky estimated ${\hat{V}}_c=[{\hat{V}}_{c,1},...{\hat{V}}_{c,l,}...{\hat{V}}_{c,N_c}],$ l=1,2,…N_c, in formulaAmass for Kronecker；

5b) it is calculated as follows the projection operator being orthogonal to clutter subspace

$P_c^{\&perp;}=I-{\hat{V}}_c{({{\hat{V}}_c}^H{\hat{V}}_c+\mathrm{\&mu;I})}^{-1}{{\hat{V}}_c}^H$

In formula, μ is positive number, to ensure matrixReversible；⊥ represents orthogonal project operator, and Η represents conjugate transposition operation, ()^-1Representing matrix is inverted, I representation unit battle array；

Steering vector when 5c) utilizing target emptyI.e. s=[s₁,s₂,…,s_k,…s_K] and the projection operator being orthogonal to clutter subspace that estimates of step 5bCalculate the filter weights vector of the space-time filtering device obtaining kth Doppler's passageK=1,2 ... K, K are Doppler's port number.Filter weights vector

Wherein $u_{\mathrm{tar}}\left(w_d\right)={\left[\begin{array}{cc}1& e^{j2\mathrm{\&pi;}{w}_d}...e^{j2\mathrm{\&pi;}(M-1)w_d}\end{array}\right]}^T$ For the time domain steering vector of target, $u_{\mathrm{tar}}\left(w_s\right)={\left[\begin{array}{cc}1& e^{j2\mathrm{\&pi;}{w}_s}...e^{j2\mathrm{\&pi;}(N-1)w_s}\end{array}\right]}^T$ For the spatial domain steering vector of target, wherein normalization Doppler frequencyNormalization spatial frequencyWhat v represented is the speed of target；θ₀What represent is the angle of pitch of main beam,What represent is the azimuth of main beam, the angle that what ψ represented is carrier aircraft velocity attitude is axial with aerial array.

The clutter suppression method of the battle array radar self adaptation subspace, airborne anon-normal side of the present invention, first pass through clutter power spectrum that matching echo data estimates with theoretical clutter curve to estimate the unknown structure parameters of clutter Two dimensional Distribution curve, the process of matching adopt LTS method improve the robustness of parameter estimation, then utilize the structure parameters that estimates to calculate clutter subspace, then the orthogonal complement space that data are corresponding to clutter subspace is carried out projection carry out clutter reduction.Its simulation experiment result shows that the method for parameter estimation of the present invention has robustness and accuracy, utilizes and estimates that the weight vector that parameter calculates can obtain good clutter recognition effect.

In reality, due to the shortage of the precision problem of inertial navigation and landform altitude knowledge, by inaccurate for the structure parameters causing clutter Two dimensional Distribution curve, now it be accomplished by utilization and be received back to wave datum to estimate these unknown structure parameters.Owing to the power spectrum of clutter is determined by clutter Two dimensional Distribution curve at the track of angle Doppler's plane, the present invention, by matching clutter power spectrum and clutter Two dimensional Distribution curve, estimates corresponding structure parameters.

In order to reduce the diffusion of the power spectrum of clutter, the present invention utilizes openness in angle Doppler's plane of clutter, adopts Corresponding Sparse Algorithm to estimate the power spectrum of echo data.Radar return data actually do not only have noise signal, it is also possible to have echo signal (main lobe target or secondary lobe target).When echo signal noise is relatively larger, after being solved by Corresponding Sparse Algorithm, the grid cell corresponding to angle Doppler's plane of target also will pass through Threshold detection.These impact points are noise spot, will have a strong impact on the Parameter Estimation Precision of least square (LS) method；In order to solve this problem, the present invention adopts minimum intercepting two to take advantage of (leasttrimmedsquares, LTS) method to replace LS method to solve.

In the present invention, the power spectrum estimated by matching echo data and theoretical clutter curve estimate the unknown structure parameters of clutter Two dimensional Distribution curve, adopt LTS method, improve accuracy and the robustness of parameter estimation in the process of matching.

In the present invention, directly utilize structure parameters and calculate the clutter subspace of unit to be detected, then the orthogonal complement space that data are corresponding to clutter subspace is carried out projection to curb clutter component.The method is without using training sample, it is to avoid the non-stationary problem of sample.

Accompanying drawing explanation

Below in conjunction with the drawings and the specific embodiments, the present invention is described in further details.

The airborne anon-normal side battle array radar clutter that Fig. 1 is the present invention suppresses flow chart.

Fig. 2 is the airborne radar geometrized structure graph of the present invention.

Fig. 3 a is that the present invention true clutter curve at noiseless takes advantage of (LTS) fitting result schematic diagram with method of least square (LS) and minimum intercepting two；Abscissa represents normalization Doppler frequency, and vertical coordinate represents normalization spatial frequency.

Fig. 3 b is that the present invention true clutter curve when there being noise spot takes advantage of (LTS) fitting result schematic diagram with method of least square (LS) and minimum intercepting two；Abscissa represents normalization Doppler frequency, and vertical coordinate represents normalization spatial frequency.

Fig. 4 is the letter miscellaneous noise ratio loss curve synoptic diagram of optimal processing, sampling covariance, the inventive method (self adaptation subspace)；Abscissa represents normalization Doppler frequency, and vertical coordinate represents letter miscellaneous noise ratio loss (dB).

Detailed description of the invention

With reference to Fig. 1, the clutter suppression method of the battle array radar self adaptation subspace, airborne anon-normal side of the present invention being described, it specifically comprises the following steps that

Step 1, sets up airborne anon-normal side battle array radar signal model, the clutter Two dimensional Distribution curve of its clutter ring

${\left(\frac{\mathrm{\&lambda;}}{d}\right)}^2{w}_s^2+{\left(\frac{f_r\mathrm{\&lambda;}}{2v}\right)}^2{w}_d^2-\left(\frac{{\mathrm{\&lambda;}}^2{f}_r\cos \mathrm{\&psi;}}{\mathrm{vd}}\right)w_s{w}_d={\cos }^2\mathrm{\&theta;}{\sin }^2\mathrm{\&psi;}$

The simplified function expression formula of above formula

C=Γ(η(λ,d,f_r,v,ψ,θ))

In formula, C is clutter curve, and Γ is the operator being calculated clutter Two dimensional Distribution curve by parameter vector, and η is the parameter vector of clutter Two dimensional Distribution curve, wherein in the parameter vector η of clutter Two dimensional Distribution curve, and λ, d, f_rFor radar system known parameters, v, ψ, θ are unknown quantity.

Concrete example illustrates as follows: wherein, as shown in Figure 2, radar operation wavelength is λ to airborne anon-normal side battle array radar arrangement, launches M pulse in coherent processing inteval, and pulse recurrence frequency is f_r, aerial array is the uniform line-array of N number of array element composition, and array element distance is d, and carrier aircraft height is H, and carrier aircraft speed is v, and the axial angle of carrier aircraft velocity attitude and aerial array is ψ, the angle of pitch of ground clutter scattering object and azimuth respectively θ andThe oblique distance of carrier aircraft and ground clutter scattering object is R.According to above signal model, echo-signal can be expressed as

$X=x_c+n=\underset{l=1}{\overset{N_c}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_l{u}_l+n$ l=1,2…N_c

In formula, n is noise component(s), x_cFor clutter component, N_cFor clutter block number independent in ground clutter ring that oblique distance is R, α_lFor the amplitude of clutter the l clutter block of ring, and separate between different clutter block, u_lSteering vector during for clutter the l clutter block of ring empty.

Without loss of generality, during clutter block empty, steering vector can be expressed as

$u(w_d,w_s)=u_t\left(w_d\right)\&CircleTimes;u_s\left(w_s\right)$

In formulaAmass for Kronecker, u_t(w_d) and u_s(w_s) respectively the time domain steering vector of clutter block and spatial domain steering vector；

$u_t\left(w_d\right)={\left[\begin{array}{cc}1& e^{j2\mathrm{\&pi;}{w}_d}...e^{j2\mathrm{\&pi;}(M-1)w_d}\end{array}\right]}^T$

$u_s\left(w_s\right)={\left[\begin{array}{cc}1& e^{j2\mathrm{\&pi;}{w}_s}...e^{j2\mathrm{\&pi;}(M-1)w_s}\end{array}\right]}^T$

In formulaFor normalization Doppler frequency,For normalization spatial frequency, T represents that transposition operates.

Utilize w_dWith w_sBetween relation carry out conversion and arrange and can obtain clutter ring at the Two dimensional Distribution curve of angle Doppler domain and be

${\left(\frac{\mathrm{\&lambda;}}{d}\right)}^2{w}_s^2+{\left(\frac{f_r\mathrm{\&lambda;}}{2v}\right)}^2{w}_d^2-\left(\frac{{\mathrm{\&lambda;}}^2{f}_r\cos \mathrm{\&psi;}}{\mathrm{vd}}\right)w_s{w}_d={\cos }^2\mathrm{\&theta;}{\sin }^2\mathrm{\&psi;}$

Now the clutter Two dimensional Distribution curve of any one clutter ring can be expressed as

C=Γ(η(λ,d,f_r,v,ψ,θ))

In formula, C is clutter curve, and Γ is the operator being calculated clutter Two dimensional Distribution curve by parameter vector, and η is the parameter vector of clutter Two dimensional Distribution curve, wherein in the parameter vector η of clutter Two dimensional Distribution curve, and λ, d, f_rFor radar system known parameters, v, ψ, θ are unknown quantity.

Step 2, carries out angle Doppler's plane discretization to echo-signal, then does weighting two-dimension fourier transform spatially and temporally and obtain fourier spectra P_b, then to fourier spectra P_bTake absolute value square, obtain correspondence clutter power spectrum；

First thresholding, the grid cell of filtering clutter Two dimensional Distribution extra curvature are set, utilize the grid cell structure search submatrix on the clutter Two dimensional Distribution curve obtainedAccording to echo-signal complex magnitude γ on angle Doppler's plane correspondence grid cell, there is sparse characteristic, pass through Corresponding Sparse Algorithm $\begin{array}{c}\min {\left|\right|\mathrm{\&gamma;}\left|\right|}_1,\\ s.t.{\left|\right|X-\widetilde{\mathrm{\&Phi;\&gamma;}}\left|\right|}_2^2\&le;\mathrm{\&epsiv;}\end{array}$ Solve the complex magnitude γ of grid cell on clutter Two dimensional Distribution curve, in formula | | | |_pRepresent L_pNorm, ε represents the limits of error that Corresponding Sparse Algorithm solves, for setting value；

According to the complex magnitude γ of grid cell on clutter Two dimensional Distribution curve, the second thresholding is set, relatively small magnitude on filtering clutter Two dimensional Distribution curve | γ_i| secondary grid cell, obtain the Power Spectrum Distribution typical grid unit in corresponding space-time two-dimensional plane distribution of reflection clutter.

Concrete sub-step illustrates as follows:

2a) by angle Doppler's plane discretization, grid cell number corresponding to Doppler frequency is N_d, grid cell number corresponding to spatial domain frequency is N_s。

Steering vector searching matrix Φ when setting up empty based on the grid cell in clutter block, concrete form $\mathrm{\&Phi;}=[u(w_{d,1},w_{s,1}),...,u(w_{d,i},w_{s,j}),...u(w_{d,N_d},w_{s,N_s})],$ Wherein u (w_d,i,w_s,j) in i=1 ... N_d;j=1,…N_sDenotation coordination be (i, the space-time two-dimensional frequency of grid cell j), namelyIn formulaAmass for Kronecker, wherein u_t(w_{D, i}) and u_s(w_s,j) respectively coordinate is (i, the time domain steering vector of grid cell j) and spatial domain steering vector in clutter block.

Then echo-signal X can be expressed as X=Φ γ,For echo-signal complex magnitude on angle Doppler's plane correspondence grid cell, wherein T represents transposition.

2b) angle Doppler's plane discretization echo-signal being done weighting two-dimension fourier transform spatially and temporally, its form isP in formula_bFor fourier spectra, t_wFor weight coefficient spatially and temporally,Amass for Hadamard, obtain fourier spectra P_b。

2c) to the fourier spectra P obtained after spatially and temporally two-dimension fourier transform_bTake absolute value square, obtain the clutter power spectrum of correspondence, and it carried out the first Threshold detection, it is assumed that detect that the power of total P1 grid cell is more than the first threshold value that the first threshold T H1, TH1 is setting, is usually set to 5-10dB.Emulation experiment takes 8dB.

The space-time two-dimensional frequency that these grid cells are corresponding is $\left[\begin{array}{ccc}u(w_{d,1}^1,w_{s,1}^1)& ...& u(w_{d,P1}^1,w_{s,P1}^1)\end{array}\right],$ Subscript 1 represents by the first thresholding TH1.

A new search submatrix is built according to these grid cells $\widetilde{\mathrm{\&Phi;}}=\left[\begin{array}{ccc}u(w_{d,1}^1,w_{s,1}^1)& ...& u(w_{d,P1}^1,w_{s,P1}^1)\end{array}\right],$ It is a subset of searching matrix Φ.

2d) owing to clutter power spectrum is along clutter Two dimensional Distribution curve distribution, and clutter component occupies main component relative to noise component(s) in echo-signal.This mean that the power spectrum of estimation only on the minority grid cell near clutter Two dimensional Distribution curve amplitude relatively big, and amplitude is only small and close to zero on the grid cell of other position.

Namely echo-signal complex magnitude γ on angle Doppler's plane correspondence grid cell has sparse characteristic, therefore can pass through Corresponding Sparse Algorithm $\begin{array}{c}\min {\left|\right|\mathrm{\&gamma;}\left|\right|}_1,\\ s.t.{\left|\right|X-\mathrm{\&Phi;\&gamma;}\left|\right|}_2^2\&le;\mathrm{\&epsiv;}\end{array}$ Solve the complex magnitude γ on grid cell, in formula | | | |_pRepresent L_pNorm, ε represents the limits of error that Corresponding Sparse Algorithm solves, for setting value.

In order to reduce operand, adopt 2c) the search submatrix that obtainsReplacing the Φ in Corresponding Sparse Algorithm, the Corresponding Sparse Algorithm namely improved is $\begin{array}{c}\min {\left|\right|\mathrm{\&gamma;}\left|\right|}_1,\\ s.t.{\left|\right|X-\widetilde{\mathrm{\&Phi;\&gamma;}}\left|\right|}_2^2\&le;\mathrm{\&epsiv;}\end{array}.$

After 2e) estimating the amplitude of complex magnitude γ on grid cell, extract big amplitude components according to crossing threshold criterion | γ_i| >=TH2, i=1 ... N_dN_s, in formula, TH2 is the second threshold value set, and is generally 10-20dB, and emulation experiment takes 16dB.

So, we are with regard to relatively small magnitude on filtering clutter Two dimensional Distribution curve | γ_i| secondary grid cell, obtain the Power Spectrum Distribution typical grid unit in corresponding space-time two-dimensional plane distribution of reflection clutter.

2f) suppose the amplitude detecting total P2 typical grid unit | γ_i| more than the second threshold T H2, the space-time two-dimensional frequency that these typical grid unit are corresponding is $\left[\begin{array}{ccc}u(w_{d,1}^2,w_{s,1}^2)& ...& u(w_{d,P2}^2,w_{s,P2}^2)\end{array}\right],$ Subscript 2 represents by the second thresholding TH2.These typical grid unit, in the distribution of space-time two-dimensional plane, reflect the Power Spectrum Distribution of clutter.

Step 3, the clutter Two dimensional Distribution curve to clutter ring

${\left(\frac{\mathrm{\&lambda;}}{d}\right)}^2{w}_s^2+{\left(\frac{f_r\mathrm{\&lambda;}}{2v}\right)}^2{w}_d^2-\left(\frac{{\mathrm{\&lambda;}}^2{f}_r\cos \mathrm{\&psi;}}{\mathrm{vd}}\right)w_s{w}_d={\cos }^2\mathrm{\&theta;}{\sin }^2\mathrm{\&psi;}$

Carrying out conversion to arrange, being write as vector form is y=a^Tβ, β=[1/v in formula²cosψ/vsin²ψcos²θ]^T, a=[(f_rλw_d/2)²-λ²f_rw_sw_d/d-1]^T,Wherein in β and η, unknown parameter v, ψ, θ are relevant, a, known parameters λ, d, f in y and η_rAnd the space-time two-dimensional frequency of clutter is relevant；

The space-time two-dimensional frequency of P2 the typical grid unit more than the second thresholding then step 2 obtained and configuration known parameter substitute into y=a^TIn a, the y of β, it is possible to obtain below equation group

$\hat{y}=\hat{A}\mathrm{\&beta;}+\mathrm{\&xi;}$

The measured value of P2 grid cell correspondence y in formula $\hat{y}={\left[\begin{array}{ccc}{\hat{y}}_1& ...& {\hat{y}}_P\end{array}\right]}^T$ For the response vector of P2 × 1 dimension, the measured value of a that P2 typical grid unit is corresponding $\hat{A}={\left[\begin{array}{ccc}{\hat{a}}_1& ...& {\hat{a}}_{N_p}\end{array}\right]}^T$ For the measurement matrix of P2 × 3 dimension, ξ is the error vector of P2 × 1 dimension, for setting value.

Step 4, adopts the equation group that minimum intercepting least square method solution procedure 3 obtainsObtain the estimated value of βAnd then solve obtain the structure parameters v of clutter Two dimensional Distribution curve, estimated value that ψ, θ are correspondingValue.

Concrete sub-step illustrates as follows:

4a) LTS is a kind of sane homing method, it is possible to be expressed as following optimization problem

$\begin{array}{c}\min \underset{i=1}{\overset{P}{\mathrm{\&Sigma;}}}{z}_i{({\hat{y}}_i-{\hat{a}}_i\mathrm{\&beta;})}^2\\ s.t.z^{T}e=Q\\ z_i\&Element;\left\{\mathrm{0,1}\right\}\end{array}$

Z in formula_iRepresent flag bit, z_i=1 represents that these data are normal point, z_i=0 represents that these data are noise spot, and e represents complete 1 vector, and Q represents the number of nonzero element in z and Q < P.

4b) optimization problem above is obtained the estimated value of βAnd then can calculate respectively and obtain estimated value corresponding to v, ψ, θIts computing formula is as follows

$\hat{v}=1/\sqrt{{\widehat{\mathrm{\&beta;}}}_{\mathrm{LTS}}\left(1\right)}$

$\widehat{\mathrm{\&psi;}}=\arccos ({\widehat{\mathrm{\&beta;}}}_{\mathrm{LTS}}\left(2\right)/\sqrt{{\widehat{\mathrm{\&beta;}}}_{\mathrm{LTS}}\left(1\right)})$

$\widehat{\mathrm{\&theta;}}=\arccos \sqrt{{\widehat{\mathrm{\&beta;}}}_{\mathrm{LTS}}\left(3\right)/(1-{({\widehat{\mathrm{\&beta;}}}_{\mathrm{LTS}}\left(2\right)/\sqrt{{\widehat{\mathrm{\&beta;}}}_{\mathrm{LTS}}\left(1\right)})}^2)}$

Step 5, utilizes the structure parameters of the clutter Two dimensional Distribution curve estimatedWith known parameters λ, d, f_r, calculate steering vector matrix during the corresponding clutter sky of clutter Two dimensional Distribution curveOrthogonal Subspaces projection operatorWith filter weights vector W.

Concrete sub-step illustrates as follows:

5a) the estimation structure parameters that will obtainWith known parameters λ, d, f_r, substitute into normalization Doppler frequencyWith normalization spatial frequencyIn, obtain the time domain steering vector of the l estimation corresponding to clutter block $u_t\left({\hat{w}}_{d,l}\right)={\left[\begin{array}{cc}1& e^{j2\mathrm{\&pi;}{\hat{w}}_{d,l}}...e^{j2\mathrm{\&pi;}(M-1){\hat{w}}_{d,l}}\end{array}\right]}^T$ With spatial domain steering vector $u_s\left({\hat{w}}_{s,l}\right)={\left[\begin{array}{cc}1& e^{j2\mathrm{\&pi;}{\hat{w}}_{s,l}}...e^{j2\mathrm{\&pi;}(N-1){\hat{w}}_{s,l}}\end{array}\right]}^T,$ During the estimation that then the l clutter block now is corresponding empty, steering vector is ${\hat{V}}_{c,l}=u_t\left({\hat{w}}_{d,l}\right)\&CircleTimes;u_s\left({\hat{w}}_{s,l}\right),$ Steering vector matrix during the clutter sky estimated ${\hat{V}}_c=[{\hat{V}}_{c,1},...{\hat{V}}_{c,l,}...{\hat{V}}_{c,N_c}],$ l=1,2,…N_c, in formulaAmass for Kronecker；

5b) it is calculated as follows the projection operator being orthogonal to clutter subspace

$P_c^{\&perp;}=I-{\hat{V}}_c{({{\hat{V}}_c}^H{\hat{V}}_c+\mathrm{\&mu;I})}^{-1}{{\hat{V}}_c}^H$

In formula, μ is positive number, to ensure matrixReversible；⊥ represents orthogonal project operator, and Η represents conjugate transposition operation, ()^-1Representing matrix is inverted, I representation unit battle array；

Steering vector when 5c) utilizing target emptyI.e. s=[s₁,s₂,…,s_k,…s_K] and the projection operator being orthogonal to clutter subspace that estimates of step 5bCalculate the filter weights vector of the space-time filtering device obtaining kth Doppler's passageK=1,2 ... K, K are Doppler's port number.Filter weights vector

Wherein $u_{\mathrm{tar}}\left(w_d\right)={\left[\begin{array}{cc}1& e^{j2\mathrm{\&pi;}{w}_d}...e^{j2\mathrm{\&pi;}(M-1)w_d}\end{array}\right]}^T$ For the time domain steering vector of target, $u_{\mathrm{tar}}\left(w_s\right)={\left[\begin{array}{cc}1& e^{j2\mathrm{\&pi;}{w}_s}...e^{j2\mathrm{\&pi;}(N-1)w_s}\end{array}\right]}^T$ For the spatial domain steering vector of target, wherein normalization Doppler frequencyNormalization spatial frequencyWhat v represented is the speed of target；θ₀What represent is the angle of pitch of main beam,What represent is the azimuth of main beam, the angle that what ψ represented is carrier aircraft velocity attitude is axial with aerial array.

Step 6, constructs corresponding space-time filtering device according to filter weights vector W, it is suppressed that the clutter in echo-signal X, obtains the final doppler spectral after current clutter block clutter recognition: Z=W^HX。

Specifically, kth Doppler channel filtering weight vector W is utilized_kConstruct corresponding space-time filtering device, it is suppressed that the clutter in echo-signal X, obtain current clutter block, kth Doppler's passage is output as:K=1,2 ... K, K are Doppler's port number, and H represents conjugate transposition operation；Obtain the final doppler spectral after current clutter block clutter recognition: Z=[Z₁,Z₂,…,Z_K]^T, T represents that transposition operates.

Below by emulation experiment, the performance of the clutter suppression method of the present invention is described in detail.

(1) emulation experiment 1 verifies the robustness of parameter estimation of the present invention

(1.1) simulation parameter

Simulation parameter is as shown in table 1, and unit oblique distance to be detected is 7km, and miscellaneous noise ratio is 50dB, add two noise spots, representing main lobe target and secondary lobe target respectively, corresponding space-time two-dimensional frequency is (-0.40.0), (0.2-0.25), and signal to noise ratio is 20dB.

Table 1 onboard radar system parameter

(1.2) emulation data processed result and analysis

The clutter curve-fitting results of least square (leastsquare, LS) and LTS both approaches when accompanying drawing 3a and Fig. 3 b sets forth noiseless point and has noise spot.' zero ' corresponding clutter point mark in figure, ' * ' corresponding Targets Dots.The clutter curve that when be can be seen that at noiseless by Fig. 3 a, two kinds of methods are estimated and true clutter curve essentially coincide, it was shown that two kinds of methods all can obtain higher Parameter Estimation Precision at noiseless.Being can be seen that the clutter curve that LS estimates when there is noise spot deviates considerably from true clutter curve by Fig. 3 b, this is because LS method solves coefficient vector by minimizing the remainder error square of all data, not considering the impact of noise spot.The clutter curve that LTS estimates still essentially coincides with true clutter curve, this is because LTS supposes there is noise spot in data, is calculated by optimum option normal data, it is possible to obtain good parameter estimation effect.

(2) emulation experiment 2 verifies the robustness of parameter estimation of the present invention

(2.1) simulation parameter

Simulation parameter, with emulation experiment 1, repeats 100 Monte Carlo experiments, and the root-mean-square error (rootmeansquareerror, RMSE) obtaining the estimation of LTS method parameter is as shown in table 2.

Table 2 parameter estimation RMSE

(2.2) emulation data processed result and analysis

The root-mean-square error that LTS method is estimated as can be seen from Table 2 is less, estimates that parameter is close to true value.This is owing to the clutter spectrum position of sparse Power estimation is comparatively accurate, Fig. 3 a and Fig. 3 b can be seen that clutter point mark is generally near true clutter spectral line, and LTS method is weighted by inhibiting the impact of noise spot.

(3) emulation experiment 3 verifies the clutter recognition performance of the present invention

(3.1) experimental technique

To believe the performance that miscellaneous noise ratio loss is criterion verification algorithm, letter miscellaneous noise ratio loss is defined as

${\mathrm{SINR}}_L=\frac{\mathrm{SINR}}{\mathrm{SNR}}=\left(\frac{w^H{R}_{s}w}{w^H\mathrm{Rw}}\right)/\left(\frac{{\mathrm{\&sigma;}}_s^2}{{\mathrm{\&sigma;}}_n^2}\mathrm{NM}\right)---\left(18\right)$

In formulaRepresent target power and the noise power of single array element pulse respectively,For the correlation matrix of target, R is theoretical clutter plus noise covariance matrix.

(3.2) emulation data processed result and analysis

Accompanying drawing 4 gives the letter miscellaneous noise ratio loss curve of optimal processing, sampling covariance, context of methods (self adaptation subspace).Sampling covariance is decreased obviously relative to optimal processing performance as seen from Figure 4, this is because sampling covariance utilizes the covariance matrix R that distance sample on average obtains_avgBiased estimation for R.The self adaptation subspace method performance of the present invention is substantially better than sampling covariance method, this is because subspace method directly utilizes structure parameters calculates clutter subspace projection operator, it is to avoid sample non-stationary problem.

Claims (4)

1. the clutter suppression method of battle array radar self adaptation subspace, an airborne anon-normal side, it is characterised in that comprise the following steps:

Step 1, sets up airborne anon-normal side battle array radar signal model, and wherein, radar operation wavelength is λ, launches M pulse in coherent processing inteval, and pulse recurrence frequency is f_r, aerial array is the uniform line-array of N number of array element composition, and array element distance is d, and carrier aircraft height is H, and carrier aircraft speed is v, and the axial angle of carrier aircraft velocity attitude and aerial array is Ψ, the angle of pitch of ground clutter scattering object and azimuth respectively θ andThe oblique distance of carrier aircraft and ground clutter scattering object is R；

Its echo-signal can be expressed as

$X=x_c+n=\underset{l=1}{\overset{N_c}{\&Sigma;}}{\mathrm{\&alpha;}}_l{u}_l+n,l=1,2...N_c$

In formula, n is noise component(s), x_cFor clutter component, N_cFor clutter block number independent in ground clutter ring that oblique distance is R, α_lFor the amplitude of clutter the l clutter block of ring, and separate between different clutter block, u_lSteering vector during for clutter the l clutter block of ring empty；

The clutter Two dimensional Distribution curve of its clutter ring is

${\left(\frac{\mathrm{\&lambda;}}{d}\right)}^2{w}_s^2+{\left(\frac{f_r\mathrm{\&lambda;}}{2v}\right)}^2{w}_d^2-\left(\frac{{\mathrm{\&lambda;}}^2{f}_{r}cos\mathrm{\&psi;}}{vd}\right)w_s{w}_d={\cos }^2{\mathrm{\&theta;sin}}^2\mathrm{\&psi;}$

The simplified function expression formula of above formula

C=Г (η (λ, d, f_r, v, Ψ, θ))

In formula, C is clutter curve, and Г is the operator being calculated clutter Two dimensional Distribution curve by parameter vector, and η is the parameter vector of clutter Two dimensional Distribution curve, wherein in the parameter vector η of clutter Two dimensional Distribution curve, and λ, d, f_rFor radar system known parameters, v, Ψ, θ are unknown quantity, w_dFor normalization Doppler frequency, w_sFor normalization spatial frequency；

Step 2, carries out angle Doppler's plane discretization to echo-signal, then does weighting two-dimension fourier transform spatially and temporally and obtain fourier spectra P_b, then to fourier spectra P_bTake absolute value square, obtain correspondence clutter power spectrum；

First thresholding, the grid cell of filtering clutter Two dimensional Distribution extra curvature are set, utilize the grid cell structure search submatrix on the clutter Two dimensional Distribution curve obtainedAccording to echo-signal complex magnitude γ on angle Doppler's plane correspondence grid cell, there is sparse characteristic, pass through Corresponding Sparse AlgorithmSolve the complex magnitude γ of grid cell on clutter Two dimensional Distribution curve, in formula | | | |_pRepresent L_pNorm, ε represents the limits of error that Corresponding Sparse Algorithm solves, for setting value；

According to the complex magnitude γ of grid cell on clutter Two dimensional Distribution curve, the second thresholding is set, relatively small magnitude on filtering clutter Two dimensional Distribution curve | γ_i| secondary grid cell, obtain the Power Spectrum Distribution typical grid unit in corresponding space-time two-dimensional plane distribution of reflection clutter；

Step 3, the clutter Two dimensional Distribution curve to clutter ring

${\left(\frac{\mathrm{\&lambda;}}{d}\right)}^2{w}_s^2+{\left(\frac{f_r\mathrm{\&lambda;}}{2v}\right)}^2{w}_d^2-\left(\frac{{\mathrm{\&lambda;}}^2{f}_{r}cos\mathrm{\&psi;}}{vd}\right)w_s{w}_d={\cos }^2{\mathrm{\&theta;sin}}^2\mathrm{\&psi;}$

Carrying out conversion to arrange, being write as vector form is y=a^Tβ, β=[1/v in formula²cosΨ/vsin²Ψcos²θ]^T, a=[(f_rλw_d/2)²-λ²f_rw_sw_d/d-1]^T,Wherein in β and η, unknown parameter v, Ψ, θ are relevant, a, known parameters λ, d, f in y and η_rAnd the space-time two-dimensional frequency of clutter is relevant；T represents that transposition operates；

The space-time two-dimensional frequency of P2 the typical grid unit more than the second thresholding then step 2 obtained and configuration known parameter substitute into y=a^TIn a, the y of β, it is possible to obtain below equation group

$\hat{y}=\hat{A}\mathrm{\&beta;}+\mathrm{\&xi;}$

The measured value of P2 grid cell correspondence y in formulaFor the response vector of P2 × 1 dimension, the measured value of a that P2 typical grid unit is correspondingFor the measurement matrix of P2 × 3 dimension, ξ is the error vector of P2 × 1 dimension, for setting value；

Step 4, adopts the equation group that minimum intercepting least square method solution procedure 3 obtainsObtain the estimated value of βAnd then solve obtain the structure parameters v of clutter Two dimensional Distribution curve, estimated value that Ψ, θ are corresponding

Step 5, utilizes the structure parameters of the clutter Two dimensional Distribution curve estimatedWith known parameters λ, d, f_r, calculate steering vector matrix during the corresponding clutter sky of clutter Two dimensional Distribution curveOrthogonal Subspaces projection operatorWith filter weights vector W；

Step 6, constructs corresponding space-time filtering device according to filter weights vector W, it is suppressed that the clutter in echo-signal X, obtains the final doppler spectral after current clutter block clutter recognition: Z=W^HX；

Specifically, kth Doppler channel filtering weight vector W is utilized_kConstruct corresponding space-time filtering device, it is suppressed that the clutter in echo-signal X, obtain current clutter block, kth Doppler's passage is output as:K is Doppler's port number, and H represents conjugate transposition operation；Obtain the final doppler spectral after current clutter block clutter recognition: Z=[Z₁, Z₂..., Z_K]^T, T represents that transposition operates.

2. the clutter suppression method of battle array radar self adaptation subspace, airborne anon-normal side according to claim 1, it is characterised in that the concrete sub-step of step 2 is:

2a) by angle Doppler's plane discretization, grid cell number corresponding to Doppler frequency is N_d, grid cell number corresponding to spatial domain frequency is N_s；

Steering vector searching matrix Ф when setting up empty based on the grid cell in clutter block, concrete formWherein u (w_{D, i}, w_{S, j}) in i=1 ... N_d；J=1 ... N_sDenotation coordination be (i, the space-time two-dimensional frequency of grid cell j), namelyIn formulaAmass for Kronecker, wherein u_t(w_{D, i}) and u_s(w_{S, j}) respectively coordinate is (i, the time domain steering vector of grid cell j) and spatial domain steering vector in clutter block；

Then echo-signal X is expressed as X=Φ γ,For echo-signal complex magnitude on angle Doppler's plane correspondence grid cell, wherein T represents that transposition operates；

2b) angle Doppler's plane discretization echo-signal being done weighting two-dimension fourier transform spatially and temporally, its form is P_b=| Φ^H(x⊙t_w) |, P in formula_bFor fourier spectra, t_wFor weight coefficient spatially and temporally, ⊙ is that Hadamard amasss, and obtains fourier spectra P_b, x is angle Doppler's plane discretization echo-signal；

2c) to the fourier spectra P obtained after weighting spatially and temporally two-dimension fourier transform_bTake absolute value square, obtain the clutter power spectrum of correspondence, and it carried out the first Threshold detection, detect power P1 the grid cell more than the first threshold T H1, the space-time two-dimensional frequency of its correspondence isSubscript 1 represents by the first thresholding TH1；And build a search submatrix It is a subset of searching matrix Ф；

2d) echo-signal complex magnitude γ on angle Doppler's plane correspondence grid cell has sparse characteristic, passes through Corresponding Sparse AlgorithmSolve the complex magnitude γ on estimation unit, in formula | | | |_pRepresent L_pNorm, ε represents the limits of error that Corresponding Sparse Algorithm solves, for setting value；

After 2e) estimating the amplitude of complex magnitude γ on grid cell, according to amplitude | γ_i| >=TH2, i=1 ... N_dN_s, in formula, TH2 is the second threshold value set, relatively small magnitude on filtering clutter Two dimensional Distribution curve | γ_i| secondary grid cell, obtain the Power Spectrum Distribution typical grid unit in corresponding space-time two-dimensional plane distribution of reflection clutter；

The space-time two-dimensional frequency that P2 typical grid unit 2f) detecting is corresponding is

Subscript 2 represents by the second thresholding TH2；Typical grid unit, in the distribution of space-time two-dimensional plane, reflects the Power Spectrum Distribution of clutter.

3. the clutter suppression method of battle array radar self adaptation subspace, airborne anon-normal side according to claim 1, it is characterised in that the concrete sub-step of step 4 is:

4a) adopt minimum intercepting least square method solving-optimizing problem

$\min \underset{i=1}{\overset{P}{\&Sigma;}}{z}_i{({\hat{y}}_i-{\hat{a}}_i\mathrm{\&beta;})}^2$

s.t.z^TE=Q

z_i∈ { 0,1}

Z in formula_iRepresent flag bit, z_i=1 represents that these data are normal point, z_i=0 represents that these data are noise spot, and e represents complete 1 vector, and Q represents the number of nonzero element in z and Q < P；

4b) optimization problem above is obtained the estimated value of β And then calculating obtains estimated value corresponding to v, Ψ, θ respectively Its computing formula is as follows

$\hat{v}=1/\sqrt{{\widehat{\mathrm{\&beta;}}}_{LTS}\left(1\right)}$

$\widehat{\mathrm{\&psi;}}=\arccos \left({\widehat{\mathrm{\&beta;}}}_{LTS}\right(2)/\sqrt{{\widehat{\mathrm{\&beta;}}}_{LTS}\left(1\right)})$

$\widehat{\mathrm{\&theta;}}=\arccos \sqrt{{\widehat{\mathrm{\&beta;}}}_{LTS}\left(3\right)/(1-{\left({\widehat{\mathrm{\&beta;}}}_{LTS}\left(2\right)/\sqrt{{\widehat{\mathrm{\&beta;}}}_{LTS}\left(1\right)}\right)}^2)}.$

4. the clutter suppression method of battle array radar self adaptation subspace, airborne anon-normal side according to claim 1, it is characterised in that the concrete sub-step of step 5 is:

5a) the estimation structure parameters that will obtainWith known parameters λ, d, f_r, substitute into normalization Doppler frequencyWith normalization spatial frequencyIn, obtain the time domain steering vector of the l estimation corresponding to clutter blockWith spatial domain steering vectorDuring the estimation that then the l clutter block now is corresponding empty, steering vector isSteering vector matrix during the clutter sky estimatedIn formulaAmass for Kronecker；

5b) it is calculated as follows the projection operator being orthogonal to clutter subspace

$P_c^{\&perp;}=I-{\hat{V}}_c{({\hat{V}}_c^H{\hat{V}}_c+\mathrm{\&mu;}I)}^{-1}{\hat{V}}_c^H$

In formula, μ is positive number, to ensure matrixReversible；⊥ represents orthogonal project operator, and H represents conjugate transposition operation, ()^-1Representing matrix is inverted, I representation unit battle array；

Steering vector when 5c) utilizing target emptyI.e. s=[s₁, s₂..., s_k... s_K] and the projection operator being orthogonal to clutter subspace that estimates of step 5bCalculate the filter weights vector of the space-time filtering device obtaining kth Doppler's passageK is Doppler's port number, filter weights vector

WhereinFor the time domain steering vector of target,For the spatial domain steering vector of target, wherein normalization Doppler frequencyNormalization spatial frequencyWhat v represented is the speed of target；θ₀What represent is the angle of pitch of main beam,What represent is the azimuth of main beam, the angle that what Ψ represented is carrier aircraft velocity attitude is axial with aerial array.