Next Article in Journal
Assessing Snow Phenology over the Large Part of Eurasia Using Satellite Observations from 2000 to 2016
Next Article in Special Issue
Digital Image Correlation of Google Earth Images for Earth’s Surface Displacement Estimation
Previous Article in Journal
Quantification of Precipitation Using Polarimetric Radar Measurements during Several Typhoon Events in Southern China
Previous Article in Special Issue
The Implications of M3C2 Projection Diameter on 3D Semi-Automated Rockfall Extraction from Sequential Terrestrial Laser Scanning Point Clouds
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

A Novel 2-D Geometry Reconstruction Approach for Space Debris via Interpolation-Free Operation under Low SNR Conditions

1
School of Physics and Optoelectronic Engineering, Xidian University, Xi’an 710071, China
2
National Key Laboratory of Science and Technology on Space Microwave, Xi’an 710100, China
3
School of Microelectronics and Communication Engineering, Chongqing University, Chongqing 400044, China
4
Chongqing Key Laboratory of Space Information Network and Intelligent Information Fusion, Chongqing University, Chongqing 400044, China
5
Chongqing Key Laboratory of Mobile Communications Technology, Chongqing University of Posts and Telecommunication, Chongqing 400065, China
*
Author to whom correspondence should be addressed.
Remote Sens. 2020, 12(12), 2059; https://0-doi-org.brum.beds.ac.uk/10.3390/rs12122059
Submission received: 5 June 2020 / Revised: 21 June 2020 / Accepted: 23 June 2020 / Published: 26 June 2020
(This article belongs to the Special Issue Remote Sensing in Engineering Geology)

Abstract

:
Two unsolved key issues in inverse synthetic aperture radar (ISAR) imaging for non-cooperative rapidly spinning targets including high computational complexity and poor imaging performance in the case of low signal-to-noise ratio (SNR) are addressed in this work. In the strip-map imaging mode of SAR, it is well known that azimuth spatial invariant characteristics exist, and inspired by this, we propose a fast ISAR imaging approach for spinning targets. Our approach involves two steps. First, a precise analytic expression in the range-Doppler (RD) domain is produced using the principle of stationary phase (POSP). Second, a novel interpolation kernel function is designed to remove two-dimensional (2-D) spatial-variant phase errors, and the corresponding fast implementation steps that only require Fourier transform and multiplications are also presented. Finally, a well-focused ISAR image is obtained by compensating the azimuth high-order terms. Compared with current imaging methods, our approach avoids multi-dimensional search and interpolation operations and exploits the 2-D coherent integrated gain; the proposed method is of low computational cost and robustness in the low SNR condition. The effectiveness of the proposed approach is confirmed by numerically simulated experiments.

1. Introduction

The rapid developments of the aerospace industry raised the issue of space debris, which pose a serious threat to spacecraft or satellites in orbit, with potentially disastrous consequences [1,2,3]. In turn, detecting, identifying, and cleaning space debris is currently an important research topic. Among the different methods to monitor space debris, inverse synthetic aperture radar (SAR) (ISAR) represents a powerful tool able to produce high-resolution images of spinning targets, and it has been widely adopted to obtain space debris images [4,5,6]. However, the trajectories of space debris are complex, and besides translational motion, they are usually characterized by rapid spinning along the main rotating axis. This makes range-Doppler (RD) [7,8] and range instantaneous Doppler (RID) algorithms for ISAR imaging [9,10] unusable and, in turn, it makes it challenging to produce well-focused high-resolution ISAR imaging for the spinning target. One of the problems is that fast rotating targets cause great range cell migration (RCM) and non-stationary Doppler frequency modulation (DFM). Moreover, the RCM and DFM exhibit the two-dimensional (2-D) spatial-variant features [11,12,13,14,15,16,17,18], which are correlated to fast and slow times dynamics, and with the positions of scattering centers. If the 2-D spatial-variant dependent phases are not properly corrected, the final RD-based images would be blurred. Another issue, arising from the spinning motion of space debris and the long imaging distance, is the strong noise, which makes the signal-to-noise ratio (SNR) of the echo signal very low and prevents a precise estimation of the motion parameters [11,17]. The overall effect is that the performance of the conventional ISAR imaging for spinning target is seriously degraded.
To overcome the above mentioned issues and produce images of rapidly spinning targets, several ISAR approaches with narrowband and wideband transmitted signals have been suggested and developed [11,12,13,14,15,16,17,18,19,20,21]. In [11] a single-range Doppler interferometry approach (SRDI) for spinning targets was developed building on the properties of the Doppler spectrum in one rotation period. SRDI detects the spinning targets by transforming the 1-D Doppler signal into a 2-D time-frequency signal, from which the energy of spinning target may be accumulated into a sharp peak. However, its performance is degraded from cross-terms, interferences, and resolution loss [11,12]. Moreover, its complexity is high, due to the need for a multi-dimensional search. To overcome the shortcomings of SRDI, a single-range matching filtering method (SRMF) was developed [13,14] to produce well-focused an ISAR image. To this aim, a match function was designed with different rotating radii and initial phases to compensate the transverse echo signal cell-by-cell. The clean technique was also exploited to migrate the high sidelobes effects. In [15], the authors developed a single-range image fusion method to produce 2-D images in the Cartesian coordinates. The polar format algorithm (PFA) [16], the modified complex valued back-projection algorithm (BPA) [17], and the tomography algorithm [18] were also used to obtain the images of rotating targets. However, all those techniques have high computational complexities due to the involved integrals, and still lack fast implementation. In [19], a real-valued generalized Radon transform (GRT) approach was presented to achieve imaging of spinning targets, which employs the sinusoidal envelop of the echo signal for spinning targets to obtain an ISAR image via non-coherent accumulation operation. As such, it cannot obtain satisfactory imaging performance in low SNR conditions. In [20], an extended Hough transform (HT) was introduced to realize ISAR imaging of rotating targets. However, the algorithm is time-consuming, due to the multi-dimensional parameter search. Upon using the sparse property of the echo signal and the matching pursuit (MP) algorithm, a micro-motion imaging approach was proposed in [21], which, however, has a high computational cost and fails when Doppler ambiguity is present. In [22], well-focused ISAR images of rapidly spinning targets were generated using the segmental pseudo Keystone transform (KT) (SPKT). However, this technique is also time-consuming, and it suffers from accuracy loss due to the need for an initial phase estimation stage and to the interpolation operation used for each scattering point of the spinning target. Overall, despite all the aforementioned algorithms having their merits, a globally effective ISAR imaging method for non-cooperative spinning targets is still lacking, especially for those situations where the SNRs are low.
In this paper, exploiting the azimuth spatial invariance of strip-map SAR imaging mode [23,24], an effective imaging approach for rapidly spinning targets is suggested, which works also with a low SNR environment. To this aim, we first formulate the RCM and time-varying DFM caused by the target’s spinning motion in terms of range and azimuth, i.e., a 2-D structure. Then, using the principle of stationary phase (POSP) [25,26], a precise analytic expression in RD domain is obtained, where the scatters located in one range cell with different azimuth positions are assigned to the same range migration trajectory to effectively eliminate the azimuth spatial variance of the initial phase. Third, a novel interpolation kernel function is designed to remove the range variance of the RCM and DFM, and an interpolation-free mapping algorithm is suggested to further reduce the computational burden. Finally, by compensating the azimuth high-order terms, a high resolution ISAR image for rapidly spinning targets is obtained. Compared to currently available imaging methods, our technique can be efficiently implemented, since it only requires Fourier transform (FT) without the need of any interpolation or multi-dimensional search operations. Furthermore, thanks to the use of a 2-D coherent integrated gain, a satisfactory imaging performance in low SNR condition may also be obtained. We also performed numerically simulated experiments to verify the correctness and the effectiveness of the proposed algorithm.
The paper is organized as follows. In Section 2, the geometry and signal model for rapidly spinning targets are provided. In Section 3, a novel efficient ISAR imaging algorithm in the RD domain using the SAR technique is described. The fast implementation of the method and some remarks about its practical applications are presented in Section 4. In Section 5, results from simulated experiments are given together with their analysis. In Section 6, we discuss the computational complexity of the method, whereas Section 7 closes the paper with some concluding remarks.

