A moving target detection method for GEO satellite-borne bistatic SAR

Abstract

The application discloses a moving target detection method for GEO satellite-borne bistatic SAR, and is applied to the technical field of radars, and aims at the problem that the prior art is not applicable to detecting moving targets in GEO satellite-borne bistatic SAR; the application firstly establishes a multi-channel echo model of GEO satellite-borne bistatic SAR; then preprocessing is carried out through channel delay and phase compensation; secondly, moving target detection is modeled as a rank-1 decomposition problem of a matrix; finally, the established optimization problem is iteratively solved through an alternating direction method of multipliers (ADMM), and two sub-problems are alternately solved through a singular value shrinkage operator of a matrix and a negative gradient descent method. The method can effectively solve the problem of performance decline of GEO satellite-borne bistatic SAR in detecting non-uniform sparse moving targets, and finally realizes detection of moving target echoes from clutter.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of radar, and particularly relates to a GEO satellite-aircraft bistatic SAR moving target detection technology. BACKGROUND

[0002] The GEO satellite-aircraft bistatic synthetic aperture radar (SAR) transmits signals by a geosynchronous orbit (GEO) satellite and receives signals by an airborne receiving station, and has the characteristics of high concealment and flexible configuration. However, when observing the earth, the echo of a moving target is usually submerged in the background of ground static clutter, and therefore how to detect the moving target is a challenging task. Moving target detection methods are mainly divided into two-channel algorithms and multi-channel algorithms. The two-channel moving target detection algorithm is mainly displaced phase center antenna (DPCA) and along track interference (ATI). The multi-channel moving target detection algorithm is mainly space time adaptive processing (STAP), and the STAP method has been widely researched and applied due to the higher system degree of freedom of the multi-channel SAR. However, due to the problem of distance non-stationarity, the method cannot be directly applied to the GEO satellite-aircraft bistatic SAR.

[0003] In recent years, researchers at Xi'an University of Electronic Science and Technology applied the robust principal component analysis (RPCA) method to SAR moving target detection in their published paper "Dong Yang, Guisheng Liao, Shengqi Zhu and Xi Yang, "RPCA based moving target detection in strong clutter background," 2015 IEEE Radar Conference (RadarCon), 2015, pp. 1487-1490." This method loads the two-dimensional echo matrix after range direction pulse compression of each channel as column vectors, and then splices the column vectors formed by each channel into a matrix. Using the low-rank characteristics of clutter and the sparse characteristics of moving targets, this separation method using RPCA completes the detection of moving targets. The disadvantage of this method is that it does not consider the influence of channel differences on the azimuth envelope in principle, in addition to which the method needs to adjust the regularization parameter to balance the influence of the original constraint condition and the sparse condition in the objective function, and the adjustment process is relatively cumbersome. More seriously, the moving target echo of the GEO satellite-aircraft bistatic SAR system also has a certain correlation between channels, which means that when the moving target also exhibits a low-rank structure, the low-rank and sparse decomposition using RPCA will also decompose the moving target into the low-rank matrix. In addition, the relative position between the target and the platform changes within the synthetic aperture time, and the energy trajectory is scattered as a migration curve, which does not meet the assumption of scattered sparse isolated strong points. Therefore, the application of RPCA to GEO satellite-aircraft bistatic SAR will cause performance degradation.

[0004] The researchers at Southeast University used a method based on RPCA for SAR moving target detection in the image domain in their published paper "G. Xu, X. Wang, Y. Huang, L. Cai and Z. Jiang, "Joint Multi-Channel Sparse Method of Robust PCA for SAR Ground Moving Target Image Indication," 2019 IEEE International Geoscience and Remote Sensing Symposium pp. 1709-1712." The method changes the position of part of the objective function and the constraint condition in the RPCA method, and then uses the Godec algorithm to solve the new constraint optimization problem. The traditional RPCA uses the nuclear norm to relax the rank constraint to obtain a convex problem, and the Godec algorithm only performs truncated rank-1 singular value decomposition when solving the low-rank matrix. However, the rank-1 matrix thus generated is not sufficient to meet the low-rank matrix structure constraint that each column is equal, and the Godec algorithm is prone to converge to a local optimal solution, so it cannot be applied to detect moving targets in GEO satellite-aircraft bistatic SAR. SUMMARY

