Next Article in Journal
Survey Assessment for Decision Support Using Self-Organizing Maps Profile Characterization with an Odds and Cluster Heat Map: Application to Children’s Perception of Urban School Environments
Next Article in Special Issue
Land-Cover Classification Using MaxEnt: Can We Trust in Model Quality Metrics for Estimating Classification Accuracy?
Previous Article in Journal
Fractional Refined Composite Multiscale Fuzzy Entropy of International Stock Indices
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Analytic Expressions for Radar Sea Clutter WSSUS Scattering Functions †

Applied Technology, Inc., King George, VA 22485, USA
This paper is an extended version of our paper published in the 2018 IEEE Radar Conference.
Current address: Applied Technology, Inc., 5200 Potomac Dr., King George, VA 22485, USA.
Entropy 2019, 21(9), 915; https://doi.org/10.3390/e21090915
Submission received: 28 August 2019 / Revised: 16 September 2019 / Accepted: 18 September 2019 / Published: 19 September 2019
(This article belongs to the Special Issue Maximum Entropy and Its Application III)

Abstract

:
Bello’s stochastic linear time-varying system theory has been widely used in the wireless communications literature to characterize multipath fading channel statistics. In the context of radar backscatter, this formulation allows for statistical characterization of distributed radar targets in range and Doppler using wide-sense stationary uncorrelated scattering (WSSUS) models. WSSUS models separate the channel from the effect of the waveform and receive filter, making it an ideal formulation for waveform design problems. Of particular interest in the radar waveform design community is the ability to suppress unwanted backscatter from the earth’s surface, known as clutter. Various methods for estimating WSSUS system functions have been studied in the literature, but to date no analytic expressions for radar surface clutter range-Doppler scattering functions exist. In this work we derive a frequency-selective generalization of the Jakes Doppler spectrum model, which is widely used in the wireless communications literature, adapt it for use in radar problems, and show how the maximum entropy method can be used to extend this model to account for internal clutter motion. Validation of the spectral and stationarity properties of the proposed model against a subset of the Australian Ingara sea clutter database is performed, and good agreement is shown.

1. Introduction

Random linear time-varying (LTV) system theory was first comprehensively described by Bello [2] and has widely been used in the wireless communications field ever since, particularly to model the multipath fading of mobile radio channels [3,4]. In particular, wide-sense stationary uncorrelated scattering (WSSUS) models, which are a subset of the category of random LTV systems, are the most common in the literature. As well, more recent work in this field has focused on nonstationary LTV communications channels [5,6,7], of which vehicle-to-vehicle (V2V) communication is a prime example [8]. A separate set of models that are distinct, but can be related to Bello’s LTV theory are the Clarke/Jakes Doppler spectrum class of models [9,10], originally only applicable for flat fading (i.e., where the symbol time is much larger than the multipath delay spread). This model is also ubiquitous in the wireless communications literature, and many extensions have been proposed, such as for varying geometries [11] and more accurate fading statistics [12]. It is also commonly coupled with Bello’s LTV theory (which is general enough to model frequency-selective fading where the multipath delay spread is much larger than the symbol duration) to model the Doppler component of multipath scattering [13]. Empirical [14] and analytic [15] expressions characterizing range-Doppler spreading of communication channels using LTV theory are commonplace.
Outside of communications applications, random LTV system models have been primarily used for sonar target detection and waveform design problems [16,17,18], although a few recent papers have applied this method to radar target detection problems as well [19,20]. Much of the literature on the subject treats the problem agnostically, treating radar and sonar as the same problem [21,22,23,24,25,26]. In both radar and sonar, the goal is to detect (usually) man-made targets of interest while suppressing thermal noise, which is internal to the receiver, and clutter, which are unwanted signal-dependent returns from the natural environment. Waveform design problems expressed in this form require characterization of the delay-Doppler “scattering functions” of the target and clutter, and the result is an optimized transmit waveform/receive filter pair.
In addition to waveform design problems, clutter Doppler spectrum characterization is useful for optimized moving target indication (MTI) and space-time adaptive processing (STAP) filtering, which optimize the receiver only. Specifically, the range-Doppler scattering function can be used to predict pulse-to-pulse correlations in the clutter return to suppress it using eigenfilter techniques [27,28]. An up-to-date summary of recent work on Doppler spectrum modeling can be found in [29]. Most papers on Doppler spectrum modeling do not use Bello’s LTV formulation; the effect of the transmit waveform, measurement system, and clutter are usually combined. This is appropriate if the goal is to generate realistic clutter samples for simulation purposes or if the goal is to design an MTI/STAP filter because the processing is done on receive. However, waveform design requires the clutter to be partitioned separately from the transmit signal, which creates a motivation to “translate” these Doppler spectrum results to the “language” of random LTV system theory so that these Doppler spectrum models can be utilized in other problem domains.
Because scattering function estimation is a fundamental component of waveform design problems, it is a topic that has been studied for decades [30,31,32,33,34,35,36,37,38,39,40,41]. However, analytic expressions for scattering functions of sea clutter for airborne radars have not been reported in the literature except for some initial work by the author [1,42,43]. In this work we derive such analytic expressions for the case when localized internal clutter motion (ICM) is small relative to radar platform motion using a frequency-selective version of the Jakes model. Extensions to this model are then proposed to connect previous work on Doppler spectrum modeling to random LTV theory. As well, we demonstrate how an estimation of spatially-local ICM spectra can be cast as a probability density estimation problem, for which solutions can be found using the Jaynes maximum entropy (MaxEnt) method [44] and directional statistics [45]. The spectrum prediction method is then validated against a subset of the Australian Ingara medium grazing angle clutter dataset [46], and good agreement is shown.

2. Results

2.1. Mathematical Background

For a monostatic radar of arbitrary polarization, the baseband backscattered signal y ( t ) will be modeled as the response of an LTV system h ( τ , t ) to a transmitted signal x ( t ) :
y ( t ) = h ( τ , t ) x ( t τ ) d τ + w ( t ) ,
where τ represents downrange delay (i.e., “fast time”), t is absolute time (“slow time”), and w ( t ) is thermal noise. In the literature, h ( τ , t ) is commonly referred to as the time-varying impulse response. In this work, h will include the transceiver, antenna, atmospheric propagation, and backscattering environment, so x ( t ) is essentially the output of the transmit D/A, and y ( t ) is the input to the receive A/D, prior to matched filtering. Polarization considerations are contained in the antenna pattern and surface scattering radar cross-section (RCS). Variations in h ( τ , t ) with respect to t represent channel-induced Doppler shifts. In the case where the entire scenario is static, h ( τ , t ) = h ( τ ) (i.e., h is LTI), and the integral in Equation (1) reduces to a convolution. In applying the LTV model to a radar system, we use the following assumptions:
  • The LTV impulse response h ( τ , t ) is only valid over a single coherent processing interval (CPI) of N transmitted pulses
  • Over a single CPI, the range to target(s) is roughly constant
  • Relative motion produces a Doppler shift on the carrier only and does not introduce time dilation of the pulse. The conditions for this assumption to be valid are that B T < 0.1 c / r ˙ , where B T is the time-bandwidth product of the waveform and r ˙ is the maximum range rate [47] (p. 382).
We will model h ( τ , t ) as a wide-sense stationary uncorrelated scattering (WSSUS) process characterized by a set of related correlation functions and show how these are a function of sensor geometry, system parameters, and environmental conditions. WSSUS processes, which are a subset of all random LTV systems, are the most widely used model in channel characterization problems [3,13]. In the following sections, we will illustrate the properties of WSSUS processes and validate the core assumptions for the radar sea clutter problem.

2.2. WSSUS Processes

