# Development of a Near-Field Bistatic Synthetic Aperture Radar for Complex Target Reconstruction

Authors: David G. Johnson; Graham M. Brooker
Journal: International Journal of Antennas and Propagation
 (2012)

---

## Abstract

This paper begins with a description of the design, construction, and characterization of a small electromagnetic anechoic chamber, developed specifically to house a bistatic ISAR system for the analysis of rock samples. Particular emphasis is given to the practicalities of construction, with the intention of assisting those in a similar position, wishing to build an anechoic chamber on a tight budget. The second part of the paper outlines efficient algorithms that may be applied to the tomographic and topographic reconstruction of complex targets within the viewing geometry of this ISAR system.


---

## Body

## 1. Introduction

This paper describes the development of a radar system designed for the 3D topographic reconstruction of unstructured surfaces. This was motivated by the need within the open-pit mining industry to determine the precise size distribution of rock fragments postblasting and precrushing. In these environments laser and camera systems are often unsuitable due to the presence of large amounts of dust. As a solution to this problem, alinear-frequency-modulated continuous-wave(LFM-CW) radar-based approach was envisioned that would be able to operate at short range, and provide sufficient range and angular resolution to compete with fair-weather optical systems. A coherent near-field multistatic methodology was chosen using synthetic-aperture-radar (SAR) techniques to minimise the number of moving parts, whilst providing the greatest number of measurement degrees of freedom and the best possible resolution [1]. A prototype bistatic system, operating in inverse-SAR(ISAR) mode within the well-behaved confines of an anechoic chamber, was therefore constructed to determine whether the required performance characteristics were achievable, prior to the development of a full-scale system on a platform, such as that shown in Figure 1. This paper describes solely the theory and measurements for this prototype system. The extension to the real world must be accompanied by further research. It should be noted that both the anechoic chamber and prototype radar system were built for a combined cost of US $35k. This compares to the quoted cost of a full professional installation of the chamber alone of in excess of US $80k.Figure 1
Artists impression of radar system on-board an electric-rope-shovel, showing the approximate radar imaging region in grey for the given bistatic antenna configuration. Further receivers could then be added along the shovel arm to extend this to a full multistatic system.Multiple papers and books have been written on the subject of anechoic chamber design for standard RCS and antenna measurements, for example, [2, 3]. However, little information can be found on the construction and installation processes, particularly when working in the near-field. The paper, therefore, includes a number of innovative cost and labour saving techniques for near-field test-facility design, test, and measurement.In the radiative near field,0.6(D3/λ)1/2≤R≤2D2/λ, as defined by the maximum extent of both the radar aperture and the targets themselves, D, relative to the wavelength,λ, plane-wave approximations become inaccurate. Spherical wave functions must then be applied to solve the system imaging problem, extending the ISAR methodology described previously to a spherical-wave form, known as SWISAR [4]. The computational burden of solving these complex functions in the Fourier domain can then be reduced by employing the method of stationary phase (MoSP). This algorithm approximates the dominant signal from each target as the region over which its phase is slow-moving or stationary, generally found around the point of closest (or furthest) approach to the radar. Whilst this is commonly described by a quadratic expansion of the phase, as in [5, 6], it can be, and for the full bistatic case, must be, more accurately described by considering higher-order polynomials.ISAR Fourier imaging is a form of convolution back-projection [7] and is, therefore, well suited to coherently processing large angular sectors, as the full azimuth spectral bandwidth must be calculated in any case. This is particularly of use at short range where a target’s translation leads to a larger rotational movement. One of the principal difficulties that then arise is the aspect dependence of dispersive target scattering. This can be dealt with either by assuming an anisotropic (i.e., non-point-like) scatterer model, or, as is described here, by applying a Fourier-domain-matched filter-bank tailored to the local curvature of the (unknown) target. In this sense, there is a strong correspondence between the current work and the use of microlocal-analysis techniques [8–12], based on geometric optics and the MoSP, which identifies the coherent scattering centres of known topography and provides an analytic solution to the subsequent imaging problem. A priori known topography is also utilised in a similar manner by [13], to identify and equate coherent regions in the Fourier domain, observed from spatially separated positions, in order to optimise a coherent change detection (CCD) process. The current work may also prove useful in terms of providing an experimental means of justifying these techniques and extending the same ideas to target detection and the mapping of unknown topography.Section2 of this paper describes the development and characterisation of the prototype SWISAR system and the near-field test chamber that was custom-built to hold it, whilst Section 3 describes measures taken to achieve a 1D range resolution of <20 mm. Sections 4 and 5 then describe algorithms to perform focussed multistatic 3D imaging of point and spherical targets, respectively, using banks of spatially-variant Fourier-domain-matched filters, extending the work of [14].


## 2. Development and Characterisation of the Near-Field Test Chamber

