Micro-Doppler Effect Analysis and Feature Extraction in

ISAR Imaging With Stepped-Frequency Chirp Signals

The MIT Faculty has made this article openly available. Please share
how this access benefits you. Your story matters.

Citation Ying Luo et al. “Micro-Doppler Effect Analysis and Feature

Extraction in ISAR Imaging With Stepped-Frequency Chirp
Signals.” Geoscience and Remote Sensing, IEEE Transactions
on 48.4 (2010): 2087-2098. © 2010 IEEE.
As Published http://dx.doi.org/10.1109/tgrs.2009.2034367
Publisher Institute of Electrical and Electronics Engineers

Version Final published version

Accessed Wed Jan 11 01:17:57 EST 2017
Citable Link http://hdl.handle.net/1721.1/61674
Terms of Use Article is made available in accordance with the publisher's policy
and may be subject to US copyright law. Please refer to the
publisher's site for terms of use.
Detailed Terms
IEEE TRANSACTIONS ON GEOSCIENCE AND REMOTE SENSING, VOL. 48, NO. 4, APRIL 2010 2087

Micro-Doppler Effect Analysis and Feature

Extraction in ISAR Imaging With
Stepped-Frequency Chirp Signals
Ying Luo, Qun Zhang, Senior Member, IEEE, Cheng-wei Qiu, Member, IEEE, Xian-jiao Liang, and Kai-ming Li

Abstract—The micro-Doppler (m-D) effect induced by the ro- helicopter with rotating rotors, a ship with scanning antennas,
tating parts or vibrations of the target provides a new approach etc.), have attracted great research attention in accurate target
for target recognition. To obtain high range resolution for the ex- identification. High-resolution time–frequency techniques have
traction of the fine m-D signatures of an inverse synthetic aperture
radar target, the stepped-frequency chirp signal (SFCS) is used been developed to extract these time-varying m-D signatures
to synthesize the ultrabroad bandwidth and reduce the require- [3], [5], and the adaptive chirplet presentation for m-D fea-
ment of sample rates. In this paper, the m-D effect in SFCS is tures from the target with rotating parts is studied in [6]. The
analyzed. The analytical expressions of the m-D signatures, which empirical-mode decomposition is utilized to estimate accurate
are extracted by an improved Hough transform method associated parameters of the rotating parts and focus the image of the
with time–frequency analysis, are deduced on the range-slow-time
plane. The implementation of the algorithm is presented, partic- main body in [7]. The image processing algorithms have also
ularly in those extreme cases of rotating (vibrating) frequencies been introduced for the separation of the m-D features, such
and radii. The simulations validate the theoretical formulation and as the Radon transform [8] and the Hough transform (HT) [9],
robustness of the proposed m-D extraction method. [10]. Furthermore, some literature has reported the application
Index Terms—Hough transform (HT), inverse synthetic aper- of the m-D effect in automatic target recognition, e.g., the
ture radar (ISAR), micro-Doppler (m-D) effect, stepped-frequency geometrical feature extraction of wheels [11], the micromotion
chirp signal (SFCS), time–frequency analysis. feature extraction of helicopters [5], the urban sensing and
indoor sensing [12], the gait and activity analyses of pedestrians
I. I NTRODUCTION
[1]–[3], [5], [6], [11], [13], etc.

T HE MICROMOTION dynamics of an object or structures