2. Imaging Model of Spinning Target and Problem Formation

The geometric configuration of ISAR imaging for spinning space debris is depicted in Figure 1, where the origin O denotes the position of the center of the rotating target. During the imaging period, the motion of the target may be decomposed into translation and rotation. The translational motion is parallel to the line-of-sight (LOS) and does not contribute to the cross-range resolution. On the other hand, rotation changes the Doppler frequency and affects the azimuth imaging. In Figure 1, we denote by P an arbitrary scattering point that rotates around the origin O with constant rotation rate ω . The spinning angular rate ω is the sum of the self-spinning angular velocity and a constant term due to the time-variant LOS angle [18,19,20,21,22]. The distance between the spinning target center and the radar is denoted by R 0 ( R 0 r p ) , where   φ p is the initial phase of the spinning scatterer P .
During the rotation, shadowing of the scatterers will cause the loss of the echo signal of the spinning target. The non-uniform rotating rate, or residual translational movement, will result in defocused ISAR images. However, the current paper is focused on designing an efficient imaging algorithm based on the presented imaging model. Here we assume that translational motion compensation (TMC) has been implemented by the available TMC algorithm [27,28] and focus on the effects of target rotation. According to Figure 1, the instantaneous slant range r ( t a ; r p ) between the radar phase center and the scatterer in P is given by
r ( t a ; r p ) = R 0 + r p sin ( ω t a + φ p )
where t a represents the azimuth slow-time variable. Assuming that the radar transmits a linear frequency-modulated (LFM) signal, the received baseband echo signal of the point target after demodulation is
S ( t r ; t a ; r p ) = σ p ω r [ t a 2 r ( t a ; r p ) c ] e x p [ j 2 π f c ( t r 2 r ( t a ; r p ) c ) ] × ω a ( t a ) e x p [ j π γ ( t r 2 r ( t a ; r p ) c ) 2 ]
where σ p is the reflection coefficient of the scatterer in   P , t r is the range fast-time variable, and ω r ( · ) and ω a ( · ) represent the envelope in the range and azimuth direction, respectively. c , f c , and γ denote the speed of light, the carrier frequency of radar signal, and the chirp rate, respectively.
Performing FT with respect to the range fast-time variable t r and removing the range frequency-modulated term, the received echo signal becomes
S ( f r ; t a ; r p ) = σ p ω a ( t a )   ω r ( f r ) e x p [ j 4 π ( f r + f c ) R 0 λ ] e x p [ j 4 π   f r λ r p sin ( ω t a + φ p ) ] × e x p [ j 4 π λ r p sin ( ω t a + φ p ) ]
where λ denotes the carrier wavelength, and f r is the range frequency variable with respect to t r .
In Equation (3), the phase consists of three terms. The first is the linear part corresponding to the range direction position, the second term is the RCM, and the third one is due to the azimuth time-varying DFM. As a matter of fact, the RCM and DFM terms contain a 2-D spatial dependence, which is due to rotation of the target and makes it difficult to obtain an effective imaging scheme for spinning targets. Indeed, unless the second and third terms are satisfactorily compensated, the final ISAR image will be blurred. Notice that the 2-D spatially dependent phase terms are different for each scatterer, and therefore, the compensation should be performed pixel by pixel. In order to face this problem while utilizing the sinusoidal envelope of the spinning targets in the range slow-time domain, existing GRT methods correct the phase via multi-dimensional searching operation along the rotating radius r p and initial phase φ p for all scatterers in the spinning target [19]. Performance is, however, poor due to the non-coherent accumulation. In the SPKT method, the initial phase estimation and the interpolation operation are required for each scatterer in the spinning target [22], thus leading to a great computational burden and a low imaging quality in the environment, inducing a low SNR. Here, by exploiting the inherent azimuth invariance of SAR, a fast ISAR imaging method in the case of spinning targets is developed to produce a well-focused ISAR image for spinning targets.
The instantaneous Doppler frequency f r is calculated by taking the derivative of the phase term in (3) with respect to t a as
f d ( t a ) = 1 2 π d [ ϕ ( t a ) ] d t a = 2 ω r p cos ( ω t a + φ p ) λ
From the Nyquist sampling theorem, the pulse repetition frequency (PRF) must satisfy PRF 2 f d m a x . Thus, the maximum of Doppler frequency is f d m a x = 2 ω r p m a x / λ , where r p m a x is the longest rotation radius of the targets. The Doppler bandwidth of scatterers is then
B d = 4 ω r p λ

