Skip to main content

MTRC compensation in high-resolution ISAR imaging via improved polar format algorithm based on ICPF

Abstract

In this paper, we present a detailed analysis on the performance degradation of inverse synthetic aperture radar (ISAR) imagery with the polar format algorithm (PFA) due to the inaccurate rotation center. And a novel algorithm is developed to estimate the rotation center for ISAR targets to overcome the degradation. In real ISAR scenarios, the real rotation center shift is usually not coincided with the gravity center of the high-resolution range profile (HRRP), due to the data-driven translational motion compensation. Because of the imprecise information of rotation center, PFA image yields model errors and severe blurring in the cross-range direction. To tackle this problem, an improved PFA based on integrated cubic phase function (ICPF) is proposed. In the method, the rotation center in the slant range is estimated firstly by ICPF, and the signal is shifted accordingly. Finally, the standard PFA algorithm can be carried out straightforwardly. With the proposed method, wide-angle ISAR imagery of non-cooperative targets can be achieved by PFA with improved focus quality. Simulation and real-data experiments confirm the effectiveness of the proposal.

1 Introduction

Inverse synthetic aperture radar (ISAR) can generate two-dimensional (2D) images for a non-cooperative moving target, which has wide military (e.g., military target classification and recognition) and civil (e.g., civil aircraft control) applications. The high range resolution is proportional to the bandwidth of the transmitted signal. And the high cross-range resolution is achieved using a synthetic antenna aperture respecting to the radar line of sight from the target rotation [1]. The range resolution ρ r is defined by ρ r  = c/(2B), where c is the velocity of light, and B is the signal bandwidth. And based on the 2D rotational model, the cross-range resolution ρ a is usually defined as ρ a  = λ/(2Θ), where λ is the signal wavelength, and Θ is the target rotation angle within the coherent processing interval (CPI).

Several useful ISAR imaging algorithms have been proposed in the literatures, such as the Range-Doppler algorithm (RDA) [2], the Keystone algorithm for precise de-coupling [36], and the polar format algorithm (PFA) [710]. The RDA is applicable when the rotation angle is small to ensure there is no significant migration through resolution cells (MTRC) and time-varying Doppler in the echoed signal. The Keystone algorithm is optimal to correct the linear MTRC caused by the rotation for all scattering centers. However, its performance is also limited by the rotation angle. In high-resolution ISAR imaging (always needs a wide rotation angle), high-order MTRC together with the linear one are present, where the Keystone algorithm would degrade significantly. The PFA was first used to compensate the MTRC in SAR imaging [7]. It was introduced in ISAR imaging later. The PFA algorithm firstly interpolates the observed signal from polar sector support region to the Cartesian rectangular region. And the fast Fourier transform (FFT) is followed to focus the ISAR image [10, 11]. Since both the range and the Doppler curve are permitted in the PFA, its applicable rotation angle is significantly wider than those of both the RDA and the Keystone algorithm.

In practical ISAR processing, PFA requires accurate information of the rotation parameters before the interpolation. However, those parameters are usually unavailable in advance. Due to the uncooperative characteristic of ISAR targets, translational motion compensation [1214] is also necessary to eliminate the translational range shift and phase errors in the echoes. This processing is usually data-driven. As a result, after translational compensation, the rotation center (RC) is commonly deviated from the measured reference range of the radar. In this case, serious blur and defocus are also present in the cross range of a PFA image.

This paper first introduces the standard PFA processing for ISAR imaging and then analyzes the performance degradation caused by the imprecise rotation center. From the analysis, one can find that the PFA image by using an inaccurate rotation center is equivalent to the response of a disturbance system with the ideal PFA image as the input. Therefore, for real applications of PFA, an accurate and robust method for rotation parameter estimation is necessary.

Many approaches have been proposed for the rotation parameter estimation in ISAR imaging. In ref. [15], the adaptive S-method is applied to reduce the effect of quadratic phase errors in each slant-range cell by calculating the conjugate correlation of signal sequences. However, it also suffers from the cross terms. In ref. [16], a method was proposed to estimate the rotation angle and center using a correlation-based function. It needs careful choosing and tracking of two scattering centers in the range-compressed phase history, which may be difficult for high-resolution ISAR. In ref [17], a method based on a concept of rotation correlation is proposed, which searches the rotation parameters by iteratively interpolating the RD image. This method may also be invalidated if the prominent scattering points are seriously interfered by their neighbors. The algorithms in ref. [18] are utilized for uniform rotation compensation in ISAR imaging, but only the rotation angle can be estimated in the processing. In ref. [19] and [20], two nonparametric methods, which utilize golden section search to estimate the rotation angle and center, are proposed. The methods use entropy as criterion and require an initial rotation angle, and the simulated and measurement data moderately show the validity. However, the concomitant risk of divergence cannot be avoided in the iteration and the computational load may become very large once the 2D map of the entropy is flat. In this paper, a method based on integrated cubic phase function (ICPF) is proposed for PFA to eliminate the effect of the disturbance system and the rotation center can be derived precisely. The ICPF extends the standard cubic phase function (CPF) and provides improved estimation performances of multi-component chirp signals at a cost of a moderate complexity increase. Accordingly, the signal after the translational compensation can be transformed into the standard polar format, and a well-focused image is finally obtained using the standard PFA algorithm. Experimental results were carried on both of the simulated and real data to demonstrate the accuracy of the analysis and the effectiveness of the improved method.

The rest of the paper is organized as follows. The signal model of ISAR imaging and the standard PFA processing are introduced in Section 2. Section 3 describes the PFA performance degradation from the imprecision of rotation center and makes a detailed analysis on the disturbance system. The improvement of PFA imaging based on ICPF is presented in Section 4. Experimental results and analysis are reported in Section 5, followed by conclusions in Section 6.

2 Signal model and standard PFA algorithm for ISAR imaging

2.1 Signal and geometry model

The geometry of the ISAR imaging is depicted in Figure 1, where the 2D rotation geometry with point scattering model is used [11]. The instantaneous distance from a scattering point p k with position (x k , y k ) to the radar can be approximated as:

Figure 1
figure 1

Geometry of the ISAR imaging for scattering point.

R k t = r a t + x k sinθ t + y k cosθ t
(1)

where r a is the range from target RC to the radar and θ(t) is the rotation angle at instantaneous time t.