A small anechoic chamber (2.9 m × 2.5 m × 2.3 m) was built at theAustralian Centre for Field Robotics(ACFR) as described in [15] and shown in Figure 2(a). The main features of this facility are as follows.(i)
A polystyrene/steel sandwich wall construction [16], with seams sealed using AT526 35 micron copper foil shielding tape and SOFT-SHIELD 2000 gaskets; lined with 300 mm pyramidal foam radar absorbing material (RAM). Some simple construction tips and a RAM mitre-joint cutting template are provided in the Appendix.(ii)
A pallet-mounted (and removable) 1.2 m diameter turntable. The turntable was constructed from 20 mmmedium density fibreboard (MDF), with the drive mechanism consisting of a geared stepper motor and a timing-belt pulley made from a bicycle wheel. An important feature of the turntable design (owing to the high torque requirement and temperature variations within the lab) was the belt-tensioners. These allowed occasional adjustment of belt tension (using a rubber mallet) to reduce belt-skip. The trucks from two skateboards provide additional load-bearing support. When operated under correct tension, the turntable was capable of rotating 300 kg rock samples with an angular accuracy of <0.01°.(iii)
A 90° sector arc, of radius 2 m, on which two independently powered trolleys could transit, each carrying a 2–18 GHz antenna linked to the radar hardware described shortly, and a 5 mega-pixel remote-operated camera (with flash). The compact cameras (Olympus FE130) were particularly useful for the formation of bistatic photographs, employing the flash from the Tx-mounted camera to illuminate the scene viewed by the Rx-mounted camera, from which bistatic scattering centres at optical frequencies could be identified and related to similar topology radar measurements.(iv)
Radar hardware installed directly above the chamber to minimise cable lengths (which pass through filtered connectors in the ceiling).(v)
Double doors provide the ability to obtain far-field measurements at MMW using the full length of the ACFR field-research-laboratory (~40 m). For an example of this functionality, see [17].(a) Cutaway schematic of anechoic chamber showing transmit and receive beams. (b) Schematic of turntable. (c) View of anechoic chamber through open doors. (d) Close-up view of antenna trolley showing camera and antenna.
(a)(b)(c)(d)At each stage in construction, measurements of reflectivity from and transmission through the chamber walls were made at H- and V-polarisations to maximise the anechoic performance of the chamber. The transmission tests involved transmitting a swept 1–19 GHz, +10 dBm signal, from a Marconi 6313 synthesiser, via 1 m of Ultraflex FM402 flexible coax, to a Q-Par Angus Broadband Horn Antenna (WBH218S) placed outside the chamber 1 m from the wall. The signal was then received by an identical antenna 1 m inside the chamber and passed via 8 m of FM402 cable to an HP8563E Spectrum Analyser. Both the Sweep Generator and Spectrum Analyser were controlled synchronously by a remote PC, thus forming a basic scalar network analyser. Reflectivity measurements were made in the same manner, with the antennas placed side by side inside the chamber. As an example, Figure3 displays the V-pol. side-wall transmission measurements (without compensation for antenna gain).Figure 3
Anechoic chamber characterisation results. Measurements were made with the transmitter outside the chamber and the receiver inside, each aligned orthogonally to the chamber wall.


## 3. Nonlinearity Correction for cm Level Resolution

Figure4 shows a schematic of the zero-IF (homodyne) LFM-CW radar and the power levels expected over the frequency of operation. The advantage of this configuration is the small number of components required, which minimises both cost and system losses. The recent development of multi-octave “ultra” broadband components for use in electromagnetic test chambers and electronic warfare devices has opened the door to radar system bandwidths in excess of 15 GHz. All radio-frequency (RF) components (including cables and connectors) are, therefore, rated from 2–18 GHz, giving a theoretical range-resolution of 9.4 mm. In reality, the amplitude and phase linearity of the components within this configuration is not perfect, leading to some degradation in resolution.Figure 4
Schematic of radar hardware showing power levels detected at various points.Although the range-resolution of LFM-CW radar is in theory only limited by its swept bandwidth, it is more realistically limited by the linearity of the frequency chirp and the spectral purity of the instantaneous transmitted signal. Before the radar data can be passed to the SWISAR algorithm for conversion to a 3D image, it must therefore be purified, that is, preprocessed to remove hardware-induced artefacts such as leakage signals, phase nonlinearities, and system losses. Numerous methods exist to linearise the frequency chirp in hardware prior to transmission. These are generally based on Phase Locked Loops (PLLs) and/or false-target calibration [18]; however, these methods are still hardware-bandwidth-limited and can only provide a certain level of improvement. Also, although many of our radar’s hardware components operate over the full 2–18 GHz band, there are no commercially available oscillators that can sweep continuously over this region, and the switching between oscillators introduces severe nonlinearities that cannot be corrected by a PLL. Indeed, over such a large swept bandwidth, nonlinearities in both the antennas and cables emerge that cannot be corrected for in hardware, therefore a software linearisation solution must be found. Other aspects of chamber calibration that should be taken into account relate to imperfections in the geometry of the measurement system. These include variations in the curvature of the steel arc on which the transmit- and receive-antenna trolleys translate, inertial sensor errors due to timing belt slippage on both of the antenna trolleys and the turntable, and the possible misalignment of the turntable beneath the centre of the steel arc.The flowchart in Figure5 describes the processing steps required to obtain a set of consistent, linearised radar samples.Figure 5
Flowchart of signal linearisation process.


### 3.1. Individual Sweep Extraction/Empty Chamber Subtraction

Despite the use of a RAM baffle installed on the receive antenna, the direct-path leakage signal between the transmit and receive antenna dominates the measured data and must be removed prior to the subsequent preprocessing steps. This is a common technique, often employed in the form of areflective power canceller, as described in [19]. Because the leakage signal remains constant irrespective of the sample within the chamber, it provides a convenient measure of both the consistency of empty-chamber-subtraction for a fixed antenna configuration, and a means of relocating a particular bistatic configuration to great accuracy. The relative positions of the transmit and receive antennas along the steel arc are measured by optical encoders on the pulley wheel that drives each trolley along the timing belt affixed to the steel-arc itself. However, due to the weight of the trolleys and the steepness of the incline in places, it is not unknown for the pulley to slip on the timing belt, which goes unnoticed by the proprioceptive sensor. To account for this, empty chamber subtraction is employed to measure the bistatic leakage distance and accurately locate the receive antenna (assuming a fixed transmit antenna position). The graph in Figure 6 is made up of range profiles corresponding to 350 receiver positions over a range of 70°, that is, steps of 0.2°, with the transmitter fixed at the end-stop of the arc. The receive antenna can then be relocated to any of these 350 positions by determining the best level of leakage-signal suppression through empty chamber subtraction. The system of antenna localisation was found to be accurate to within ±0.1° of arc as described. Whilst having cm-level two-target resolution, it was found that after phase-linearisation the system could follow the phase progression of a single dominant target to sub-mm accuracy.Figure 6
Empty chamber measurements for varying bistatic angles. Colour represents signal intensity in dBm.