on the object may induce frequency modulations on the re-
turned signal, which generate sidebands of the target’s Doppler
The m-D signatures describe the fine micromotion features of
a radar target. To more effectively exploit the m-D signatures of
an ISAR target, the broad-bandwidth transmitted radar signal is
frequency shift, which is called the micro-Doppler (m-D) effect necessary to obtain high range resolution. The chirp train with
[1]–[4]. Micro-Doppler features can be regarded as a unique stepped carrier frequencies, i.e., the stepped-frequency chirp
signature of an object with movements, providing additional in- signal (SFCS), is one kind of effective signals for obtaining
formation for the classification, recognition, and identification high-resolution downrange profiles of targets. The main advan-
of the object. On the other hand, the existence of m-D effect tage of this signal is that the actual instantaneous bandwidth
may contaminate the inverse synthetic aperture radar (ISAR) of the radar is quite small, compared to the total processing
image of the main body and bring great difficulties to the image bandwidth, which allows the transmission of waveforms with
interpretation [1]–[9]. Therefore, it is necessary to separate and extremely wide overall bandwidth without expensive hardware
extract the m-D signals from the returns of the main body. to support the wide instantaneous bandwidth [14]. Thus, this
Since Chen introduced it into the radar imaging in 2000 [1], kind of signals can be utilized to introduce imaging capability
the m-D effects, which widely exist in the radar targets (e.g., a to existing narrow-bandwidth radar [15]. In contrast to the
stepped-frequency signal, the SFCS possesses longer detection
Manuscript received March 27, 2009; revised June 4, 2009. First published distance and lower grating lobes in the range response by means
December 4, 2009; current version published March 24, 2010. This work of replacing the fix-frequency subpulses with chirp subpulses
was supported in part by the National Natural Science Foundation of China
under Grant 60971100 and Grant 60672032 and in part by the Natural Science
[16]. However, due to the frequency-stepped processing, it is
Foundation Research Program of Shaanxi Province under Grant 2007F28. sensitive to the radial velocity of the target. When a target con-
Y. Luo, X. Liang, and K. Li are with the Institute of Telecommunication tains micromotion structures, the movements of these structures
Engineering, Air Force Engineering University, Xi’an 710077, China.
Q. Zhang is with the Institute of Telecommunication Engineering, Air Force will modulate the high-resolution range profile (HRRP) and
Engineering University, Xi’an 710077, China, and also with the Key Labora- exhibit different m-D characteristics from the traditional chirp
tory of Wave Scattering and Remote Sensing Information (Ministry of Edu- signal. Most of the present studies focus on the m-D effects
cation), Fudan University, Shanghai 200433, China (e-mail: zhangqunnus@
gmail.com). with single-frequency signals and chirp signals [1], [6]–[9].
C. Qiu is with the Research Laboratory of Electronics, Massachusetts Although the m-D effects with SFCS are briefly discussed in
Institute of Technology, Cambridge, MA 02139-4307 USA. our previous work [17], detailed analyses of the m-D effects
Color versions of one or more of the figures in this paper are available online
at http://ieeexplore.ieee.org. and the algorithm for m-D signature extraction are not studied,
Digital Object Identifier 10.1109/TGRS.2009.2034367 which motivates this paper.

0196-2892/$26.00 © 2009 IEEE

2088 IEEE TRANSACTIONS ON GEOSCIENCE AND REMOTE SENSING, VOL. 48, NO. 4, APRIL 2010

In this paper, we first analyze the m-D signatures and in-

vestigate the range profile’s migration and wrapping in SFCS
in Section II. Particularly, we discuss the m-D effects induced
by rotation and vibration in detail and deduce the analytical
expressions of the m-D signatures on the range-slow-time plane
since the rotating and vibration dynamics are the most common
micromotions in radar targets. Then, we propose an algorithm
based on an improved HT and time–frequency analysis for
the m-D signature extraction in Section III. Due to the range Fig. 1. Frequency variety in a burst of SFCS.
profile’s migration and wrapping, the m-D signal of a rotating
or vibrating scatterer on the target no longer satisfies the respectively; f0 + iΔf is the carrier frequency; Tr and N
sinusoidal function on the range-slow-time plane. As a result, denote the subpulse repetition interval and number of subpulses
the HT method for sinusoid detection and estimation in [9] in a burst, respectively; and θi is the initial phase. B1 = μT1
is not suitable for the m-D signature extraction in SFCS. In stands for the bandwidth of the subpulse, which should be equal
this section, we construct new HT equations and associate with to or a bit larger than the stepped frequency value Δf to avoid
the time–frequency analysis to extract the m-D signatures of the interspaces in the synthetic spectrum.
the rotating or vibrating parts with different moving periods Assuming a target of a point scatterer with unit scattering
and radii. Different from [9], which focuses on the imaging coefficient, the echo of the ith subpulse from the target is
for targets with large-radius and low-frequency rotating parts  
such as helicopters, this paper mainly deals with the m-D t−iTr −2R/c  
s(t, i) = rect · exp jπμ(t−iTr −2R/c)2
signature extraction for the rotating or vibrating parts with T1
extreme frequencies (very high or very low) and radii that · exp (j2π(f0 + iΔf )(t − iTr − 2R/c) + jθi ) (2)
are larger/smaller than the range resolution. Simulations are
presented in Section IV, and more discussions are presented in where R is the distance between the target and the radar, and c is
Section V. the wave propagation velocity. The “deramping” is performed
In addition, it should be mentioned that the “deramping” by multiplying the echo signal with the complex conjugate of
method [18] is utilized to synthesize the HRRP for simplicity in the reference signal [18], and the reference point can be chosen
this paper. In fact, the analyses of m-D effect and the signature among the scatterers with the strongest scattering coefficients
extraction method can easily be extended to imaging systems on target. The reference signal can be written as
utilizing other HRRP synthesis methods, such as the time-  
domain synthesis [19] and the frequency-domain synthesis [20]. t−iTr −2R0 /c  
s0 (t, i) = rect · exp jπμ(t−iTr −2R0 /c)2
Furthermore, considering the small size of general ISAR target Tref
and the narrow bandwidth of the chirp subpulses, we assume
· exp (j2π(f0 + iΔf )(t − iTr − 2R0 /c) + jθi ) (3)
that the length of target in line of sight (LOS) is shorter than the
coarse range resolution of SFCS. Thus, complicated extraction
where R0 is the distance between the reference point and the
processing of HRRP from redundant data can be avoided, and
radar, Tref is the duration of the reference signal, which is a bit
the analyses of m-D effect are simplified.
larger than T1 . After the “deramping” processing, it yields
II. m-D E FFECT A NALYSIS
sc (t, i) = s(t, i)s∗0 (t, i)
In this section, we first introduce the HRRP synthesis tech-  
t − iTr − 2R/c
nique with SFCS and then analyze the m-D effects induced = rect
T1
by the micromotions, particularly the rotating motion and the    
vibrating motion. 4πμ 2R0
· exp −j t − iTr − RΔ
c c
A. HRRP Synthesis With SFCS    
4π 4πμ 2
The SFCS transmitted by the radar is made up of a series of · exp −j (f0 + iΔf )RΔ · exp j 2 RΔ
c c
bursts, each of which consists of a sequence of chirp subpulses
of stepped carrier frequency. Fig. 1 shows the variety of the (4)
signal frequency in a burst.
where RΔ = R − R0 . Replacing (t − iTr − 2R0 /c) by t , ap-
Suppose that the ith subpulse of a burst transmitted by the
plying the Fourier transform to (4) with respect to t , and
radar is
removing the Residual Video Phase (RVP), one can obtain
Ui (t) = u(t − iTr ) · exp (j2π(f0 + iΔf )(t − iTr ) + jθi ) ,   
4πμ
0 ≤ i ≤ N − 1 (1) Sc (ω, i) = T1 sinc T1 ω + RΔ
c
 
where u(t) = rect(t/T1 ) · exp(jπμt2 ) is a chirp with μ and 4π
· exp −j (f0 + iΔf ) RΔ . (5)
T1 being its frequency modulation slope and time duration, c
LUO et al.: MICRO-DOPPLER EFFECT ANALYSIS AND FEATURE EXTRACTION 2089

Thus, the HRRP peaks at

4π 4π 8π
k = Φ (i) = − Δf · RΔ0 − f0 Tr v − Δf iTr v (9)
c c c
Fig. 2. Nonrigid target with micromotion. where Φ(i) = −4π(f0 + iΔf )(RΔ0 + iTr v)/c is the phase
term in (8). From (9), one can see that the first term on the right
side denotes the initial location of the point during the burst
It can be found that the peak value of |Sc (ω, i)| appears at time, the second term demonstrates that the location of Q in the
ω = −4πμRΔ /c, and in fact, it is the coarse range profile of
HRRP has an offset proportional to v, and the last term is the
the point target created by the ith subpulse exactly. N subpulses coupling term between the velocity and the stepped frequency,
will create N coarse range profiles. Assuming the length of which is negligibly small, compared with the first two terms
target in LOS to be shorter than the coarse range resolution (e.g., when f0 = 35 GHz, Tr = 1 μs, i = 1, Δf = 5 MHz,
c/(2B1 ), the range profile of the target obtained by chirp sub- v = 20 m/s, and RΔ0 = 5 m, the value of this term is only
pulse will have only one peak at ω = −4πμRΔ /c. Therefore, −8.3776 × 10−6 , whereas the values of the first and second
keeping ω = −4πμRΔ /c and taking the Fourier transform of terms are −1.0472 and −0.0293, respectively). However, this
these coarse range profiles with respect to i, one can obtain term may make the peaks of the HRRP seriously expanded,
    particularly when v is quite large.
4πΔf 4πf0
S(k) = Cn · sinc k+ RΔ · exp −j RΔ (6) Ignoring the last term, (9) is rewritten as
c c
4π 4π
k=− Δf · RΔ0 − f0 Tr v (10)
where Cn is a constant. From (6), it can be found that the peak c c
value of |S(k)| appears at k = −4πΔf · RΔ /c, i.e., the HRRP and (10) denotes the peak location of Q in the HRRP. The
is obtained. Divide k by −4πΔf /c, and then, the location of location of Q in range can be obtained by dividing (10) by
the scatterer in the range can be scaled. Each burst creates −4πΔf /c, i.e.,
an HRRP. By transmitting a number of bursts and taking the
Fourier transform to the HRRPs with respect to slow time, the r = RΔ0 + f0 Tr v/Δf. (11)
ISAR image of the target can be achieved [17].
Because the HRRP is obtained by N -point Fourier transform,
the value range of k is within [−π, π]. If v is too large, k in (10)
B. Micromotion Signatures in SFCS may exceed this scope, and the peak location of the HRRP may
wrap and appear on the other side of the reference point. That
In many cases, the radar target cannot be regarded as a rigid is to say, to avoid the wrapping of the HRRP, k must satisfy the
object, e.g., a helicopter with rotating rotors and a bird with following condition:
wings beating. When a target contains micromotion parts, the
existence of micromotions will affect the synthesis of the HRRP |k| < π (12)
and bring in m-D effect.
As shown in Fig. 2, point P is the reference point, and point i.e.,
Q is a micromotion scatterer with the velocity v relative to P c c
along the LOS. Since the duration of a burst is very short, v can − <r< . (13)
4Δf 4Δf
be considered as a constant in burst time approximately. The
distance between P and Q at the initial time of a burst is RΔ0 . Equation (13) demonstrates that the peak location of the mi-
Due to the movement of Q, in a burst time, RΔ in (4) should be cromotion point is limited to the unambiguous range window.
rewritten as Thus, v must also satisfy
1 c  1 c 
RΔ = RΔ0 + iTr v. (7) − +Δf · RΔ0 < v < −Δf · RΔ0 . (14)
f0 Tr 4 f0 Tr 4
When the movement of the micromotion point does not satisfy
Generally, the displacement of Q in a burst cannot exceed the
the conditions of (12)–(14), due to the wrapping of the HRRP,
coarse range resolution. Therefore, (5) can be rewritten as
the peak location of the micromotion point in the HRRP will be
   determined by
4πμ
Sc (ω, i) = T1 sinc T1 ω + (RΔ0 + iTr v)  
c c c c
  r̂ = mod r + , − (15)
4π 4Δf 2Δf 4Δf
· exp −j (f0 + iΔf )(RΔ0 + iTr v)
c where mod(a, b) returns the remainder of a divided by b.
  
4πμ
≈ T1 sinc T1 ω + RΔ0
c C. Micro-Doppler Signatures of Rotation and Vibration
 
4π In order to more profoundly comprehend the m-D effect in
· exp −j (f0 + iΔf )(RΔ0 + iTr v) . (8)
c SFCS, we take further analyses of the m-D signatures induced
2090 IEEE TRANSACTIONS ON GEOSCIENCE AND REMOTE SENSING, VOL. 48, NO. 4, APRIL 2010

be a sinusoidal curve on the range-slow-time plane. The

frequency of the sinusoid is equal to the rotating or
vibrating rate of the scatterer. Divide the amplitude of the
sinusoid by 1 + ω02 f02 Tr2 /Δf 2 , and then, the radius of
Fig. 3. Target with rotating or vibration scatterers.
the rotation or vibration ρ is obtained.
2) Since r(tm ) in (19) cannot reach its maximum when
by rotating and vibrating parts in a radar target. The rotation and R(tm ) in (16) reaches its maximum, owing to the phase
vibration are the most common kinds of micromotions in the shift φ, it means that, at the moment of the sinusoid
radar targets, e.g., a helicopter with rotating rotors, a ship with maximum, the rotating/vibrating scatterer does not reach
scanning antennas, and a truck with vibrating surface induced the maximum distance from the reference point.
by a running engine. 3) When the micromotion parameters of the rotating/
As shown in Fig. 3, when the micromotion scatterer Q is vibrating scatterer fail to satisfy the condition in (12), the
taking the rotation or vibration motion around C with the m-D signature of the rotating/vibrating scatterer appears
rate ω0 = 2πΩ and radius ρ, the distance between Q and the to be a wrapped sinusoidal curve on the range-slow-time
reference point P in LOS is expressed as plane. The frequency of the curve is equal to the rotating
or vibrating rate of the scatterer.
R(tm ) = RQ + ρ cos(ω0 tm + θ0 ) (16)
These characteristics of the m-D effect are different from that
where RQ indicates the distance between the rotating center C in chirp signals. In chirp signals, the sinusoid on the range-
and the reference point, θ0 is the initial phase, and tm is the slow-time plane would never be wrapped, and the amplitude
slow time, i.e., the transmitting time of each burst. The velocity of the sinusoid is equal to the radius of the rotating or vibrating
of the point relative to the radar is scatterer [9]. In addition, the sinusoid crest also means that the
rotating or vibrating scatterer is reaching the maximum distance
v(tm ) = R (tm ) = −ρω0 sin(ω0 tm + θ0 ). (17) from the reference point at the moment. These differences be-
tween the SFCS and chirp signal result from the range profile’s
According to (10), the HRRP of Q peaks at migration in (10) and the wrapping of the HRRP in (15).
4π 4π
k(tm ) = − Δf · R(tm ) − f0 Tr v(tm )
c c
  III. M ICRO -D OPPLER F EATURE E XTRACTION