In the ISAR imaging system, a linear frequency modulated (LFM) waveform is usually utilized and it can be formulated as:

s R t ^ , t =rect t ^ / T p exp - j 2 π f c t + 1 2 γ t ^ 2
(2)

where

rect u = 1 , u 1 / 2 0 , u > 1 / 2
(3)

and f c , T p , and γ are carrier frequency, pulse width, and chirp rate of the transmitted signal, respectively. The full time t= t ^ + t m is expressed as the form of the fast time t ^ and the slow time t m  = m PRI, where m are integers and PRI is the pulse repetition interval. Let us suppose that the target consists of K point scatterers, and the returned signals can be expressed by:

s R t ^ , t = k = 1 K σ k rect t ^ - 2 R k t / c / T p × exp - j 2 π f c t ^ - 2 R k t / c + 1 2 γ t ^ - 2 R k t / c 2 ,
(4)

where σ k , c, and f c are the backward scattering coefficient, light speed, and carrier frequency, respectively. For the target without high velocity, the ISAR imaging process satisfies the go-stop model [11], thus the full time t in (4) can be replaced by the slow time t m . In order to reduce the receiver bandwidth, the echo is often dechirped on reception [21]. Suppose that the reference point (RP) is o and its range to radar is r0(t m ). Then after the dechirping, the received signal is written as:

s R 2 t ^ , t m = k = 1 K σ k rect t ^ - 2 R k t m - r 0 t m / c / T p × exp - j 4 π c f c + γ t ^ - 2 r 0 t m / c R ps _ k t m exp j ψ k t m ,
(5)

where

R ps _ k t m = R k t m - r 0 t m , ψ t m = 4 π γ c 2 R ps _ k 2 t m ,
(6)

and Rps _k(t m ) is the instantaneous range from P k to o and ψ(t m ) is the residual video phase (RVP). Generally, the RVP can be removed easily [22], when dealing with small targets.

Let’s denote f=γ t ^ - 2 r 0 t m / c , and by substituting f into Equation (5), one can obtain:

s R 3 f , t m = k = 1 K σ k rect f / γ - 2 R k t m - r 0 t m / c T p × exp - j 4 π c f c + f R ps _ k t m .
(7)

Assuming that the RP are located at the RC, i.e., r0(t m ) = r a (t m ), then one can obtain Rps _k from (1):

R ps _ k t m = r a t m + x k sin θ t m + y k cos θ t m - r 0 t m = x k sin θ t m + y k cos θ t m
(8)

Since R k (t m ) ≈ r0(t m ) for small targets, the dechirped echo signals can be finally expressed as:

s R 4 f , θ = k = 1 K σ k rect f B × exp - j 4 π c f c + f x k sin θ t m + y k cos θ t m ,
(9)

where B = γT p denotes the signal bandwidth.

2.2 Standard polar format algorithm model for ISAR imaging

In order to obtain 2D ISAR image, FFT is usually carried on the slow time t m to yield the range-Doppler (RD) image of a target. However, in high-resolution ISAR imaging, RD algorithm will be invalid because of the presence of the MTRC and time-varying Doppler. The Keystone algorithm can correct the linear MTRC, but it is invalid for the high-order MTRC and phase errors. After the correction of linear MTRC, the ISAR image with large rotation angle would be still blurred. Since both the range and Doppler curve are permitted in the PFA, both range-azimuth coupling and time-varying Doppler can be compensated well and yields a well-focused image. In this subsection, the standard PFA algorithm will be introduced and the maximum rotation angles of PFA will be analyzed with comparison to RDA and Keystone algorithm. The efficiency of PFA in MTRC compensation will be also illustrated.

Let K R  = -(4π/c) (f c  + f), then the phase term of (7) can be expressed as:

Φ k θ , K R = K R x k sin θ + y k cos θ = K x x k + K y y k ,
(10)

where (θ, K R ) is the polar format of (K x , K y ) in a Cartesian coordinate system. If the total target rotation angle Θ is known, the target image can be reconstructed from 2D Fourier transform (2D-FT). However, the FT algorithm requires a rectangular raster distribution of the data points, and the signal model can acquire a sector support region. Therefore, the interpolation of the observed signal on a polar sector support region to a Cartesian rectangular system is required.

As it is shown in Figure 2, the radial wavenumber K R and the rotation angle θ satisfies:

Figure 2
figure 2

Wavenumber spectrum after de-chirping.

K R K R 1 , K R 2 = 4 π c f c - B 2 , f c + B 2 θ - Θ 2 , Θ 2
(11)

The azimuth wavenumber is (4π/λ)Θ.

An illustration of the PFA algorithm is shown in Figure 3. If KR 2 cos Θ ≤ KR 1, the interpolation error will rise with Θ increasing, which results in a large vacant space in the wavenumber spectrum filled with stripe in Figure 3. To avoid this phenomenon, the maximum rotation angle for PFA imaging is constrained as:

Figure 3
figure 3

Wavenumber spectrum for K R 2 cos ΘK R 1.

Θ P F A =2 cos - 1 K R 1 K R 2 =2 cos - 1 1 - 0.5 μ 1 + 0.5 μ ,
(12)

where μ = B/f c . According to (10), ΘPFA is merely a function of μ independent of the target size and the range resolution. For comparison, the maximum rotation angles for RD and Keystone are given below. In order to avoid the MTRC, the maximum rotation angle for the RD algorithm in ISAR imaging should satisfy [2]:

1 ρ r L X 2 Θ < 1 1 ρ a L Y 2 Θ < 1 ,
(13)

where L X and L Y denote the target size on the cross range and slant range, respectively. For a fixed size target, the maximum rotation angle for RD algorithm is:

Θ R D =min 2 ρ r L X , λ L Y .
(14)

When the rotation theta Θ in CPI is larger than ΘRD, the MTRC and time-varying Doppler will occur failing RDA to generate optimal images.

The linear MTRC can be compensated using the Keystone algorithm [5, 6]. Given the target size L X and L Y , the Keystone algorithm restricts the maximum rotation angle by:

Θ Keystone =min 4 ρ r L Y , λ L Y .
(15)