[0005] To solve the above technical problems, the present application proposes a moving target detection method for GEO satellite-aircraft bistatic SAR.

[0006] The technical scheme adopted by the present application is as follows: a moving target detection method for GEO satellite-aircraft bistatic SAR, based on a radar system including a GEO satellite and an airborne platform, the signal transmitted by the GEO satellite includes:

[0007] S1, the signal transmitted by the GEO satellite is received by the airborne platform, thereby obtaining a multi-channel two-dimensional time domain signal;

[0008] S2, compensating for the delay and phase of the multi-channel two-dimensional time domain signal of step S1 to obtain a processed echo matrix;

[0009] S3, modeling the echo matrix obtained in step S2 as a matrix rank-1 recovery problem;

[0010] S4, solving the matrix rank-1 recovery problem in step S3 by iterative optimization through ADMM to realize moving target detection.

[0011] The beneficial effects of the present application: the present application proposes a method based on matrix rank 1 decomposition to detect moving targets, the objective function in the optimization design is the weighted sum of the nuclear norm of the auxiliary low-rank matrix and the 1-norm of the sparse matrix, the multi-baseline design problem is modeled as an optimization problem of maximizing the correlation coefficient, and the low-rank auxiliary matrix and the rank 1 coefficient vector are solved alternately by ADMM algorithm. The beneficial effects of the present application are that the traditional RPCA method overcomes the requirement of uniform sparse mode of sparse matrix in low-rank sparse decomposition. The echo detection performance of the present application in non-uniform sparse mode is better than that of the traditional RPCA method. BRIEF DESCRIPTION OF DRAWINGS

[0012] Figure 1 is a flow chart of detecting moving target signals of the present application;

[0013] Figure 2 is a GEO satellite moving target detection geometry configuration diagram;

[0014] Wherein, (a) is the transmission process, (b) is the receiving process;

[0015] Figure 3 is a signal superposition diagram given by signal clutter 0dB;

[0016] Wherein, (a) is the superposition signal, (b) is the moving target signal, and (c) is the clutter signal;

[0017] Figure 4 is the detection result of the present method given by signal clutter ratio-26dB;

[0018] Wherein, (a) is the signal before processing, (b) is the signal after processing. DETAILED DESCRIPTION

[0019] In order to facilitate the description of the content of the present application, first of all, the following belongs to is explained:

[0020] Term 1: GEO satellite-aircraft bistatic SAR

[0021] GEO satellite-aircraft bistatic SAR refers to a bistatic synthetic aperture radar with a geosynchronous orbit (Geosynchronous, GEO) satellite as the transmission source and an airborne platform as the receiving station.

[0022] The present application provides a moving target detection method for GEO satellite-aircraft bistatic SAR, which comprises:

[0023] Step S1: Establishing a radar echo model of GEO satellite-aircraft bistatic SAR system

[0024] The multi-channel two-dimensional time domain signal transmitted by the GEO satellite and received by the airborne platform is:

[0025]

[0026] where τ denotes the range dimension, L is the number of range samples, t denotes the azimuth dimension, and M is the number of pulses. c denotes the speed of electromagnetic wave, f c denotes the carrier frequency of the system. n = 0, 1,..., N - 1 denotes the channel number, and there are N total receive channels; V denotes the velocity of the airborne receiving station along the track direction, x is the initial position of the ground target along the track direction, X n is the initial position of the channel n along the track direction;

[0027] In equation (1), R T (t) denotes the transmitted range, which is represented by the position vector r s (t) of the GEO satellite and the position vector r T (t) of the moving target. Because the synthetic aperture time is short, the first-order motion parameters of the satellite and the moving target are used to approximate the motion of the satellite and the moving target, which means that R T (t) = ||Δr T + Δv T t||, where Δr T = (r s0 - r t0 - v t0 τ nsgs ), denotes the equivalent relative displacement vector of the high-orbit transmitting station and the ground target under the "non-stop-and-go" assumption. A coordinate system OXYZ is established with the ground reference point as the center to describe the motion of the ground moving point, the airborne receiving station, and the satellite. In this coordinate system, r s0 = (x s0 , y s0 , z s0 ) T denotes the initial position vector of the GEO satellite relative to the reference point, where x s0 , y s0 , z s0 are the components of the initial position vector of the satellite in the coordinate system OXYZ, and the superscript T denotes the transpose operation. r t0 = (x t0 , y t0 , z t0 ) T , v t0 = (v xt , v yt , v zt ) T denote the initial position vector and the velocity vector of the ground target relative to the reference point, respectively, where x t0 , y t0 , z t0 are the components of the initial position vector of the ground target in the coordinate system OXYZ, and vxt ,v yt ,v zt are the components of the velocity vector of the ground target in the coordinate system OXYZ. nsgs denotes the signal propagation time under the "no-stop-go" assumption. T = v s0 - v t0 denotes the velocity difference vector between the GEO satellite and the ground target, where v s0 = (v xs , v ys , v zs ) T denotes the velocity vector of the GEO satellite, v xs , v ys , v zs are the components of the velocity vector of the satellite in the coordinate system OXYZ.

[0028] Using the Taylor expansion approximation for the launch distance, we obtain:

[0029]

[0030] The Taylor expansion coefficients for the launch distance are:

[0031]

[0032] r r (t) denotes the position vector of the airborne receiving station with respect to channel 0, R R (t, n) = ||A r R - n d + A v R t|| denotes the reception distance of channel n. d = (dx, dy, dz) T denotes the interval vector using a uniform array. A r R = r t0 - r r0 + (v t0 - v r0 ) t nsgs denotes the equivalent relative displacement vector of the reference channel and the ground target under the "no-stop-go" assumption. r r0 denotes the initial position vector of the reference channel. v r0 = (v xr , v yr , v zr ) T denotes the velocity vector of the airborne receiving station. A v R = v t0 - v r0 denotes the relative velocity vector between the target and the airborne platform.

[0033] Step S2: echo signal preprocessing

[0034] The specific implementation method of this step is time delay and phase compensation operation on the input signal. To achieve this operation, first transform the signal to the frequency domain, and get

[0035]

[0036] In the above formula, f τ represents the distance to frequency, and B represents the signal bandwidth.

[0037] For clutter, the receiving distance between adjacent channels approximately satisfies where d is the channel spacing, so Therefore, consider the channel delay operation, delay the channel n echo Get:

[0038]

[0039] where, represents the double base distance of channel 0;

[0040] Next, use the reference point to construct the following phase compensation function:

[0041]

[0042] Next, perform a phase compensation operation on each channel by the following formula, then:

[0043]

[0044] Back to the two-dimensional time domain by the inverse distance Fourier transform operation:

[0045]

[0046] Step S3: Model the moving target detection as a matrix rank 1 recovery problem

[0047] Since the echo matrix is a complex matrix, the real part and the imaginary part are detected and processed respectively. Since the methods are the same, take the real part as an example. Let the real part of the preprocessed echo matrix s3(t,τ,n) be denote the real set. Each column of it corresponds to a slow-time signal, and it is loaded as a column vector That is, x n =vec(X n ), where vec(·) performs column loading operation. The obtained column vector is assembled into a new matrix That is, X = [x0,…,x N ];

[0048] Since the echo data is superimposed by clutter and moving targets, i.e. X = L + S, where L and S are the corresponding clutter matrix and moving target matrix, respectively. Since the clutter signal s3(t, τ, n) is approximately equal among channels, L is a rank-1 matrix, i.e. rank(L) = 1, and can be written as L = p1 T ;

[0049] where p is the detected column vector, and 1 is the all-1 column vector. Since the number of moving targets is small, S is modeled as a sparse matrix. Therefore, moving targets can be detected by the reconstruction of rank-1 matrix and sparse matrix, which means solving the following optimization problem:

[0050]

[0051] In the above formula, ‖·‖ * represents the kernel norm of the matrix, i.e. the sum of all singular values. ‖·‖1 is the 1-norm of the matrix, λ is the penalty parameter, and M is the auxiliary variable.

[0052] Step S4: Iterative solution of optimization problem

[0053] The rank-1 matrix recovery algorithm iteratively optimizes the designed optimization problem by ADMM. And the two sub-problems of auxiliary matrix M and vector p need to be alternately optimized by the singular value shrinkage operator and the negative gradient descent method;

[0054] In step 4, the two sub-problems of auxiliary matrix M and vector p are alternately optimized by inexact ADMM to effectively solve, and the specific steps are as follows:

[0055] Step S4-1: Initialize iteration parameters

[0056] (1) Reload the formed matrix as X, and the penalty factor as λ;

[0057] (2) Initialize the Lagrange multiplier matrix as Y, and the augmented Lagrange penalty factor as μ;

[0058] Step S4-2: If the number of iterations is less than the maximum number of iterations, execute step S4-3, otherwise, go to S4-7; In this embodiment, the maximum number of iterations is 20.

[0059] Step S4-3: At the m+1th iteration, first fix p m and Y m , solve M m+1

[0060] (1) This step solves the following optimization sub-problem:

[0061]

[0062]

[0063] In the above formula, <A,B> = tr(A T B) represents the inner product of the trace of the matrix, ||·||F F represents the Frobenius norm of the matrix, which is numerically equal to the 2-norm of the matrix, whose square is equal to the inner product of the trace of the matrix and itself;

[0064] (2) p m , Y m respectively represent the vector p and the Lagrange multiplier Y obtained at the mth iteration, since both are fixed when solving the (m+1)th time, the optimization problem is converted to:

[0065]

[0066] After solving, the solution is Here D(A,a) represents the singular value shrinkage operator of the matrix, D(A,a) = US(Σ,a)V T , where A = UΣV T represents the singular value decomposition of the matrix A, and S(B,b) represents the soft thresholding operator of the matrix;

[0067] Step S4-4: At the m+1th iteration, fix M m+1 and Y m Solve p m+1

[0068] (1) At this step, solve the following optimization subproblem:

[0069]

[0070]

[0071] (2) Let That is, optimize the following problem:

[0072]

[0073] Solving the above unconstrained optimization problem using the negative gradient descent method can obtain p m+1 , the gradient vector of the objective function is:

[0074]

[0075] In the above formula, sgn(·) is the sign function in the sense of matrix elements, therefore, the gradient descent subproblem is:

[0076] Step S4-4-1 initializes the optimization variable p 0 and the iteration step size ρ 0

[0077] Step S4-4-2 If the iteration number is less than the maximum iteration number, step S4-4-3 is performed, otherwise, go to S4-4-5

[0078] Step S4-4-3 Gradient descent and step length update are performed

[0079]

[0080] p j+1 = p j * q (18)

[0081] In the above formula, p j represents the vector p at the jth iteration in the gradient descent subproblem, 0 < q < 1 represents that the iteration step length gradually decreases during the iteration process;

[0082] Step S4-4-4 j <- j + 1, indicating that the iteration number is incremented by 1, and go to step S4-4-2

[0083] Step S4-4-5 The subproblem iteration optimization is completed, that is, the solution is obtained

[0084] Step S4-5: Update the augmented Lagrange multiplier matrix Y m+1

[0085] Y m+1 = Y m + u (M - p1 T ) (19)

[0086] Step S4-6: m <- m + 1, indicating that the iteration number is incremented by 1, and go to step S4-2

[0087] Step S4-7: The iteration optimization is completed, and the sparse matrix is obtained

[0088] The above steps are respectively decomposed for the real part and the imaginary part of the echo matrix s3(t,tau,n), one column of the sparse matrix is expanded into a matrix by column, and the matrices obtained by respectively decomposing the real part and the imaginary part are envelope synthesized to obtain the moving target indication matrix. After the above steps are performed, the indication result of the echo of the clutter submerged moving target in the signal can be obtained. Then, the traditional moving target parameter estimation method can be used to further estimate the motion parameters of the moving target.