−4πΔf
= RQ + ρ 1 + ω02 f02 Tr2 /Δf 2
c Since the curves on the range-slow-time plane contain the
 micromotion information of the target, we can extract the m-D
× cos(ω0 tm + θ0 + φ) (18) signatures by detecting the parameters of these curves. The
HT has been recognized as one of the most popular methods
for the detection of line segments. It was first introduced by
where φ = arccos(1/ 1 + ω02 f02 Tr2 /Δf 2 ). When k(tm ) sat-
Hough [21] and was later improved to detect analytical curves
isfies the condition in (12), the location of the point in range is
[22], [23]. The major advantage of the HT is its robustness
r(tm ) = R(tm ) + f0 Tr v(tm )/Δf to noise and discontinuities in the pattern to be detected. In
 general, the HT can be used to extract the curves’ feature with
= RQ + ρ 1 + ω02 f02 Tr2 /Δf 2 cos(ω0 tm + θ0 + φ). the following basic principle: the parameters (a1 , a2 , . . . , aN )
are used to represent the curve, and each parameter accumulator
(19) simply counts the number of points (xi , yi ) fitted into the curve
equation
If k(tm ) does not satisfy the condition in (12), the HRRP will
be wrapped, and according to (15), the location of the point in
range is f (a1 , a2 , . . . , aN , xi , yi ) = 0. (21)
 
r̂(tm ) = mod RQ + ρ 1 + ω02 f02 Tr2 /Δf 2 The accumulator, which contains the local maxima in the HT
space, may correspond to the presence of a curve.

c c c The m-D effects induced by different forms of micromotions
× cos(ω0 tm + θ0 + φ) + , − . (20) appear as different curves on the range-slow-time plane. If the
4Δf 2Δf 4Δf
analytical expressions of the curves are obtained, the parame-
It can easily be found that (19) is just a subset of (20) under ters of the m-D features can be extracted by HT by constructing
the condition in (12). Therefore, (20) describes the m-D effect corresponding transform equations. In the following, we take
in SFCS induced by the rotation and vibration. From (19) and the rotation and vibration as examples to discuss the m-D
(20), three conclusions can be drawn. extraction algorithm.
1) When the micromotion parameters of the rotating/ According to (20), the parameters (RQ , ρ, ω0 , θ0 ) have deter-
vibrating scatterer satisfy the condition in (12), the m-D mined the curve’s equation on the range-slow-time plane. We
signatures of the rotating/vibrating scatterer appears to construct the HT equation as follows.
LUO et al.: MICRO-DOPPLER EFFECT ANALYSIS AND FEATURE EXTRACTION 2091

First, (20) is rewritten in the equivalent form The rotating or vibrating center can be identified from the
c c parameter RQ . The radius ρ can also be estimated using the
r̂(tm ) + +M · following:
4Δf 2Δf
 
= RQ + ρ 1 + ω02 f02 Tr2 /Δf 2 cos(ω0 tm ρ = A/ 1 + ω02 f02 Tr2 /Δf 2 . (28)
c
+ θ0 + φ) + , M = 0, ±1, ±2, . . . (22)
4Δf The accuracy of RQ and A is mainly dependent on the range
resolution of the HRRP, and the accuracy of Ω and ζ is
i.e., determined by the set interval between the adjacent discrete
c elements in [Ωmin , Ωmax ] and [ζmin , ζmax ], respectively. The
RQ = r̂(tm ) + M ·
2Δf more accurate the Ω and ζ, the more the computer storage and
 computational time required.