To present the analysis in an explicit manner, ISAR images of a simulated aircraft with the three different algorithms are presented in Figure 4. The aircraft is composed of 49 point scatterers with equal reflectivity. The size of the aircraft shown in Figure 4a is 20 m × 20 m. The radar transmits LFM signals with a carried frequency of 15 GHz and a bandwidth of 2 GHz. The range resolution is ρ r  = 0.075. From Equations (12), (14), and (15), the maximum rotation angles of RD, Keystone, and PFA for the target in Figure 4 are ΘRD = 0.0075 rad = 0.43°, ΘKeystone = 0.032 rad = 1.82°, and ΘPFA = 1.01 rad = 57.91°, respectively. Obviously, ΘPFA is much larger than ΘRD and ΘKeystone. It must be declared that, it is difficult to realize a uniform rotation for an azimuth angle of ΘPFA = 57.91°. In practical ISAR processing, a rotation angle 6°  10° is adaptive to obtain a relevant cross-range resolution. Even though, the value of the rotation angle is over the limitations of the RD algorithm and Keystone transform. When the radar observes the target flying stably, the rotation velocity can be approximate as constant in the range of 6°  10°. Therefore, the PFA is more efficient and significant than RD algorithm and Keystone transform in high-resolution ISAR imaging.

Figure 4
figure 4

Images of an aircraft with the three algorithms. (a) Target model. (b) RD image. (c) Image with Keystone transform. (d) Image with PFA.

In the above analysis, an underlying assumption is that the translational motion is compensated ideally and only the rotational motion is involved. With the development of the translational motion compensation algorithm, the assumption can be satisfied in much case. In order to obtain high cross-range resolution matching with the range resolution, the target coherent rotation angle is set to Θ = 0.2 rad and the rotational velocity is set to ω = 0.01 rad/s. In the dataset, the pulse repetition frequency (PRF) is 50 Hz, i.e., 1,024 pulses within dwell time [-10, 10] (s) are used. Meanwhile, the number of range bins in high-resolution range profile (HRRP) is N = 1,000. The RC and RP are both located at (0 m, 0 m) in the target coordinate system.

Figure 4b shows the ISAR image achieved by RDA. Since the rotation angle Θ is larger than ΘRD, the blurring effects caused by MTRC and time-varying Doppler can be found in both of the range cell and Doppler cell. The obtained ISAR image using the Keystone algorithm is shown in Figure 4c. It can be seen that the blurring in range cells caused by MTRC is compensated, but the blur in the Doppler dimension still remains. Figure 4d presents the ISAR image obtained by the PFA algorithm, where all target scatterers are optimally focused.

According to the algorithm analysis and simulations, the PFA has advantage on MTRC compensation with the resolution of radar increasing.

3 Analysis of the disturbance system in PFA

3.1 The performance degradation of PFA when Δr ≠ 0

According to Equations (8 to 10), one notes that when the RC locates at the cell of the RP, the phase term of the dechirped signals can be expressed in the polar format. Thus we can use the standard PFA algorithm to obtain a focused ISAR image, as shown in Figure 4. However, in practical applications, since the target is non-cooperative and the error of range measurement is unavoidable in the radar system, after translational motion compensation, RC is usually not coincided with RP. Assuming that the range difference between RC and RP is Δr = r a  - r0, Equations (8) and (9) can be expressed as:

R ps _ k =Δr+ x k sinθ t m + y k cosθ t m ,
(16)
s R 4 f , θ = k = 1 K σ k rect f B × exp - j 4 π c f c + f Δ r + x k sin θ t m + y k cos θ t m .
(17)

Taking scattering point p k for example, its echo can be formulated as:

s R 4 _ k f , θ = σ k rect f B × exp - j 4 π c f c + f Δ r + x k sin θ t m + y k cos θ t m
(18)

Then, the phase term of Equation (18) is expressed by:

Φ k θ , K R = K R Δ r + x k sin θ t m + y k cos θ t m = K R Δ r + x k K R sin θ t m + y k K R cos θ t m .
(19)

Equation (19) indicates that the disturbed phase term K R Δr can cause the performance degradation when the PFA algorithm is applied. For clarity, we provide some simulation results in Figure 5 to indicate the degradation of PFA image from the imprecise RC. The target is the one in Figure 4. The new given dechirped signals in Figure 5 have the RC highlighted by the red dashed line, which locates at (0 m, 10 m) in the target coordinate system. The RP is highlighted by the white solid line in Figure 5a, locating at (0 m, 0 m). The ISAR image achieved by the RDA is shown in Figure 5b and the imaging result achieved by PFA is given in Figure 5d. Compared to the image in Figure 4d, the blur in Figure 5b was not fully overcome by PFA.

Figure 5
figure 5

Analysis of the PFA algorithm. (a) HRRPs of the simulated aircraft. (b) Image by RDA. (c) Image after Keystone transform. (d) Image by PFA.

3.2 Analysis and simulation of the disturbance system in PFA

According to the simulation above, the blur in the cross-range dimension of the image by PFA occurs when the RC is deviated from the RP. In this subsection, the influence of the disturbed phase term K R Δr on PFA will be analyzed in detail, by adopting the equivalence of a disturbance system. From Equation (10), since:

K x = K R sin θ K y = K R cos θ ,
(20)

the echo phase term of p k with respect to (K x , K y ) can be formulated as:

Φ k K x , K y = K x 2 + K y 2 Δr+ K x x k + K y y k .
(21)

Supposing that the PFA interpolation has been carried on the echo in Equation (18), the echo s R 4 _ k K x , K y of (K x , K y ) in a Cartesian coordinate system can be obtained. Taking 2D FFT on the s R 4 _ k K x , K y , one can obtain the 2D spectrum of scattering point p k :

S R 4 _ k X , Y = s R 4 _ k K x , K y e - j 2 π X K x e - j 2 π Y K y d K x d K y = σ k rect f B exp j Φ k K x , K y e - j 2 π X K x e - j 2 π Y K y d K x d K y = σ k rect f B exp j K x 2 + K y 2 Δ r + K x x k + K y y k e - j 2 π X K x e - j 2 π Y K y d K x d K y
(22)

Obviously, when Δr = 0, the 2D spectrum in Equation (22) has the form of the 2D Sinc function:

S R 4 _ k X , Y Δ r = 0 =sinc 2 Θ λ X - x k , 2 B c Y - y k .
(23)