### 3.2. Bandpass Filtering

To further improve the subsequent nonlinearity-measurement and normalisation steps, the “chamber-subtracted” signal is bandpass filtered to isolate the region of interest defined by the range-extent of the turntable.


### 3.3. Normalisation

Following an analysis of cable and antenna losses between 2 and 18 GHz, as described in [15], a log-linear power gain-correction factor of 0–25 dB is applied to the filtered signal at this point. Quantisation noise then reduces the effective-number-of-bits (ENOB) of the ADC to ~8.5 bits at 18 GHz; however, this is still adequate for the subsequent linearisation stage.


### 3.4. Interoscillator Phase-Shift Correction

The Marconi 6313 sweep generator is made up of a number ofyttrium-iron-garnet (YIG) tuned oscillators [20], which although inherently linear, experience a band switching delay of 0.5 μsec. This switching between oscillator sections leads to a significant, although consistent, non-linearity in phase that can initially be dealt with by applying a constant phase shift to data within each of the upper (8–12.4 GHz and 12.4–18 GHz) oscillator bands.A technique has recently been put forward in the SAR literature by Meta et al. [21], claiming to solve the non-linearity problem entirely in software, although it does require precise knowledge of the original non-linearity. For the case of the anechoic environment, this can be extracted from an isolated test target (such as a sphere). The signal-phase from this target is then unwrapped (following the pre-processing steps described above), which is a simple matter for a strong complex signal, and a linear fit of the result removed to produce a residual phase characteristic as seen in Figure 7. Because this signal includes the effects of target location and antenna configuration, this process is repeated for a number of movements and an average non-linearity is derived.Determination of LFM-CW system non-linearity using phase fill-in and unwrapping applied to an isolated target.
(a)(b)(c)(d)After theinteroscillator phase-shift correction process, there then remains a period of irregular phase caused by variations in amplitude as the oscillators are switched, which is too nonlinear to be dealt with by Meta’s algorithm. However, by applying a phase fill-in procedure, as depicted in Figure 7, that is, replacing the irregular phase section by a signal with a linear phase variation based on the well-behaved portion of the signal.The Meta technique then expresses the nonlinear beat signal of a target received after a time delay,τ, as the product of the ideal linear beat signal, sif,lin(t), and the non-linear error function, sϵ(t), applied both on transmission and reception:(1)sif(t)=ej2π(fcτ+αtτ+(1/2)ατ2+ϵ(t)-ϵ(t-τ))=ej2π(fcτ+αtτ+(1/2)ατ2)ej2πϵ(t)e-j2πϵ(t-τ)=sif,lin⁡(t)sϵ(t)sϵ(t-τ)*,
where (i) fc; The central carrier frequency of the transmitted signal. (ii) α; The slope of the frequency chirp (Hz/s). (iii) ϵ(t); The time dependent non-linearity applied to the voltage ramp.Assuming the error function,sϵ(t), is well known, the ideal linear beat signal can be recovered by a 3-step process:(1)
multiplying the non-linear beat signal by the conjugate of the error function,sϵ(t)*, to remove the effects of the transmitted non-linearity;(2)
applying a range-dependent time shift to the received signal to remove the range/time dependence of the received non-linearity. This is performed by aresidual video phase (RVP) removal filter in Fourier space;(3)
multiplying the filtered signal by the conjugate of the error function after it too has been passed through the RVP filter.The process described in the flowchart of Figure5 is now shown for the example of two metallic spheres placed in isolation within the anechoic chamber. In Figure 8, each step of the phase-linearisation process is shown with the final result that the two spheres, with a bistatic range separation of 100 mm, can be seen clearly 15 dB above the next spurious peak. Here, values in dBm relate to the direct measurement made into a 50 Ω load, and changes in antenna gain and cable losses with frequency are subsumed into the empty chamber subtraction step.Figure 8
Signal linearisation of two targets closely separated in range.The two spheres shown in Figure9(a) have then moved to produce a closer bistatic range separation of 20 mm and still have clearly separated peaks as seen in Figure 9(b), although there is a 3 dB difference in intensity caused by the variation in gain of both the transmitter and receiver [22].(a) Two 150 mm spherical targets within an otherwise empty chamber. (b) Linearised results of both targets in isolation, and together.
(a)(b)It should be noted that a side-effect of the Meta linearisation technique is to add ~2 dBm to all dominant targets, as energy from the error function is essentially added in to the linearised result. This can be seen in both Figure8 (inset) and Figure 9(b).


## 3.1. Individual Sweep Extraction/Empty Chamber Subtraction

Despite the use of a RAM baffle installed on the receive antenna, the direct-path leakage signal between the transmit and receive antenna dominates the measured data and must be removed prior to the subsequent preprocessing steps. This is a common technique, often employed in the form of areflective power canceller, as described in [19]. Because the leakage signal remains constant irrespective of the sample within the chamber, it provides a convenient measure of both the consistency of empty-chamber-subtraction for a fixed antenna configuration, and a means of relocating a particular bistatic configuration to great accuracy. The relative positions of the transmit and receive antennas along the steel arc are measured by optical encoders on the pulley wheel that drives each trolley along the timing belt affixed to the steel-arc itself. However, due to the weight of the trolleys and the steepness of the incline in places, it is not unknown for the pulley to slip on the timing belt, which goes unnoticed by the proprioceptive sensor. To account for this, empty chamber subtraction is employed to measure the bistatic leakage distance and accurately locate the receive antenna (assuming a fixed transmit antenna position). The graph in Figure 6 is made up of range profiles corresponding to 350 receiver positions over a range of 70°, that is, steps of 0.2°, with the transmitter fixed at the end-stop of the arc. The receive antenna can then be relocated to any of these 350 positions by determining the best level of leakage-signal suppression through empty chamber subtraction. The system of antenna localisation was found to be accurate to within ±0.1° of arc as described. Whilst having cm-level two-target resolution, it was found that after phase-linearisation the system could follow the phase progression of a single dominant target to sub-mm accuracy.Figure 6
Empty chamber measurements for varying bistatic angles. Colour represents signal intensity in dBm.