To apply the WSSUS assumption to a radar scattering environment, we will express the impulse response in the following form:
h ( τ , t ) = F ρ 1 [ η ( τ , ρ ) ] = η ( τ , ρ ) e j 2 π ρ t d ρ
where F 1 is the inverse Fourier operator, ρ is Doppler shift, and η ( τ , ρ ) is known as the “delay-Doppler spread function” [13,31]. This function is a stochastic description of the range-Doppler map of the targets. Substituting Equation (2) into Equation (1) yields:
y ( t ) = η ( τ , ρ ) x ( t τ ) e j 2 π ρ t d ρ d τ + w ( t ) .
In Equation (3) it is clear to see that y ( t ) can be expressed as the superposition of delayed (by τ ), frequency-shifted (by ρ ), and scaled (by the complex gain η ( τ , ρ ) d τ d ρ ) copies of x ( t ) . It should be noted that this definition is general enough to include all scatterers in the radar’s field of view—clutter as well as useful targets. In this case the overall delay-Doppler spreading function is:
η total = η targets + η clutter ,
due to the linearity of the system.
To characterize the behavior of h as a stochastic process, we will assume it is zero mean and compute its second-order statistics:
R h ( τ , τ , t , t + Δ t ) = E h * ( τ , t ) h ( τ , t + Δ t ) .
Here will we apply the two fundamental assumptions of a WSSUS system: (1) The system is wide-sense stationary (WSS), implying that the autocorrelation depends only on Δ t along the t axis, and (2) the scattering at different lags τ are uncorrelated (US):
E h * ( τ , t ) h ( τ , t + Δ t ) = R h ( τ , τ , Δ t ) = A h ( τ , Δ t ) δ ( τ τ ) .
The function A h ( τ , Δ t ) is known as the system correlation function [3,48]. By taking the Fourier transform of A h ( τ , Δ t ) with respect to Δ t , we can view the Doppler spectrum of the return as a function of delay, referred to in the literature as the “scattering function” S h ( τ , ρ ) :
S h ( τ , ρ ) = F Δ t A h ( τ , Δ t ) .
The scattering function can be thought of as the true range-Doppler map of the scattering environment independent of the waveform used to probe it. It can also be shown that the autocorrelation function of the delay-Doppler spread function η for a WSSUS system is:
E η * ( τ , ρ ) η ( τ , ρ ) = S h ( τ , ρ ) δ ( τ τ ) δ ( ρ ρ ) .
It is sometimes easier to work with this form when deriving new expressions. This result shows that under a WSSUS model, the scattering not only at different ranges, but also different Doppler shifts are uncorrelated.
The full set of WSSUS system functions are shown in Figure 1, all of which are related by Fourier transforms of the different temporal variables. Knowledge of any one of these functions is a sufficient description of the second-order statistics of the system. These functions are standard tools in wireless communication channel modeling, and the diagram shown in Figure 1 can be found in virtually any textbook on the subject [3].
Note that nowhere in this discussion have we specified a distribution for the samples of h ( τ , t ) . If we assume that h is a circularly-symmetric complex Gaussian random process, then its distribution is completely specified by the WSSUS system functions described previously, and the magnitude envelope | y | of the return will be Rayleigh distributed [3]. However, in many cases a distribution with heavier tails, such as the Weibull or K distribution, is more appropriate, particularly for low grazing angles or when the radar can resolve individual sea spikes [49]. The results derived in the remainder of this paper, however, are valid regardless of the distribution of h. Note, however, that if h is non-Gaussian, higher-order statistics (third-order and above) are required to uniquely specify all of its properties.

2.3. Simulation Geometry

The simulation geometry is shown in Figure 2, where the airborne radar is at an altitude H above the surface traveling with velocity v = [ v X , v Y , v Z ] T . We will assume a coordinate system fixed to the phase center of the radar antenna, such that the antenna is at the origin, the X-axis is pointed parallel to the surface and in the plane of symmetry of the aircraft, the Y-axis is pointed out the right side of the aircraft and parallel to the surface, and the Z-axis is pointing down. In this work, we will assume that v always lies in the X-Z plane (i.e., no crabbing). Because we are seeking to model the scattering characteristics as a function of t and τ , we will need to express all spatial quantities in terms of these temporal variables.
We will model the radar signal using an approach similar to Barrick [50,51], i.e., as a spherical wave emanating from the source located at the origin as is shown in Figure 3. Each small segment of the wave reflects off of successive isorange rings on the sea surface; the total return is therefore the superposition of the returns from each isorange ring. The incremental power gain as a function of delay can be obtained using the radar range equation:
d P R P T = G 2 ( ϕ ( τ , ρ ) , θ ( τ ) ) λ 2 σ 0 ( α ( τ ) ) d A ( τ ) ( 4 π ) 3 1 2 c τ 4
where:
  • P T = transmit power
  • d P R = incremental received power from isorange ring at delay τ
  • ϕ , θ = azimuth and depression angles relative to platform
  • α = grazing angle, which equals θ in a flat earth model
    α = θ = sin 1 ( H / ( 1 2 c τ ) )
  • G ( ϕ , θ ) = one-way antenna power gain
  • λ = carrier wavelength
  • σ 0 = σ 0 ( α , λ , sea state , ) = surface normalized radar cross section (NRCS)
  • d A = π 2 c 2 τ d τ = incremental area of isorange ring at delay τ .

2.4. WSSUS System Function Derivations