If Δr ≠ 0, the 2D spectrum of scattering point p k is given by:

S R 4 _ k X , Y = σ k rect f B exp j K x x k + K y y k exp j K y 2 + K x 2 Δ r e - j 2 π X K x e - j 2 π Y K y d K x d K y = s R 4 _ k K x , K y exp j K x 2 + K y 2 Δ r e - j 2 π X K x e - j 2 π Y K y d K x d K y .
(24)

Let ψ K x , K y =exp j K y 2 + K x 2 Δ r and its 2D FFT spectrum is Ψ(X, Y), then S R 4 _ k X , Y satisfies:

S R 4 _ k X , Y = s R 4 _ k K x , K y ψ K x , K y e - j 2 π X K x e - j 2 π Y K y d K x d K y = S R 4 _ k X , Y Δ r = 0 Ψ X , Y ,
(25)

where operation symbol  denotes the operator of the 2D convolution. The characteristics of Ψ(X, Y) are the key point of the analysis of the influence on S R 4 _ k X , Y . We define Ψ(X, Y) as:

Ψ X , Y = ψ K x , K y e - j 2 π X K x e - j 2 π Y K y d K x d K y = K y 1 , K y 2 K x 1 , K x 2 ψ K x , K y e - j 2 π X K x d K x e - j 2 π Y K y d K y = K y 1 , K y 2 K x 1 , K x 2 exp j K x 2 + K y 2 Δ r e - j 2 π X K x d K x e - j 2 π Y K y d K y .
(26)

From Equation (26), the spectrum of Ψ(X, Y) is independent on the location (x k , y k ) of the scattering point. As a result, when Δr ≠ 0, the degraded PFA is equivalent to the output of the standard PFA model passing through a disturbance system denoted by Ψ(X, Y). The schematic diagram of the disturbance system Ψ(X, Y) is shown in Figure 6.

Figure 6
figure 6

The schematic diagram of the disturbance system Ψ ( X , Y ).

According to Figure 3, the value zone of K x and K y after the PFA interpolation is:

K x K x 1 , K x 2 = - K R 1 tan Θ / 2 , K R 1 tan Θ / 2 K y K y 1 , K y 2 = K R 1 , K R 2 cos Θ / 2 .
(27)

Since ψ(K x , K y ) have the symmetrical format with respect to K x and K y , we only need to analyze the FFT distribution of either one. The formula g(x) with respect to x is given by:

g x =exp j Δ r a 2 + x 2 ,
(28)

where a and Δr are constant. When Δr > 0, the FFT of g(x) is given by:

G X = x 1 x 2 exp j Δ r a 2 + x 2 e - j 2 π X x d x .
(29)

Figure 7a,b shows the curves of the g(x) phase term under the conditions x2 > x1 0 and x2 = - x1, respectively. By comparing between Figure 7a and b, one can observe that the phase term is approximately linear to x when x2 > x1 0, while it has quadratic phase if x2 = -x1. Furthermore, the smaller x is, the more distinct the quadratic feature of g(x) phase term is. For a limited time signal with linear phase and rectangular envelope, its spectrum performs as a Sinc function [23], since the signal can be approximated by a point frequency signal. The signal with quadratic phase is somewhat like a LFM signal, and its spectrum will be widened around the center frequency. In a nutshell, the spectrum of G(X) is related with the integral region [x1, x2]. Simulations in Figure 8 are given to analyze of the disturbance system Ψ(X, Y). Take the signal model in Section 3.1 for example, the interval of K x and K y can be calculated as:

Figure 7
figure 7

The curve of g ( x ) under different conditions. (a) The curve of g(x) when x2 > x1 0. (b) The curve of g(x) when x2 = -x1.

Figure 8
figure 8

The simulation results of scattering point p k (13, 9). (a) The location of scattering point p k (13, 9). (b) The ISAR image of RDA. (c) The ISAR image under standard PFA model. (d) The 2D distribution of Ψ(X, Y). (e) The 2D spectrum of the FFT results on the K x dimension. (f) The 1D slice of the 2D spectrum in (e). (g) The 2D distribution of the FFT transform result on K y dimension. (h) The 1D slice of the 2D spectrum in (g). (i) The ISAR image after PFA interpolating with Δr ≠ 0.

K x - K R 1 tan Θ / 2 , K R 1 tan Θ / 2 = - 58.84 , 58.84 K y K R 1 , K R 2 cos Θ / 2 = 586.43 , 666.85 .
(30)

For generality, we take the scattering point p k at (x k , y k ) = (13, 9) for example, as shown in Figure 8a. By applying RD and PFA algorithms to the simulated standard PFA signal, the obtained ISAR images of p k are shown in Figure 8b,c, respectively. The blurring in range cells and Doppler cells caused by the MTRC in the RD model is compensated. The 2D spectrum of disturbance system Ψ(X, Y) is shown in Figure 8d. Figure 8e,f,g,h illustrates the one-dimensional (1D) distributions of Ψ(X, Y) along cross-range axes X and range axes Y. Since Ky 2 > Ky 1 0 and Kx 2 = -Kx 1, one can see from Figure 8d, e, f, g, h that the spectra in the cross range are much wider than the ones in the slant range. According to the definition and characteristics of 2D convolution [24, 25], the output S R 4 _ k X , Y Δ r = 0 of the ideal PFA image will be widened in the cross range caused by the disturbance system Ψ(X, Y). As a result, the final output is blurred in the cross range and focused in slant range, (Figure 8i).

Furthermore, the influence of Δr and Θ values on the disturbance system Ψ(X, Y) is analyzed. When Θ = 0.2 rad is fixed, the distributions of Ψ(X) and Ψ(Y) under different values of Δr are shown in Figure 9, where Δr {2, 6, 10, 14, 18} m. With the distance of RC and RP increasing, Ψ(X, Y) is widened in the cross range and shifted in the slant range. The distributions of Ψ(X) and Ψ(Y) under different values of rotation angle Θ is shown in Figure 10, where Θ {0.05, 0.1, 0.15, 0.2, 0.25, 0.3} rad and Δr = 10 m. Similarly to the results in Figure 9, Ψ(X, Y) is widened in the cross range and shifted in the slant range with Θ increasing. Particularly, when the rotation angle Θ is small enough, the spectrum of Ψ(X, Y) in the cross range trends to the impulse function. At this level, one can obtain the formulas

Figure 9
figure 9

The influence of the value of Δr on Ψ ( X , Y ). (a) Distributions of Ψ(X) under different values of Δr. (b) Distributions of Ψ(Y) under different values of Δr.

Figure 10
figure 10

The influence of the value of Θ on Ψ ( X , Y ). (a) Distributions of Ψ(X) under different values of Θ. (b) Distributions of Ψ(Y) under different values of Θ.

K x = K R sin θ K R θ K y = K R cos θ K R 1 - θ 2 / 2 K x 2 + K y 2 K y .
(31)

Then, the phase term of P k echo is approximately equal to:

Φ k θ = K R Δ r + x k sin θ + y k cos θ = K R Δ r + K x x k + K y y k K y Δ r + K x x k + K y y k .
(32)

Apparently, Δr only produces translation in the cross range without the blurring in the slant range. Actually, the signal model becomes the RD model when Θ is small.

From the analysis and simulations above, we can note that, because of the range measurement errors and the target's non-cooperative motion, the disturbed phase term remained after data-driven translational compensation and the blur in the cross-range dimension is induced. Furthermore, with the increasing of Δr and Θ, the blur in the cross-range dimension becomes serious. As a conclusion, the accurate rotation center is essentially significant for PFA imaging for ISAR targets and only accurate estimation of RC can overcome the blur after PFA.

4 PFA revision based on integrated cubic phase function

According to the theoretical analysis of the PFA algorithm, if the RC on the slant range is not overlapped with the RP, the blur in cross-range cells remains after. In order to solve this problem, the dechirped signals after the translational compensation are shifted in the slant range, and the RC is located at the center of the HRRP. As a result, the effect of the disturbance system caused by Δr can be eliminated.

It is well known that the RP is located at the center of the HRRP. Therefore, the location of RC in HRRP is r a t =Δr= r a t - r 0 t . When the estimation value Δ r ^ of Δr is obtained, s R 4 f , θ is shifted in the slant range, and the new dechirped signals sR 5(f, θ) can be written as:

s R 5 f , θ = s R 4 f , θ exp - j 4 π c f c + f Δ r ^ = k = 1 K σ k rect f B × exp j 4 π c f c + f Δ r - Δ r ^ + x k sin θ t m + y k cos θ t m .
(33)

If the value of Δ r ^ is equal to Δr, then Δr- Δ r ^ =0. The shifted signals in Equation (33) become:

s R 5 f , θ = s R 4 f , θ exp - j 4 π c f c + f Δ r ^ = k = 1 K σ k rect f B exp j 4 π c f c + f x k sin θ t m + y k cos θ t m .
(34)

From formula (34), the shifted signals have the same format as Equation (9). s R 4 f , θ is converted into a standard PFA formation. When the PFA algorithm is performed on the revised signal sR 5(f, θ), the MTRC of s R 4 f , θ can be compensated accurately.

Note that, in Equation (33), the RC should be estimated firstly. Herein, we propose a method based on ICPF to estimate the rotation motion information. The ICPF extends the CPF and provides improved estimation of multi-component chirp signals [26].

By applying the FFT to the dechirped signal, the HRRP is obtained.

S Y , t m = sinc 2 B c Y - y k × exp - j 4 π λ x k sin θ t m + y k cos θ t m
(35)

During the CPI, if the rotation angle of the target is not too large (e.g., Θ < 20°), cos θ(t) and sin θ(t) can be approximated by the second-order Taylor expansion. Only accounting the second-order phase terms, the signal expression can be rewritten as:

S Y , t m = sinc 2 B c Y - y k × exp - j 4 π λ y k + x k θ t m + y k θ t m 2 .
(36)

Based on the uniform rotation assumption, the rotation angle θ(t) can be expressed as θ(t) = ωt, where ω is the rotation velocity . Then, S(Y, t) can be given as:

S Y , t m = sin c 2 B c Y - y k × exp - j 4 π λ y k + x k ω t m + y k ω 2 t m 2 .
(37)

It can be noted that the signal in a certain range cell is a chirp signal in azimuth, with chirp rate Γ = 2y k ω2/λ. For the rotation center, the chirp rate of the signal equals to zero, due to y k  = 0. When RC is distant from RP, the location of RC in ISAR image can be marked as p RC (x RC , Δr). Therefore, RC can be estimated by searching the range cell with the minimum chirp rate. The cubic phase function is a 2D bilinear mapping operator to estimate the parameters of a chirp signal [14]. For a single-component cubic phase signal,

x t = b 0 e j ϕ t = b 0 exp j a 0 + a 1 t + a 2 t 2 + a 3 t 3 ,
(38)

where ϕ(t) is the expression of a cubic phase. The instantaneous chirp rate is given by:

ICR= d 2 ϕ t d t 2 =2 a 2 + 3 a 3 t .
(39)

The cubic phase function is defined as a 2D bilinear transform of the signal as:

C P F t , Ω = 0 + x t + τ x t - τ e - j Ω τ 2 d τ = b 0 2 exp j 2 a 0 + a 1 t + a 2 t 2 + a 3 t 3 0 + exp j 2 a 2 + a 3 t - Ω τ 2 d τ
(40)

where Ω represents the instantaneous chirp rate index for the spectrum of the cubic phase function. Obviously, CPF (t, Ω) reaches the peaks along the curve Ω = 2(a2 + 3a2t). Thus, the peaks of the CPF can be used for the ICR estimation. The instantaneous chirp rate estimate is:

I C ^ R= arg max Ω C P F t , Ω .
(41)

For a chirp signal, the energy of the CPF is concentrated along straight lines Ω = 2a2 in the t - Ω domain. The CPF has been verified to be highly efficient in analyzing the mono-component signal even at low signal-to-noise ratio. However, when the algorithm is applied to multi-component signals, the distinct cross terms and spurious peaks will occur. The auto terms of the CPF are distributed around the straight lines which are parallel to the time axis. The interference including cross terms and spurious peaks disperse over the time and frequency rate domain. For multi-component chirp signals, the number of cross terms is K2 - K, and that will merge into (K2 - K)/2 spurious peaks in the time and chirp rate domain, where K is the number of the signal components.

To enhance the auto terms and suppress the cross terms, ICPF is proposed [26, 27]. For either single- or multi-component linear frequency modulated signals, ICPF provides improved estimation accuracy and better performance in cross-term rejection than the standard CPF-based approach. In this paper, the ICPF is exploited to estimate the chirp rate of the multi-component linear frequency modulated signals. ICPF integrates the energy of auto terms, yielding distinct peaks in the time and frequency rate domain, by which, the chirp rate can be easily determined by finding the highest peak.

ICPF Ω = 0 + | C P F t , Ω | 2 d t
(42)
I C ^ R= arg max Ω ICPF Ω
(43)

The procedure of the RC estimation based on ICPF is summarized as follows. First, conventional translational motion compensation, including range alignment and the phase compensation [2831], should be performed. Second, the Keystone transform is used to correct the linear MRTC as much as possible. It should be emphasized that, in the wide-angle ISAR imaging, the Keystone usually cannot correct the distinctive high-order MRTCs, but it still plays a significant role for a better estimation of the chirp rate with ICPF. Then, the estimate of the quadratic phase for the multi-component linear frequency modulated signals in different range cells can be obtained by:

Γ ^ n = 1 2 I C ^ R
(44)

where n [1, N] denotes the serial number of the range cell and N is the number of range cells. Therefore, the position of RC can be estimated by:

Δr= arg min n Γ ^ n - N / 2 ρ r .
(45)

Once Γ ^ n is obtained, the rotation angle Θ during CPI can be calculated by:

Θ n = ω T p = 1 2 Γ n λ / y n T p y n = n - N / 2 ρ r - Δ r .
(46)

It should be noted that in some slant-range cells there are no presence of dominant scatterers, where significant estimate error would involve. To improve the accuracy of the estimated Θ, the least mean square error (LMSE) method can be applied for fitting a linear function of the chirp rate, and an improved Θ ^ is achieved.

Θ ^ = E LMSE Θ n
(47)

Figure 11 shows the estimated coefficient of quadratic phase terms along the slant range with the ICPF. By the proposed LMSE fitting, the RC is determined at the range gate 10.1 m and the estimated rotation angle is Θ ^ = 11.2 .As it is shown in Figure 12a, the RC was shifted to the center of the HRRPS according to the estimated RC in Figure 11. Taking PFA on the shifted signals, the compensated ISAR image is obtained in Figure 12b. Since the disturbance of the echoes in Section 3 is revised, the MTRC was fully compensated and the PFA image was well focused.

Figure 11
figure 11

Coefficient of quadratic phase estimated by ICPF method.

Figure 12
figure 12

The results of the improved PFA in Figure5. (a) HRRPs after range shift. (b) Image with improved PFA.

5 Experiment and performance analysis

In the first experiment, simulated data of aircraft (MIG-25) [1] were used. The aircraft is composed of 120 point scatterers with equal reflectivity. The radar transmits LFM signals with a carried frequency of 9 GHz and a bandwidth of 512 MHz. It was assumed that the aircraft only has rotational motion.

Figure 13a shows the image obtained by using RDA without MTRC and time-varying Doppler compensation, where the blurring effects are evident. The RC is located at -4.1 m, which was estimated by ICPF in Figure 13g and highlighted by a dashed line in Figure 13c. The estimated rotation angle is Θ ^ = 18.5 . Figure 13b presents the ISAR image after MTRC compensation with the Keystone transform, where the blur in the cross-range direction still exists. As shown in Figure 13e, since the RC was not located at the same position as the RP, the image obtained by straight forward application of PFA was blurred. With the proposed ICPF estimation of RC, the HRRPs of shifted signals are shown in Figure 13d and the image with PFA is shown in Figure 13f, where compensation of both rotation induced MTRC and time-varying Doppler is achieved in an optimal manner. And well-focused ISAR imagery of the aircraft is obtained. By these result, the effectiveness of the proposed RC estimation for PFA imaging is justified clearly.

Figure 13
figure 13

Simulated experiment with MIG-25 data. (a) RD image. (b) ISAR image after the Keystone algorithm. (c) Original HRRPs. (d) HRRPs after range shift. (e) ISAR image after PFA. (f) ISAR image of shifted signals after PFA. (g) Coefficient of quadratic phase estimated by ICPF method. (h) Entropy map of the ISAR image after MTRC compensation with different rotation angles and centers. (i) Compensated image by the method in ref. [20]. (j) Estimated values under different SNRs. (k) Estimated errors under different SNRs.

For comparison, another method using PFA presented in ref. [20] was used to compensate the simulated data. The method presented in [20] is a non-parametric algorithm via golden section search (GSS), which applies the PFA on the translational compensated data iteratively to search the optimization solution based on minimum entropy criterion. The accurate rotation angle and center are obtained when the entropy of the ISAR image is minimized. Meanwhile, the MTRC can be compensated and the ISAR image will be best focused. Figure 13h depicts the ISAR image entropy after MTRC compensation with different rotation angles and centers. The minimum entropy is achieved at the position (18.3°, -4.2 m, 6.24), where 18.3° and -4.2 m are the rotation angle and center, respectively. After 2D iteratively searching, the estimated rotation center and angle are -4.2 m and 18.3° and the compensated image is shown in Figure 13i. Comparing Figure 13f with Figure 13i, we can note that the two methods have similar compensational performance. However, in the aspect of computational efficiency, as the interpolation in PFA is time-costly, the method in ref. [20] is much less efficient than the proposed one. The consuming CPU time of the method in ref. [20] and proposed one is 35.39 s and 6.23 s, respectively. Additionally, we plot the estimated results and errors corresponding to different SNRs in Figure 13j, k, in order to show the robustness of the proposed method. The true values of the rotation angle and center are 18.5° and -4.1 m, respectively. The estimated error is denoted by the difference between the estimated value and the true value. From Figure 13j,k, we note that effective rotation parameter estimation is achievable even if SNR decreases down to -10 dB and the estimate accuracy degrades slightly with the increase of strong noise. Contrarily, the estimation by the method in ref. [20] is sensitive to the strong noise, leading to serious errors in the result.

In order to show the importance of precise RC for PFA imaging, we performed the following test. Figure 14 depicts the entropy and contrast of the ISAR images from polar format interpolating with different shift ranges. The more focused the image is, the smaller the entropy is and the larger the contrast is [14, 31, 32]. From Figure 14, the minimized entropy and the maximized contrast are achieved when the shift range is -4.1 m. It illustrates that, the closer the shift range and RC are, the more focused the compensated ISAR image is and consequently the more efficient the PFA is.