## 3.2. Bandpass Filtering

To further improve the subsequent nonlinearity-measurement and normalisation steps, the “chamber-subtracted” signal is bandpass filtered to isolate the region of interest defined by the range-extent of the turntable.


## 3.3. Normalisation

Following an analysis of cable and antenna losses between 2 and 18 GHz, as described in [15], a log-linear power gain-correction factor of 0–25 dB is applied to the filtered signal at this point. Quantisation noise then reduces the effective-number-of-bits (ENOB) of the ADC to ~8.5 bits at 18 GHz; however, this is still adequate for the subsequent linearisation stage.


## 3.4. Interoscillator Phase-Shift Correction

The Marconi 6313 sweep generator is made up of a number ofyttrium-iron-garnet (YIG) tuned oscillators [20], which although inherently linear, experience a band switching delay of 0.5 μsec. This switching between oscillator sections leads to a significant, although consistent, non-linearity in phase that can initially be dealt with by applying a constant phase shift to data within each of the upper (8–12.4 GHz and 12.4–18 GHz) oscillator bands.A technique has recently been put forward in the SAR literature by Meta et al. [21], claiming to solve the non-linearity problem entirely in software, although it does require precise knowledge of the original non-linearity. For the case of the anechoic environment, this can be extracted from an isolated test target (such as a sphere). The signal-phase from this target is then unwrapped (following the pre-processing steps described above), which is a simple matter for a strong complex signal, and a linear fit of the result removed to produce a residual phase characteristic as seen in Figure 7. Because this signal includes the effects of target location and antenna configuration, this process is repeated for a number of movements and an average non-linearity is derived.Determination of LFM-CW system non-linearity using phase fill-in and unwrapping applied to an isolated target.
(a)(b)(c)(d)After theinteroscillator phase-shift correction process, there then remains a period of irregular phase caused by variations in amplitude as the oscillators are switched, which is too nonlinear to be dealt with by Meta’s algorithm. However, by applying a phase fill-in procedure, as depicted in Figure 7, that is, replacing the irregular phase section by a signal with a linear phase variation based on the well-behaved portion of the signal.The Meta technique then expresses the nonlinear beat signal of a target received after a time delay,τ, as the product of the ideal linear beat signal, sif,lin(t), and the non-linear error function, sϵ(t), applied both on transmission and reception:(1)sif(t)=ej2π(fcτ+αtτ+(1/2)ατ2+ϵ(t)-ϵ(t-τ))=ej2π(fcτ+αtτ+(1/2)ατ2)ej2πϵ(t)e-j2πϵ(t-τ)=sif,lin⁡(t)sϵ(t)sϵ(t-τ)*,
where (i) fc; The central carrier frequency of the transmitted signal. (ii) α; The slope of the frequency chirp (Hz/s). (iii) ϵ(t); The time dependent non-linearity applied to the voltage ramp.Assuming the error function,sϵ(t), is well known, the ideal linear beat signal can be recovered by a 3-step process:(1)
multiplying the non-linear beat signal by the conjugate of the error function,sϵ(t)*, to remove the effects of the transmitted non-linearity;(2)
applying a range-dependent time shift to the received signal to remove the range/time dependence of the received non-linearity. This is performed by aresidual video phase (RVP) removal filter in Fourier space;(3)
multiplying the filtered signal by the conjugate of the error function after it too has been passed through the RVP filter.The process described in the flowchart of Figure5 is now shown for the example of two metallic spheres placed in isolation within the anechoic chamber. In Figure 8, each step of the phase-linearisation process is shown with the final result that the two spheres, with a bistatic range separation of 100 mm, can be seen clearly 15 dB above the next spurious peak. Here, values in dBm relate to the direct measurement made into a 50 Ω load, and changes in antenna gain and cable losses with frequency are subsumed into the empty chamber subtraction step.Figure 8
Signal linearisation of two targets closely separated in range.The two spheres shown in Figure9(a) have then moved to produce a closer bistatic range separation of 20 mm and still have clearly separated peaks as seen in Figure 9(b), although there is a 3 dB difference in intensity caused by the variation in gain of both the transmitter and receiver [22].(a) Two 150 mm spherical targets within an otherwise empty chamber. (b) Linearised results of both targets in isolation, and together.
(a)(b)It should be noted that a side-effect of the Meta linearisation technique is to add ~2 dBm to all dominant targets, as energy from the error function is essentially added in to the linearised result. This can be seen in both Figure8 (inset) and Figure 9(b).


## 4. A Coherent Multistatic Near-Field Imaging Algorithm for Point Targets

This section provides a brief overview of the bistatic SWISAR imaging algorithm for point targets, with an accompanying metric of performance for the case of a set of disparate sized spheres. Further details of the monostatic point-target algorithm can be found in [5], whilst the 3D multi/bistatic imaging algorithm applied to a selection of rocks has been previously presented in [14, 23].Figure10 depicts a bistatic (BS) SWISAR geometry, but can equally be employed to visualise the monostatic (MS), or quasi-monostatic (QMS) in which case a single, or pair of, antenna(s) centred on “Tx” will simply use R=Rt in the formulae for the spatially variant matched filter (SVMF), F, for each voxel in the cylindrical-polar imaging coordinate system, following Fortuny [4]:(2)I(ρ,ϕ,z)=8c3∫ff2∫θsin⁡⁡θ×∫ϕ′Es(f,ϕ′,θ)F(ρ,ϕ-ϕ′,z;f,θ)dfdθdϕ′,
where(3)F(ρ,ϕ,z;f,θ)=(RRa)2exp⁡[jkf(R-Ra)],(4)R=Ra2+ρ2+z2-2Ra[ρsin⁡⁡θcos⁡⁡ϕ+zcos⁡⁡θ].
Equation, (2) describes a circular convolution and can therefore be calculated by a 1D azimuth discrete-Fourier-transform (DFT) and matched-filter multiply. Fortuny observed that the azimuth Fourier transform of (3) can be rewritten (for fixed f and θ) in the form of a 1D Fourier-integral-operator (FIO) focussing function, ℱ(·):(5)F(ρ,kϕ,z;f,θ)=∫-ππ(R(ϕ)Ra)2exp⁡[jp(ϕ)]dϕ,
where, for the monostatic case, R(ϕ)=RMS(ϕ)and p(ϕ)=pMS(ϕ)(6)pMS(ϕ)=A[RMS(ϕ)-Ra]-kϕϕ,RMS(ϕ)=B-Ccos⁡⁡ϕ,(7)A=kf=4fc,B=Ra2+ρ2+z2-2Razcos⁡⁡θMS,C=2Raρsin⁡⁡θMS.
The MoSP can then be used to determine analytically the asymptotic approximation to this azimuth DFT, as for a rapidly oscillating FIO, such that describing the phase history of a point target on a rotating turntable, the major contribution to the integral comes from the regions around the stationary phase points.Figure 10
Bistatic SWISAR geometry (for which the monostatic geometry is a subset).For the case when two 1st-order stationary phase points are present atϕ- and ϕ+,(8)pMS′(ϕ-)=∂pMS(ϕ)∂ϕ|ϕ-=ACsin⁡⁡ϕ-2B-Ccos⁡⁡ϕ--kϕ=0,pMS′(ϕ+)=∂pMS(ϕ)∂ϕ|ϕ+=ACsin⁡⁡ϕ+2B-Ccos⁡⁡ϕ+-kϕ=0.
This can be solved algebraically to give(9)ϕ±=arccos[2α2±β2-4α2(1-α2)β],
where (10)α=kϕAB,β=CB.
Above a certain threshold of kϕ, the two stationary phase points given in (9) coalesce to a single point ϕ0 given by:(11)ϕ0=arccos[1-1-β2β].
These expressions can then be resubstituted when equating the higher derivatives of pMS(ϕ)at the stationary phase points, as they relate to the asymptotic approximations for the focussing template of (5):(12)FMS(ρ,kϕ,z;f,θ)≍j2πRa2(RMS2(ϕ-)pMS′′(ϕ-)ejpMS(ϕ-)+RMS2(ϕ+)pMS′′(ϕ+)ejpMS(ϕ+)),FMS(kϕ)≍2πejp(ϕ0)(RMS(ϕ)2Ra2)(2|pMS′′′(ϕ0)|)1/3×Ai(-[2|pMS′′′(ϕ0)|]1/3pMS′(ϕ0)),
where Ai(·)denotes the Airy function [24]. An extended description of this derivation can be found in the Appendix A of [15]. Because these expressions involve only A, B, and C, they are independent of any target properties and need only be regenerated if the sensing geometry is altered.The meaning of the two regions ofkϕin (12) can be explained with reference to Figures 11–13 below.Four representations of data corresponding to a single simulated point target from a single bistatic antenna pair: (a) the “Fourier” domain, (b) the “raw-data” domain, (d) the “time” domain, and (c) shown for completeness but not used anywhere.
(a)(b)(c)(d)Figure 12
Geometric description of boundaries ofk-space.Time domain solutions for a single-simulated point target corresponding to the use of (a):IFFTkϕ,kfof first-order (Kϕ≤K̂ϕ) Fourier domain region only, (b):IFFTkϕ,kfof second-order (Kϕ>K̂ϕ) Fourier domain region only and, (c): IFFTkϕ,kfof full kϕFourier domain regions. (The units, although not important to the figure, are nanoseconds on the x-axis and degrees on the y-axis.)
(a)(b)(c)Figure11 depicts the four time/frequency domain combinations available for SAR processing. The top-left “Fourier domain” matched-filter template corresponds to a target at a particular position as seen in the bottom-right. Also obvious, is a triangular region, that is itself a subset of a rectangular region (which for the most part is some 50 dB weaker).Targets positioned at points A1 and A2 in Figure12 will fall along the line kϕ=0, whilst the points B1 and B2 refer to an area within the triangular region itself.The boundary of the triangular region in Figure11(a) relates to the point kϕ≤k̂ϕ. If the antenna was in the planez=0, that is, θ=90°, then this would correspond to the maximum spatial bandwidth of the geometry and the point C in Figure 12. Outside of this triangular region,kϕ>k̂ϕ would then correspond to an evanescent-wave region and hence provide negligible response for ranges much greater than the wavelength of the radar. Because we are dealing with a slant-plane circular-SAR geometry, however, this region does have an affect on the subsequent time-domain solution, as seen in Figure 13.When the transmitter and receiver are separated along the arc, as seen in Figure10, applying the monostatic to bistatic equivalence theorem (MBET) works well for small bistatic angles, that is, ≲6° [25]. that is, by assuming that the bistatic-bisector angle, θb, can be used to provide a quasimonostatic (QMS) estimation of θ in (5). For larger bistatic separations, this assumption can no longer be made.Figure14 and Table 1 describe two equivalent experimental setups, which will be used to explain and in the subsequent text demonstrate a full bistatic SWISAR formulism.Table 1
Parameters for simulated and experimental SWISAR imaging.
ParameterSymbolValueNumber of MS/BS antenna configurationsN12Experiment
1: (Quasi-) Monostatic configurationMSn Tx antenna angular range (n=1:12)θt,MSn39° : −2° : 17°MSn Rx antenna angular range (n=1 : 12)θr,MSn33° : −2° : 11°QMS1n-equivalent angle =(θt,MSn+θr,MSn)/2θQMS1,n36° : −2° : 14°Experiment
2: Bistatic configurationBS0 Tx antenna angle (fixed)θt,BS039° (fixed)BSn Rx antenna angular range (n=1 : 12)θt,BSn33° : −4° : −11°QMS2n-equivalent angle =(θt,BS0+θr,BSn)/2θQMS2,n36° : −2° : 14°Other parameters:Frequency rangef2–18 GHzTurntable angular rangeϕ0° : 0.5° : 359.5°Radial coordinate rangeρ0 : 0.01 : 0.7 mHeight coordinate rangez−0.2 : 0.01 : 0.4 mCartesian coordinate rangex,y−1 : 0.01 : 1 mEquivalent QMS paths shown for (a) a pair of closely located antennas moving in sync across the arc and (b) a stationary transmitter and a moving bistatic receiver.
(a)(b)In the bistatic SWISAR context, to produce a dataset equivalent to the QMS antenna-pair measurements shown in Figure14(a), for which a pair of antennas are moved synchronously across the arc, requires the experimental geometry of Figure 14(b), using a fixed transmitter and a receiver moved at twice the angular-increments. In the most extreme case, with values of θt,BS0=+39° and θr,BS12=-11°(with θ=0 at the zenith), the bistatic pair of antennas span the centre of the arc with a bistatic-angle of 50°. This is a fairly reasonable requirement if an eventual multi-static radar is required to obtain measurements over a wide θ-range, operating with multiple static receivers rather than the single dynamic unit used by the prototype. The validity of applying a point-target approximation over this large bistatic angle is discussed later.To form a bistatic asymptotic expansion for the focussing function,ℱBS, similar to that given in (12) for the (quasi-) monostatic case,ℱMS, requires that the stationary phase points (ϕ± or ϕ0, depending on the size of kϕ) be obtained in an accurate and efficient manner. For the bistatic case, (5) still applies; however, (6) must be adjusted to take the full bistatic geometry into account:(13)pBS(ϕ)=A[Rt(ϕ)2+Rr(ϕ)2-Ra]-kϕϕ,Rt,r(ϕ)=Bt,r-Ct,rcos⁡⁡ϕ,
where(14)Bt,r=Ra2+ρ2+z2-2Razcos⁡⁡θt,r,Ct,r=2Raρsin⁡⁡θt,r.
Splitting p(ϕ)into the contributions from each antenna has the effect of complicating the solution of the stationary-phase points in the bistatic case, as the first-order derivative then hassinϕ terms for both transmitter and receiver simultaneously. In fact, within the region of Kϕ-space corresponding to two unique stationary points, these must be solved by a 6th-order polynomial equation in ϕ, whose full expansion has been shown in [15].Clearly the calculation of the polynomial coefficients adds to the complexity of the bistatic SWISAR algorithm compared to the monostatic case, particularly as the roots must then be solved approximately byeigenvalue-analysis [26] and then verified for consistency by substitution into the first-order derivative equation:(15)pBS′=A2[Ctsin⁡⁡ϕ±2Bt-Ctcos⁡⁡ϕ±+Crsin⁡⁡ϕ±2Br-Crcos⁡⁡ϕ±]-kϕ<δp′,
where δp′is an arbitrarily small threshold. The result of this procedure is, however, considerably more accurate than the QMS solution, as seen in Figure 15, and is significantly faster than computing the full DFT, particularly if the resolution in the ϕ domain is high, as is required by a large radial dimension.SVMF phase (and derivatives) for a large bistatic angle.
(a)(b)(c)(d)In the instance where no results forϕ±<δp′ can be found from the sixth-order polynomial, the second-order derivative of the “phase” term tends to zero at a single root as before, which for the bistatic geometry now becomes(16)pBS(ϕ)′′=A(-Ct2sin⁡2⁡ϕ8Rt(ϕ)3+Ctcos⁡⁡ϕ4Rt(ϕ)-Cr2sin⁡2⁡ϕ8Rr(ϕ)3+Crcos⁡⁡ϕ4Rr(ϕ))=0.
Following the same procedure as before, (16) has a seventh-order polynomial-root solution, allowing the second- and third-order phase derivatives to be calculated explicitly at these points.The algorithm then proceeds as before, with the appropriate bistatic terms substituting those for the monostatic case in (12) to form the focussing template,ℱBSfor each frequency, f, and look-angle, θb (which does still make use of the bistatic bisector angle, without any noticeable side-effects). The reflectivity of the target volume in cylindrical-polar coordinates can then be recovered as follows:(17)IBS(ρ,ϕ,z)=8c3IFFTkϕ[∑ff2∑θbsin⁡⁡θbFFTϕ[EBS(f,ϕ,θb)]×[FBS]kϕ∑ff2∑θbsin⁡⁡θbFFTϕ[EBS(f,ϕ,θb)]]ΔfΔθbΔϕ,
where[ℱBS]kϕ=ℱBS(ρ,kϕ,z;f,θt,θr)for Kϕ=-Nϕ/2,-Nϕ/2+1,…,Nϕ/2-1.To demonstrate the improvement in applying thebistatic (BS) approach to determineℱBS, over the existing ℱMSapproach used in [4], the example shown in Figure 15 uses the extreme Tx/Rx positions listed previously.The uppermost phase plot in Figure15 is formed with a full expansion of the DFT described by (5), (13), and (14). The blue lines in the remaining plots are then simple derivatives of this term, while the red line in the 2nd plot shows the value of the phase derivative assuming that the bistatic-bisector angle, θb, is used to provide a quasimonostatic approximation to the geometry. Admittedly this is far outside the expected region of validity for the bistatic angle of 6° given by the MBET. For the particular values of kϕ=25, and f=10GHz, two 1st-order stationary points can be seen, requiring the solution and validation of 6th-order polynomial roots. The dashed black line in the 2nd plot then shows a 6th order polynomial expansion of (15), two roots of which, shown by the blue markers, correspond exactly with the DFT solution. The variation in position of these stationary-phase points with the (red) QMS phase derivative is obvious from the markers.As stated previously, in order to test the volumetric imaging performance of the full Bistatic-SWISAR algorithm, the prototype ISAR system was configured with two geometries as described by Table1 and Figure 14. A set of spherical targets was then placed within the anechoic chamber as seen in Figure 16 and described by Table 2. Simulation and experimental results relating to these two geometries and common target set are then shown in Figures 17, 18, 19.Table 2
Parameters of spherical targets used in point target simulations. based on the graphically estimated true target positions seen in Figure16.
Target no.X, mY, mZ, mDiameter, m1−0.2−0.110.150.1530.110.090.1250.1080.065−0.1650.10.05Figure 16
Experimental setup used to test monostatic (MS), quasimonostatic (QMS), and bistatic (BS) SWISAR algorithms, and the basis for equivalent simulated results.Simulated point-spread-function results for a single target at [X,Y,Z]=[0,0.5,0.05] cut through each plane. (a) (Q)MS antennas, (Q)MS algorithm. (b) BS antennas, (Q)MS algorithm. (c): BS antennas, BS algorithm.
(a)(b)(c)Simulated tomographic reconstruction of the spherical targets listed in Table 2. (a) (Q)MS antennas, (Q)MS algorithm. (b): BS antennas, (Q)MS algorithm. (c) BS antennas, BS algorithm.
(a)(b)(c)Experimental tomographic reconstruction of the targets seen in Figure16. (a) (Q)MS antennas, (Q)MS algorithm. (b) BS antennas, (Q)MS algorithm. (c): BS antennas, BS algorithm.
(a)(b)(c)A summary of the experiments and a brief analysis of their main features are shown in Table3.Table 3
Description for point-target SWISAR experiments.
Simulated point target at[X,Y,Z]=[0,0.5,0.05]Simulated spherical targetsExperimental spherical targetsMono Static (MS)Figure17(a)Figure18(a)Figure19(a)θt= 39° : −2° : 17° θr= 33° : −2° : 11° θMS= 36° : −2° : 14° Computation time = 34 minsFocussed to a single point. Effects of limited aperture inθclearly visible.Clear correspondence to upper (visible) surface of spheresSome residual energy in annuliClear correspondence to upper (visible) surface of spheresQuasi Mono Static (QMS)Figure17(b)Figure18(b)Figure19(b)θt= 39° : −2° : 17° θr= 33° : −2° : 11° θMS= 36° : −2° : 14° Computation time = 34 minsPoorly focussedSmeared target responseGreater residual energy in annuli.Peak energy density away from target surfaceBi Static (BS)Figure17(c)Figure18(c)Figure19(c)θt= 39° (fixed)  θr= 33° : −4° : −11°  θBS= 36° : −2° : 14°Computation time = 7.5 hrsResult “as good as” that of MSResult “as good as” that of MSResult “as good as” that of MSThe (white) “dead-pixels” seen both at theX-Y origin and for large values of both X and Y in all images are due to the postprocessing interpolation from cylindrical-polar to Cartesian coordinates. The limited variation in θ also leads to strong sidelobes forming annuli above and below the target.When the point target SWISAR algorithm is applied to nonpoint targets, the energy is distributed around the surface of that target (corresponding to the specular region of that target at each rotation angle).To allow these algorithms to be compared more accurately, a metric is required that determines how closely the imaged data conforms to “ground-truth.” Such a metric is given by the average power density over concentric spheres whose origins correspond to the targets of the multi-ball scenario, as seen in Figure 16.The results of applying this metric are seen in Figure20. The simulated MS and BS results correctly place the peak intensity at the desired result (within measurement error). The simulated QMS results meanwhile show that the classified spherical radius, as determined by the peak intensity, is overestimated by 5–10 mm for each target. The peak-responses for the MS and BS results are also significantly narrower than the equivalent QMS results (by ~50%). The hump in intensity for a spherical radius of approximately 0.25 m, for all MS and BS simulated targets, is believed to be due to the annuli artefacts earlier discussed, and not the other targets in the image.Spherical-surface classification performance of point-target SWISAR algorithms. Left: simulated and, Right: experimental data. (a, b) (Q)MS antennas, (Q)MS algorithm. (c, d) BS antennas, (Q)MS algorithm. (e, f) BS antennas, BS algorithm.
(a)(b)(c)(d)(e)(f)Simulated data for targets 1, 3 and 8 (on the left) are compared to the experimental results for all 8 targets (on the right). These results correspond to concentric spheres whose radii increase in increments of 5 mm. The small (high-frequency) oscillations seen in the plots are, therefore, due to the discrete nature of the imaged data on a 10 mm voxel 3D grid.


## 5. A Coherent Multistatic Near-Field Imaging Algorithm for Spherical Targets

For a fixed geometry, determining focussing templates “on the fly” is a far more time consuming process than retrieving precomputed templates from a memory store (by up to three orders of magnitude). This being the case, the templates may as well be created by any tractable means. Taking this into consideration, the simulator used to test the SWISAR algorithm in the previous section was extended to allow the formation of templates corresponding to arbitrary shapes placed anywhere within the target volume. A spherical target, made up of distributed point-elements, is simulated at each required rotation angle around the turntable, and the (geometric optics) phase contributions calculated from the individual elements located within the frequency-dependent cophase region, that is, those points on the sphere, at each angle and frequency, whose bistatic range fall within a half-wavelength of the shortest path. This, therefore, forms a numerical approximation to the “wavefront-set” of the target, as recognised by microlocal-analysis, and can be easily calculated for a large number of scatterers. It can also be computed for any antenna configuration, and if necessary, arbitrary trajectories. A bank of generic SVMFs based on spherical false-targets can then be built-up to determine locally coherent scattering regions on complex experimental targets. The Fourier-domain matched filter template, equivalent to that produced by (2), is then formed by a 1D azimuth FFT of the simulated frequency/angle data, which can then be applied to detect the location of any target within the cylindrical-polar coordinate space. The simulation process is then repeated for targets of arbitrary size, under all required antenna-pair configurations. An example of the focussing template, that is produced for a 50 mm diameter sphere located at [X,Y,Z]=[0,0.5,0.05] m (or [ρ,ϕ,z]=[0.5,0,0.05]), is shown in Figure 21 for the bistatic antenna pair: θTx=39°,θRx=13°.Four representations of the false-target-generated focussing template revolution data for a single antenna pair (BS6): (a) the “raw-data” domain, (b) the “time” domain, (c) the “Fourier” domain, and (d) shown for completeness but not otherwise used.
(a)(b)(c)(d)One additional point to note is that although a focussed intensity peak is obtained at the centre of a sphere, most targets are not perfect spheres and will, therefore, have a variable local curvature, which will lead to defocussing of the peak intensity. Whilst windowing of the focussing templates may be employed to highlight angular subregions of certain curvature, this remains an area of future work.Focus templates were generated for a range of spherical radii from 0 to 100 mm in steps of 5 mm. As an example, the maximum intensity values for the mid-range value ofrs=50 mm are shown in Figure 22. While the false-target template generation method is also amenable to antenna gain-correction, this was not applied in keeping with previous results.Figure 22
Maximum intensity inX,Y,Z dimensions for the 3 spherical targets described in Table 2, when imaged by the MS spherical-target SWISAR algorithm.These results show that when the data are correctly focussed at the centre of the 100 mm diameter “Target 3”, the intensity remains off-centre for both the smaller and larger targets; indeed the main effect seems to be merely to push the intensity peak downwards inz as rsis increased. This can be expected due to the limited extent in θ.Also, whilst the intensity relating to each sphere increases to a peak as the correct radius is reached, the total energy remains distributed as for the original data, that is, the peak intensity for any search radius corresponds to the signal from the largest sphere, while the signal from the smallest remains at all times the weakest. This is most clearly summarised in the intensity distribution plots shown in Figure23.Intensity distribution per voxel within the proximity of simulated targets 1, 3 and 8.
(a)(b)(c)Here, the distribution of power within a subset of the target volume (for each target, a histogram of voxel intensity is applied over a subvolume defined by ground truth position ±75 mm inXandY and over the full range in Z) is plotted against spherical radius to highlight the fact that (i)
the energy related to each target is focussed towards a peak at the correct radius, which typically occurs within a single voxel of the target volume, and is approximately 10 dB greater than the energy in the surrounding cells;(ii)
the energy from large targets is significantly greater than for smaller targets, even when the radius is far from its true value, as the SWISAR algorithm continues to focus data as a mismatched filter;(iii)
the noise floor of large targets is significantly higher than the peak of smaller targets due to the imaging artefacts previously discussed.For the experimental results, the choice of ground-truth origin for each sphere plays an important part in determining the algorithm’s performance and in practice must be determined from the data as the target structure is “unknown.” In this instance, the spherical origin was chosen by clustering regions of high intensity and choosing the peak intensity of each cluster. As can be seen from the results in Figure24, this is an area for improvement.Variation in intensity versus spherical search radius for each target. The origin of each target was determined by peak intensity of each target cluster. Actual target radii are indicated by the dashed/dotted lines.
(a)(b)(c)(d)(e)(f)(g)(h)Reconstruction may then be performed by rebuilding complex targets as the sum of spheres of different sizes. This is a time-consuming process but has the benefit of identifying both large specular regions, as well as smaller textural features. This work has been reported in [23].