−ρ 1 + ω02 f02 Tr2 /Δf 2 cos(ω0 tm + θ0 + φ). (23) Generally, the algorithm previously presented can extract
the m-D signatures from most of the targets with rotating or
According to the assumption that the length of the target in vibrating parts, such as the horizontal rotors on a helicopter and
LOS is shorter than the coarse range resolution c/(2B1 ), i.e., the rotating antenna on a ship. However, two special situations
−c/(4B1 ) < RQ < c/(4B1 ), and B1 ≥ Δf , it yields must be taken into consideration.
c c
− < RQ < . (24)
4Δf 4Δf
A. When A < c/(4N Δf )
Therefore, (23) can be rewritten as
  In this situation, the amplitude of the curve is less than half
RQ = mod r̂(tm ) − ρ 1 + ω02 f02 Tr2 /Δf 2 of the range resolution. This situation also widely exists in
radar target, e.g., a truck with vibrating surface induced by
 a running engine or a person with chest fluctuation induced
c c c
× cos(ω0 tm + θ0 + φ) + , − . (25) by breathing. When A < c/(4N Δf ), the curve will appear
4Δf 2Δf 4Δf
as a straight line on the range-slow-time plane, and the m-D

Let A = ρ 1 + ω02 f02 Tr2 /Δf 2 , Ω = ω0 /(2π), and ζ = θ0 + signatures cannot directly be extracted by HT. Supplementary
φ, and we have processing is needed before the implementation of the HT.
 Reconsider (8); when letting ω = −4πμRΔ0 /c, the HRRP
RQ = mod r̂(tm ) − A cos(2πΩtm + ζ) can be synthesized by applying fast Fourier transform (FFT) to
the right side of (8) with respect to i. Ignoring the influence
 of the coupling term between v and i, the result of FFT is
c c c
+ , − . (26) approximately obtained, i.e.,
4Δf 2Δf 4Δf
 
Thus, the detection of curves on the range-slow-time plane 4πΔf 4π
S(k) = Cn · sinc k + RΔ0 + f0 Tr v
is transformed to the peak value detection in the parameters c c
(RQ , A, Ω, ζ) space by utilizing the HT according to (26). The  
4π
implementation of the HT is given as follows: · exp −j f0 RΔ0 (29)
c
Step 1) Set up a discrete accumulator matrix
Ac (RQ , A, Ω, ζ), where RQ ∈ [RQ min , RQ max ],
where Cn is a constant. According to (16) and (17), it can be
A ∈ [Amin , Amax ], Ω ∈ [Ωmin , Ωmax ], and ζ ∈
obtained as
[ζmin , ζmax ] are within the prospective value scales.
Zeroize matrix Ac .  
4πΔf 4π
Step 2) Apply the HT to all the points on the range- S(k, tm ) = Cn · sinc k + R(tm ) + f0 Tr v(tm )
c c
slow-time plane whose module values are greater  
than a certain threshold, i.e., calculate the curves 4π
(R̂Q , Â, Ω̂, ζ̂) in the HT space corresponding to · exp −j f0 R(tm ) . (30)
c
the point (r̂(t̂m ), t̂m ) on the range-slow-time plane
according to (26), and then apply Because A < c/(4N Δf ), |S(k, tm )| appears as a straight line
on the range-slow-time plane; hence, we have the peak location
Ac (R̂Q , Â, Ω̂, ζ̂) = Ac (R̂Q , Â, Ω̂, ζ̂) + 1. (27)
of |S(k, tm )| at
Step 3) Search the local maxima in the HT space, and their
coordinates Ac (RQ , A, Ω, ζ) are the parameters of k = −4πΔf R(tm )/c − 4πf0 Tr v(tm )/c ≈ C n (31)
the curves on the range-slow-time plane.
The rotating or vibrating rate of the micromotion scatterer where Cn is a constant. Let k = Cn , and take the derivative to
ω0 = 2πΩ can be obtained from the parameters detected by HT. the phase term of the right side of (30) in terms of tm . Then, the
2092 IEEE TRANSACTIONS ON GEOSCIENCE AND REMOTE SENSING, VOL. 48, NO. 4, APRIL 2010

variance of the instantaneous frequency of S(k, tm )|k=Cn with

respect to tm is obtained as
1  2
l(tm ) = ψ (tm ) = f0 ρω0 sin(ω0 tm + θ0 ). (32)
2π c
It can be found that l(tm ) is a sinusoidal function. There-
fore, the m-D signatures can be extracted by analyzing
S(k, tm )|k=Cn in the joint time–frequency domain first and
then utilizing the HT to detect the curve on the time–frequency
plane. Similarly, l(tm ) will also be wrapped if |l(tm )| >
BRF/2, where the BRF is the burst repetition frequency, i.e., Fig. 4. Range profile’s migration and wrapping.
the sample frequency in slow-time domain. Therefore, the curve According to the analysis in Section II, the HRRP of the
on the time–frequency plane is depicted as micromotion scatterer is synthesized by taking FFT to (8) with

respect to i. When the value of v is quite large, the coupling
ˆl(tm ) = mod 2 f0 ρω0 sin(ω0 tm + θ0 )
c term between v and i will more significantly affect the HRRP,
 and the Fourier transform is not suitable to describe the m-D
BRF BRF
+ , BRF − . (33) signatures. From (9), the linear relationship between k and i can
2 2
be observed when assuming v to be a constant approximately in
Then, the HT equation can be deduced accordingly as a burst time. Considering that k is the instantaneous frequency
 
of Sc (ω, i) in (8) at ω = −4πμRΔ /c, Sc (ω, i)|ω=−4πμRΔ /c
mod ˆl(tm ) − A sin(2πΩtm + ζ), BRF = 0 (34)
is an approximate linear frequency modulation (LFM) signal.
where A = 2f0 ρω0 /c, Ω = ω0 /(2π), and ζ = θ0 . Thus, the Taking time–frequency analysis to Sc (ω, i)|ω=−4πμRΔ /c , it will
m-D signature extraction is equivalent to the peak value detec- appear as a segment of straight line on the time–frequency
tion in the parameters (A , Ω, ζ) space. plane. Since it is a simple LFM signal, the STFT is competent
Since S(k, tm )|k=Cn is a sine frequency modulation sig- to analyze it without cross-term interferences.
nal, considering the deficient time–frequency concentration After taking the STFT to Sc (ω, i)|ω=−4πμRΔ /c of each burst,
of short-time Fourier transform (STFT) and the cross- the corresponding time–frequency planes can be obtained, and a
term interference in quadratic time–frequency distribution new frequency-slow-time plane can be constructed by stitching
such as Wigner-Ville distribution (WVD) and Smoothing up these planes in order. On that new plane, the m-D signatures
Pseudo WVD, we choose the Gabor transform to analyze of rotating or vibrating scatterers also appear as the curves
S(k, tm )|k=Cn . Gabor transform has been recognized as one of depicted by (20), but the slow-time sampling number in a cycle
the best methods to represent pattern signals, owing to its per- has increased to N times. Therefore, the HT can then be carried
fect time resolution, frequency resolution, and the absence of out to accurately extract the m-D signatures according to (26).
cross-term interference [24]. Hence, analyzing S(k, tm )|k=Cn The aforementioned algorithm is just for the m-D signature
with Gabor transform could give many advantages to the fol- extraction of rotation and vibration. Nevertheless, if other mi-
lowing HT processing. cromotions can be depicted by certain analytical expressions,
one can also extract their m-D signatures by constructing their
respective HT equations.
B. When the Rotating or Vibrating Rate Is
Extraordinarily High
IV. S IMULATIONS
When BRF is given in a radar system, the higher the rotating
In this section, we verify the theoretical derivation and the
or vibrating rate, the less the slow-time-domain samplings of
proposed m-D extraction method using simulated data. The
the curve in a cycle. If the number of samplings in a cycle is
simulations are composed of three parts: 1) simulations to
too low, the detection capability of the HT will be depressed,
validate the range profile’s migration and wrapping; 2) simu-
and the parameters of the curve cannot accurately be extracted.
lations to validate the m-D signatures induced by rotation; and
Simulation results have demonstrated that the number of sam-
3) simulations of m-D signature extraction. Simulation results
plings in a cycle should be larger than 10–15 to ensure precise
of vibration are not presented, because it is very close to the
detection ability of the HT. For example, in a radar system with
case with rotation.
BRF = 200 Hz, the parameters of the curve are difficult to
directly be detected from the range-slow-time plane when the
A. Validation of the Range Profile’s Migration and Wrapping
rotating or vibrating rate is higher than 20 Hz. (Please refer
to the simulations in Section IV-C3.) Furthermore, when the The parameters of the transmitted radar signal are f0 =
value of v increases, the value of the third term in (9) increases 35 GHz, Tr = 78.125 μs, Δf = 4.6875 MHz, N = 64, and
accordingly and more seriously expands the peaks in the HRRP. the total bandwidth B = 300 MHz. Four micromotion scat-
Therefore, some supplementary processing are also necessary terers are at the same location with the distance from the
when extracting the m-D signatures of micromotion parts with reference point RΔ0 = 5 m but with different instantaneous
high rotating or vibration rates, such as the tail rotors on a velocities. According to (14), it can easily be calculated that v ∈
helicopter or missile. (−36 m/s, 18.8571 m/s) to avoid the range profile’s wrapping.
LUO et al.: MICRO-DOPPLER EFFECT ANALYSIS AND FEATURE EXTRACTION 2093

Fig. 5. Micro-Doppler signatures induced by rotation. (a) Curves of the scatterers (2 Hz, 2 m) and (5 Hz, 2 m) on the range-slow-time plane. (b) Curves of
the scatterers (2 Hz, 2 m) and (8 Hz, 2 m) on the range-slow-time plane. The “◦”-dotted and “”-dotted dashed curves are determined by (20) when the motion
parameters of the scatterer are (2 Hz, 2 m) and (8 Hz, 2 m), respectively.

As shown in Fig. 4, the peak locations of the scatterers in the TABLE I

C URVE PARAMETERS D IRECTLY O BTAINED VIA HT
HRRP are different from each other because of their different F ROM THE R ANGE -S LOW-T IME P LANE
instantaneous velocities. In particular, the profile’s wrappings
are validated in Fig. 4(b). When v is 19 m/s close to the critical
point 18.8571 m/s, the peak of the scatterer in the HRRP is split
and appears at both sides of the HRRP. When v exceeds the
limited value of (14), such as −60 m/s, the corresponding peak
of the scatterer is wrapped and appears at the right side of the
peak of the scatterer with v = −15 m/s. It can also be observed
that the peaks are more and more seriously expanded versus the
horizontal axis with the increase in |v|.
curves depicted by the equation on the range-slow-time plane
using the “◦”-dotted and “”-dotted dashed lines. It can be
B. Validation of m-D Signatures Induced by Rotation observed that the curves determined by (20) coincide well with
the simulation results.
The parameters of the transmitted radar signals are given as
follows: f0 = 10 GHz, Tr = 78.125 μs, Δf = 4.6875 MHz,
C. Validation of m-D Signature Extraction Algorithm
N = 64, BRF = 200 Hz, and the total bandwidth B =
300 MHz. The range resolution of HRRP is ΔR = 0.5 m. Two The parameters of the transmitted radar signals are the same
scatterers both rotate around one center with the distance from as those in Section IV-B. In this section, the effectiveness of the
the reference point RQ = −3 m. The radii are both 2 m, and m-D signature extraction algorithm is validated in Sections IV-
the initial phase θ0 are both π/2. As shown in Fig. 5(a), when C1–IV-C3, and the robustness of the algorithm is discussed in
the rotating rate is 2 and 5 Hz, respectively, it can be calculated Section IV-C4.
that |k|max = 2.6824 < π according to (18), which satisfies the 1) Detect Curves Using HT on the Range-Slow-Time Plane
conditions limited by (12). The magnifying coefficient of the si- Directly Under the Most Common Conditions: Assume that
nusoid’s amplitude 1 + ω02 f02 Tr2 /Δf 2 is 2.3209 and 5.3306, two scatterers are both rotating around one center with the dis-
respectively. Hence, the amplitude of the sinusoid in Fig. 5(a) tance from the reference point RQ = −3 m. The radii are both 2
is the rotating radius multiplied by 1 + ω02 f02 Tr2 /Δf 2 , which m, and the rotating rates are 2 and 8 Hz, respectively. The initial
is 4.6418 and 10.6613 m, respectively. phase θ0 are both π/2. The corresponding curves on the range-
In the simulation, the scatterer with a rotating rate of 5 Hz slow-time plane have been shown in Fig. 5(b). The parameter
was set to reach the critical position with the largest distance in pairs (RQ , A, Ω, ζ) detected by HT are shown in Table I. It
LOS from the reference point at tm = 0.05 s, 0.25 s, 0.45 s, . . .. can be found that the rotating rates of the two scatterers and
Taking tm = 0.25 s into consideration, although the scatterer the value of RQ are accurately obtained. The average value of
has moved to the furthest position in LOS from the reference parameter A is about 4.75 and 17 m, respectively. According to
point, the corresponding sinusoid on the range-slow-time plane (28), the radii of rotating scatterers are 2.0466 and 2.0149 m,
has not reached the crest at the moment, such as point A shown which agree quite well with the set value of 2 m.
in Fig. 5(a). The phase shift φ can be calculated as 79.1875◦ 2) Detect Curves When A < c/(4N Δf ): Assume that two
according to (18). As a result, the sinusoid crests at the position scatterers rotate around the center with the distance from the
indicated by point B in Fig. 5(a). reference point RQ = −3 m. The radii are both 0.05 m, and the
When the rotating rate of the scatterer (5 Hz, 2 m) increases rotating rate is 2 and 4 Hz, respectively. The initial phase θ0 are
up to 8 Hz, it can be calculated that |k|max = 3.9023 > π. both 0 rad. It can be calculated that A1 = 0.1160, A2 = 0.2153,
Therefore, the sinusoid will be wrapped as shown in Fig. 5(b). and c/(4N Δf ) = 0.25, which satisfy A1 < c/(4N Δf ) and
For the sake of validating the correctness of (20), we draw the A2 < c/(4N Δf ). Because the amplitudes are smaller than
2094 IEEE TRANSACTIONS ON GEOSCIENCE AND REMOTE SENSING, VOL. 48, NO. 4, APRIL 2010