3. Proposed Imaging Approach Description

In this section, we focus on the 2-D spatial-variant RCM and DFM fast correction under low SNR conditions. A fast imaging approach for non-cooperative spinning targets is developed. Before going into the details of the imaging algorithm, let us summarize the principles of azimuth invariant SAR strip-map imaging mode. A schematic representation is given in Figure 2; panel 2(a) shows the geometric model of multiple targets located in different range and azimuth cells at a given time, whereas panels 2(b) and 2(c) show the situation after range compression. For targets located in the same range cell at different azimuth times, the migration trajectories are overlapped in the RD domain, and this is the well-known azimuth invariance of SAR strip-map imaging mode [23,24,25,26].
In our proposed method, in the Doppler domain and based on the above-mentioned invariance, echo signals of the same distance are only processed once, thus reducing the computational complexity.

3.1. Analytical Expression Derivation in RD Domain

Scatterers with the same rotating radius have overlapped trajectories in the range Doppler domain, and thus the RCM correction located in same range gate with different azimuth time instant can be realized just once, thus improving imaging efficiency. In order to obtain an expression in the range Doppler domain after range compression, let us perform FT of the slow-time variable t a in Equation (2), thus arriving at
S ( r ; f a ; r p ) = σ p ω a ( t a ) s i n c [ 2 π B r c ( r r ( t a ; r p ) ) ] e x p [ j 4 π r ( t a ; r p ) λ ] × e x p ( j 2 π f a t a )   d t a
where B r is the transmitted LFM signal bandwidth, and s i n c denotes the range envelope. The integral in (6) cannot be evaluated analytically due to the sin function term in the slant range formula r ( t a ; r p ) , and to obtain an analytic expression we resort to POSP as in the following. The phase expression in (6) is
Φ ( t a ) = 4 π r p sin ( ω t a + φ p ) λ 2 π f a t a
According to the POSP, the value in (7) is calculated by the stationary point, that is,
d Φ ( t a ) d t a = 4 π [ ω r p cos ( ω t a + φ p ) ] λ 2 π f a = 0
By solving the Equation (8), one obtains the stationary point
t a * = { 1 ω arc cos ( f a λ 2 r p ω ) φ p ω , ω t a + φ p [ 0 , π ] 1 ω arc cos ( f a λ 2 r p ω ) φ p ω , ω t a + φ p [ π , 0 ]
which is then substituted in Equation (1), leading to the instantaneous slant distance expression in the RD domain
r ( f a ) = { r p sin [ arc cos ( f a λ 2 r p ω ) ] , ω t a + φ p [ 0 , π ] r p sin [ arc cos ( f a λ 2 r p ω ) ] , ω t a + φ p [ π , 0 ]
S ( r ; f a ; r p ) = { σ p ω a ( f a ) sin c [ 2 π γ T p c ( r r ( f a ) ) ] e x p ( j 2 π f a 1 ω a r c cos ( f a λ 2 r p ω ) ) , e x p ( j 2 π f a φ p ω ) e x p [ j 4 π r ( f a ) λ ] , ω t a + φ p [ 0 , π ] σ p ω a ( f a ) sin c [ 2 π γ T p c ( r r ( f a ) ) ] e x p ( j 2 π f a 1 ω a r c cos ( f a λ 2 r p ω ) ) e x p ( j 2 π f a φ p ω ) e x p [ j 4 π r ( f a ) λ ] , ω t a + φ p [ π , 0 ]
Equation (10) reveals that r ( f a ) includes a trigonometric function sin , which is generated by the fast rotation of the targets. Finally, the expressions of the echo signal in the range-Doppler in Equation (11) is obtained, where the second exponential term contains the initial phase φ p , and the first and the third exponentials represent high order phase terms.

3.2. The Proposed ISAR Imaging Method

In Equation (11), in the high order phase term or the sin c envelop, the initial phase φ p is absent. The range sin c envelop is, in turn, only related to the rotating radius, showing that scatterers in the same range cell and different azimuth instantaneous times have the same trajectories in the RD domain. Let us now describe in some detail the steps needed to derive the proposed ISAR imaging method.
In Equation (11), the range envelope term sin c [ · ] depends on the radial range r ( f a ) and the azimuth frequency f a , which cause the RCM. Because of the RCM, the energy of the same target is scattered in different range cells in the original ( r f a ) plane. In order to efficiently remove the RCM, we introduce a novel interpolation kernel function [29,30], given by
Δ r ( f a ) = { r p sin [ a r c cos ( f a λ 2 r p ω ) ] r p , ω t a + φ p [ 0 , π ] r p sin [ a r c cos ( f a λ 2 r p ω ) ] r p , ω t a + φ p [ π , 0 ]
After applying the interpolation, the signal in the RD domain is given in (13), where RCM has been removed. The phenomenon is illustrated in Figure 3. The original discrete signal in the ( r f a ) plane is shown in Figure 3a, where the sample points of the Doppler frequency are equally spaced. The signal after traditional KT in the ( r f a ) plane is illustrated in Figure 3b, where Doppler frequency intervals are varied with range sampling points, but the intervals between range sampling points are identical. The signal of the proposed method in the ( r f a ) plane is demonstrated in Figure 3c, where different Doppler frequencies have different range sampling intervals.
S R C M ( r ; f a ; r p ) = { σ p ω a ( f a ) sin c [ 2 π γ T p c ( r r p ) ] e x p ( j 2 π f a 1 ω a r c cos ( f a λ 2 r p ω ) ) e x p ( j 2 π f a φ p ω ) e x p [ j 4 π r ( f a ) λ ] , ω t a + φ p [ 0 , π ] σ p ω a ( f a ) sin c [ 2 π γ T p c ( r r p ) ] e x p ( j 2 π f a 1 ω a r c cos ( f a λ 2 r p ω ) ) e x p ( j 2 π f a φ p ω ) e x p [ j 4 π r ( f a ) λ ] , ω t a + φ p [ π , 0 ]
Upon comparison of the original RD domain expression in (11) to the RD domain expression after transformation in (13), we see that along the sin c trajectory, the energies of the scatterers are transformed into one range cell r p , with the sin c [ · ] term depending only on the rotating radius r p . The azimuth frequency f a has been removed from the range envelope, which means that the RCM has been properly corrected (i.e., compensated).
Let us now address the azimuth processing procedure. The compensation function H ( r , f a ) is designed to compensate the high order phase terms in (13) as
H ( r , f a ) = { e x p ( j 2 π f a 1 ω a r c   c o s ( f a λ 2 r p ω ) ) e x p [ j 4 π r ( f a ) λ ] , ω t a + φ p [ 0 , π ] e x p ( j 2 π f a 1 ω a r c   c o s ( f a λ 2 r p ω ) ) e x p [ j 4 π r ( f a ) λ ] , ω t a + φ p [ π , 0 ]
Multiplying (13) by (14), the signal is given by
S H C ( r , f a ) = σ p ω a ( f a ) sin c [ 2 π γ T p c ( r r ( f a ) ) ] e x p ( j 2 π f a φ p ω )
such that the energy of a single scatterer P is still in the same range cell r p . Finally, by applying an azimuth inverse fast Fourier transform (IFFT), the ISAR image expression becomes
S ( r , t a ) = σ p sin c [ 2 π γ T p c ( r r p ) ] sin c [ π B a ω ( ω t a + φ p ) ]
From (16), one obtains the coordinates ( r p ,   φ p ) of the scatterer in   P and, in turn, an accurate imaging of the target.

4. A Fast Implementation of the Proposed Method and Some Remarks

4.1. Fast Implementation

In general, the novel interpolation operation in the RD domain proposed in Section 3 may increase the computational complexity and also introduce loss of accuracy. On the other hand, the data for non-sampling points may be obtained through a zero-padding in another domain, where the interpolation operation is replaced by FFTs and phase multiplication [31,32,33]. As a consequence, the drawback mentioned above may be eliminated. In fact, the interpolation function may be reformulated as
Δ r ( f a ) = r p · α f a ( r p   )
where the nonlinear shifting coefficient α f a ( r p   ) is given by
α f a ( r p   ) = { sin [ a r c cos ( f a λ 2 r p ω ) ] 1 , ω t a + φ p [ 0 , π ] sin [ a r c cos ( f a λ 2 r p ω ) ] 1 , ω t a + φ p [ π , 0 ]
In order to obtain a fast implementation of the interpolation-free mapping method, the shift length Δ r ( f a ) for the azimuth bin is decomposed into two terms, given by
Δ r ( i ) = { n i · Δ r inter bin + δ i · Δ r sub bin }
where n i is an integer, δ i ( 0 , 1 ) , and Δ r is the range cell size. Then, the interpolation operation may be accomplished through the following steps:
(1)
According to (7), one extracts the range cell data for each azimuth bin in the RD domain, i.e.,
S r = S ( r , f a i )
and calculates n i and δ i according to (10) and (12);
(2)
Zero-padding S r i to S r of length 4 to achieve
S r = [ 0   0    S r i   0   ]
and indexing
I i = [ i 2   i 1   i   i + 1   ] n i
circularly in [ 1 , ,   N ] , where N is original data length;
(3)
FFT of S r to obtain S f r ;
(4)
Phase multiplication
S f r = S f r · e x p { j δ i Δ r · f r }
where f r is the new range frequency;
(5)
IFFT of S f r to complete the sub-bin shift of S r ;
(6)
Coherent summation
F [ I i ] = F [ I i ] + S r
where F is the interpolation mapped data;
(7)
Iterating step (1) to step (6) for all bins to achieve interpolation-free mapping of the full image.
The interpolation-free mapping method replaces the RD domain interpolation operation through zero-padding in the wavenumber domain, which reduces the computational complexity and ISAR image errors. The flowchart of the algorithm and its fast implementation is shown in Figure 4.

4.2. Some Remarks about the Practical Application of the Proposed Method

(1) Range and Azimuth Resolution Analysis:
The range resolution ρ r of ISAR imaging is determined by the transmitted signal bandwidth B r
ρ r = c 2 B r
whereas the azimuth resolution ρ a depends on the azimuth Doppler frequency resolution
ρ a = λ Δ f d   2 ω
where Δ f d is Doppler bandwidth, i.e., the inverse of the coherent imaging accumulation time Δ f d = 1 / T . The azimuth resolution may be thus rewritten as
ρ a = λ   2 ω T = λ   2 Δ θ
where Δ θ is the rotation angle of the spinning target relative to the radar platform during the coherent accumulation time. In summary, the range resolution of ISAR is determined by the radar transmitting signal bandwidth, whereas the azimuth resolution is obtained by the target’s rotation angle relative to the radar’s line of sight during the coherent imaging accumulation.
(2) Effects of Estimation of ω on ISAR Images:
As mentioned before, the proposed method is based on the knowledge of the angular velocity ω . In practice, however, the estimation of the angular velocity may not completely accurate. In order to assess the effects of this uncertainty on the final ISAR images, let us remember that the imaging results are obtained by coherent accumulation through the FT, which means the final phase error Δ ψ should be restricted in | Δ ψ | < π 2 . Assuming that the uncertainty of ω is δ ω , the phase error of the spinning target during the coherent time Δ t ( 0 < Δ t < π / ω ) may be expressed by
Δ ψ Δ t 4 π λ r m a x { cos [ ( ω + δ ω ) Δ t + φ p ] cos ( ω Δ t + φ p ) } = S π λ   r m a x sin ( δ ω Δ t 2 ) sin ( ω Δ t + φ p + δ ω Δ t 2 )
where r m a x is maximum rotating radius for the spinning target. Assuming δ ω Δ t / 2 1 , we obtain sin ( δ ω Δ t / 2 ) δ ω Δ t / 2 and sin ( ω Δ t + φ p + δ ω Δ t / 2 ) 1 . The constrained condition in (28) can be written as
|   δ ω ω π | < λ   2 r m a x
Equation (29) provides a sufficient and unnecessary condition for obtaining well-focused ISAR images for the spinning target and should be taken into account in selecting the system parameters.
(3) Removing the effects of the sidelobes:
It may happen that the sidelobes of a strong target hide the mainlobe of a weak target. Therefore, clean technology [34] is employed to suppress the sidelobes. The process of clean method proceeds as follows:
(a)
Estimate the initial phase φ p and rotation radius r p of all scatterers. The foregoing imaging algorithm may have already estimated these parameters.
(b)
According to the estimated initial phase φ p and rotation radius r p , a point spread function (PSF) is designed. Assuming that the coordinate of the maximum peak is ( r p ξ , φ p ξ ) , the PSF is given by
X ( r , ψ ) = sin c [ 2 π γ T p c ( r r p ξ ) ] sin c [ 2 π B a ( ψ φ p ξ ) ]
(c)
Using a minimum norm criterion, the scattering coefficient at ( r p ξ , φ p ξ ) is estimated by
min σ ^ [ I ( σ ^ ) ] = min σ ^ s σ ^ X = min σ ^ r , ψ | s ( r , ψ ) σ ^ X ( r , ψ ) | 2
where σ ^ is the estimated value of the scattering coefficient. Performing the derivative of I ( σ ^ ) with respect to σ ^ leads to
I ( σ ^ )   σ ^   = 2 r , ψ [ s ( r , ψ ) σ ^ X ( r , ψ ) ] X * ( r , ψ )
where X * ( r , ψ ) represents the conjugate function of X ( r , ψ ) . When I / σ = 0 , the reflection coefficient of scatterer in coordinate ( r p ξ , φ p ξ ) is obtained. Under this condition, the estimated value of the scattering coefficient is
σ ^ = s ( r , ψ ) X * ( r , ψ ) | X ( r , ψ ) | 2
(d)
Steps (a) to (c) may be biased because of the effects of sidelobes, which could affect the usage of the clean approach. The coordinate estimation in Step (c) initializes the following fine search.
min σ , r p ξ ,   φ p ξ I ( σ ^ ) = min σ , r p ξ ,   φ p ξ s ( r , ψ ) σ ^ X ( r , ψ ) = min σ , r p ξ , φ p ξ r , φ | s ( r , ψ ) σ ^ X ( r p ξ ,   φ p ξ , ψ ) | 2
Using this technique, the resulting r p ξ ,   φ p ξ ,   σ ^ are accurate estimates of the coordinates and the reflection coefficient. Subtract the contribution of scatterer ( r p ξ ,   φ p ξ ) from the original signal, thus arriving at
s c l e a n ( r , ψ ) = s σ ^ X
(e)
Use s c l e a n ( r , ψ ) as the new input to step (a), and repeat until convergence.

5. Simulation Results and Discussion

In this section, results from numerically simulated experiments are provided to confirm the effectiveness and the robustness of the proposed ISAR imaging algorithm for rapidly spinning targets. The simulation settings are given in Table 1.

5.1. Imaging Results for a Single Scatterer Target

Results for the case of a single scatterer are shown in Figure 4. The rotation radius of the target was 3   m , and the initial phase was 0   r a d . The distance between the rotating center and the antenna phase center was set to 20   k m .
Figure 5a shows the scatterer distribution in Cartesian coordinates. Figure 5b provides the 2-D echo signal data. The profile after range compressed for echo signal is shown in Figure 5c, where a trigonometric function may be clearly recognized. The envelop of signal in the RD domain is shown in Figure 5d. Obviously, the energy of a single spinning target in the RD domain spanned multiple range cells. We now separated the energies in the RD domain into ψ [ 0 , π ] and ψ [ π , 0 ] parts and performed an interpolation-free mapping operation according to (12) for each part. The results are depicted in Figure 5e and Figure 5f, respectively. It can be seen that the main energy of a single target was concentrated into one range cell, which shows that the RCM was well corrected using the proposed algorithm. The ISAR image after high order phase compensation and combining the ψ [ 0 , π ] and ψ [ π , 0 ] parts of energy are presented in Figure 5g. The imaging result after clean post-processing is shown in Figure 5h. Finally, the interpolation-free mapping is shown in Figure 5i.

5.2. Imaging Results for Surface Targets

The simulation results for a set of targets on a surface are depicted in Figure 6. The surface was composed of eight scatterers, whose rotation radii were 2 m , 2 m , 2 m , 2 m , 3 m , 3 m , 3 m , 3 m , respectively, and the initial phases were 0 ,   π / 2 , π / 2 , π ,   π / 4 , π / 4 ,   3 π / 4 , 3 π / 4 , respectively. The reflection coefficients were set to one, and the distance between the antenna phase center and the rotating center was set to 20   k m .
Figure 6a illustrates the scatterer distribution under the Cartesian coordinates. The profile after compressing range for echo signal is shown in Figure 6b, which resembles a set of trigonometric functions with all the trajectories for the scattering points intersecting with each other. In the GRT algorithm, the energy of the target is accumulated along the profile, as shown in Figure 6b. Thus, searching different initial phases and rotating radii was needed here. The ISAR image after high order phase compensation is shown in Figure 6c, and after clean processing in Figure 6e. Meanwhile, to demonstrate the effectiveness of the proposed algorithm, the comparisons with the GRT imaging result and GRT- clean imaging result are presented in Figure 6d and Figure 6f, respectively.
As it is apparent from the figures, both of the images obtained after clean processing may be used to reconstruct the structure of the surface target. However, the proposed method is computationally efficient since the 2-D search is not needed here.

5.3. Imaging Results under Different SNRs Condition

To further show the robustness of different methods, Gaussian white noise was added. Figure 7 shows the imaging results obtained by using the GRT method of [19], the SPKT method of [22], and our proposed method for two values of the SNR: 0 dB and –35 dB, respectively. All the GRT, SPKT, and the proposed method have coherent integration in the range direction.
However, thanks to the coherent integration in the azimuth direction, our method and the SPKT showed a higher SNR gain than the GRT method due to the incoherent integration in the azimuth direction. The imaging results under different SNRs were consistent with the output SNR gain analysis, which verifies the correctness and effectiveness of the proposed method. However, compared with the SPKT method [22], our proposed method has great advantages in computational cost, since it avoids the initial phase estimation and interpolation operation used for each scatterer of the spinning target.

5.4. Computational Complexity

The running times of GRT, SPKT, and the proposed algorithm are plotted versus imaging scene sizes. The results are illustrated in Figure 8, where MATLAB running on an Intel quad-core processor, CPU clocked frequency at 3.2 GHz, memory 8 GB, and Windows 10 were used. As it is apparent from Figure 8, our method was faster than GRT and SPKT methods. Moreover, as the size increased, the difference became more striking.

6. Discussion

In this section, we compare our method with the existing GRT [19] and SPKT [22] methods in terms of their computational complexity and denote   by   N r and N a the number of samples in the range and the azimuth dimension, respectively. In GRT, the computational complexity of range compression consisting of N a -times N r -point FT and complex multiplications is O ( N a N r log 2 ( N r ) + N a N r ) . Let Q and M denote, respectively, the search numbers for rotating radius and initial phase. The computational complexity of energy accumulation is O ( M Q N a ) . Therefore, the overall computation complexity of GRT is
C GRT = O ( N a   N r log 2 ( N r ) + N a   N r + M Q N a )
The SPKT method is based on the segmental pseudo Keystone transform, which is designed to realize migration through resolution cell correction. The whole data is divided into P segments ( P 1 , and N a / P is integer) and the number of scatterers is H . By ignoring the low computational cost of a few steps, the computational cost of the SPKT method is approximately given by
C SPKT = O (   N r H ( N a 2 P 2 N a P log 2 P + ( 1 + 1 P ) N a log 2 N a ) )
In our proposed algorithm, implementation mainly includes range compression,   N r -times azimuth FFT, N a N r point interpolation operation, complex multiplication, and azimuth IFFT. Among them, the range compression step consists of a range FFT with complexity O ( N a N r log 2 ( N r ) ) , a complex multiplication with complexity O ( N a N r ) , and a range IFFT with complexity O ( N a N r log 2 ( N r ) ) , then followed by N r -times azimuth FFT with a computational complexity of O ( N r N a log 2 ( N a ) ) . The computational complexity of the interpolation operation is 2 ( 2 N k e r 1 ) N a N r , where N k e r is the length of the interpolation kernel. This step could be time consuming, and therefore, the interpolation-free mapping method is considered here. The interpolation-free mapping method involves zero-padding, FFT, phase multiplication, and IFFT. Since only a few zeros are added to the data, the total computational cost of this step is O ( 2 N a N r log 2 ( N r ) + N a N r ) . After this, the high order phase is compensated through a complex multiplication with complexity O ( N a N r ) . Finally, N r -times azimuth IFFT with complexity O ( N r N a log 2 ( N a ) ) is conducted to obtain the ISAR image. Therefore, the overall computation complexity of the proposed method is
C proposed = O ( 4 N a N r log 2 ( N r ) + 2 N r N a log 2 ( N a ) + 3 N a N r )
The above analysis confirms that, in general, our proposed algorithm is convenient in terms of computational complexity. In addition, compared to GRT and SPKT methods, the proposed algorithm does not need to estimate the initial phase.

7. Conclusions

Rapidly spinning space debris generates echo signals with great RCM during a short coherent processing interval. Upon exploiting the inherent azimuth invariance of SAR, we suggested and discussed in detail an efficient ISAR imaging algorithm. At first, a precise analytic expression in the RD domain was developed based on the POSP principle. Then, a novel interpolation kernel function was designed to remove 2-D space dependent phase errors. At the same time, an efficient interpolation-free mapping algorithm was exploited to reduce the computational complexity, and the energy spanning multiple range cells were transformed into one range cell. Finally, a well-focused ISAR image was produced by compensating the azimuth high-order modulating terms. Simulation results prove that the proposed ISAR imaging algorithm is robust, effective, and time efficient compared to GRT and SPKT methods.

Author Contributions

X.L., L.G., and D.L. proposed the method, conceived and designed the experiments, and analyzed the data; H.L. and M.Q. performed the experiments and wrote the paper; and D.L. revised the paper. All authors have read and agreed to the published version of the manuscript.

Funding

This work was funded by the National Key Laboratory Foundation of China (6142411193202), the National Natural Science Foundation of China (61971075, 61627901, 61431010), the Foundation for Innovative Research Groups of the National Natural Science Foundation of China (61621005), the National Stability Support Foundation of China (2018SSFNKLSMT-10), the Chongqing Research Program of Basic Research and Frontier Technology (cstc2018jcyjAX0351), the Key Project of Application and Development of Chongqing (cstc2019jscx-fxyd0354), the Pre-Research Fund Project (61404130114, 61404130219), and the Fundamental Research Funds for the Central Universities (2019CDQYTX012).

Acknowledgments

The authors are very grateful to the editor and reviewers for their constructive comments that have an important role in further improving this work.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Pate-Cornell, E.; Sachon, M. Risks of particle hits during space walks in low Earth orbit. IEEE Trans. Aerosp. Electron. Syst. 2001, 37, 134–146. [Google Scholar] [CrossRef]
  2. Luo, Y.; Zhang, Q.; Yuan, N. Three-dimensional precession feature extraction of space targets. IEEE Trans. Aerosp. Electron. Syst. 2014, 50, 1313–1329. [Google Scholar] [CrossRef]
  3. Garrett, H.B.; Pike, C.P. Collision frequency of artificial satellites: Creation of a debris belt. J. Geophys. Res. Space Phys. 1978, 83, 2637–2646. [Google Scholar]
  4. Li, D.; Liu, H.Q.; Zhan, M.Y.; Zhang, X.Z. ISAR Imaging of Nonuniformly Rotating Target Based on the Multi-Component CPS Model under Low SNR Environment. IEEE Trans. Aerosp. Electron. Syst. 2017, 53, 1119–1135. [Google Scholar] [CrossRef]
  5. Bai, X.R.; Xing, M.; Zhou, F.; Bao, Z. High-resolution three-dimensional imaging of spinning space debris. IEEE Trans. Geosci. Remote Sens. 2009, 47, 2352–2362. [Google Scholar]
  6. Berizzi, F.; Mese, E.D.; Diani, M.; Martorella, M. High-resolution ISAR imaging of maneuvering targets by means of the range instantaneous Doppler technique: Modeling and performance analysis. IEEE Trans. Image Process. 2001, 10, 1880–1890. [Google Scholar] [CrossRef]
  7. Walker, J.L. Range-Doppler imaging of rotating objects. IEEE Trans. Aerosp. Electron. Syst. 1980, AES-16, 23–52. [Google Scholar] [CrossRef]
  8. Lakkis, I.; Wang, Y.; el Assad, S.; Saillard, J. Range-Doppler imaging of rotating objects using state space model. In Proceedings of the Progress in Electromagnetics Research Symposium, Noordwijk, The Netherlands, 22 July 1994. [Google Scholar]
  9. Wang, Y.; Jiang, Y.C. ISAR imaging of a ship target using product high-order matched-phase transform. IEEE Geosci. Remote Sens. Lett. 2009, 6, 658–661. [Google Scholar] [CrossRef]
  10. Bao, Z.; Wang, G.; Luo, L. Range-instantaneous Doppler imaging in ISAR. Acta Electron. Sin. 1998, 24, 34–36. [Google Scholar]
  11. Sato, T. Shape estimation of space debris using single-range Doppler interferometry. IEEE Trans. Geosci. Remote Sens. 1999, 37, 1000–1005. [Google Scholar] [CrossRef] [Green Version]
  12. Li, J.; Qiu, C.W.; Zhang, L.; Xing, M.D. Time-frequency imaging algorithm for high speed spinning targets in two dimensions. IET Radar Sonar Navig. 2010, 4, 806–817. [Google Scholar] [CrossRef]
  13. Wang, Q.; Xing, M.D.; Lu, G.Y.; Bao, Z. Single Range Matching Filtering for Space Debris Radar Imaging. IEEE. Geosci. Remote Sens. Lett. 2007, 4, 576–580. [Google Scholar] [CrossRef]
  14. Wang, Q.; Xing, M.D.; Lu, G.Y.; Bao, Z. SRMF-CLEAN imaging algorithm for space debris. IEEE Trans. Antennas Propag. 2007, 55, 3524–3533. [Google Scholar] [CrossRef]
  15. Wang, H.X.; Quan, Y.H.; Xing, M.D.; Zhang, S.H. Single-range image fusion for spinning space debris radar imaging. IEEE Geosci. Remote Sens. Lett. 2010, 7, 626–630. [Google Scholar] [CrossRef]
  16. Yuan, Y.; Sun, J.; Mao, S. PFA algorithm for airborne spotlight SAR imaging with non-ideal motions. Proc. Inst. Elect. Eng.—Radar Sonar Navig. 2002, 149, 174–182. [Google Scholar] [CrossRef]
  17. Bai, X.R.; Sun, G.C. Narrow-band radar imaging of spinning targets. Sci. China Inf. Sci. 2011, 54, 873–883. [Google Scholar] [CrossRef]
  18. Ding, X.F. Narrowband imaging for spatial rotating targets based on tomography algorithm. Signal Process. 2010, 26, 648–653. [Google Scholar]
  19. Hansen, K.V.; Toft, P.A. Fast curve estimation using preconditioned generalized radon transform. IEEE Trans. Image Process. 1996, 5, 1651–1661. [Google Scholar] [CrossRef] [Green Version]
  20. Zhang, Q.; Yeo, T.S.; Tan, H.; Luo, Y. Imaging of a moving target with rotating parts based on the Hough transform. IEEE Trans. Geosci. Remote Sens. 2008, 46, 291–299. [Google Scholar] [CrossRef]
  21. Gong, J.; Teng, L.; Li, Z.L.; Li, C. Micro-motion ISAR imaging method based on MP sparse decomposition for fast moving target. Adv. Mater. Res. 2013, 760, 1482–1485. [Google Scholar] [CrossRef]
  22. Huo, K.; Liu, Y.X.; Hu, J.M.; Jiang, W.D.; Li, X. A novel imaging method for fast rotating targets based on the segmental pseudo Keystone transform. IEEE Trans. Geosci. Remote Sens. 2011, 49, 1464–1472. [Google Scholar] [CrossRef]
  23. Tang, S.Y.; Zhang, L.R.; Guo, P.; Liu, G.; Sun, G.-C. Acceleration model analyses and imaging algorithm for highly squinted airborne spotlight-mode SAR with maneuvers. IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens. 2015, 8, 1120–1131. [Google Scholar] [CrossRef]
  24. Li, D.; Lin, H.; Liu, H.Q.; Liao, G.; Tan, X. Focus improvement for high-resolution highly squinted SAR imaging based on 2-D spatial-variant linear and quadratic RCMs correction and azimuth-dependent Doppler equalization. IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens. 2017, 10, 168–183. [Google Scholar] [CrossRef]
  25. Huang, P.H.; Liao, G.S.; Yang, Z.W.; Xia, X.-G.; Ma, J.-T.; Zhang, X.P. A fast SAR imaging method for ground moving target using a second-order WVD transform. IEEE Trans. Geosci. Remote Sens. 2016, 54, 1940–1956. [Google Scholar] [CrossRef]
  26. Wong, F.H.; Cumming, I.G.; Neo, Y.L. Focusing bistatic SAR data using the nonlinear chirp scaling algorithm. IEEE Trans. Geosci. Remote Sens. 2008, 46, 2493–2505. [Google Scholar] [CrossRef] [Green Version]
  27. Li, D.; Zhan, M.Y.; Liu, H.Q.; Liao, G.S. A robust translational motion compensation method for ISAR imaging based on keystone transform and fractional Fourier transform under low SNR environment. IEEE Trans. Aerosp. Electron. Syst. 2017, 53, 2140–2156. [Google Scholar] [CrossRef]
  28. Liu, L.; Zhou, F.; Tao, M.; Sun, P.; Zhang, Z. Adaptive translational motion compensation method for ISAR imaging under low SNR based on particle swarm optimization. IEEE J. Sel. Topics Appl. Earth Obs. Remote Sens. 2015, 8, 5146–5157. [Google Scholar] [CrossRef]
  29. Ruan, H.; Wu, Y.; Jia, X.; Ye, W. Novel ISAR imaging algorithm for maneuvering targets based on a modified keystone transform. IEEE Trans. Geosci. Remote Sens. Lett. 2014, 11, 128–132. [Google Scholar] [CrossRef]
  30. Davidson, G.W.; Cumming, I.G. Signal properties of spaceborne squint-mode SAR. IEEE Trans. Geosci. Remote Sens. 1997, 35, 611–617. [Google Scholar] [CrossRef] [Green Version]
  31. Jehanzeb, B.; Christopher, F.B. Interpolation-free algorithm for SAR 2D aperture synthesis. In Proceedings of the SPIE 6237, Algorithms for Synthetic Aperture Radar Imagery XIII, Orlando (Kissimmee), FL, USA, 10 May 2006. [Google Scholar]
  32. Li, Z.; Wang, J.; Wu, J.J.; Liu, H.Q. A fast radial scanned near-field 3-D SAR imaging system and the reconstruction method. IEEE Trans. Geosci. Remote Sens. 2015, 53, 1355–1363. [Google Scholar] [CrossRef]
  33. Li, Z.; Wang, J.; Liu, H.Q. Interpolation-Free Stolt Mapping for SAR Imaging. IEEE Trans. Geosci. Remote Sens. 2014, 11, 926–929. [Google Scholar] [CrossRef]
  34. Tsao, J.; Steinberg, B.D. Reduction of sidelobe and speckle artifacts in microwave imaging: The CLEAN technique. IEEE Trans. Antennas Propag. 1988, 36, 543–556. [Google Scholar] [CrossRef]
Figure 1. Geometry of ISAR imaging for spinning space debris.
Figure 1. Geometry of ISAR imaging for spinning space debris.
Remotesensing 12 02059 g001
Figure 2. Illustration of the azimuth invariance of strip-map SAR.
Figure 2. Illustration of the azimuth invariance of strip-map SAR.
Remotesensing 12 02059 g002
Figure 3. Data format comparison. (a) Original ( r f a ) plane. (b) KT method ( r f a ) plane. (c) Proposed method ( r f a ) plane.
Figure 3. Data format comparison. (a) Original ( r f a ) plane. (b) KT method ( r f a ) plane. (c) Proposed method ( r f a ) plane.
Remotesensing 12 02059 g003
Figure 4. The flowchart of the proposed ISAR imaging method (left) and details of its fast implementation (right).
Figure 4. The flowchart of the proposed ISAR imaging method (left) and details of its fast implementation (right).
Remotesensing 12 02059 g004
Figure 5. Results of numerically simulated reconstruction of a single target by our proposed method. (a) 2-D distribution of scatterers. (b) 2-D echo signal. (c) The profile of echo signal after range compression. (d) The signal envelope in the RD domain. (e) Interpolation for ψ [ 0 , π ] . (f) Interpolation for ψ [ π , 0 ] . (g) 2-D imaging result. (h) Imaging result after clean processing. (i) Imaging result with interpolation-free mapping after clean processing.
Figure 5. Results of numerically simulated reconstruction of a single target by our proposed method. (a) 2-D distribution of scatterers. (b) 2-D echo signal. (c) The profile of echo signal after range compression. (d) The signal envelope in the RD domain. (e) Interpolation for ψ [ 0 , π ] . (f) Interpolation for ψ [ π , 0 ] . (g) 2-D imaging result. (h) Imaging result after clean processing. (i) Imaging result with interpolation-free mapping after clean processing.
Remotesensing 12 02059 g005
Figure 6. Results of numerically simulated reconstruction of rapidly spinning scatterers. (a) 2-D distribution of scatterers. (b) The profiles after range compression. (c) Result of the proposed method without clean. (d) Result after GRT. (e) Result of the proposed method with clean. (f) Result of GRT- clean.
Figure 6. Results of numerically simulated reconstruction of rapidly spinning scatterers. (a) 2-D distribution of scatterers. (b) The profiles after range compression. (c) Result of the proposed method without clean. (d) Result after GRT. (e) Result of the proposed method with clean. (f) Result of GRT- clean.
Remotesensing 12 02059 g006
Figure 7. Results of numerically simulated reconstruction with different SNRs. (a) GRT with SNR= 0 dB. (b) SPKT with SNR = 0 dB. (c) The proposed method with SNR= 0 dB. (d) GRT with SNR= –35 dB. (e) SPKT with SNR = –35 dB. (f) The proposed method with SNR = –35 dB.
Figure 7. Results of numerically simulated reconstruction with different SNRs. (a) GRT with SNR= 0 dB. (b) SPKT with SNR = 0 dB. (c) The proposed method with SNR= 0 dB. (d) GRT with SNR= –35 dB. (e) SPKT with SNR = –35 dB. (f) The proposed method with SNR = –35 dB.
Remotesensing 12 02059 g007
Figure 8. Comparison of computational complexity.
Figure 8. Comparison of computational complexity.
Remotesensing 12 02059 g008
Table 1. Parameters for simulated ISAR Imaging.
Table 1. Parameters for simulated ISAR Imaging.
Parameter Value
Carrier frequency8 GHz
Sample frequency6 GHz
Transmit bandwidth4 GHz
Pulse width1 us
Pulse repetition frequency2400 Hz
Rotating angle velocity6.28 rad/s

Share and Cite

MDPI and ACS Style

Luo, X.; Guo, L.; Li, D.; Liu, H.; Qin, M. A Novel 2-D Geometry Reconstruction Approach for Space Debris via Interpolation-Free Operation under Low SNR Conditions. Remote Sens. 2020, 12, 2059. https://0-doi-org.brum.beds.ac.uk/10.3390/rs12122059

AMA Style

Luo X, Guo L, Li D, Liu H, Qin M. A Novel 2-D Geometry Reconstruction Approach for Space Debris via Interpolation-Free Operation under Low SNR Conditions. Remote Sensing. 2020; 12(12):2059. https://0-doi-org.brum.beds.ac.uk/10.3390/rs12122059

Chicago/Turabian Style

Luo, Xi, Lixin Guo, Dong Li, Hongqing Liu, and Mengyi Qin. 2020. "A Novel 2-D Geometry Reconstruction Approach for Space Debris via Interpolation-Free Operation under Low SNR Conditions" Remote Sensing 12, no. 12: 2059. https://0-doi-org.brum.beds.ac.uk/10.3390/rs12122059

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop