SAR Signal Formation and Image Reconstruction of a Moving Sea Target

Abstract

Maritime application of synthetic aperture radar (SAR) technology for sea-target surveillance and imaging is considered in this study. A SAR scenario, including the kinematics of a SAR satellite and a ship moving on the sea, along with the geometry of the target, are analytically described. A linear frequency modulation (LFM) waveform is applied for the target’s illumination. Based on the target’s geometry, SAR and target kinematics and the LFM waveform, a SAR signal model is synthesized. It is proven that the process of signal formation is a transformation of the three-dimensional (3D) image into a two-dimensional (2D) signal, whereas the target’s 2D imaging is an inverse transformation of the 2D signal into the target’s 2D image. SAR signal components, linear Fourier terms and higher-order phase terms are analytically derived and discussed in detail. Moreover, it is proven that SAR image reconstruction is a motion-compensation procedure, i.e., it removes all phases induced by first- and higher-order motion. Based on the SAR signal analysis, an illustrative iterative image-reconstruction algorithm is derived. The quality of the imaging is evaluated by an entropy cost function. Simulation experiments are carried out to verify the correctness of the theoretical statements in respect of SAR signal formation and image reconstruction.

1. Introduction

Aperture synthesis is the registration of an object’s complex signature as a SAR complex signal. The radar systems using the principle of aperture synthesis are synthetic aperture radar (SAR) and inverse synthetic aperture radar (ISAR). The reflectivity of the object modulating the signal is a SAR/ISAR image and can be retrieved from the complex signals by applying adequate image-reconstruction techniques.

The trajectory of a moving target or a SAR system limited by the radar’s antenna pattern or observation time is the synthetic aperture. SAR/ISAR is a coherent instrument for imaging of objects, using highly informative wideband pulses to realize the range resolution and gradient of the Doppler frequency induced by mutual displacement of the observed object and the system of observation to achieve azimuth or cross-range resolution. By definition, SAR is a moving observation system with a stationary object, whereas ISAR is a stationary observation system with a moving object. If the SAR antenna pattern tracks the object’s mass center, the system can be regarded as inverse SAR based on the apparent object’s rotation with respect to its mass center.

SAR systems are used for remote sensing of electromagnetic waves from the Earth’s surface and objects on it. SAR provides high-resolution surveillance in all meteorological and all day conditions. Automatic ship target detection is a common problem in maritime surveillance. Traditional and deep learning detection methods have been proposed in [1,2].

Standard SAR signal processing and image-reconstruction techniques are applied for static observed regions and objects. In the case of moving objects, their locations in a SAR image are moved in a cross-range direction which causes defocusing. In [3], methods for moving objects’ detection, imaging and position estimation are discussed.

A high-resolution method for moving target imaging is presented in [4]. A mathematical model of a moving object in a SAR scene is proposed in [5], where the object exhibits second- and third-order motions.

A SAR super-resolution multiple-signal classification (MUSIC) algorithm is discussed in [6]. A focusing method for SAR ship target imaging applying an ISAR approach is presented in [7]. Based on the sparsity of ship targets in a SAR image, an image-reconstruction method for multiple ship targets is proposed in [8]. A regression approach combined with convolutional neural networks to estimate radial velocity and slant range components of moving objects is considered in [9]. A target-detection algorithm in SAR images by applying the fusion of multiscale superpixel segmentations is proposed in [10].

Due to the weakness of object reflection and limited pixel numbers, it is difficult to detect small ships with a classical SAR system. In contrast, due to its abundant target-scattering information, a polarimetry SAR has been utilized for small ship detection [11]. Algorithms for ship target detection and effective exclusion from background clutter in SAR images with the application of information instruments are considered in [12,13,14].

To realize better target detection in maritime contexts, a SAR system has been combined with an Automatic Identification System and data fusion has been performed [15].

In case the signal-to-noise ratio (SNR) is low, standard imaging algorithms are not effective for motion compensation. To cope with this constraint, an original approach based on the generalized Radon–Fourier transform and gradient descent optimization is proposed in [16,17]. An important processing step in ISAR ship imaging is the selection of a suitable time window. A selection algorithm using a time window to extract Doppler frequencies of multiple scattering points for SAR ship imaging is developed in [18].

Based on dynamic azimuthal constraints used in SAR image reconstruction, a SAR target classification method is developed in [19]. SAR systems generate images of two and three dimensions, which gives rise to problems of range alignment. An algorithm for coarse and precise alignment in SAR imaging is presented in [20]. Phase correction or autofocusing (known as translational motion compensation) is a significant step in ISAR imaging. A non-parametric autofocusing algorithm based on contrast maximization is proposed in [21]. Due to the target’s translational displacement during aperture synthesis, moving objects in SAR images are blurred. To achieve a focused SAR image, a motion-compensation algorithm is developed in [22]. Based on features of Fisher vectors, robust ship detection algorithms are developed in [23]. Passive SAR/ISAR are effective tools for target detection and recognition. An imaging algorithm for a multi-static passive ISAR system is presented in [24].