[0089] The present application is mainly verified by simulation experiment method, and all steps and conclusions are verified correct on Matlab2018. The present application will be further described in detail in combination with the simulation parameters shown in Table 1.

[0090] Table 1 Simulation parameters

[0091] Parameter Value Parameter Value Orbital semi-major axis 42241 Km Signal bandwidth 50 MHz Orbital eccentricity 0.012 Carrier frequency 1.25 GHz Right ascension of ascending node 120° Antenna separation 2m Orbital inclination 50° Pulse repetition frequency 100 Hz Argument of perigee 20° Range sampling frequency 100 MHz On-board receiving station altitude 10000m On-board platform velocity 200 m / s Number of channels 3 Synthetic aperture time 2s

[0092] Step one: generate GEO bistatic SAR multi-channel signals for clutter and moving target respectively, system parameters are as follows Figure 2 , wherein the clutter and moving target are both located at the non-scene center point.

[0093] First, set the signal-to-clutter ratio to 0dB, Figure 3 The superposition result of the signals is shown. Figure 3 (a) shows the two-dimensional time domain distribution of the superimposed signal formed by mixing the clutter and moving target signals with equal amplitudes. Figure 3 (b) and Figure 3 (c) are the moving target signal and the clutter signal respectively. Since the amplitudes of the clutter and moving target signals are set to be equal here, the clutter in (a) does not drown out the moving target signal. Next, adjust the signal amplitude and set the signal-to-clutter ratio to -26dB. Figure 3

[0094] Step two: perform delay and phase compensation preprocessing operations.

[0095] Step three: complete the rank-1 recovery of the matrix by iteratively solving the optimization problem through ADMM. And the two sub-problems need to be alternately optimized and solved through the singular value shrinkage operator and the negative gradient descent method respectively;

[0096] Figure 4 The signal processed by the method of the application is shown. Figure 4 (a) shows the two-dimensional time domain distribution of the superimposed signal formed by mixing the clutter and moving target signals. It can be clearly seen that since the amplitude of the clutter is much higher than that of the moving target signal, the moving target signal is completely drowned in the clutter and cannot be accurately detected.

[0097] Figure 4 (b) shows the detected moving target result after the method of the application is processed. Comparing the moving target envelope of Figure 3 (b) and Figure 4 the detection result of the moving target of Figure 4 (b), it can be seen from (b) that the method can detect the moving target and display the time domain range of the moving target signal.

[0098] Figure 4 The detected moving target result in (b) is clear, indicating that the method of the application can detect the moving target from the clutter after range compression.

[0099] Compared with the case where the moving and stationary targets are mixed together after separation by the traditional method and the separation is not clean, the separation performance of the method of the application is obviously better than that of the traditional method.

[0100]

[0100] Those skilled in the art will appreciate that the embodiments described herein are presented for purposes of illustration and understanding of the principles of the application and should not be construed as limiting the scope of the application to such specifically enumerated embodiments. Various modifications and changes can be made to the application by those skilled in the art which will be apparent from this disclosure without departing from the spirit and principles of the application. Any modifications, equivalent substitutions, improvements, etc. made to the application should be included within the scope of the application as defined in the following claims.

Claims