To create a model for the scattering function S h ( τ , ρ ) , we need to characterize the signal return as a function of delay and Doppler. To do so, we will consider the Doppler spectrum generated from the backscatter from a single isorange ring and then apply the principle of superposition to obtain the total response from all ranges. The approach taken in this section can be considered a generalization of the Clarke model for flat fading [9,12], applied to modeling radar surface clutter.
We will assume that the impulse response from this ring is the superposition of returns from N equiangular patches in azimuth with random amplitudes and phases. We will also assume that the scatterers on this ring are located at delay τ ˜ and write the incremental impulse response of this thin ring as follows:
d h ( τ , t ; τ ˜ ) = δ ( τ τ ˜ ) n = 0 N 1 d a ( ϕ n ) e j ( 2 π ρ ( ϕ n ) t + γ n )
where ϕ n = 2 π n / N is the azimuth angle to patch n, d a ( ϕ n ) is the infinitesimal amplitude gain at angle ϕ n , ρ ( ϕ n ) is the Doppler shift at angle ϕ n , and γ n is a random phase shift. The incremental autocorrelation function of the impulse response at τ ˜ is then given by:
d R h ( τ , τ , t , t + Δ t ; τ ˜ ) = E d h * ( τ , t ; τ ˜ ) d h ( τ , t + Δ t ; τ ˜ ) = δ ( τ τ ˜ ) δ ( τ τ ˜ ) m , n E d a ( ϕ m ) d a ( ϕ n ) e j [ 2 π ( ρ ( ϕ n ) ( t + Δ t ) ρ ( ϕ m ) t ) + γ n γ m ] .
We will also assume that the scatterer amplitudes d a ( ϕ n ) are mutually independent, and we will assume that the amplitude d a ( ϕ n ) is independent of the Doppler shift ρ ( ϕ n ) , therefore:
d R h ( τ , τ , t , t + Δ t ; τ ˜ ) = δ ( τ τ ˜ ) δ ( τ τ ˜ ) n = 0 N 1 E | d a ( ϕ n ) | 2 E e j 2 π ρ ( ϕ n ) Δ t = d R h ( τ , τ , Δ t ; τ ˜ )
Note that E | d a ( ϕ n ) | 2 is simply the backscattered power gain from scatterer n, therefore by using Equation (9) and scaling it to account for the fact that the surface area of each patch is smaller by a factor of 1 / N , we obtain:
E | d a ( ϕ n ) | 2 = 1 N d P R P T = λ 2 σ 0 d τ ˜ 8 π 2 c 2 τ ˜ 3 1 N G 2 ( n Δ ϕ , θ ) .
Note that N = 2 π / Δ ϕ , where Δ ϕ is the angular spacing between patches, which upon substituting in Equation (13) yields:
E | d a ( ϕ n ) | 2 = λ 2 σ 0 d τ ˜ 16 π 3 c 2 τ ˜ 3 G 2 ( n Δ ϕ , θ ) Δ ϕ .
Combining this with Equation (12) and taking the limit as N yields
d R h ( τ , τ , Δ t ; τ ˜ ) = δ ( τ τ ˜ ) δ ( τ τ ˜ ) λ 2 σ 0 d τ ˜ 16 π 3 c 2 τ ˜ 3 · lim N n = 0 N 1 G 2 ( n Δ ϕ , θ ) E e j 2 π ρ ( ϕ n ) Δ t Δ ϕ = δ ( τ τ ˜ ) δ ( τ τ ˜ ) λ 2 σ 0 d τ ˜ 16 π 3 c 2 τ ˜ 3 2 π G 2 ( ϕ , θ ) E e j 2 π ρ ϕ , θ Δ t d ϕ .
Note that the expression E e j 2 π ρ ϕ , θ Δ t is the characteristic function of the random angular frequency 2 π ρ ϕ , θ [52]. In this work, we will define the characteristic function k ( Δ t | ϕ , θ ) as:
k ( Δ t | ϕ , θ ) = p ( ρ | ϕ , θ ) e j 2 π ρ Δ t d ρ = F ρ 1 p ( ρ | ϕ , θ ) ,
where p ( ρ | ϕ , θ ) is the probability density function (pdf) of the random frequency ρ ϕ , θ . Substituting this expression in Equation (15) yields:
d R h ( τ , τ , Δ t ; τ ˜ ) = δ ( τ τ ˜ ) δ ( τ τ ˜ ) λ 2 σ 0 d τ ˜ 16 π 3 c 2 τ ˜ 3 · lim N n = 0 N 1 G 2 ( n Δ ϕ , θ ) E e j 2 π ρ n Δ t Δ ϕ = δ ( τ τ ˜ ) δ ( τ τ ˜ ) λ 2 σ 0 d τ ˜ 16 π 3 c 2 τ ˜ 3 · 2 π G 2 ( ϕ , θ ) k ( Δ t | ϕ , θ ) d ϕ
One interpretation of the function p ( ρ | ϕ , θ ) that is used in the wireless communications literature is that the pdf of the random Doppler shift can be thought of as a normalized power spectral density (PSD) [3,8]. Thus if one has some model of the local Doppler spectrum, this information is accounted for in Equation (16).
To obtain the autocorrelation function for all delays τ we will use the US property that the return from each isorange ring τ ˜ is uncorrelated and thus we can integrate over τ ˜ to obtain the total response:
R h ( τ , τ , Δ t ) = d R h ( τ , τ , Δ t ; τ ˜ ) = δ ( τ τ ) λ 2 σ 0 16 π 3 c 2 τ 3 2 π G 2 ( ϕ , θ ) k ( Δ t | ϕ , θ ) d ϕ .
It is clear from Equation (18) that A h ( τ , Δ t ) is therefore:
A h ( τ , Δ t ) = λ 2 σ 0 16 π 3 c 2 τ 3 2 π G 2 ( ϕ , θ ) k ( Δ t | ϕ , θ ) d ϕ .
To find the clutter scattering function S h ( τ , ρ ) , we take the Fourier transform of A h ( τ , Δ t ) with respect to Δ t :
S h ( τ , ρ ) = F Δ t A h ( τ , Δ t ) = λ 2 σ 0 16 π 3 c 2 τ 3 2 π G 2 ( ϕ , θ ) F Δ t F ρ 1 p ( ρ | ϕ , θ ) d ϕ = λ 2 σ 0 16 π 3 c 2 τ 3 2 π G 2 ( ϕ , θ ) p ( ρ | ϕ , θ ) d ϕ
The expression in Equation (20) is significant because it gives an analytical expression for the full clutter spectrum, not just the mainlobe clutter, for any antenna pattern and provides a mechanism for supplying a priori information about the local Doppler spectra as a function of look angle. In this integral it can be seen that the antenna pattern has the effect of performing a weighted average of the local Doppler spectra over azimuth.

2.5. Important Special Cases

In this section, we will use Equation (20) to derive analytic expressions for the scattering function for several useful cases.

2.5.1. No ICM

In the degenerate case where there is no ICM, i.e., the surface motion is small relative to the platform motion, then the localized Doppler shifts ρ ϕ , θ are deterministic and only depend on the look angle relative to platform motion. This means that the local Doppler spectra are each just an impulse, and the scattering function Equation (20) reduces to [42]:
S h ( no ICM ) ( τ , ρ ) = λ 2 σ 0 16 π 3 c 2 τ 3 2 π G 2 ( ϕ , θ ) δ ( ρ β ( ϕ , θ ) ) d ϕ = λ 2 σ 0 ( α ( τ ) ) 16 π 3 c 2 τ 3 G 2 ( ϕ , θ ) + G 2 ( ϕ , θ ) ρ X ( τ ) 2 ρ ρ Z ( τ ) 2 , | ρ ρ Z ( τ ) | < ρ X ( τ ) 0 , else ,
where the nominal Doppler shift β ( ϕ , θ ) due to platform motion is given as:
β ( ϕ , θ ) = ρ X cos θ cos ϕ + ρ Z sin θ ,
where ρ X = 2 v X / λ and ρ Z = 2 v Z / λ are the Doppler contributions due to motion in the X Z plane, ρ X ( τ ) = ρ X cos ( θ ( τ ) ) and ρ Z ( τ ) = ρ Z sin ( θ ( τ ) ) are the effective Doppler contributions in the observation direction, and ϕ ( τ , ρ ) = arccos ( ( ρ ρ Z ( τ ) ) / ρ X ( τ ) ) is the nominal azimuth angle. The result in Equation (21) is a frequency-selective generalization of the Clarke/Jakes spectrum for airborne radar geometries [9,10].
A constant- τ cut of Equation (21) without the antenna pattern is plotted in Figure 4, with G = 1 and ρ Z = 0 . The singularity of the scattering function at ρ ρ Z ( τ ) ρ X ( τ ) is extremely critical for modeling nose-aspect clutter, as it serves to narrow the Doppler spread of the antenna pattern. For side-looking antennas, however, the singularities are irrelevant as the Jakes spectrum is flat near | ρ ρ Z ( τ ) | 0 .

2.5.2. Side-Looking Antenna, Level Flight Path