Based on the original vector geometrical description of a SAR scenario with a sea target, satellite SAR and LFM waveform, the goal of the present study is the analytical derivation and physical interpretation of all steps and components of SAR signal formation and SAR image reconstruction, proving their Fourier and phase structures to be of high order. In contrast to the spectral analysis of SAR signal formation well known in the ISAR/SAR literature, in the present work, SAR signal formation is considered in the spatial–temporal domain. Thus, the kinematics of the SAR carrier and sea target are depicted in a 3D coordinate system of observation, whereas the target’s geometry is drawn using a 3D separate coordinate system. The target is illuminated with highly informative wideband LFM waveforms. To illustrate the solution of the defined problem, a range-Doppler image-reconstruction algorithm is built which includes phase compensation of a high order, inverse Fourier transformation and entropy evaluation of image quality.

In accordance with the defined problem, the rest of the paper is organized as follows. Section 2 describes the SAR scenario with geometrical and kinematical equations. Section 3 illustrates a process of SAR signal formation in a 3D space. Section 4 presents the algorithm for SAR imaging combined with phase correction of a high order and image entropy minimization. Section 5 illustrates the results of the numerical experiment. Section 6 presents a discussion of the theoretical and experimental results. Section 7 contains some concluding remarks.

2. Analytical Description of the SAR Scenario

Assume that a SAR carrier is moving in a 3D space of observation defined by a coordinate system $Oxyz$ (Figure 1). A sea target depicted in a three-dimensional (3D) coordinate system $O’XYZ$ is moving rectilinearly. The distance vector $R_{Sijk}(p)$ measured from the point S where the SAR is placed to the object’s scattering point can be written as

$$R_{Sijk}(p)=R_{S0’}(p)+A·R_{ijk}$$

(1)

where

$R_{S0’}(p)=R_S(p)+R_{0’}(p)$ denotes the distance vector SAR target’s mass center, $R_{ijk}={\left[X_{ijk},Y_{ijk},Z_{ijk}\right]}^T$ denotes the ijk-th scattering point’s distance vector in the coordinate system $O’XYZ$, where $X_{ijk}=i(\Delta X)$, $Y_{ijk}=j(\Delta Y)$ and $Z_{ijk}=k(\Delta Z)$ denote scattering point coordinates and $\Delta X$, $\Delta Y$ and $\Delta Z$ denote the dimensions of the resolution cell.

$p=\overline{0,N-1}$ denotes the transmitted waveform’s index; N denotes the number of waveforms illuminating the object during the aperture synthesis.

$$R_S(p)=R_S(0)+V_S·p·T_p,$$

(2)

denotes the SAR distance vector.

$V_S$ denotes the SAR carrier vector velocity with guiding cosines $\cos {\mathsf{\alpha }}_S$, $\cos {\mathsf{\beta }}_S$ and $\cos {\mathsf{\gamma }}_S=\sqrt{1-{\cos }^2{\mathsf{\alpha }}_S-{\cos }^2{\mathsf{\beta }}_S}$, and a module $V_S$

$$R_{0’}(p)=R_{0’}(0)+V_{0’}·p·T_p$$

(3)

denotes the distance vector of the target’s mass center, while $V_{O’}$ denotes the velocity vector of the target’s mass center.

$V_{0’}$ denotes the object’s velocity vector with guiding cosines $\cos {\mathsf{\alpha }}_{0’}$, $\cos {\mathsf{\beta }}_{0’}$ and $\cos {\mathsf{\gamma }}_{0’}=\mathsf{\pm }\sqrt{1-{\cos }^2{\mathsf{\alpha }}_{0’}-{\cos }^2{\mathsf{\beta }}_{0’}}$ and a module $V_{O’}$.

$R_S(0)$ denotes the SAR distance vector measured at the moment $p=0$ (the moment of imaging).

$R_{0’}(0)$ denotes the distance vector of the object’s mass center measured at the moment $p=0$.

$A$ denotes the Euler rotation matrix in the 3D coordinate space, the elements of which are defined by the equations

$$\begin{array}{l}{A}_{11}=\cos \mathsf{\psi }\cos \mathsf{\varphi }-\sin \mathsf{\psi }\cos \mathsf{\theta }\sin \mathsf{\varphi };\\ A_{12}=-\cos \mathsf{\psi }\sin \mathsf{\varphi }-\sin \mathsf{\psi }\cos \mathsf{\theta }\cos \mathsf{\varphi };\\ A_{13}=\sin \mathsf{\psi }\sin \mathsf{\theta };\\ A_{21}=\sin \mathsf{\psi }\cos \mathsf{\varphi }+\cos \mathsf{\psi }\cos \mathsf{\theta }\sin \mathsf{\varphi };\\ A_{22}=-\sin \mathsf{\psi }\sin \mathsf{\varphi }+\cos \mathsf{\psi }\cos \mathsf{\theta }\cos \mathsf{\varphi };\\ A_{23}=-\cos \mathsf{\psi }\sin \mathsf{\theta };\\ A_{31}=\sin \mathsf{\theta }\sin \mathsf{\varphi };\\ A_{32}=\sin \mathsf{\theta }\cos \mathsf{\varphi };\\ A_{33}=\cos \mathsf{\theta }.\end{array}$$

(4)