1. A method for moving target detection for GEO satellite-airborne bistatic SAR, based on a radar system comprising a GEO satellite and an airborne platform, the signal transmitted by the GEO satellite, characterized in that, The method comprises the following steps: S1, receiving a signal transmitted by a GEO satellite through an airborne platform, thereby obtaining a multi-channel two-dimensional time domain signal; S2, performing a delay and phase compensation operation on the multi-channel two-dimensional time domain signal of step S1, thereby obtaining a processed echo matrix; S3, modeling the echo matrix obtained in step S2 as a matrix rank 1 recovery problem; step S3 is specifically: denoting the real part of the echo matrix obtained in step S2 as denoting the real number set, as a pulse number, each column of which corresponds to a slow-time signal, and loading it as a column vector , that is, , wherein performing the column loading operation; assembling the obtained column vector into a new matrix , that is, ;​​ The echo data is superimposed by clutter and moving targets, that is where the corresponding clutter matrix and moving target matrix are denoted as and respectively; since the clutter signals are approximately equal among the channels, the clutter matrix is a rank-1 matrix, that is , and is written as ; wherein is a column vector of detections, is a column vector of full is a column vector of detections; Since the number of moving targets is small, the matrix is modeled as sparse matrix; Therefore, the moving target is detected through reconstruction of the rank 1 matrix and the sparse matrix, and an optimal problem expression is obtained as follows: ； In the above formula, denotes the nuclear norm of a matrix; denotes the 1-norm of a matrix, is a penalty parameter, is an auxiliary variable; S4, performing an iterative optimization on the matrix rank 1 recovery problem in step S3 through ADMM, thereby achieving moving target detection; step S4 specifically comprises: S41, initialize iteration parameters, including: reload the formed matrix as , the penalty factor is ; initialize the Lagrange multiplier matrix as , the augmented Lagrange penalty factor is ; S42, if the current iteration number is less than the set maximum iteration number, step S43 is executed, otherwise, step S46 is executed; S43, record the current iteration number as m+1, first fix and , solve : ； ； wherein denotes the trace inner product of matrices, denotes the Frobenius norm of a matrix, which is numerically equal to the 2-norm of the matrix, whose square is equal to the trace inner product of the matrix with itself; denote the vector obtained at the th iteration and the Lagrange multiplier Since both are fixed at the th solution, the optimization problem is transformed into: ； After solving, the solution is ; here denotes the singular value shrinkage operator of a matrix, where denotes the singular value decomposition of a matrix and denotes the soft thresholding operator of a matrix, the superscript T denotes the transpose operation; Then fix And Solve : ； ； Recall i.e. to optimize the problem of ； The solution is obtained using the negative gradient descent method. The gradient vector of the objective function is: ； wherein is the sign function in the sense of a matrix element; S44, updating the augmented Lagrange multiplier matrix : ； S45、 the number of iterations plus 1, and the process goes to step S42. S46, iteration optimization ends, get sparse matrix .

2. The method according to claim 1, wherein, The multi-channel two-dimensional time domain signal expression of step S1 is as follows: ； wherein, denotes the time in distance direction, and let the number of sampling points in distance direction be , denotes the time in azimuth direction, and B denotes the signal bandwidth, denotes the wave velocity of electromagnetic wave, denotes the transmitting distance, denotes the receiving distance of channel , denotes the carrier frequency of the system, denotes the channel number, is the total number of receiving channels; denotes the velocity of airborne receiving station along the track direction, is the initial position of ground target along the track direction, is the initial position of channel along the track direction.

3. The method according to claim 2, wherein, Step S2 specifically comprises: Firstly, the multi-channel two-dimensional time domain signal is transformed into a range frequency domain, thereby obtaining: ； wherein represents the distance to frequency, B represents the signal bandwidth; Considering the channel delay operation, the channel echo delay , we get: ； wherein is the channel spacing, , denotes the double base distance of the channel and, ; Next, a phase compensation function is constructed using a reference point: ； The phase compensation operation is performed on each channel through the following formula, thereby obtaining: ； Then, the inverse Fourier transform operation is performed on the distance vector to return to the two-dimensional time domain, thereby obtaining: 。 4. The method according to claim 3, wherein, The solution process is as follows: A1, initialize optimization variable and iteration step ; A2, if the iteration number is less than the maximum iteration number, step A3 is executed, otherwise, step A5 is executed; A3, performing a negative gradient descent and step length update: ； ； wherein denotes the vector at the iteration of the gradient descent subproblem , denotes that the iteration step size is gradually reduced during the iteration process; A4、 , indicating the number of iterations plus 1, go to step A2; A5, solve for .

Images