In many airborne radars, the antenna broadside is pointed at ϕ = π / 2 . In this case, we will create a new angular variable ψ = ϕ + π / 2 to represent angular deviation from broadside. The antenna pattern and local Doppler spectrum with respect to this coordinate will be denoted G ˜ and p ˜ , respectively. If we assume a level flight path (i.e., ρ Z = 0 ), then:
β ( ϕ , θ ) = ρ X cos θ cos ϕ = ρ X cos θ sin ψ ( ρ X cos θ ) ψ = ρ X ψ ,
which, for an antenna with beamwidth of less than 20°, is a very good approximation. Substituting into Equation (20) to compute the scattering function yields:
S h ( side ) ( τ , ρ ) = λ 2 σ 0 16 π 3 c 2 τ 3 π π G ˜ 2 ( ψ , θ ) p ˜ ( ρ | ψ , θ ) d ψ λ 2 σ 0 16 π 3 c 2 τ 3 G ˜ 2 ( ψ , θ ) p ˜ ( ρ | ψ , θ ) d ψ = λ 2 σ 0 16 π 3 c 2 τ 3 1 ρ X G ˜ 2 ( β / ρ X , θ ) p ˜ ( ρ | β / ρ X , θ ) d β ,
where the second equality is due to the fact that the antenna gain is nearly zero near the edges of the angular limit, so extending the limits of the integral to infinity does not change the result, and the third equality converts the integral from the angular domain to the Doppler domain. If we simplify the integrals by defining G ̲ ( β | τ ) = G ˜ ( β / ρ X ( τ ) , θ ( τ ) ) and p ̲ ( ρ | β , τ ) = p ˜ ( ρ | β / ρ X ( τ ) , θ ( τ ) ) , which are just the antenna pattern and local Doppler spectrum projected into delay-Doppler space, we can see that the scattering function is just a linear transformation of the antenna pattern with the local Doppler spectrum acting as the kernel function:
S h ( side ) ( τ , ρ ) = λ 2 σ 0 16 π 3 c 2 τ 3 1 ρ X G ̲ 2 ( β | τ ) p ̲ ( ρ | β , τ ) d β .
Note that p ̲ ( ρ | β , τ ) is the Doppler spectra observed by a moving platform, meaning that its center frequency is being modulated by the platform Doppler shift β ( ϕ , θ ) . If we assume that p ̲ does not change much relative to β , which is a reasonable assumption due to the fact that radar antennas usually have a very small azimuthal beamwidth, we can express p ̲ in terms of a “baseband” Doppler spectrum, i.e., the Doppler spectrum caused only by ICM that would be seen by a stationary radar with an infinitesimal beamwidth. The ICM spectrum is usually what is discussed in papers on Doppler spectrum modeling (e.g., [28,29,53]). We will denote this quantity as b ̲ ( ρ | τ ) and note that under the narrow beamwidth assumption, p ̲ ( ρ | β , τ ) b ̲ ( ρ β | τ ) , therefore:
S h ( side ) ( τ , ρ ) = λ 2 σ 0 16 π 3 c 2 τ 3 1 ρ X G ̲ 2 ( β | τ ) b ̲ ( ρ β | τ ) d β = λ 2 σ 0 ( τ ) 16 π 3 c 2 τ 3 1 ρ X ( τ ) G ̲ 2 ( ρ | τ ) * ρ b ̲ ( ρ | τ ) .
If we note that the coefficients outside the brackets only depend on τ , we can simplify Equation (27):
S h ( side ) ( τ , ρ ) = f ( τ ) G ̲ 2 ( ρ | τ ) * ρ b ̲ ( ρ | τ ) .
The relationship in Equation (27) is a common model of the relationship between ICM spectrum, antenna pattern, and observed Doppler spectrum used in side-looking airborne radars (e.g., [54]); this derivation explicitly shows the conditions that are required for this model to be accurate.
We can also define a baseband analogue to the characteristic function of Equation (16) as:
k ̲ ( Δ t | τ ) = F ρ 1 b ̲ ( ρ | τ ) ,
and note that for a side-looking radar:
A h ( side ) ( τ , Δ t ) = F ρ 1 S h ( side ) = f ( τ ) F ρ 1 G ̲ 2 ( ρ | τ ) * ρ b ̲ ( ρ | τ ) = f ( τ ) F ρ 1 G ̲ 2 ( ρ | τ ) F ρ 1 b ̲ ( ρ | τ ) = A h ( platform ) ( τ , Δ t ) · k ̲ ( Δ t | τ ) ,
where A h ( platform ) is the autocorrelation due to platform motion only. Because k ̲ serves the function of windowing the correlation function (which in this case is also a covariance because the random process is zero mean) in Δ t , it is referred to in the STAP literature as a covariance matrix taper (CMT) [27,28].

2.5.3. Arbitrary Orientation

In general, for a narrow-azimuthal beamwidth antenna, the scattering function will be:
S h ( τ , ρ ) S h ( no ICM ) ( τ , ρ ) * ρ b ̲ ( ρ | τ ) ,
where S h ( no ICM ) is the scattering function defined in Equation (21), which includes the Jakes spectrum peaks, and b ̲ is the ICM spectrum.

2.6. Output Time-Frequency Power Distribution

We will assume the matched-filter output is the correlation of N pulses against an infinite pulse train input, represented by the periodic ambiguity function (PAF) | χ N T | , given by [55]:
| χ N T ( τ , ρ ) | = 1 N T r 0 N T r x ( t ) x * ( t + τ ) e j 2 π ρ t d t ,
where x ( t ) is the normalized unit-energy input pulse train and T r is the pulse repetition interval (PRI). The pulse repetition frequency (PRF) F r = 1 / T r . The PAF has the following properties that are relevant to our application:
| χ N T ( τ , ρ ) | = | χ N T ( τ + T r , ρ ) | ,
and:
| χ N T ( τ , ρ ) | = | χ 1 T ( τ , ρ ) | sin N π ρ T r N sin π ρ T r .
One important result from Green’s work [56] was that the time-frequency power distribution is proportional to the convolution of the ambiguity function with the WSSUS system scattering function:
P ( τ , ρ ) = E x S h ( τ , ρ ) * * | χ N T ( τ , ρ ) | 2 = E x S h ( τ , ρ ) | χ N T ( τ τ , ρ ρ ) | 2 d τ d ρ ,
where E x is the energy of the pulse train. Note that since | χ N T | 2 is periodic along the τ axis, the result of the convolution integral of Equation (34) is periodic along this axis as well. If we assume the support of S h is limited to [ 0 , Q T r ] for some positive integer Q along the τ axis, then the output power distribution can be written as a circular convolution:
P ( τ , ρ ) = E x 0 Q T r S h ( τ , ρ ) | χ N T ( τ τ , ρ ρ ) | 2 d ρ d τ = E x S h ( τ , ρ ) τ * ρ | χ N T ( τ , ρ ) | 2 ,
where ⊛ represents circular convolution and * represents ordinary convolution. The problem was cast into this form so that discrete implementations can use the Fast Fourier Transform (FFT) to efficiently compute the periodic output power distribution along the τ axis without zero padding.

2.7. ICM Spectrum Characterization

2.7.1. Maximum Entropy Prior

As it has been noted in Section 2.4, b ̲ ( ρ | τ ) can be interpreted simultaneously as either (a) a PSD corresponding to a small patch of sea surface, or (b) a pdf of a random Doppler shift introduced by surface motion. Thus the spectrum characterization problem is reduced to a problem in prior probability density estimation.
Choosing a prior distribution is an important component of Bayesian estimation, and thus, there is a wide selection of literature available on the subject. If one has prior knowledge of the Doppler distribution from measurements or physical calculations, they can immediately compute the appropriate clutter taper using Equation (16). However, this is not usually the case, as the Doppler spectra for any given operating wavelength λ may depend on a multitude of factors, such as sea state, wind speed, temperature, salinity, etc.
In the absence of such detailed information, one tool for selecting an appropriate prior is the principle of maximum entropy (MaxEnt) [44]. It is based on the objectivist Bayesian philosophy that probabilities represent a state of knowledge rather than a degree of belief, and thus the selection of a prior distribution should be based on objective criteria, such as known expectations or known support of the random variable, in such a way that the amount of “assumed” information is minimized.
Formally, this is expressed as choosing some distribution function b ̲ ( x ) with support ( a , b ) of the random variable X that maximizes the differential entropy H:
H ( X ) = a b b ̲ ( x ) log b ̲ ( x ) d x ,
subject to the constraints imposed by the known “testable information” F i :
E f i ( X ) = a b b ̲ ( x ) f i ( x ) d x = F i , i = 1 , , n .
The solution to this optimization problem is found using the calculus of variations:
b ̲ ( x ) = c 0 m ( x ) exp i λ i f i ( x ) ,
where c 0 is a normalization constant, m ( x ) is a partition function that is constant in the support region x ( a , b ) and zero elsewhere, and λ i are the Lagrange multipliers. The shape of the resulting distribution depends on the amount of known “testable information”.

2.7.2. Known Mean, Unknown Variance

We often know the wind speed v w and direction ϕ w but no other information about the shape of the Doppler prior. If we start with the following assumptions:
  • The mean wave speed (and hence mean Doppler shift) is proportional to the wind speed. For example, it is commonly assumed that the wave speed is 1/8th the wind speed [57], and
  • The wind velocity vector has no Z component, and
  • The waves move in the same direction as the wind,
then the following mathematical statements can then be constructed to express this testable information:
  • μ ρ v w cos θ cos ( ϕ ϕ w ) ,
  • ρ sgn ( μ ρ ) [ 0 , ) .