Fig. 6. Micro-Doppler signature extraction when A < c/(4N Δf ). (a) Curves on the range-slow-time plane appear to approximately be a straight line.
(b) Sinusoidal curve obtained by applying the Gabor transform to the row of straight lines.

TABLE II
C URVE PARAMETERS O BTAINED VIA HT W HEN A < C/(4N Δf )

Fig. 7. Micro-Doppler signature extraction when the rotating rate is extraordinarily high. (a) Blurry curves on the range-slow-time plane. (b) STFT result
of Sc (ω, i) in (8) at ω = −4πμRΔ /c. It appears to be segments of straight lines on the time–frequency plane. (c) New frequency (range)-slow-time plane
constructed by stitching up the planes of all the bursts in order.

half-range resolution, the curves appear to be a straight line around one center with the distance from the reference point
approximately, as shown in Fig. 6(a). In this condition, the HT RQ = −3 m. The radii are both ρ = 1 m, and the rotating rates
is unable to detect the parameters from the range-slow-time are the same at 20 Hz. The initial phase θ0 is 0 and π/2, respec-
plane. Extracting the row of straight line, taking interpolations tively.
The magnifying coefficient of the sinusoid’s amplitude is
for four times, and then implementing Gabor transform to the 1 + ω02 f02 Tr2 /Δf 2 = 20.9678, i.e., A = 20.9678 according
row, the sinusoidal curves are obtained as shown in Fig. 6(b). to (28).
The synthesis window in the Gabor transform is a Gaussian The curves on the range-slow-time plane are shown in
window. The number of Gabor coefficients in time is 100, Fig. 7(a). The curves are very blurry because of the high rotat-
and the degree of oversampling is 20. Then, the HT can be ing rates; as a result, the parameters of the curves cannot exactly
implemented to obtain the parameters (A , Ω, ζ), and the results be detected by the HT from the range-slow-time plane, as
are shown in Table II. shown in Table III. It can be found that the detected parameters
It can be observed in Table II that the rotating rates are fail to fit to the set values.
accurately detected as Ω = 2 and 4 Hz. The average values of Fig. 7(b) shows the STFT result of Sc (ω, i) in (8) at ω =
parameter A are about 43.33 and 84.75, respectively. Accord- −4πμRΔ /c. It can be seen that Sc (ω, i)|ω=−4πμRΔ /c is an
ing to (34), the rotating radii are calculated as ρ = A c/(2f0 · approximate LFM signal on the time–frequency plane, i.e.,
2πΩ) = 0.0517 m and 0.0506 m, which are quite close to the the i−k plane. Stitching up the planes of all the bursts in
set value of 0.05 m. In fact, the errors are mainly dependent on order, a new frequency-slow-time plane can be constructed,
the frequency resolution of the Gabor transform. as shown in Fig. 7(c). For clarity of the figure and ease
3) Detect Curves When the Rotating or Vibrating Rate Is of comprehension, Fig. 7(c) presents just a part of the new
Extraordinarily High: Assume that two scatterers both rotate plane from 0 to 0.1 s in the slow time, and the frequency
LUO et al.: MICRO-DOPPLER EFFECT ANALYSIS AND FEATURE EXTRACTION 2095

TABLE III As shown in Fig. 8(a), the simulation conditions are the same
C URVE PARAMETERS O BTAINED VIA HT FROM THE R ANGE -S LOW-T IME
P LANE W HEN THE ROTATING R ATE I S E XTRAORDINARILY H IGH as those in Section IV-C1, and the body’s returns are added
with SNR = −20 dB. It can be found that the m-D curves are
relatively very weak and that the HT may have errors when
detecting the parameter pairs of curves. The result of the noise
suppression processing is shown in Fig. 8(b). Applying the HT
to this figure, the parameter pairs can accurately be detected.
However, similar with many other algorithms, the m-D ex-
traction algorithm presented herewith cannot effectively work
when the noise is overwhelming. As shown in Fig. 8(c), when
the SNR decreases to −30 dB, the noise suppression does not
TABLE IV effectively work, so that the curve with a frequency of 8 Hz
C URVE PARAMETERS O BTAINED VIA STFT AND HT W HEN
THE ROTATING R ATE I S E XTRAORDINARILY H IGH can no longer be detected by HT. In this case, it is necessary to
improve the noise suppression algorithm in future works.
When A < c/(4N Δf ) or the rotating rate is extraordinarily
high, the Gabor transform or the STFT needs to be implemented
before the HT. In these cases, the body’s returns will bring in
adverse effects because of the limitation of the time–frequency
resolution. Therefore, the antinoise ability of the algorithm
will be degraded to some extent. In particular, for the STFT,
the multicomponent signal of the body’s returns will seriously
contaminate the m-D curves on the frequency-slow-time plane.
axis is transformed to the range axis for consistency with Fig. 9 shows the robustness of the proposed algorithm when
Fig. 7(a). A < c/(4N Δf ). The simulation conditions are the same as
Thus, the curves in Fig. 7(c) are more obvious, and the those in Section IV-C2, and SNR = −16 dB. The results of
HT can be utilized to detect the parameters exactly. The de- Gabor transform and the noise suppression processing are
tected parameters (RQ , A, Ω, ζ) of the two curves are given in shown in Fig. 9(a) and (b), respectively, and the HT can
Table IV and fit well to the set values. accurately detect the parameters of the curves.
4) Robustness Discussion: In actual situations, the micro- Fig. 10 shows the robustness of the algorithm when the
motion parts are always masked by the target’s body. As a rotating rate is extraordinarily high. The simulation conditions
result, the induced m-D signals are contaminated by the target are the same as those in Section IV-C3, and SNR = −10 dB.
body’s returns. In the processing of m-D feature extraction, Applying the HT to Fig. 10(c), the m-D features can success-
these body returns must be considered as “noise,” because fully be extracted.
they will bring in mistakes. Meanwhile, for simplicity, we Time-varying complex scattering coefficients: It is highly
have assumed the target consists of point scatterers with unit possible that the scattering coefficients are time-varying com-
scattering coefficients. However, for a rotating target with a plex values if the rotating scatterers have large radii or high
large radius or high frequency in real-world situations, it is frequencies. When the rotating scatterers have large radii, we
highly possible that the backscattering coefficients have time- directly detect the curves from the range-slow-time plane;
varying complex values. In the following, we further study the therefore, the time-varying complex scattering coefficients just
performance of the proposed algorithm in these situations. lead to the varying moduli of the curves. When the rotating rates
SNR: After sufficiently accurate translation compensa- of the scatterers are extraordinarily high, we detect the curves
tion, the body scatterers can be considered as micromotion on the frequency (range)-slow-time plane after taking the STFT
scatterers with zero rotating frequency approximately; there- to each burst and stitching up the time–frequency planes in
fore, the corresponding curves on the range-slow-time plane order, and in each burst time, the scattering coefficients can be
or the frequency-slow-time plane appear to be straight lines. considered as constants approximately. In such case, the time-
According to the difference between the m-D curves and the varying complex scattering coefficients also just lead to the
body’s straight lines, a primary and simple noise suppression varying moduli of the curves on the frequency (range)-slow-
method is designed as follows: assuming that the data matrix time plane. Due to the robustness of the HT to the discontinu-
on the range-slow-time plane or the frequency-slow-time plane ities of the pattern, the adverse effects of the curves’ varying
is M with X × Y points, the straight-line-removing processing moduli are very limited.
can be expressed as Fig. 11 shows the curves on the range-slow-time plane,