## 6. Conclusions

A prototype 3D imaging radar system has been developed utilising software linearisation and the highest available component bandwidth of 2–18 GHz. This bandwidth requires multiple YIG oscillator sections, which have been linearised to one another to reduce phase noise and achieve close to theoretical range resolutions. Two-target range discrimination performance of better than 20 mm has been obtained for a single FMCW measurement within a custom anechoic test environment.Bistatic antennas in an ISAR configuration have been used to demonstrate the concept of a multistatic spherical SAR system. The Fourier domain point-target focussing templates for this bistatic antenna configuration have been derived using the MoSP. Due to the differing trigonometric terms for two antennas, the phase derivative terms must be solved using polynomial root equations, which although quicker than a full DFT expansion are considerably slower than the equivalent QMS calculation. Results, however, are considerably improved, particularly when the bistatic antenna separation crosses the spherical-polar coordinate system’s zenith.When the point target SWISAR algorithm is applied to nonpoint targets, the imaged intensity is distributed around the surface of that target, but the shallow drop-off in power away from the surface makes size discrimination difficult. A MS spherical-target SWISAR algorithm has, therefore, been developed, acting as a matched filter for targets of a particular radius.Applying the MoSP to the detection of spherical targets in a bistatic antenna configuration is intractable due to the high number of polynomial root equations that must be solved. Alternatively, the Fourier-domain focussing templates for a spherical target can be precomputed using a geometric-optics false-target model, based on numerical rather than asymptotic-analytic integration of elements within the wavefront-set of the target.Experimental results show that while the MS and BS methods provide positive identification of spherical targets of varying diameter, applying the MS algorithm to the BS antenna dataset (i.e., QMS-mode) has significantly poorer results.While further work remains to achieve the original aim of producing size-distributions of blasted rock fragments, this work has then been extended to the 3-D reconstruction of more regular rocks as reported in [14, 23]. It also remains to incorporate the frequency-dependent antenna gain variation, inherent with multioctave measurements, which might otherwise lead to suboptimal image fidelity.



---

*Source: 761816-2012-06-11.xml*