There are cases when these statements are not true, such as when there is a sudden change of wind direction, but for a fully developed wind-wave these are not particularly controversial statements [49]. The MaxEnt prior that satisfies these criteria can be shown to be:
b ̲ ( ρ | τ ) = 1 | μ ρ ( τ ) | exp ρ μ ρ ( τ ) u ρ sgn μ ρ ( τ ) ,
where u ( · ) is the Heaviside step function. This is an exponential distribution with respect to the wind direction. The intuition behind this is that, while the mean Doppler shift will be proportional to the wind speed, there will be some wave components that move much faster, but none that move in the opposite direction of the wind.

2.7.3. Known Mean and Variance

If the mean E ρ = μ ρ and variance E ( ρ E ρ ) 2 = σ ρ 2 of the Doppler frequency are known, then the MaxEnt prior is a truncated Gaussian:
b ̲ ( ρ | τ ) = c 0 m ( ρ ) exp ( ρ μ ρ ( τ ) ) 2 2 σ ρ 2 ( τ ) .
In general μ ρ μ ρ and σ ρ 2 σ ρ 2 due to the truncation of the tails of the distribution by m ( ρ ) . The Doppler spectra model given in §3.8.2 of [49] can be considered a version of this with a and b , in which case μ ρ = μ ρ and σ ρ = σ ρ .

2.7.4. Distribution Comparison

The claimed exponential shape of the MaxEnt Doppler spectrum may seem puzzling to experienced readers who have experience with Doppler spectra being commonly modeled as having a Gaussian shape. The reader is reminded that Equation (39) is not the spectrum the radar “sees”; the shape of the spectrum a (side-looking) radar will see is determined by the convolution (along ρ ) of the Doppler spectra b ̲ ( ρ | τ ) with the antenna pattern G ̲ 2 ( ρ | τ ) , which is then convolved (along τ and ρ ) with the waveform ambiguity squared | χ N T ( τ , ρ ) | 2 . Each of these processing steps alters the resolution of the radar.
As an illustrative example, we will consider the Doppler spectra modeled in §3.8.2 of [49]. Since the radar is stationary, G ̲ 2 ( ρ | τ ) δ ( ρ ) . The predicted spectrum in this scenario is described using the following parameters:
  • F r = 2 kHz
  • N = 32 pulses
  • Dolph-Chebyshev window with −55 dB sidelobes
  • μ ρ = 62.5 Hz
  • σ ρ 2 = 20 Hz2
  • Overall clutter to noise ratio (CNR) = 20 dB
The computed MaxEnt spectra for both the unknown and known variance cases are shown in Figure 5. The predicted power distribution P ( ρ ) = | χ N T ( 0 , ρ ) | 2 * b ̲ ( ρ ) is shown in Figure 6. It can be seen in Figure 5 that the specification of the variance narrows the Doppler pdf b ̲ ( ρ ) significantly, but if too short of a CPI is used, as is the case in Figure 6, the observed spectrum after matched filtering will not change significantly.

2.8. Experimental Validation

2.8.1. Doppler Spectrum Modeling

To validate our framework, we will compare the modeled clutter spectra to empirically measured spectra. For our comparison, we will use the Ingara dataset, which is a medium grazing angle dataset containing measured returns from an airborne radar flying in a circular path at a speed of approximately 200 knots, operating in spotlight mode (i.e., illuminating the patch of sea at the center of the circle.) Data were collected at grazing angles from 15 to 45 degrees. The radar operates at a center frequency of 10.1 GHz using a linear frequency modulated (LFM) waveform with a pulse width of 20 μs and a bandwidth 200 MHz, leading to a range resolution of about 0.75 m.
The data in the flight that we will use for comparison was collected from an altitude of 0.5 nmi with the plane flying in a circular path of radius 1.9 nmi, leading to a grazing angle of approximately 15 degrees. This run contains the measured radar IQ returns over a 305 second interval, or about 1.4 revolutions around the circular path. A plot of the downrange video at each slow time is shown in Figure 7. Close inspection of this picture shows that the individual wave crests and troughs are clearly visible.
To test the ability of the WSSUS framework to predict Doppler spectra, we will compare the predicted spectra to the Doppler spectra observed while the radar is looking downwind. The observed range-Doppler map seen when the radar is looking downwind is shown in Figure 8. It can be seen that the spectrum is biased in the negative Doppler direction, as expected. Significant variation in the mean Doppler frequency in each range bin can be seen; this is likely due to the high range resolution of the radar. Because these minute fluctuations are essentially random and difficult to predict, we will take the ensemble average of the normalized spectra over τ . In this dataset, the elevation beampattern was removed as a preprocessing step, so it is expected that the Doppler spectra in each range bin will have similar statistics with the exception of path loss-induced amplitude decay. The ensemble average spectrum is shown in Figure 9, where it is labeled “empirical”.
Because the Ingara radar is side-looking and is flying at a constant altitude, we will use the simplifications of the WSSUS model in Section 2.5.2 to compute the scattering function S h , which we will convolve with the ambiguity to obtain the delay-Doppler map P ( τ , ρ ) . The azimuth antenna pattern and waveform parameterization were supplied with the dataset, so the only unknown in our prediction efforts is the ICM Doppler prior b ̲ ( ρ | τ ) .
We will use the MaxEnt procedure of Section 2.7.1 to estimate the ICM prior. For these models, it is necessary to predict mean and variance. Because the spectrum is periodic, we estimated the mean μ ρ Doppler shift using the circular mean [45]:
μ ρ F r 2 π arg E P e j 2 π ρ / F r ,
where E P denotes expectation is taken with respect to the empirically observed spectrum P ( ρ ) , normalized to unit area. From the data in Figure 9, μ ρ = −47.9 Hz, which corresponds to a radial wave speed of −0.71 m/s. We will assume that μ ρ accounts for both the mean ICM Doppler shift as well as antenna pointing errors.
For estimating the variance σ ρ 2 of the prior, we need to account for the fact that, according to Equations (27) and (35), the spectral width is due to the convolution of the prior with the antenna with the waveform ambiguity, therefore each step “broadens” the spectrum. In the elementary case when all three functions are Gaussian, the observed spectral variance σ obs 2 is the sum of the component variances:
σ obs 2 = σ ρ 2 + σ a 2 + σ w 2 ,
where σ a and σ w are the RMS Doppler widths of the antenna pattern and waveform, respectively. This can be used to estimate the variance of the clutter prior:
σ ρ σ obs 2 σ a 2 σ w 2 .
However, because the observed distribution P ( ρ ) and the waveform ambiguity cut | χ N T ( 0 , ρ ) | 2 are both periodic, the appropriate way to measure σ obs and σ w are using the circular standard deviation [45]:
σ obs = 2 ln E P e j 2 π ρ / F r σ w = 2 ln E | χ | 2 e j 2 π ρ / F r ,
where E | χ | 2 is expectation with respect to the normalized ambiguity squared.
The predicted spectrum using the WSSUS framework with the MaxEnt priors, labeled “WSSUS u.v.” for the unknown variance case, and “WSSUS k.v.” for the known variance case, as well a plot of the “standard” Gaussian Doppler spectrum model for P ( ρ ) , are shown in Figure 9. It can be seen that both predicted spectra capture the shape of the empirical spectrum; the specification of the variance produces a spectrum that is within approximately 1 dB of the empirical spectrum at all frequencies. The Gaussian model for P ( ρ ) produces excelled agreement for in the main lobe but does not predict the presence of the sidelobe clutter floor of approximately –29 dB. This could lead to excessively optimistic predictions of radar performance in the sidelobe clutter region.
It is expected that for a side-looking radar in level flight the clutter Doppler spectrum will be centered at 0 Hz. The negative frequency bias of the observed Doppler spectrum in Figure 9 could be due to:
  • Wave motion away from the radar, as described in Section 2.7.2, of approximately 0.71 m/s,
  • Aircraft crabbing (i.e., translational motion in the Y-direction),
  • Antenna pointing errors; using Equation (23) yields a pointing error of approximately 0.4 degrees,
  • Some combination of the above effects.

2.8.2. Validation of WSSUS Assumption

Another validation test that was performed using the Ingara dataset was to determine the validity of the WSSUS assumption. If the clutter impulse response h ( τ , t ) is non-WSS in t, then the autocorrelation A h ( τ , Δ t ) would not just depend on the time-difference Δ t , but also the absolute time t.
To characterize the change in the autocorrelation over time, we computed the autocorrelation of the data over a sliding window of width 1 s over the entire record length of 305 s in range bin 1. The start time of each window was advanced by 10 ms, leading to a 99% overlap between windows. We will denote this function A t ( Δ t ) . To characterize the autocorrelation at each time instant t, we measure the decorrelation time T c , defined implicitly by:
A t ( T c ) A t ( 0 ) = e 1 ,
i.e., T c is the time it takes for the autocorrelation centered at time t to decay to 1 / e from its peak value. The measured values of T c versus time are plotted in Figure 10.
We will use the stability of T c over time as a proxy to for the stationarity properties of the clutter. It is clear that the clutter is nonstationary over long periods of time, but for short instants of time it may be approximated as being locally WSS. Ironically, even researchers that study the nonstationary nature of clutter spectra must implicitly assume the clutter signal is locally WSS for them to be able to apply the Wiener-Khinchin theorem to estimate the clutter power spectrum in the first place! So the important question is not whether or not clutter is stationary, but how long does it remain locally stationary? Our modeling so far has assumed that the stationarity duration is larger than a CPI.
To test this assumption, we created a statistic called the “stability duration” of T c , which is defined as the region of time over which T c changes less than X%. (We used 10 percent, but this is admittedly arbitrary.) Because T c = T c ( t ) is itself a function of time, the stability duration was computed relative to the value of T c ( t ) at each timestep. A histogram of the stability times over the entire record length is shown in Figure 11. It can be seen that the median stability time is about 1.5 s, with the mean being even higher. Thus for any waveform with a PRI on the order of one millisecond or less (i.e., PRF > 1 kHz, which is very common) and integrating less than 1000 pulses, this is a perfectly reasonable assumption. This assumption may be violated for low-PRF waveforms with long integration times, such as those used in surveillance radars.

3. Discussion and Conclusions

Effective suppression of unwanted sea clutter returns requires an accurate characterization of their spatio-temporal statistics. We have created a WSSUS model for sea clutter that allows for the prediction of the clutter range-Doppler spectrum and provided a mechanism by which internal clutter motion may be accounted for and quickly estimated. Validation against the Ingara medium grazing angle dataset shows agreement between the WSSUS model and the ensemble average spectrum to within 1 dB at all frequencies.
Future work could apply this model to developing signal processing techniques such as adaptive transmit waveform design and improved MTI filtering to optimize signal detection in the presence of heterogeneous clutter.

4. Materials and Methods

The Ingara dataset, which was used for all experiments in this paper, is not publicly available; access is controlled by the Australian Defence Science and Technology (DST) Group [46]. Data analysis methods are described fully in the paper.
References

Funding

This work was performed in support of US NRL Contract N00173-17-F-2010.

Acknowledgments

Adam Anderson (Tennessee Technological University, Cookeville, TN, USA, email: [email protected]) has assisted in proofreading and has provided feedback in the capacity of an academic advisor.

Conflicts of Interest

The author declares no conflict of interest. The sponsor of this work had no role in the design, execution, interpretation, or writing of the study.

References

  1. Cooke, C.D. Scattering function approach for modeling time-varying sea clutter returns. In Proceedings of the 2018 IEEE Radar Conference (RadarConf), Oklahoma City, OK, USA, 23–27 April 2018. [Google Scholar]
  2. Bello, P. Characterization of randomly time-variant linear channels. IEEE Trans. Commun. Syst. 1963, 11, 360–393. [Google Scholar] [CrossRef]
  3. Goldsmith, A. Wireless Communications; Cambridge University Press: Cambridge, UK, 2005. [Google Scholar]
  4. Matz, G.; Bolcskei, H.; Hlawatsch, F. Time-frequency foundations of communications: Concepts and tools. IEEE Signal Process. Mag. 2013, 30, 87–96. [Google Scholar] [CrossRef]
  5. Matz, G. On non-WSSUS wireless fading channels. IEEE Trans. Wirel. Commun. 2005, 4, 2465–2478. [Google Scholar] [CrossRef]
  6. Matolak, D.W. Channel modeling for vehicle-to-vehicle communications. IEEE Commun. Mag. 2008, 46, 76–83. [Google Scholar] [CrossRef]
  7. Paetzold, M.; Gutierrez, C.A. Definition and analysis of quasi-stationary intervals of mobile radio channels–invited paper. In Proceedings of the 2018 IEEE 87th Vehicular Technology Conference (VTC Spring), Porto, Portugal, 3–6 June 2018; pp. 1–6. [Google Scholar] [CrossRef]
  8. Walter, M.; Shutin, D.; Dammann, A. Time-variant Doppler PDFs and characteristic functions for the vehicle-to-vehicle channel. IEEE Trans. Veh. Technol. 2017, 66, 10748–10763. [Google Scholar] [CrossRef]
  9. Clarke, R.H. A statistical theory of mobile-radio reception. Bell Syst. Tech. J. 1968, 47, 957–1000. [Google Scholar] [CrossRef]
  10. Jakes, W.C. Microwave Mobile Communications; Wiley-IEEE Press: Hoboken, NJ, USA, 1974. [Google Scholar]
  11. Aulin, T. A modified model for the fading signal at a mobile radio channel. IEEE Trans. Veh. Technol. 1979, 28, 182–203. [Google Scholar] [CrossRef]
  12. Feng, T.; Field, T.R. Statistical analysis of mobile radio reception: An extension of Clarke’s model. IEEE Trans. Commun. 2008, 56, 2007–2012. [Google Scholar] [CrossRef]
  13. Matz, G. Time-frequency characterization of random time-varying channels. In Time-Frequency Signal Analysis and Processing: A Comprehensive Reference; Boashash, B., Ed.; Academic Press: Cambridge, MA, USA, 2016; Chapter 9; pp. 554–561. [Google Scholar]
  14. Kivinen, J.; Zhao, X.; Vainikainen, P. Empirical characterization of wideband indoor radio channel at 5.3 GHz. IEEE Trans. Antennas Propag. 2001, 49, 1192–1203. [Google Scholar] [CrossRef]
  15. Sadowsky, J.S.; Kafedziski, V. On the correlation and scattering functions of the WSSUS channel for mobile communications. IEEE Trans. Veh. Technol. 1998, 47, 270–282. [Google Scholar] [CrossRef]
  16. Ziomek, L.J. A Scattering Function Approach to Underwater Acoustic Detection and Signal Design. Ph.D. Thesis, The Pennsylvania State University, Centre County, PA, USA, 1981. [Google Scholar]
  17. Ziomek, L.J.; Sibul, L.H. Broadband and narrow-band signal-to-interference ratio expressions for a doubly spread target. J. Acoust. Soc. Am. 1982, 72, 804–819. [Google Scholar] [CrossRef]
  18. Ricker, D.; Cutezo, A. Estimation of coherent detection performance for spread scattering in reverberation-noise mixtures. J. Acoust. Soc. Am. 2000, 107, 1978–1986. [Google Scholar] [CrossRef] [PubMed]
  19. Zhu, Z.; Kay, S.; Raghavan, R.S. Waveform design for doubly spread targets detection. In Proceedings of the 2018 IEEE Radar Conference (RadarConf18), Oklahoma City, OK, USA, 23–27 April 2018; pp. 49–54. [Google Scholar] [CrossRef]
  20. Zhu, Z.; Kay, S.; Raghavan, R.S. Locally optimal radar waveform design for detecting doubly spread targets in colored noise. IEEE Signal Process. Lett. 2018, 25, 833–837. [Google Scholar] [CrossRef]
  21. Van Trees, H.L. Detection, Estimation, and Modulation Theory; John Wiley & Sons: Hoboken, NJ, USA, 1971; Volume 3. [Google Scholar]
  22. Sibul, L.E. Optimum detection and estimation in stochastic transmission media. In Applied Stochastic Processes; Academic Press: Cambridge, MA, USA, 1980; pp. 247–267. [Google Scholar]
  23. Sibul, L.H.; Titlebaum, E.L. Signal design for detection of targets in clutter. Proc. IEEE 1981, 69, 481–482. [Google Scholar] [CrossRef]
  24. Ziomek, L.J.; Sibul, L.H. Maximation of the signal-to-interference ratio for a doubly spread target: Problems in nonlinear programming. Signal Process. 1983, 5, 355–368. [Google Scholar] [CrossRef]
  25. Sibul, L.H.; Fowler, M.L.; Harbison, G.J. Optimum array processing and representation of nonstationary random scattering. In Proceedings of the Twenty-Second Asilomar Conference on Signals, Systems and Computers, Pacific Grove, CA, USA, 31 October–2 November 1988; Volume 1, pp. 431–435. [Google Scholar] [CrossRef]
  26. Sibul, L.H.; Roan, M.J.; Hillsley, K.L. Wavelet transform techniques for time-varying propagation and scattering characterization. In Proceedings of the Conference Record of Thirty-Second Asilomar Conference on Signals, Systems and Computers (Cat. No.98CH36284), Pacific Grove, CA, USA, 1–4 November 1998; Volume 2, pp. 1644–1649. [Google Scholar] [CrossRef]
  27. Guerci, J.R. Theory and application of covariance matrix tapers for robust adaptive beamforming. IEEE Trans. Signal Process. 1999, 47, 977–985. [Google Scholar] [CrossRef]
  28. Guerci, J.R.; Bergin, J.S. Principal components, covariance matrix tapers, and the subspace leakage problem. IEEE Trans. Aerosp. Electron. Syst. 2002, 38, 152–162. [Google Scholar] [CrossRef]
  29. Rosenberg, L. Comparison of bi-modal coherent sea clutter simulation techniques. IET Radar Sonar Navig. 2019, 13, 1519–1529. [Google Scholar] [CrossRef]
  30. Gaarder, N. Scattering function estimation. IEEE Trans. Inf. Theory 1968, 14, 684–693. [Google Scholar] [CrossRef]
  31. Nguyen, L.T.; Senadji, B.; Boashash, B. Scattering function and time-frequency signal processing. In Proceedings of the 2001 IEEE International Conference on Acoustics, Speech, and Signal Processing. Proceedings (Cat. No.01CH37221), Salt Lake City, UT, USA, 7–11 May 2001; Volume 6, pp. 3597–3600. [Google Scholar] [CrossRef]
  32. Kay, S.M.; Doyle, S.B. Rapid estimation of the range-Doppler scattering function. IEEE Trans. Signal Process. 2003, 51, 255–268. [Google Scholar] [CrossRef]
  33. Artes, H.; Matz, G.; Hlawatsch, F. Unbiased scattering function estimators for underspread channels and extension to data-driven operation. IEEE Trans. Signal Process. 2004, 52, 1387–1402. [Google Scholar] [CrossRef]
  34. Allemand, V. Estimation of the Scattering Function of Fading Channels For Acoustic Communications in Shallow Waters; Florida Atlantic University: Boca Raton, FL, USA, 2005. [Google Scholar]
  35. Van Walree, P. Channel Sounding for Acoustic Communications: Techniques and Shallow-Water Examples; Technical Report FFI-rapport 2011/00007; Norwegian Defence Research Establishment (FFI): Horten, Norway, 2011. [Google Scholar]
  36. Yin, X.; Zuo, Q.; Zhong, Z.; Lu, S.X. Delay-Doppler frequency power spectrum estimation for vehicular propagation channels. In Proceedings of the 5th European Conference on Antennas and Propagation (EUCAP), Rome, Italy, 11–15 April 2011; pp. 3434–3438. [Google Scholar]
  37. Boayue, A. Characterization of Underwater Acoustic Communication Channels: Statistical Characteristics of the Underwater Multipath Channnels. Master’s Thesis, Institutt for Elektronikk og Telekommunikasjon, Oslo, Norway, 2013. [Google Scholar]
  38. Oktay, O.; Pfander, G.E.; Zheltov, P. Reconstruction of the scattering function of overspread radar targets. IET Signal Process. 2014, 8, 1018–1024. [Google Scholar] [CrossRef] [Green Version]
  39. Pfander, G.E.; Zheltov, P. Identification of stochastic operators. Appl. Comput. Harmon. Anal. 2014, 36, 256–279. [Google Scholar] [CrossRef]
  40. Pfander, G.E.; Zheltov, P. Sampling of stochastic operators. IEEE Trans. Inf. Theory 2014, 60, 2359–2372. [Google Scholar] [CrossRef]
  41. Pfander, G.E.; Zheltov, P. Estimation of overspread scattering functions. IEEE Trans. Signal Process. 2015, 63, 2451–2463. [Google Scholar] [CrossRef]
  42. Jarvis, D.E.; Cooke, C.D. Application of Statistical Linear Time-Varying System Theory to Modeling of High Grazing Angle Sea Clutter; Technical Report NRL/MR/5750–17-9747; U.S. Naval Research Laboratory: Washington, DC, USA, 2017. [Google Scholar]
  43. Cooke, C.D.; Jarvis, D.E. Application of a Statistical Linear Time-Varying System Model of High Grazing Angle Sea Clutter for Computing Interference Power; Technical Report NRL/MR/5753–17-9744; U.S. Naval Research Laboratory: Washington, DC, USA, 2017. [Google Scholar]
  44. Jaynes, E.T. Probability Theory: The Logic of Science; Cambridge University Press: Cambridge, UK, 2003. [Google Scholar]
  45. Mardia, K.V.; Jupp, P.E. Directional Statistics; John Wiley & Sons: Hoboken, NJ, USA, 2000. [Google Scholar]
  46. Rosenberg, L.; Watts, S.; Bocquet, S.; Ritchie, M. Characterisation of the Ingara HGA dataset. In Proceedings of the 2015 IEEE Radar Conference (RadarCon), Arlington, VA, USA, 10–15 May 2015; pp. 27–32. [Google Scholar] [CrossRef]
  47. Rihaczek, A.W. Principles of High Resolution Radar; McGraw-Hill: New York, NY, USA, 1969. [Google Scholar]
  48. Proakis, J.; Salehi, M. Digital Communications, 5th ed.; McGraw-Hill: New York, NY, USA, 2008. [Google Scholar]
  49. Ward, K.; Tough, R.; Watts, S. Sea Clutter: Scattering, the K Distribution and Radar Performance; Institution of Engineering and Technology: Stevenage, UK, 2013. [Google Scholar]
  50. Barrick, D.E. The interaction of HF/VHF radio waves with the sea surface and its implications. In Proceedings of the AGARD Conference Proceedings on Electromagnetics of the Sea, Paris, France, 22–25 June 1970; Volume 77, pp. 1–25. [Google Scholar]
  51. Barrick, D.E. Remote sensing of sea state by radar. In Proceedings of the Ocean 72—IEEE International Conference on Engineering in the Ocean Environment, Newport, RI, USA, 13–15 September 1972; pp. 186–192. [Google Scholar] [CrossRef]
  52. Leon-Garcia, A. Probability, Statistics, and Random Processes for Electrical Engineering; Pearson Prentice Hall: Upper Saddle River, NJ, USA, 2008. [Google Scholar]
  53. Billingsley, B.; Farina, A.; Gini, F.; Greco, M.; Lombardo, P. Impact of experimentally measured Doppler spectrum of ground clutter on MTI and STAP. In Radar Systems (RADAR 97); Institution of Engineering and Technology: Stevenage, UK, 1997. [Google Scholar]
  54. Rosenberg, L. Characterization of high grazing angle X-band sea-clutter Doppler spectra. IEEE Trans. Aerosp. Electron. Syst. 2014, 50, 406–417. [Google Scholar] [CrossRef]
  55. Levanon, N.; Mozeson, E. Radar Signals; IEEE Press/Wiley-Interscience: Hoboken, NJ, USA, 2004. [Google Scholar]
  56. Green, P.E, Jr. Radar Astronomy Measurement Techniques; Technical Report 282; MIT Lincoln Laboratory: Lexington, MA, USA, 1962. [Google Scholar]
  57. Watts, S.; Rosenberg, L. Radar Clutter Modelling and Exploitation. In Proceedings of the 2018 IEEE Radar Conference Tutorial Sessions, Oklahoma City, OK, USA, 23–27 April 2018. [Google Scholar]