Y  where the simulation conditions are the same as those in
 1  Section IV-C1. The moduli and phases of scattering coefficients
M (x, y) = |M(x, y)| − |M(x, y)| · MY (35)
Y y=1 randomly vary within [0, 1] and [0, 2π], respectively. Applying
the HT to Fig. 11, the parameter pairs of the curves can
where MY is a 1 × Y vector with all elements equal to one. accurately be detected.
Then, the HT can be implemented to M (x, y) for m-D feature Fig. 12 shows the curves on the frequency (range)-slow-time
extraction. plane, where the simulation conditions are identical to those
2096 IEEE TRANSACTIONS ON GEOSCIENCE AND REMOTE SENSING, VOL. 48, NO. 4, APRIL 2010

Fig. 8. Robustness validation of the proposed algorithm under different SNRs when detecting curves on the range-slow-time plane. (a) Curves and straight lines
correspond to the micromotion scatterers and the body scatterers on the range-slow-time plane, respectively. The straight lines are considered as “noise” when
extracting the m-D features. The SNR is −20 dB. (b) Result of the noise suppression processing when SNR = −20 dB. (c) Result of the noise suppression
processing when SNR = −30 dB.

Fig. 9. Robustness validation when A < c/(4N Δf ) and SNR = −16 dB. (a) Result of Gabor transform. (b) Result of the noise suppression processing.

Fig. 10. Robustness validation when the rotating rate is extraordinarily high and SNR = −10 dB. (a) Blurry curves correspond to the micromotion scatterers,
and straight lines correspond to the body scatterers on the range-slow-time plane. (b) New frequency (range)-slow-time plane constructed by stitching up the STFT
results of all the bursts in order. (c) Result of the noise suppression processing.

in Section IV-C3. The moduli and phases of scattering coeffi-

cients also randomly vary within [0, 1] and [0, 2π], respectively.
Applying the HT to Fig. 12, the parameter pairs of the curves
can also accurately be detected.
From these analyses and simulations, it can be found that the
proposed algorithm possesses good antinoise performance and
robustness to the time-varying complex scattering coefficients.

V. D ISCUSSION AND C ONCLUSION

The m-D information has been verified to be very useful
for radar target identification and recognition. In this paper,
Fig. 11. Curves with time-varying complex scattering coefficients on the the m-D effects in SFCS radar have been investigated, and
range-slow-time plane. an algorithm for m-D signature extraction has been proposed.
LUO et al.: MICRO-DOPPLER EFFECT ANALYSIS AND FEATURE EXTRACTION 2097

using simulated data at present. In our future work, we will

introduce the algorithm to the real data processing.

R EFERENCES
[1] V. C. Chen, “Analysis of radar micro-Doppler signature with
time–frequency transform,” in Proc. Stat. Signal Array Process., Pocono
Manor, PA, 2000, pp. 463–466.
[2] V. C. Chen, “Micro-Doppler effect of micromotion dynamics: A review,”
Proc. SPIE, vol. 5102, pp. 240–249, Apr. 2003.
[3] V. C. Chen, F. Li, S. S. Ho, and H. Wechsler, “Micro-Doppler effect in
radar: Phenomenon, model, and simulation study,” IEEE Trans. Aerosp.
Electron. Syst., vol. 42, no. 1, pp. 2–21, Jan. 2006.
[4] T. Sparr and B. Krane, “Micro-Doppler analysis of vibrating targets in
SAR,” Proc. Inst. Elect. Eng.—Radar Sonar Navig., vol. 150, no. 4,
Fig. 12. Curves with time-varying complex scattering coefficients on the pp. 277–283, Aug. 2003.
frequency (range)-slow-time plane after taking the STFT. [5] T. Thayaparan, S. Abrol, E. Riseborough, L. Stankovic, D. Lamothe, and
G. Duff, “Analysis of radar micro-Doppler signatures from experimental
helicopter and human data,” IET Radar Sonar Navig., vol. 1, no. 4,
The proposed m-D signature extraction method is mainly based pp. 289–299, Aug. 2007.
on the HT to detect the parameters of the curves on the range- [6] J. Li and H. Ling, “Application of adaptive chirplet representation for
ISAR feature extraction from targets with rotating parts,” Proc. Inst. Elect.
slow-time plane or the frequency-slow-time plane. It is known Eng.—Radar Sonar Navig., vol. 150, no. 4, pp. 284–291, Aug. 2003.
that, when using the HT to detect patterns, the computer storage [7] X. Bai, M. Xing, F. Zhou, G. Lu, and Z. Bao, “Imaging of micromo-
tion targets with rotating parts based on empirical-mode decomposition,”
and the computational time exponentially grow with the number IEEE Trans. Geosci. Remote Sens., vol. 46, no. 11, pp. 3514–3523,
of dimensions in the parameter space [22]. The HT that we Nov. 2008.
utilized here is a four-parameter mapping (when using the HT [8] L. Stankovic, I. Djurovic, and T. Thayaparan, “Separation of target rigid
body and micro-Doppler effects in ISAR imaging,” IEEE Trans. Aerosp.
on the range-slow-time plane directly under the most common Electron. Syst., vol. 42, no. 4, pp. 1496–1506, Oct. 2006.
conditions or when the rotating rate is extraordinarily high) or a [9] Q. Zhang, T. S. Yeo, H. S. Tan, and Y. Luo, “Imaging of a moving target
three-parameter mapping [when A < c/(4N Δf )]. In our sim- with rotating parts based on the Hough transform,” IEEE Trans. Geosci.
Remote Sens., vol. 46, no. 1, pp. 291–299, Jan. 2008.
ulation, the computation time for obtaining the parameter pairs [10] S. He, Y. F. Zhu, H. Z. Zhao, J. X. Zhou, and F. Qiang, “Analysis of
(A , Ω, ζ) in Section IV-C2 was approximately 175 s on a per- rotating structures for stepped frequency radar,” in Proc. Int. Conf. Radar,
sonal computer with an Intel Pentium Dual T2330 (1.60 GHz) 2008, pp. 386–390.
[11] A. Ghaleb, L. Vignaud, and J. M. Nicolas, “Micro-Doppler analysis of
processor, whereas that for obtaining the parameter pairs wheels and pedestrians in ISAR imaging,” IET Signal Process., vol. 2,
(RQ , A, Ω, ζ) from Fig. 7(c) was approximately 1311 s. In no. 3, pp. 301–311, Sep. 2008.
order to reduce the computation cost for real-time applications, [12] P. Setlur, M. Amin, and F. Ahmad, “Urban target classifications using
time–frequency micro-Doppler signatures,” in Proc. 9th Int. Symp. Signal
some supplementary processing for speeding up the computa- Process. Appl., 2007, pp. 1–4.
tion of the algorithm can be used, e.g., the fast method presented [13] Y. Kim and H. Ling, “Human activity classification based on micro-
Doppler signatures using a support vector machine,” IEEE Trans. Geosci.
in [25] to detect fixed-period sinusoids based on the HT by Remote Sens., vol. 47, no. 5, pp. 1328–1337, May 2009.
decomposing the parameters into two groups and estimating [14] Q. Zhang and Y. Q. Jin, “Aspects of radar imaging using frequency-
them in turn can be extended to detect curves in this paper. The stepped chirp signals,” EURASIP J. Appl. Signal Process., vol. 2006,
no. 13, pp. 85 823-1–85 823-8, 2006.
parallel processing technique can also be used to accelerate the [15] A. Freedman, R. Bose, and B. D. Steinberg, “Thinned stepped frequency
computation of the HT. waveforms to furnish existing radars with imaging capability,” IEEE
It is mentioned that the main purpose of this paper is Aerosp. Electron. Syst. Mag., vol. 11, no. 11, pp. 39–43, Nov. 1996.
[16] N. Levanon and E. Mozeson, “Nullifying ACF grating lobes in stepped-
the m-D effect analysis and signature extraction in SFCS frequency train of LFM pulses,” IEEE Trans. Aerosp. Electron. Syst.,
radars. In fact, based on the m-D signature extraction, the vol. 39, no. 2, pp. 694–703, Apr. 2003.
final ISAR image of the target body can be improved after [17] Q. Zhang, Y. S. Zeng, Y. Q. He, and Y. Luo, “Avian detection and iden-
tification with high-resolution radar,” in Proc. IEEE Radar Conf., Rome,
eliminating the curves from the range-slow-time plane. The Italy, May 26–30, 2008, pp. 2194–2199.
elimination algorithm for chirp radars has been discussed in [18] H. Schimpf, A. Wahlen, and H. Essen, “High range resolution by means
our previous work [9], and it can also be extended to focus a of synthetic bandwidth generated by frequency-stepped chirps,” Electron.
Lett., vol. 39, no. 18, pp. 1346–1348, Sep. 2003.
clear body image in SFCS radars. Of course, considering the [19] R. T. Lord and M. R. Inggs, “High resolution SAR processing using
m-D characteristic differences between the chirp radar and the stepped-frequencies,” in Proc. IEEE IGARSS, Singapore, Aug. 1997,
SFCS radar, the algorithm should be modified accordingly. vol. 1, pp. 490–492.
[20] D. J. Rabideau, “Nonlinear synthetic wideband waveforms,” in Proc.
For example, in this paper, the Gabor transform and STFT IEEE Radar Conf., Long Beach, CA, May 2002, pp. 212–219.
are adopted before the HT, and therefore, how the signals in [21] P. V. C. Hough, “Methods and means for recognizing complex patterns,”
U.S. Patent 3 069 654, Dec. 18, 1962.
range-slow-time domain could be reconstructed after curves [22] D. H. Ballard, “Generalization the Hough transform to detect arbitrary
eliminated from the time–frequency domain will be a subject shapes,” Pattern Recognit., vol. 13, no. 2, pp. 111–122, 1981.
for further research. The ISAR imaging algorithm based on the [23] D. S. Mckenzie and S. R. Protheroe, “Curve description using the in-
verse Hough transform,” Pattern Recognit., vol. 23, no. 3/4, pp. 283–290,
m-D signature extraction will be considered in a separate paper. Mar. 1990.
The simulations in Section IV have validated the accuracy [24] S. Qian and D. Chen, “Discrete Gabor transform,” IEEE Trans. Signal
of m-D effect analysis and the effectiveness of the proposed Process., vol. 41, no. 7, pp. 2429–2438, Jul. 1993.
[25] Z. Changchun and S. Ge, “A Hough transform-based method for fast
signature extraction algorithm. Due to the lack of the required detection of fixed period sinusoidal curves in images,” in Proc. 6th Int.
real data, we can just examine the performance of the algorithm Conf. Signal Process., 2002, pp. 909–912.
2098 IEEE TRANSACTIONS ON GEOSCIENCE AND REMOTE SENSING, VOL. 48, NO. 4, APRIL 2010