The rotation angles $\mathsf{\psi },\mathsf{\theta }$, $\mathsf{\varphi }$ determine the target’s spatial orientation and are defined by kinematic characteristics of the SAR carrier and target.

3. SAR Signal Formation Based on the 3D Object’s Geometry

3.1. SAR LFM Waveform

The SAR illuminates the object by a sequence of LFM waveforms modeled by

$$S(t)={\displaystyle \sum _{p=0}^{N-1}rect\frac{t}{T}·exp\left\{-j\left[\mathsf{\omega }\left(t-p{T}_p\right)+b{(t-p{T}_p)}^2\right]\right\}}$$

(5)

where $rect\frac{t}{T}=\{\begin{array}{l}1\mathrm{if}0\le \frac{t}{T}<1\\ 0\mathrm{otherwise}\end{array}$ is the rectangular signal time limiting function, $\mathsf{\omega }=2\mathsf{\pi }c/\mathsf{\lambda }$ denotes the angular frequency, $c=3\times {10}^8$ m/s denotes the light velocity, $b=\mathsf{\pi }\Delta F/T$ denotes the rate of the linear frequency modulation, $T$ denotes the timewidth of the LFM waveform and $\Delta F$ denotes the bandwidth of the linear frequency modulation waveform providing the resolution on the range direction $\Delta R=c/2\Delta F$.

Consider $\tilde{t}=t-p·T_p$ as a slow time, t as a fast time. It can be written $t=\tilde{t}\mathrm{mod}{T}_p$. The slow time determines the structure of the SAR signal for the azimuth direction, whereas the fast time determines the structure of the SAR signal for the range direction. In other words, the slow time and fast time define a 2D structure of the SAR signal.

3.2. SAR Signal Formation

The 2D SAR signal reemitted from the surface of the 3D object and received by the SAR antenna at the p-th slow time moment measured for the azimuth direction and the k-th fast time moment measured for the range direction can be expressed as

$$S(k,p)={\displaystyle \sum _{k=0}^{K-1}{\displaystyle \sum _{j=0}^{J-1}{\displaystyle \sum _{i=0}^{I-1}rect\frac{t-t_{Sijk}}{T}·a_{ijk}·}}}\exp \left\{j\left[\mathsf{\omega }\left(t-t_{Sijk}\right)+b·{\left(t-t_{Sijk}\right)}^2\right]\right\}$$

(6)

where $a_{ijk}$ is the intensity of the ijk-th scattering point and $I,J$ and $K$ as upper bounds of the sums determine the number of discrete coordinates on the $X,Y$ and $Z$ coordinate axes, respectively.

$t_{Sijk}=\frac{2R_{Sijk}(p)}{c}$ denotes the time delay of the SAR signal reemitted by the ijk-th scattering point; $t=t_{Sijk\min }+(k-1)\Delta T$ denotes the discrete fast time.

$k=\overline{1,K+\widehat{L}}$ denotes the fast time or range sample’s index, K in this case denotes the total number of fast time or range samples, $\Delta T$ denotes the timewidth of the waveform’s sample or the range sample, $t_{ijk\min }=\frac{2R_{ijk\min }(p)}{c}$ denotes the minimum signal time delay from the object, $t_{Sijk\max }=\frac{2R_{Sijk\max }(p)}{c}$ denotes the maximum time delay of the signal from the object, $\widehat{L}=\lceil \frac{t_{Sijk\max }-t_{Sijk\min }}{\Delta T}\rceil$ denotes the time-dependent relative object’s dimension,

$$rect\frac{t-t_{Sijk}}{T}=\{\begin{array}{l}1\mathrm{if}0\le \frac{t-t_{Sijk}}{T}<1\\ 0\mathrm{otherwise}\end{array}$$

(7)

For each p, demodulation of the received SAR signal is performed. This involves multiplication of the received signal with a complex conjugated emitted LFM waveform $rect\frac{t}{T}·exp\left[j\left(\mathsf{\omega }t+b{t}^2\right)\right]$. Apply $t_{Sijk}=\frac{2R_{Sijk}(p)}{c}$, then the demodulated SAR signal can be expressed by

$$\widehat{S}(k,p)={\displaystyle \sum _{k=0}^{K-1}{\displaystyle \sum _{j=0}^{J-1}{\displaystyle \sum _{i=0}^{I-1}rect\frac{t-\frac{R_{Sijk}(p)}{c}}{T}·a_{ijk}\exp \left[-j\left(2{\mathsf{\omega }}_k\frac{R_{Sijk}(p)}{c}-b{\left(\frac{2R_{Sijk}(p)}{c}\right)}^2\right)\right]}}}$$

(8)

where ${\mathsf{\omega }}_k=\mathsf{\omega }+2b(k\Delta T)$ denotes the current angular frequency.

Expression (8) can be regarded as a spatial-time transform, and the 2D SAR signal can be interpreted as a 2D complex spatial spectrum for the 3D object image.

3.3. Analysis and Physical Interpretation of a SAR Signal’s Structure

Expression (8) is used to model a demodulated LFM SAR signal reflected by all scattering points from the entire 3D object’s space. The distance $R_{Sijk}(p)$ and the image $a_{ijk}$ in (8) are defined in a three-dimensional space, whereas the synthetic aperture radar signal is described in a 2D signal plane $(k,p)$. Hence, the multiplication of $a_{ijk}$ with an exponential term limited by rectangular function and summation on 3D discrete coordinates in Equation (8) performs a spatial transformation of the 3D object’s image function $a_{ijk}$ into a 2D SAR signal function $\widehat{S}(k,p)$. Analytically, the transformation from 3D to 2D coordinates can be revealed as well by Taylor expansion of the argument

$$\widehat{\mathsf{\Phi }}(k,p)=2{\mathsf{\omega }}_k\frac{R_{Sijk}(p)}{c}-b{\left(\frac{2R_{Sijk}(p)}{c}\right)}^2$$

(9)

in the projection exponential term for the azimuth and range directions, i.e., for slow time $\tilde{t}\approx p·T_p$ and fast time t at the moment of imaging, $p=0$ and $k=0$, i.e.,

$$\begin{array}{l}\widehat{\mathsf{\Phi }}(k,p)=a_0+a_1·(p{T}_p)+a_2·{(p{T}_p)}^2+\dots +a_m{(p{T}_p)}^m\\ +b_1·(k\Delta T)+b_2·{(k\Delta T)}^2+\dots +b_n{(k\Delta T)}^n,\end{array}$$

(10)

where $a_m=\frac{{\widehat{\mathsf{\Phi }}}^{(m)}(0,0)}{m!}$ denotes the m-th derivative of the phase term in respect of slow time measured for the cross-range (azimuth) direction for $p=0$, and $k=0$.

$b_n=\frac{{\widehat{\mathsf{\Phi }}}^{(n)}(0,0)}{n!}$ denotes the n-th derivative of the phase term in respect of fast time measured for the range direction for $p=0$, and $k=0$.

A term $a_0=\frac{{\widehat{\mathsf{\Phi }}}^{(0)}(0,0)}{0!}$ denotes the constant phase term at the moment of imaging which can be omitted.

Based on the structure of $\widehat{\mathsf{\Phi }}(k,p)$, the coefficients $b_n=0$ for $n\ge 2$. The coefficients $a_1$ and $b_1$ can be determined as follows.

In case m = 1,

$$a_1=\frac{{\widehat{\mathsf{\Phi }}}^{(1)}(0,0)}{1!}=\frac{d\widehat{\mathsf{\Phi }}(0,0)}{d\tilde{t}}=\left(2\mathsf{\pi }\frac{2}{\mathsf{\lambda }}-\frac{8·b}{c^2}{R}_{Sijk}(0)\right)·\frac{d{R}_{Sijk}(0)}{d\tilde{t}}$$

(11)

The term $\frac{8·b}{c^2}{R}_{Sijk}(0)$ is much less than $2\mathsf{\pi }\frac{2}{\mathsf{\lambda }}$ and can be omitted. Then,

$$\frac{d\widehat{\mathsf{\Phi }}(0,0)}{d\tilde{t}}=2\mathsf{\pi }\frac{2}{\mathsf{\lambda }}\frac{d{R}_{Sijk}(0)}{d\tilde{t}}=2\mathsf{\pi }\frac{2}{\mathsf{\lambda }}{v}_{Sijk}(0)=2\mathsf{\pi }·f_{Sijk}(0)$$

(12)

where $v_{Sijk}(0)=\frac{d{R}_{Sijk}(0)}{d\tilde{t}}$ denotes the scattering point’s radial velocity at the moment of imaging, $f_{Sijk}(0)=\frac{2}{\mathsf{\lambda }}{v}_{Sijk}(0)$ denotes the scattering point’s Doppler frequency at the moment of imaging.

The first-order term $a_1·(p{T}_p)$ in the Taylor expansion can be modified as $a_1·(p{T}_p)=2\mathsf{\pi }·f_{Sijk}(0)·p,T_p$. Define the time parameter $T_p$ as $T_p=\frac{\tilde{T}}{N}$, where $\tilde{T}$ is the total observation time, i.e., the coherent processing period. Define the Doppler resolution as $\Delta f_D=\frac{1}{\tilde{T}}$. Then, the linear term can be written as

$$a_1·(p{T}_p)=\frac{2\mathsf{\pi }}{N}·\frac{f_{Sijk}(0)}{\Delta f_D}·p$$

(13)

Denote the index of the scattering point’s Doppler frequency as $\widehat{p}=\frac{f_{Sijk}(0)}{\Delta f_D}$, and the linear term (13) can be rewritten as

$$a_1·(p{T}_p)=\frac{2\mathsf{\pi }}{N}·p·\widehat{p}$$

(14)

In case $n=1$, the coefficient $b_1$ of the linear term $b_1·(k\Delta T)$ is modified as

$$b_1=\frac{{\widehat{\mathsf{\Phi }}}^{(1)}(0,0)}{1!}=\frac{d\widehat{\mathsf{\Phi }}(0,0)}{dt}=2·b·\frac{2R_{Sijk}}{c}=2·\frac{\mathsf{\pi }·\Delta F}{T}·\frac{2R_{Sijk}}{c}$$

(15)

Assume the time resolution of the LFM emitted waveform equal to the time sample to be $\Delta T={\left(\Delta F\right)}^{-1}$, then

$$\frac{d\widehat{\mathsf{\Phi }}(0,0)}{dt}=2·\frac{\mathsf{\pi }·}{T}·\frac{2R_{Sijk}(0)}{c·\Delta T}$$

(16)

The linear term $b_1·(k\Delta T)$ can be modified as $b_1·(k\Delta T)=2·\frac{\mathsf{\pi }·\Delta F}{T}·\frac{2R_{Sijk}(0)}{c}(k·\Delta T)$, i.e.,

$$b_1·(k\Delta T)=\frac{2·\mathsf{\pi }·}{\frac{T}{\Delta T}}·\frac{2R_{Sijk}(0)}{c·\Delta T}·k$$

(17)

Denote the total number of the time samples as $K=\frac{T}{\Delta T}$, the range resolution as $\Delta R=\frac{c·\Delta T}{2}$, then the range index of the scattering point is defined as $\widehat{k}=\frac{R_{Sijk}(0)}{\Delta R}$ and the linear term $b_1·(k\Delta T)$ can be rewritten as

$$b_1·(k\Delta T)=\frac{2·\mathsf{\pi }}{K}·k·\widehat{k}$$

(18)

Hence, the 2D Taylor expansion of $\widehat{\mathsf{\Phi }}(k,p)$ can be rewritten as

$$\widehat{\mathsf{\Phi }}(k,p)=\mathsf{\Phi }(0,0)+\frac{2\mathsf{\pi }·p·\widehat{p}}{N}+\frac{2\mathsf{\pi }·k·\widehat{k}}{K}+\mathsf{\Phi }(k,p)$$

(19)

where $\frac{2\mathsf{\pi }·p·\widehat{p}}{N}$ and $\frac{2\mathsf{\pi }·k·\widehat{k}}{K}$ are Fourier terms, the linear phases for the cross-range and range directions, respectively, and $\mathsf{\Phi }(k,p)$ denotes the higher-order phase term approximated by a higher-order polynomial.

The 3D image function $a_{ijk}$ is transformed into a 2D image function defined by discrete coordinates $\widehat{p}$ and $\overleftarrow{k}$ and can be expressed as $a_{\widehat{p}\overleftarrow{k}}$. Equation (8) can be rewritten as

$$\widehat{S}(k,p)={\displaystyle \sum _{\widehat{k}=0}^{\widehat{K}-1}{\displaystyle \sum _{\widehat{p}=0}^{\widehat{P}-1}{a}_{\widehat{p}\ddot{k}}rect\frac{t-t_{\widehat{p}\ddot{k}}(p)}{T_0}}}·\exp \left[-j\left(\frac{2\mathsf{\pi }·p·\widehat{p}}{N}+\frac{2\mathsf{\pi }·k·\widehat{k}}{K}+\mathsf{\Phi }(k,p)\right)\right]$$

(20)

for each $p=\overline{0,N-1}$, and $k=\overline{1,K}$.

In matrix form, Equation (20) can be written as $S=P·A·K^T{\mathsf{\Phi }}^T$, where $S=S[K\times N]$ is the SAR signal matrix, $P=P[N\times \widehat{P}]$ is the Fourier transformation matrix for the cross-range direction, $K=K[N\times \widehat{K}]$ is the Fourier transformation matrix for the range direction and $\mathsf{\Phi }=\mathsf{\Phi }[N\times K]$ is the higher-order phase matrix.

The 2D SAR signal can be regarded as a direct projection of the 2D image A onto the 2D SAR signal plane implemented by 2D Fourier transformation of the image A multiplied with a higher-order phase matrix $\mathsf{\Phi }$.

4. SAR Image Reconstruction and Signal Phase Correction-Based Entropy Minimization

4.1. 2D SAR Image Reconstruction

Based on Equation (20), an image-reconstruction procedure can be implemented by a 2D inverse Fourier transformation of the phase-corrected demodulated SAR signal, i.e.,

$$a_{\widehat{k}\widehat{p}}=\frac{1}{N·K}{\displaystyle \sum _{p=0}^{N-1}{\displaystyle \sum _{k=0}^{K-1}\widehat{S}(k,p)}}·\exp \left[j\left(\frac{2\mathsf{\pi }·p·\widehat{p}}{N}+\frac{2\mathsf{\pi }·k·\widehat{k}}{K}+\mathsf{\Phi }(k,p)\right)\right]$$

(21)

for each $\widehat{p}=\overline{0,N-1}$, and $\widehat{k}=\overline{0,K-1}$.

Equation (20) performs an inverse transformation of the two-dimensional SAR signal into a two-dimensional SAR image and in a matrix form is presented as

$$A=P^T·S·\mathsf{\Phi }K$$

(22)

The higher-order phase-correction matrix is unknown and at the beginning stage of the SAR image reconstruction this matrix is accepted as identity one. Hence, the demodulated SAR signal first undergoes a two-dimensional Fourier transformation, i.e.,

$$a_{\widehat{k}\widehat{p}}=\frac{1}{N·K}{\displaystyle \sum _{p=0}^{N-1}{\displaystyle \sum _{k=0}^{K-1}\widehat{S}(k,p)}}·\exp \left[j\left(\frac{2\mathsf{\pi }·p·\widehat{p}}{N}+\frac{2\mathsf{\pi }·k·\widehat{k}}{K}\right)\right]$$

(23)

for each $\widehat{p}=\overline{0,N-1}$, and $\widehat{k}=\overline{0,K-1}$.

Equation (23) is regarded as a correlation procedure searching for coordinates $\widehat{k}$ and $\widehat{p}$ and maximizing intensities of scattering points $a_{\widehat{k}\widehat{p}}$ in the SAR imaging plane, with a compensation of first-order motion removing the phases defined by scattering points’ radial velocities, except for phases at the moment of imaging that define the complex nature of the SAR image.

4.2. Higher-Order Phase Correction and Image Entropy Minimization

In case the SAR image retrieved by (23) is not focused, a higher-order phase correction of $\widehat{S}(k,p)$ is needed. Taking into account that for $n\ge 2$ the coefficients $b_n=0$, the higher-order term of the Taylor expansion is presented as a one-dimensional matrix function $\mathsf{\Phi }(p)$. This means that each row of the $\widehat{S}(k,p)$ matrix is to be multiplied with a higher-order phase correction term $\exp \left[j\mathsf{\Phi }(p)\right]$, i.e.,

$$\tilde{S}(k,p)=\widehat{S}(k,p)·\exp \left[j\mathsf{\Phi }(p)\right]$$

(24)

where

$$\mathsf{\Phi }(p)=a_2·{(p{T}_p)}^2+\dots +a_m{(p{T}_p)}^m$$

(25)

Coefficients $a_m$ are computed iteratively. First, $a_2$ is computed until the image entropy function has achieved a minimal value. Second, $a_3$ is calculated, and so on. The coefficient $a_2$ most influences the phase correction of a higher order.

Let ${\mathsf{\Phi }}_s(p)$ be the phase-correction function, defined at the s-th stage.

Implement the higher-order phase correction

$${\tilde{S}}_s(k,p)=\widehat{S}(k,p)\exp \left[j{\mathsf{\Phi }}_s(p)\right]$$

(26)

Implement an image reconstruction by range and azimuth compression

$$a_{\widehat{k}\widehat{p}}=\frac{1}{N·K}{\displaystyle \sum _{p=0}^{N-1}{\displaystyle \sum _{k=0}^{K-1}{\tilde{S}}_s(k,p)}}·\exp \left[j\left(\frac{2\mathsf{\pi }·p·\widehat{p}}{N}+\frac{2\mathsf{\pi }·k·\widehat{k}}{K}\right)\right]$$

(27)

Calculate the normalized intensities of the image

$$I_s(\widehat{k},\widehat{p})=\frac{{\left|a_{\widehat{k}\widehat{p},s}\right|}^2}{{\displaystyle \sum _{p=0}^{N-1}{\displaystyle \sum _{k=0}^{K-1}{\left|a_{\widehat{k}\widehat{p},s}\right|}^2}}}$$

(28)

Calculate the entropy function of the image with normalized intensities

$$H_s=-{\displaystyle \sum _{p=0}^{N-1}{\displaystyle \sum _{k=0}^{K-1}{I}_s(\widehat{k},\widehat{p})\ln [I_s(\widehat{k},\widehat{p})]}}$$

(29)

Increase $s=s+1$. If $H_{s+1}<H_s$, continue the calculation. The estimate of the polynomial coefficient $a_2$ corresponds to the minimum of the minimal value of the image entropy

$$a_2=\arg \underset{\mathsf{\Phi }}{\min }\left\{H\right[I_s(\widehat{k},\widehat{p})\left]\right\}$$

(30)

The phase-correction procedure lasts until the global minimum of the image entropy $H_s$ is achieved.

Based on Expression (21), the following image-reconstruction algorithm can be derived.

4.3. SAR Imaging Algorithm with Higher-Order Motion Compensation Based on Entropy Minimization

Let $\widehat{S}(k,p)$ be the demodulated 2D SAR signal obtained during aperture synthesis and registered in a 2D signal matrix with dimension $(K+\widehat{L})\times N$.

• Initialization of the second order coefficient $a_2=0$.
• Image reconstruction by a two-dimensional Fourier transformation for the range and azimuth (cross-range) directions

  $${\widehat{a}}_{\widehat{k}\widehat{p}}=\frac{1}{N·K}{\displaystyle \sum _{p=0}^{N-1}{\displaystyle \sum _{k=0}^{K-1}\tilde{S}(k,p)\exp \left[j\left(\frac{2\mathsf{\pi }·p·\widehat{p}}{N}+\frac{2\mathsf{\pi }·k·\widehat{k}}{K}\right)\right]}},$$

Range compression of the SAR signal by the first Fourier transformation

$$S(\widehat{k},p)=\frac{1}{K}{\displaystyle \sum _{k=0}^{K-1}\tilde{S}(k,p)\exp \left[j\left(\frac{2\mathsf{\pi }·k·\widehat{k}}{K}\right)\right]}.$$

The azimuth compression of the range-compressed SAR signal by the second Fourier transformation

$${\widehat{a}}_{\widehat{k}\widehat{p}}=\frac{1}{N}{\displaystyle \sum _{p=0}^{N-1}S(\widehat{k},p)\exp \left[j\left(\frac{2\mathsf{\pi }·p·\widehat{p}}{N}\right)\right]}.$$

3.

Normalization of the SAR image amplitude needed for calculation of the entropy cost function

$$I_{\widehat{k},\widehat{p}}=\frac{{\left|{\widehat{a}}_{\widehat{k},\widehat{p}}\right|}^2}{{\displaystyle \sum _{\widehat{p}=0}^{\widehat{P}-1}{\displaystyle \sum _{\widehat{k}=0}^{\widehat{K}-1}{\left|{\widehat{a}}_{\widehat{k},\widehat{p}}\right|}^2}}}$$

4.

Calculation of the entropy cost function

$$H=-{\displaystyle \sum _{\widehat{p}=0}^{\widehat{P}-1}{\displaystyle \sum _{\widehat{k}=0}^{\widehat{K}-1}{I}_{\widehat{k},\widehat{p}}\ln (I_{\widehat{k},\widehat{p}})}}.$$

5.

If H accepts a minimal value, stop the procedure.

6.

Otherwise, increment the coefficients $a_2=a_2+1$ and calculate $\mathsf{\Phi }(p)$.

7.

Apply a higher-order phase correction of the SAR signal by multiplication

$$\tilde{S}(k,p)·\exp \left[j\mathsf{\Phi }(p)\right],$$

(31)

and go to step 3.

5. Simulation Experiment

To prove the precision of the 3D SAR scenario’s analytical and geometrical description, signal modeling and the developed SAR imaging algorithm, including range compression by a Fourier transformation, azimuth compression by a Fourier transform, higher-order phase correction based on image entropy minimization, the results of a numerical experiment are provided. Assume that a sea object is moving rectilinearly with vertical fluctuations in a Cartesian coordinate system $Oxyz$. The object’s geometry is described in a Cartesian coordinate system $O’XYZ$ with a 3D resolution cell’s dimensions $\Delta X=\Delta Y=\Delta Z=1$ m (Figure 2). The numbers of the reference points on the axes $O’X$, $O’Y$ and $O’Z$ are as follows: $I=100$, $J=100$ and $K=10$. The mathematical expectation of the normalized intensities of the scattering points placed on the ship target is $a_{ijk}=0.01$.

The reflection from the water surface superimposes correlated noise on the SAR signal. Assume the noise is whited and removed by application of a standard Gaussian filter. This is the reason why, in the experiment, the reflection of the water surface can be ignored.

The ship target’s trajectory parameters are as follows:

• Vector velocity modulus $V=20$ m/s;
• Angle between ship’s vector velocity and $Ox$ axis ${\mathsf{\alpha }}_{0’}=\mathsf{\pi }/4$;
• Angle between ship’s vector velocity and $Oy$ axis ${\mathsf{\beta }}_{0’}=\mathsf{\pi }/4$;
• Angle between ship’s vector velocity and $Oz$ axis ${\mathsf{\gamma }}_{0’}=\mathsf{\pi }/2$;
• Vertical movement of the ship along the $Oz$ axis with angular velocity 2 rad/s and an amplitude of 2 m.

The SAR carrier’s trajectory parameters are as follows:

• Modulus of the SAR carrier’s vector velocity $V=1500$ m/s;
• Angle between SAR carrier vector velocity and $Ox$ axis ${\mathsf{\alpha }}_S=\mathsf{\pi }/2$;
• Angle between SAR carrier vector velocity and $Oy$ axis ${\mathsf{\beta }}_S=0$;
• Angle between SAR carrier vector velocity and $Oz$ axis ${\mathsf{\gamma }}_S=\mathsf{\pi }/2$;
• SAR carrier’s initial coordinates $x_S(0)=40$ km, $y_S(0)=30$ km, $z_S(0)=100$ km;
• Ship target’s initial coordinates $x_{0’}(0)=10$ km, $y_{0’}(0)=50$ km, $z_{0’}(0)=0$ km;
• SAR LFM waveform’s parameters:
• Wavelength is $\mathsf{\lambda }=3$ cm;
• LFM pulse waveform’s repetition period $T_p=0.005$ s;
• LFM pulse waveform’s timewidth $T={10}^{-6}$ s;
• LFM pulse waveforms’ number during aperture synthesis $N=256$;
• LFM signal bandwidth $\Delta F=150$ MHz;
• LFM rate $b=4.71\times {10}^{14}$;
• Number of samples of LFM transmitted signal and range samples $K=256$;
• Carrier frequency $f=10$ GHz;
• Dimension of the LFM pulse sample $\Delta T=T/K=7.8$ μs.

The quadrature complex components of the complex SAR signal, the real part (a) and imaginary part (b), are presented in Figure 3.

The relative time signal length is greater than 256 due to the relative timewidth of the target. The structure of the complex SAR signal can be interpreted as a 2D complex spatial spectrum of a 3D image projected onto a 2D SAR imaging plane.

Multiplication of the SAR signal with a complex conjugate LFM waveform generates a demodulated SAR signal. The quadrature components of the demodulated complex SAR signal, the real part (a) and the imaginary part (b), are presented in Figure 4.

The quadrature components, the real part (a) and imaginary part (b), of the range-compressed SAR signal by applying an inverse Fourier transform on range direction are presented in Figure 5. Displacement of the spectrum in the upper part of the picture is due to the fact that a frequency shift operation after the first Fourier transformation was not applied.

The quadrature components, the real part (a) and the imaginary part (b), of the azimuth-compressed SAR signal derived by applying an inverse Fourier transform for the cross-range (azimuth) direction after range compression of the SAR signal, are presented in Figure 6.

The real and imaginary parts of the azimuth-compressed SAR have wavelike interferential properties and define the complex character of the SAR image.

The SAR images of the ship target before (a) and after (b) application of the higher-order phase correction, presented as normalized intensities of the scattering points, are depicted in Figure 7. An optimal higher-order phase correction is achieved for a coefficient $a_2=10$. For this coefficient, the image entropy accepts a minimal value H = 7.97.

In the experiment, the phase correction polynomial of a higher order is limited with a second-order term, i.e., $\mathsf{\Phi }(p)=a_2·{(p{T}_p)}^2$. The variation of the image entropy function in accordance with a polynomial coefficient $a_2$ during the higher-order phase correction is illustrated in

Table 1. Image entropy H and a corresponding polynomial coefficient $a_2$.

a2	0	5	10	15	20
H	8.01	7.98	7.97	7.99	8.02

. Only a limited number of entropy values around the minimal entropy that differ in the second decimal place are presented in Table 1 to evaluate the dynamics of higher-order phase correction.

The small differences in the entropy values are due to a comparatively clear target image after range compression and azimuth compression of the complex SAR signal. The effect of the phase-correction procedure can be noticed in the stern and an upper and lower board of the ship. In Figure 7a, the unfocused SAR image of the ship is presented, whereas, in Figure 7b, the focused SAR image of the ship is depicted.

6. Discussion

In the present work, SAR kinematics, geometry and the process of signal formation and image extraction are analytically described and physically interpreted. The SAR scenario, including a moving SAR carrier and ship target, is depicted in a 3D coordinate system. The target is described in a 3D coordinate grid, in nodes of which scattering points are placed in accordance with the target’s geometry. Kinematic equations defining the positions of the SAR carrier and the target are derived. Based on the LFM waveform and 3D target geometry, a 2D SAR signal model is synthesized. The process of the signal formation is regarded as a spatial transformation of the 3D image into a 2D SAR signal. Through 2D Taylor expansion of the projection operator’s argument, i.e., the argument of the exponential term in Expression (8), the 3D coordinates of the target’s scattering points are transformed into 2D coordinates in the image plane. Dimensionless cross-range and range coordinates of scattering points are analytically derived and physically interpreted. An original interpretation of SAR signal formation as a 2D Fourier transformation of the image function with higher-order phase correction is given. The image reconstruction is an inverse projection procedure realized by a 2D inverse Fourier transformation of a higher-order phase-corrected SAR signal in the image plane defined by unknown dimensionless range and cross-range coordinates. The inverse Fourier transformation for the cross-range direction is a correlation procedure searching for and detecting unknown coordinates of the target’s scattering points. Moreover, the inverse Fourier transformation performs a motion compensation of the first order. It compensates linear phases induced by the radial displacement of scattering points during aperture synthesis, except phases in point of imaging that define the complex properties of the SAR image. The results of the simulation experiment illustrate the capabilities of the geometric, kinematic, SAR signal models, the imaging algorithm including a 2D Fourier transformation and the higher-order motion compensation implemented by a phase correction based on the minimization of image entropy. The simulation experiment was carried out with data close to real data. The ship target exhibits rectilinear movement with a linear vector velocity and vertical fluctuation (displacement), with an angular velocity of 2 rad/s and an amplitude of 2 m.

Different to the phase gradient autofocusing methods that require dominant scattering points in the observed scene, phase retrieval and nonparametric autofocusing based on contrast maximization [25,26,27,28], in the present work, the higher-order phase correction of the SAR signal is performed iteratively by applying entropy as an image quality evaluation function, avoiding measurement of dominant scattering points’ phase patterns. The suggested iterative image-focusing algorithm runs at high speed and distinguishes with fast convergence.

Due to the comparably clear image of the target after 2D inverse Fourier transformation of a demodulated signal, the variation of the polynomial coefficient of second order and respective image entropy is in a narrow interval.

7. Conclusions

In the present study, an original geometrical and analytical approach has been applied in the analysis and synthesis of a SAR signal and the structure of the imaging algorithm in a SAR scenario with a moving sea target. The main contributions are the mathematical modeling and the physical interpretation of the SAR signal components and SAR signal formation, the kinematic SAR and the geometric target parameters, along with the creation of the logical chain of target image extraction. An iterative image-reconstruction algorithm has been created based on the use of entropy as an image-quality evaluation function. The simulation experimental results prove the correctness of the theoretical statements, mathematical models and analytical expressions.