Figure 1. Wide-sense stationary uncorrelated scattering (WSSUS) system function relationships. Knowledge of any one of these functions grants complete knowledge of the second-order statistics of the system.
Figure 1. Wide-sense stationary uncorrelated scattering (WSSUS) system function relationships. Knowledge of any one of these functions grants complete knowledge of the second-order statistics of the system.
Entropy 21 00915 g001
Figure 2. Flat earth geometry in the X-Z plane. Because the goal is to ultimately produce a range-Doppler map, we need to express all quantities in terms of downrange delay τ and absolute time t.
Figure 2. Flat earth geometry in the X-Z plane. Because the goal is to ultimately produce a range-Doppler map, we need to express all quantities in terms of downrange delay τ and absolute time t.
Entropy 21 00915 g002
Figure 3. Return from a single infinitesimal isorange ring on the sea surface. The autocorrelation function of the impulse response of this ring is computed at delay τ ˜ , then the result is integrated over all τ ˜ to get the total autocorrelation as a function of delay τ and time offset Δ t .
Figure 3. Return from a single infinitesimal isorange ring on the sea surface. The autocorrelation function of the impulse response of this ring is computed at delay τ ˜ , then the result is integrated over all τ ˜ to get the total autocorrelation as a function of delay τ and time offset Δ t .
Entropy 21 00915 g003
Figure 4. Constant- τ cut of the no-internal clutter motion (ICM) scattering function Equation (21) with the antenna pattern removed (i.e., G = 1 ) and no vertical motion ( ρ Z = 0 ). The true scattering function will be windowed by the antenna pattern to focus on a specific region of ρ —a side-looking antenna will be focused near ρ = 0 , whereas a nose-aspect antenna will be focused near ρ = ρ X .
Figure 4. Constant- τ cut of the no-internal clutter motion (ICM) scattering function Equation (21) with the antenna pattern removed (i.e., G = 1 ) and no vertical motion ( ρ Z = 0 ). The true scattering function will be windowed by the antenna pattern to focus on a specific region of ρ —a side-looking antenna will be focused near ρ = 0 , whereas a nose-aspect antenna will be focused near ρ = ρ X .
Entropy 21 00915 g004
Figure 5. Maximum entropy (MaxEnt) Doppler priors for unknown variance (u.v.) and known variance (k.v.) cases. For the u.v. case, the MaxEnt prior is an exponential distribution, whereas for the k.v. case, the MaxEnt prior is a Gaussian. Note that additional knowledge of the variance significantly reduces the predicted spectrum width.
Figure 5. Maximum entropy (MaxEnt) Doppler priors for unknown variance (u.v.) and known variance (k.v.) cases. For the u.v. case, the MaxEnt prior is an exponential distribution, whereas for the k.v. case, the MaxEnt prior is a Gaussian. Note that additional knowledge of the variance significantly reduces the predicted spectrum width.
Entropy 21 00915 g005
Figure 6. Output power distributions for unknown variance (u.v.) and known variance (k.v.) cases. Note that in this case the predicted spectrum when the variance is known becomes waveform-limited because σ ρ F r / N .
Figure 6. Output power distributions for unknown variance (u.v.) and known variance (k.v.) cases. Note that in this case the predicted spectrum when the variance is known becomes waveform-limited because σ ρ F r / N .
Entropy 21 00915 g006
Figure 7. Ingara data from trial SCT04, flight F42, run 34,877. The Ingara radar contains 1024 range gates covering a range span of approximately 767 m, captured from an airborne radar operating in spotlight mode flying in a circular path at an altitude of 0.5 nmi and observing the ocean at a grazing angle of 15°. In this run, the plane completes approximately 1.4 revolutions of its circular flight path.
Figure 7. Ingara data from trial SCT04, flight F42, run 34,877. The Ingara radar contains 1024 range gates covering a range span of approximately 767 m, captured from an airborne radar operating in spotlight mode flying in a circular path at an altitude of 0.5 nmi and observing the ocean at a grazing angle of 15°. In this run, the plane completes approximately 1.4 revolutions of its circular flight path.
Entropy 21 00915 g007
Figure 8. Range-Doppler map when the radar is looking in the downwind direction, obtained over a coherent processing interval (CPI) of 128 pulses (0.52 s). The wind speed was 8.5 m/s. Note that the Doppler shifts are biased towards negative frequencies. This bias could be due to wave motion away from the radar as well as antenna pointing errors.
Figure 8. Range-Doppler map when the radar is looking in the downwind direction, obtained over a coherent processing interval (CPI) of 128 pulses (0.52 s). The wind speed was 8.5 m/s. Note that the Doppler shifts are biased towards negative frequencies. This bias could be due to wave motion away from the radar as well as antenna pointing errors.
Entropy 21 00915 g008
Figure 9. Empirical ensemble average spectrum plotted versus WSSUS predicted spectra for the unknown variance (u.v.) and known variance (k.v.) cases. As well, the “standard” Gaussian clutter spectrum when the variance is known is plotted for reference. It can be seen that the Gaussian spectrum model vastly underestimates the clutter floor caused by antenna and waveform sidelobes.
Figure 9. Empirical ensemble average spectrum plotted versus WSSUS predicted spectra for the unknown variance (u.v.) and known variance (k.v.) cases. As well, the “standard” Gaussian clutter spectrum when the variance is known is plotted for reference. It can be seen that the Gaussian spectrum model vastly underestimates the clutter floor caused by antenna and waveform sidelobes.
Entropy 21 00915 g009
Figure 10. Clutter decorrelation times measured over time in range bin 1. Note that the decorrelation time remains roughly constant over several intervals much longer than a CPI. This provides support for the assumption that the clutter statistics are locally WSS for sufficiently long instants in time for the assumptions in Section 2.4 to apply.
Figure 10. Clutter decorrelation times measured over time in range bin 1. Note that the decorrelation time remains roughly constant over several intervals much longer than a CPI. This provides support for the assumption that the clutter statistics are locally WSS for sufficiently long instants in time for the assumptions in Section 2.4 to apply.
Entropy 21 00915 g010
Figure 11. Histogram of stability durations. It can be seen that the decorrelation time T c is roughly constant over intervals much greater than one CPI for medium- and high-PRF waveforms.
Figure 11. Histogram of stability durations. It can be seen that the decorrelation time T c is roughly constant over intervals much greater than one CPI for medium- and high-PRF waveforms.
Entropy 21 00915 g011

Share and Cite

MDPI and ACS Style

Cooke, C. Analytic Expressions for Radar Sea Clutter WSSUS Scattering Functions. Entropy 2019, 21, 915. https://doi.org/10.3390/e21090915

AMA Style

Cooke C. Analytic Expressions for Radar Sea Clutter WSSUS Scattering Functions. Entropy. 2019; 21(9):915. https://doi.org/10.3390/e21090915

Chicago/Turabian Style

Cooke, Corey. 2019. "Analytic Expressions for Radar Sea Clutter WSSUS Scattering Functions" Entropy 21, no. 9: 915. https://doi.org/10.3390/e21090915

APA Style

Cooke, C. (2019). Analytic Expressions for Radar Sea Clutter WSSUS Scattering Functions. Entropy, 21(9), 915. https://doi.org/10.3390/e21090915

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

Article Metrics

Back to TopTop