Ying Luo was born in Hunan, China, in 1984. He Cheng-wei Qiu (M’06) was born in Zhejiang, China,
received the M.S. degree in electrical engineering, in on March 9, 1981. He received the B.Eng. degree
2008, from the Institute of Telecommunication Engi- from the University of Science and Technology of
neering, Air Force Engineering University (AFEU), China, Hefei, China, in 2003 and the Ph.D. de-
Xi’an, China, where he is currently working toward gree from the Joint Ph.D. Program between the
the Ph.D. degree in electrical engineering. He is cur- National University of Singapore, Singapore and
rently with the Radar and Signal Processing Labo- the L’École Supérieure d’Électricité (SUPELEC),
ratory, Institute of Telecommunication Engineering, Gif-sur-Yvette, France, in 2007.
AFEU. His research interests include signal process- Since 2009, he has been a Postdoctoral Fel-
ing and ATR in SAR and ISAR. low with the Research Laboratory of Electronic,
Massachusetts Institute of Technology, Cambridge,
MA. He has authored more than 40 journal papers and one book chapter, and
has given a few invited talks at conferences. His research interests include the
electromagnetic wave theory of complex media (e.g., chiral, anisotropic, and
bianisotropic materials), invisibility cloaks, metamaterials, and imaging.
Dr. Qiu was the recipient of the SUMMA Graduate Fellowship in Advanced
Electromagnetics in 2005; the IEEE AP-S Graduate Research Award in 2006;
the Young Scientist Travel Grant in Japan in 2007; the Travel Grant for
Metamaterial ’07 in Rome, Italy; and the URSI Young Scientist Award in 2008.

Qun Zhang (M’02–SM’07) received the M.S. de- Xian-jiao Liang was born in Hunan, China, in 1986.
gree in mathematics from Shaanxi Normal Univer- She is currently working toward the M.S. degree in
sity, Xi’an, China, in 1988 and the Ph.D. degree electrical engineering in the Institute of Telecom-
in electrical engineering from Xidian University, munication Engineering, Air Force Engineering Uni-
Xi’an, in 2001. He was a Research Engineer from versity, Xi’an, China. Her research interests include
2001 to 2003 and a Research Fellow from 2005 to signal processing in SAR and ISAR.
2006 with the Department of Electrical and Com-
puter Engineering, National University of Singapore,
Singapore. He is currently a Professor with the
Telecommunication Engineering Institute, Air Force
Engineering University, Xi’an, and an Adjunct Pro-
fessor with the Key Laboratory of Wave Scattering and Remote Sensing
Information (Ministry of Education), Fudan University, Shanghai, China. He
has authored more than 100 papers on journals and conference proceedings. Kai-ming Li received the M.S. degree in electrical
His research interests include signal processing, clutter suppression, and its engineering from the Institute of Telecommunica-
application in SAR and ISAR. tion Engineering, Air Force Engineering University
Dr. Zhang is a Senior Member of the Chinese Institute of Electronics (CIE); (AFEU), Xi’an, China, in 2009. He is currently with
a Committee Member of the Radiolocation Techniques Branch, CIE; and a the Radar and Signal Processing Laboratory, Insti-
member of the Signal Processing Council of Shaanxi, China. He was also a tute of Telecommunication Engineering, AFEU. His
Vice Co-Chair of the Technical Program Committee of Workshop for Space, research interests include signal processing and ATR
Aeronautical and Navigational Electronics in 2009. He was the recipient of the in SAR and ISAR.
First-Grade and the Second-Grade Prizes of Shaanxi Natural Science Excellent
Academic Paper in 2006 and 2008, respectively.