Figure 14
figure 14

The entropy and contrast of the image after PFA with different shifted ranges.

In the second experiment, the real data provided by ATR National Defense Science and Technology Key Laboratory were used to analyze the performance of the proposed method. The data were collected by a ground-based imaging radar. The target is a stably flying airplane of Boeing 737, whose model is shown in Figure 15a. The carrier frequency is 10 GHz and the bandwidth is 1 GHz, providing a range resolution of 0.15 m. In order to obtain high cross-range resolution matching with the range resolution, the rotation angle in CPI should exceed 5.7°. Thus 512 pulses are used in the dataset and Θ is about 8°.

Figure 15
figure 15

Real data experiment of Boeing 737. (a) Model of Boeing 737. (b) Aligned profiles. (c) HRRPs after range shift. (d) Image after RDA. (e) Image after Keystone algorithm. (f) ISAR image after PFA on the data in (b). (g) Shifted HRRPs after PFA interpolation. (h) ISAR image after improved PFA. (i) Region 1 in (d). (j) Region 1 in (e). (k) Region 1 in (f). (l) Region 2 in (d). (m) Region 2 in (e). (n) Region 1 in (h). (o) Region 2 in (h).

Herein, the translational compensation was firstly carried on the echoes. The aligned profiles are shown in Figure 15b. By using the ICPF algorithm, the estimated rotation angle is Θ ^ = 8.03 and the RC locates at range -15 m, which is highlighted in Figure 15b. As it is shown in Figure 15c, the RC was shifted to the center of the HRRPS according to the estimated RC. For comparison, the RDA, Keystone algorithm, and traditional PFA were carried on the data in Figure 15b, and the images are shown in Figure 15d, e, f. From Figure 15d, the blurring effects caused by MTRC in both of the range and Doppler cell can be found in the RD image. The image using the Keystone algorithm in Figure 15e illustrates that the MTRC in the slant-range dimension is compensated but the blur in the Doppler dimension still remains. As shown in Figure 15g, when one utilizes PFA on the dechirped signals directly, the image is blurred in the cross-range direction. The aligned profiles and 2D image after PFA interpolation on the shifted data are shown in Figure 15g, h. In order to be clearer, the detailed exhibitions of the local regions are given in Figure 15i, j, k, l, m, n, o. Comparing Figure 15i, j, k , l, m with Figure 15n, o, one can observe that more focused image was obtained by the proposed method.From Figure 15n, o, two dominant scatterers are chosen in range 29.3 m and 11.5 m. We plot the dominant scatterers corresponding to different algorithms respectively along the cross-range and range dimensions in Figure 16, in order to show the improvement on the point spread function (PSF) by the use of the improved PFA. From Figure 16, we note that the blur caused by MTRC both in range and cross-range dimensions was compensated by the proposed method.

Figure 16
figure 16

PSF images corresponding to different algorithms of two dominant scattering centers. (a) Images of scatterer 1 along cross-range dimension. (b) Images of scatterer 2 along range dimension.

What is more, we use image entropy and contrast value to evaluate the focused images in Figure 15, which are presented in Table 1. It is notable in Table 1 that a well focused ISAR image was obtained by using the improved PFA.

Table 1 Entropy and contrast evaluation

6 Conclusions

In this paper, the limitations of the PFA algorithm for high-resolution ISAR imaging are analyzed. The difference between RC and RP produces a disturbance system and results in the distribution of the 2D scatterer mismatch in the polar format. As a result, the image is blurred in the cross range from PFA. To solve this problem, RC estimation via ICPF is proposed before PFA interpolation. Then the dechirped signals after translational compensation are shifted along the slant range in order to make the rotation center be located at the center of the HRRP. Consequently, well-focused wide-angle ISAR images can be obtained with PFA. Experiments with both simulated and real data show that the proposed method provides a promising performance. However, this model restricts the rotational motion model, leading to limitations on specific scenarios, such as the space target and the civil aviaton aircraft moving stably. The complex motion of maneuvering target, such as battleplane and missile, could include 3D rotation and deviation. In this case, the MTRC compensation in high-resolution ISAR imaging is still open. Thereby, extending the analysis and compensation to the 3D plane or non-uniform rotation is a significant and continual work in the future.

References

  1. Chen CC, Andrews HC: Target motion induced radar imaging. IEEE Trans. Aerosp. Electron. Syst. 1980, 16: 2-14.

    Article  Google Scholar 

  2. Zhou F, Bai X, Xing M, Bao Z: Analysis of wide-angle radar imaging. IET Radar. Sonar. Navig. 2011, 5: 449-457. 10.1049/iet-rsn.2010.0076

    Article  Google Scholar 

  3. Xing M, Wu R, Lan J, Bao Z: Migration through resolution cell compensation in ISAR imaging. IEEE Geosci. Remote Sens. Lett. 2004, 1: 141-144. 10.1109/LGRS.2004.824766

    Article  Google Scholar 

  4. Wang Y, Jiang YC: A novel algorithm for estimating the rotation angle in ISAR imaging. IEEE Geosci. Remote Sens. Lett. 2008, 5: 608-609.

    Article  Google Scholar 

  5. Ruan H, Wu YH, Jia X: Novel ISAR imaging algorithm for maneuvering targets based on a modified keystone transform. IEEE Geosci. Remote Sens. Lett. 2014, 11: 128-132.

    Article  Google Scholar 

  6. Zheng JB, Zhu WT, Liu QH: ISAR imaging of targets with complex motions based on the Keystone Time-Chirp rate distribution. IEEE Geosci. Remote Sens. Lett. 2014, 11: 1275-1279.

    Article  Google Scholar 

  7. Carrara WG, Goodman RS: Spotlight Synthetic Aperture Radar Signal Processing Algorithms. Artech House Inc, MA, Boston; 1995:152-154.

    MATH  Google Scholar 

  8. LeCaillec JM: SAR remote sensing analysis of the sea surface by polynomial filtering. IEEE Trans. Signal Process. 2007, 24: 105-107.

    Article  Google Scholar 

  9. Yu T, Xing MD, Zhen B: The polar format imaging algorithm based on double chirp-Z transforms. IEEE Geosci. Remote Sens. Lett. 2008, 5: 610-614.

    Article  Google Scholar 

  10. Yang ZM, Sun GC, Wu YF, Xing MD: A new fast back projection algorithm using polar format algorithm. J. Electron. Inf. Technol. 2013, 36: 373-376.

    Google Scholar 

  11. Mensa DL: High Resolution Radar Imaging. Artech House Inc, Dedham, MA; 1981:87-89.

    Google Scholar 

  12. Wu HQ, Grenier D, Delisle GY: Translational motion compensation in ISAR image processing. IEEE Trans. Image. Proc. 1995, 4: 1561-1571. 10.1109/83.469937

    Article  Google Scholar 

  13. Zhang L, Sheng JL, Duan J: Translational motion compensation for ISAR imaging under low SNR by minimum entropy. EURASIP J. Adv. Signal. Process. 2013, 33: 1-17.

    Google Scholar 

  14. Martorella M, Berizzi F, Bruscoli S: Use of genetic algorithms for contrast and entropy optimization in ISAR autofocusing. EURASIP J. Adv. Signal. Process 2006, 87298: 11.

    Google Scholar 

  15. Stanković L, Thayaparan T, Popovi V, Djurović I, Daković M: Adaptive s-method for SAR/ISAR imaging. EURASIP J. Adv. Signal. Process 2008, 593216: 10.

    MATH  Google Scholar 

  16. Prodi F: ISAR cross-range scaling using a correlation based functional. In Radar Conference, 2008. RADAR '08. IEEE, Rome; 2008:1-6.

    Chapter  Google Scholar 

  17. Yeh M, Xu J, Peng YN, Wang XT: Cross-range scaling for ISAR based on image rotation correlation. IEEE Geosci. Remote Sens. Lett. 2009, 6: 597-601.

    Article  Google Scholar 

  18. Hu JM, Zhou W, Fu YW: Uniform rotational motion compensation for ISAR based on phase cancellation. IEEE Geosci. Remote Sens. Lett. 2011, 8: 636-641.

    Article  Google Scholar 

  19. Liu Y, Zou JW, Xu SY, Chen ZP: Nonparametric motion compensation technique for high-resolution ISAR imaging via golden section search. Prog Electromagnetics Res M. 2014, 36: 67-76.

    Article  Google Scholar 

  20. Liu Y, Li H, Li N, Xu SY, Chen ZP: MTRC compensation in high-resolution ISAR imaging via improved polar format algorithm. In 2014 Security & Defence Conference. Proceedings of SPIE, Amsterdam, Netherlands; 2014. doi:10.1117/12.2066089

    Google Scholar 

  21. de Witt JJ, Nel WAJ: Range Doppler dynamic range considerations for dechirp on receive radar. In Proceeding of 5th European Radar Conference. Amsterdam, Netherlands; 2008:136-139.

    Google Scholar 

  22. Zhu DY: A novel approach to residual video phase removal in spotlight SAR image formation. In Radar Systems, 2007 IET International Conference. IET, Edinburgh, UK; 2007:15-18.

    Google Scholar 

  23. Lathi BP: Signal Processing and Linear Systems. Oxford University Press, New York; 1998:256-265.

    Google Scholar 

  24. Ran XN, Ray Liu KJ: Fast algorithms for 2-D circular convolutionsand number theoretic transforms based on polynomial transforms over finite rings. IEEE Trans. Signal Process. 1995, 43: 569-578. 10.1109/78.370613

    Article  Google Scholar 

  25. Radhakrishnan C, Jenkins WK: The 2-D modulated disctrete Fourier transform for 2-D fast convolution and digital filtering. In Proc. Int. Symp. Circ. Syst. IEEE, Rio de Janeiro, Brazil; 2011:1508-1511.

    Google Scholar 

  26. Wang P, Li HB, Himed B: Parameter estimation of linear frequency-modulated signals using integrated cubic phase function. In 42nd Asilomar Conference on Signals, Systems and Computers. IEEE, Pacific Grove, California; 2008:487-491.

    Chapter  Google Scholar 

  27. O'Shea P: A new technique for estimating instantaneous frequency rate. IEEE Signal. Process. Lett. 2002, 9: 251-252.

    Article  Google Scholar 

  28. Wang JF, Liu XZ: Improved global range alignment for ISAR. IEEE Trans. Aerosp. Electron. Syst. 2007, 43: 1070-1075.

    Article  Google Scholar 

  29. Zhu DY, Wang L, Yu YS, Tao QN, Zhu ZD: Robust ISAR range alignment via minimizing the entropy of the average range profile. IEEE Geoscience. Remote. Sensing Lett. 2009, 6: 204-208.

    Article  Google Scholar 

  30. Berizzi F, Corsini G: Autofocusing of inverse synthetic aperture radar images using contrast optimization. IEEE Trans. Aerosp. Electron. Syst. 1996, 32: 1185-1191.

    Article  Google Scholar 

  31. Wang J, Liu X, Zhou Z: Minimum-entropy phase adjustment for ISAR. IEE Proc. Radar Sonar Navig. 2004, 151: 203-209. 10.1049/ip-rsn:20040692

    Article  Google Scholar 

  32. Wang JF, Kasilingam D: Global range alignment for ISAR. IEEE Trans. Aerosp. Electron. Syst. 2003, 39: 351-357. 10.1109/TAES.2003.1188917

    Article  Google Scholar 

Download references

Acknowledgements

This work was supported in part by the National Natural Science Foundation of China (No. 61002025). The authors would also like to thank the editors and the anonymous reviewers for their helpful comments and suggestions.

Author information

Authors and Affiliations

Authors

Corresponding author

Correspondence to Yang Liu.

Additional information

Competing interests

The authors declare that they have no competing interests.

Authors’ original submitted files for images

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (https://creativecommons.org/licenses/by/4.0), which permits use, duplication, adaptation, distribution, and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Liu, Y., Xu, S., Chen, Z. et al. MTRC compensation in high-resolution ISAR imaging via improved polar format algorithm based on ICPF. EURASIP J. Adv. Signal Process. 2014, 179 (2014). https://doi.org/10.1186/1687-6180-2014-179

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/1687-6180-2014-179

Keywords