Journal of Global Positioning Systems (2009)

Vol.8, No 1: 88-100

Details of an Imaging System Based on Reflected GPS Signals and

Utilizing SAR Techniques
M Usman, D W Armitage
School of Electrical and Electronic Engineering, The University of Manchester, UK

Abstract GPS as a bistatic radar illuminator for the purposes of air

This research endeavor reveals a method for utilizing target detection. The reflected GPS signals can be used
reflected Global Positioning System (GPS) signals to over the ocean to calculate mean sea height, wind speed,
form an image of targets within a region of interest. The wind direction and significant wave height [5]. Over land
principle is based upon a type of bistatic synthetic to measure soil moisture content, biomass and bistatic
aperture radar (SAR) in which a matched filter technique imaging [6] and over ice to ascertain ice age, thickness
is employed to perform the image reconstruction. A major and surface ice density [7].
challenge was the appalling signal to noise ratio
associated with the received reflected GPS signals. One of the pioneers in the effort of utilizing reflected GPS
Moreover, the reconstruction method resulted in an signals is J. Garrison of Purdue University, USA. His
undesirable point spread function (PSF) which degraded research primarily involved in finding and developing
the reconstructed image. The entire GPS signal generation new applications of the Global Positioning System (GPS)
and image reconstruction process was simulated as and using the GPS signal itself as a remote sensing
faithfully as possible and it has been demonstrated that a instrument [8].
spatial resolution of the order of the GPS L1 frequency
(wavelength of 19 cm) is possible. The smeared image
obtained due to poor quality PSF was improved by means
The idea of 3D multi-static SAR imaging
of a Wiener filter based deconvolution method. An system, which utilized reflected GPS signals
imaging system based on a stationary receiver has been from objects on the Earth's surface was
practically demonstrated with the successful identification presented by Chris Rizos et al. [9]. The output
of a 0.5 m2 spherical target. of the matching filter was mentioned as: -
Ts

Keywords: GPS, Imaging system, Bistatic synthetic 2

aperture radar, Correlation. Sij ( ) =  Ts

SRij (t +  ) SRij  (t) dt

2 (1)
____________________________________________ where Ts is the observation time, τ being the delay and
SRij(t) is the received signal and considered as an
1. Introduction approximately linear FM signal. In fact, bistatic radar
GPS was developed as a military navigation system with systems have been studied and built since the earliest days
potential commercial applications [1]. Direct GPS signals of radar. As an example, the Germans used the British
are utilized for the measurement of position and other Chain Home radars as illuminators for their Klein
navigation parameters with the help of dedicated Heidelberg bistatic system during world war II.
hardware. However, reflected GPS signals, which are Fortunately, the advent of GPS solved many of the
considered a nuisance factor in navigation, have the synchronization and timing problems that have previously
potential to be used for various remote-sensing limited the performance of bistatic radar systems [10].
applications.
M Cherniakov and co-authors have reported much work
The possibility of using GPS signals reflected from the in the field of bistatic radars in general and utilizing
Earth’s surface as a new remote-sensing opportunity was different systems as transmitters of opportunity in
first described in 1993 by the European Space Agency particular. In an initial account they studied a Bistatic
[2]. Similarly [3] and [4] explored the idea of using the
Journal of Global Positioning Systems (2009)
Vol.8, No 1: 88-100
system with the aim to analyze the possibility of air target detection using
bistatic radar and specified the associated problems and The authors are pleased to demonstrate that a spatial
challenges [11]. The research utilized the parameters of resolution of the order of GPS L1 wavelength (19 cm) is
Iridium LEO satellite system and concentrated on the possible provided that sufficient coherent correlator
detection of air targets against a background of white integration time is allowed (in case of static receiver for at
noise, clutter and direct satellite signal. M Cherniakov et least two hours) and that the direct signals can be
al. also examined the power budget analysis of signal sufficiently suppressed for extended imaging range.
detection in bistatic synthetic aperture radar, with global However, the generated image exhibited significant
navigation satellite systems (GNSS) acting as non- blurring, this smeared image obtained due to poor quality
cooperative transmitters (NCT) [12]. PSF was improved by means of a Wiener filter based
deconvolution method. In this context the deconvolution
One of the key parameters around which a remote sensing method is rather simplistic and this issue can form the
system is designed is the spatial resolution. In the SAR basis of further research.
system analysis of the point spread function (PSF) is used
The author’s research endeavor describes a bistatic SAR
to characterize spatial resolution [13]. G P Cardillo used
imaging system, which utilizes reflected GPS signals
vector gradients to derive equations for range resolution, from targets and objects on the Earth's surface to form an
Doppler resolution and cross range resolution in meters image of the area of interest. The GPS satellites, a
and implemented them in computer code. The computer modified GPS receiver and its signal detection
simulations for 1000 km altitude transmitter, with radar components, form the bistatic SAR system. Since the GPS
bandwidth of 500 MHz and integration time of 10 satellites and / or receiver platform are in motion during
seconds showed the ground region for which the range the integration time, the signal obtained at the receiver
and cross range resolution are each less than one meter position will be a Frequency Modulated (FM) GPS signal
[14]. also termed as a chirp signal owing to the constantly
changing Doppler shift.
However, M Cherniakov et al. argued that these equations
can only apply to the rectangular spectrum shape of the The auto-correlation function of a chirp signal exhibits a
ranging signal [15]. They based the resolution analysis narrow pulse property. We know from the basic radar
with an introduction of BSAR spatial (Generalized concept that range resolution is an inverse ratio of the
Ambiguity Function) GAFs and first represented it in the bandwidth of the signal. Chirp signals are the preferred
delay-Doppler domain, and later on expanded to the transmitted signals used to obtain wide frequency band
spatial (coordinates) domain. Using this spatial within a relatively wide time width. The linear FM pulse
presentation of GAF the range and azimuth resolution, as is called chirp signal and can be represented by the
well as a resolution cell area are considered relevant to following simplified equation for a linear increase of
different directions and coordinate systems. Taking the frequency with time: -
idea further they presented equations of resolution of
BSAR [16]. With the assumptions that the target area and f(t) = cos(ωot + πF/T t²) for -T/2 < t < T/2 (3)
synthetic apertures are small, the transmitted signal is
narrow band signal and used the stop-and-go where F = f2 – f1 is the sweep bandwidth with f1 and f2
approximation [13]. The equation derived for range being the initial and final instantaneous frequencies and T
resolution is: - is the time duration [17]. It is possible for the bistatic
SAR imaging system to obtain enough range and
c resolution simultaneously, since the frequency modulated
 range = GPS signal has superior pulse compression properties. A
2bcos( /2) (2) matched filter is employed in the received signal

Transmitter (Tx)
Where b is the bandwidth and β is the bistatic angle. The
equation reveals that the highest resolution will he GPS Satellite
achieved for a monostatic case when β = 0 and range
Direct Signal
resolution is limited by the bandwidth of the transmitted
Imaging System
signal. The theoretical range and cross-range resolution with Modified
Reflected Signal
for the SS-BSAR has been mentioned as 5 meters and 10 GPS Receiver
(Rx)
meters respectively.
Land or Sea Targets
Fig. 1. Scenario - Airborne receiver
Journal of Global Positioning Systems (2009)
Vol.8, No 1: 88-100 Transmitter

processing to improve system resolution. The matched ( GPS Satellite)

filter is an optimum filter for signal detection in a white Direct Signal

noise environment [18]. Each reflected GPS signal is

expected to have a unique chirp or signature which can be
detected
and coherently correlated with the direct signal provided
that the geometry of the satellite is optimum. The location
of the target will have a strong correlation value and can Reflected Signals Static Imaging

thus be identified in the area of interest. Detailed System

simulations and later on practical experiments have Valley

revealed that the bistatic SAR system has the potential to

develop acceptable quality and low cost images of a Few km

localised area. Fig. 2. Scenario – Stationary receiver

The concept of generating an image by the correlation of and relevant problems, shortcomings and challenges have
direct and reflected GPS signals is based on the SAR been recorded for all three scenarios that can prove
principle in which the receiver or bistatic radar is required beneficial for later research efforts.
to sustain its movement in spatial domain. There can be
three different types of imaging scenarios based on the The main objective of the author’s research was to model
source of synthetic aperture, receiver kept airborne, these imaging arrangements and scenarios in a software
ground moving receiver and static receiver. environment with the aim to remove major existing
bottlenecks, optimize the design and later on, attempt to
The transmitter, target and receiver geometry in case of an assemble a proof of concept hardware of the imaging
airborne receiver can be understood with the help of device in order to demonstrate a practical system.
figure (1), as per the traditional SAR concept the receiver
is kept airborne and can detect the presence of ships or The purpose of this research paper is to elaborate how this
other terrestrial targets. The second option is to mount the goal was achieved by careful antenna design, suitable
imaging system on a ground vehicle and drive along a electronic circuit, appropriate GPS signal acquisition and
target thus simplifying the geometry and eliminating the tracking algorithms and efficient reconstruction code. The
need for airborne platform. However, the requirement for document explains the efforts put forward to generate an
precise receiver position and velocity and associated image by the correlation of direct and reflected GPS
logistics issue rendered this option non-practical. signals and based on bistatic SAR techniques.

In order to simplify the experimental set up, it is possible There are many advantages of this type of system (a) GPS
to generate an image with a static receiver, at least target detection and imaging has the attraction that the
theoretically, if sampling is performed over longer period user can take advantage of the expensive GPS
of time, such a scenario is depicted in figure (2). With the infrastructure maintained for navigation purposes by
receiver being static the change is geometry is provided employing it as a ‘illuminator of opportunity’ thus there is
by the moving GPS satellite, though the change in no requirement for a dedicated transmitter. (b) There are
geometry is a relatively slow process and results in very several measurement opportunities, one for each GPS
long data acquisition times of the order of hours. One of satellite in view, thus forming a multi-static radar system.
the numerous possible applications with this type of set (c) The imaging system has the additional advantage of
up are remote sensing for monitoring of landslides, long- being very cost effectiveness, as it can be built for a
term seismic studies and similar applications. fraction of the cost of traditional radars, space-borne
equipment and other sensors. (d) The operation will be
The possibility of using this technique for imaging covert as no signal will be transmitted locally unlike
purposes has been mentioned in various research papers mono-static radar. (e) Bistatic radar may also have a
as a number of institutions study bistatic SAR with counter stealth capability, since target shaping to reduce
asymmetric configuration, but there is little information the mono-static RCS of a target will in general not reduce
regarding detailed simulations to confirm technical its bistatic RCS [10].
feasibility or experimental results. In this context the
efforts carried out to simulate and practically validate the We know that power spectral density of any given GPS
results with the help of static or stationary GPS receiver C/A code signal is below the power spectral density of
are a novel achievement. Nonetheless, simulation details noise. It is converted to a usable (signal to noise ratio)
Journal of Global Positioning Systems (2009)
Vol.8, No 1: 88-100
- x (meters) + x (meters)
SNR by correlation of a locally generated C/A code
sequence that provides an effective processing gain to the +y (meters)
SNR ratio. The theoretical processing gain for the C/A
code is 1023 or 30.1 dB [1]. High autocorrelation peak
and low cross-correlation values provide a wide dynamic
range for signal acquisition. Hence, on the basis of
correlation theory the code is able to be detected, although
it is buried in noise. This important property of the GPS
signal plays a vital part in our simulation to detect the
target with the help of reflected GPS signals [19].

When the GPS signal arrives at a target, it is reflected in a -y (meters)

Fig. 4. Results for integration time of 0.1 sec.
manner that depends on the surface roughness and the
(0, 0, 6378000) One pixel or element
dielectric properties of the target. The reflected signals North Pole
that are of most use for imaging purposes are those Fig. 3. Candidate position for the target
reflected in a non-specular manner. Such a ‘reflection’
results from the signal being scattered over a wide angle, signal paths are calculated based on satellite and receiver
and will thus be subjected to further attenuation as a trajectories.
function of distance traveled. We could consider the
surface of a complex target to consist of a number of The reconstruction engine divides the region of interest
‘point reflectors’. In reality ‘point reflectors’ must have into an array of square bins and attempts to coherently
finite size, and for our purposes it would seem expedient correlate the received signal with an individual template
to model these ‘point reflectors’ as spheres with a chirp unique to each bin. Only where the received signal
diameter comparable to the resolution of the imaging contains a matching chirp does a given bin achieve a high
method. score. In order to generate the matching templates, the
reconstruction engine independently calculates the signal
2. Simulation transit delays and PRN codes in the same way as the
The simulation work was carried out in a Matlab ® signal synthesizer, but for all candidate target locations.
environment and during initial phase of the research a
simple code was generated in order to verify the imaging The search area in which the targets can be located is
principle and concept. After successful initial simulations, shown in figure (3). The position (0, 0, 6378000) in ECEF
the code was modified and expanded to facilitate image coordinate system was fixed in the middle of the search
generation under practical scenario and efforts were put in area. These are actually the coordinates of North Pole,
to improve the computational efficiency. employed to simplify coordinate transformations.

The simulation was deliberately split into two distinct The matched filter processing plays a vital part in the
parts: the first part being responsible for synthesizing the simulation. In the Matlab environment this is performed
signal emanating from the analogue to digital converter of by multiplying two signals sample by sample and than
the GPS front end, the second part being responsible for accumulating. The process was carried out for each
the image reconstruction. Care was taken to ensure that sample and the relevant value was stored in the form of a
the only information available to the reconstruction two dimensional matrix. The matrix was scanned for each
engine was the digitized GPS signal, and ‘globally sample in time, rendering the procedure equivalent to a
available’ data such as the satellite ephemerides and the matched filter process. A different value of correlation
receiver location. The signal synthesizer determines the was obtained for each instant in time. The location of a
satellite orbits together with the individual PRN code target is expected to have a maximum value of coherent
offsets and combines the direct signals from all satellites correlation. A grey scale image was constructed using
in view with the reflected signals from the targets together above-mentioned information. Thus, the image produced
with the receiver noise. Transit delays associated with the directly pin points the target at the desired coordinate in
the candidate position.
3. Results and Analysis of 0.1 seconds was selected and the results were
The Matlab® simulations were carried out under comparable to that presented in [9].
different conditions and scenarios by varying
various parameters. Initially an integration time As depicted in the figure (4), the two targets can easily be
recognized with acceptable spatial resolution. For the next
Journal of Global Positioning Systems (2009)
Vol.8, No 1: 88-100
case the integration time was increased to 0.2 seconds, 1000
Simulation performed with Integration Time of 21160 seconds

although the simulation took many hours longer to 800

execute but resulted in improved resolution. The longer
600
integration time produced larger number of samples and
400
improved the PSNR (peak signal to noise ratio) in the
reconstructed image. 200

Metres Y
0

As mentioned earlier the other configuration for this type -200

of system would involve a static GPS receiver. Figure (5) -400

reveals that the resolution was limited only by the image -600

pixel size of 20 m in case of such a scenario, when sampling -800

is performed every few seconds over a longer period of -1000

-1000 -800 -600 -400 -200 0 200 400 600 800 1000
time, for example six hours and later on integrated. The Metres X

targets are clearly visible with no distortion. In fact the Fig. 5. Results for integration time of 21160 sec (6 hours)
longer integration time has resulted in larger number of
samples and cancellation of uncorrelated noise. Simulations were performed to explore the limits within
which the targets can be resolved in terms of the
3.1. Direct Signal Suppression magnitude of the direct and reflected signals. The direct
Under practical scenario the strength of GPS signal is signal was attenuated to a significant amount in an effort
very low and buried in noise and only with the help of to find out the threshold of target detection, the indirect
correlation theory and deterministic nature of the GPS signal was also attenuated in incremental form. It was
C/A codes the signal is detected. The reflected signal is deduced that in order to bring both signals comparable in
even weaker but it is expected that the targets can be terms of magnitude and make sure that the autocorrelation
resolved by integration over longer periods of time. tails of direct signal are smaller that the reflected signal it
However, in case the reflected signal is further reduced, is imperative to attenuate and scale down the direct signal.
due to the type of reflecting surface and distance covered, It is assumed that this large suppression figure could be
it can become equal or lower than the tails of auto achieved by a combination of careful antenna design and
correlation function of direct signal and can be signal subtraction within the signal processing, since the
permanently lost. No matter how long we integrate the receiver channel for the reflected GPS signals may be
reflected signal cannot be separated and the targets are no locked to the direct signal. Other options including the
longer resolved. This phenomenon can be understood utilizing of multi-beam antenna with main lobe for the
with the help of the figure (6). reflected signal and side lobe for the direct signal can be
explored.
In numeric terms, the tails of autocorrelation of direct
signals are 24dB lower than the direct signal peak. However, there are some simple ways to improve the
Theoretically, the reflected signal must not be lower than SNR of the reflected signal, for example, carry out
24dB as compared to the correlation peak of direct correlation for longer period of time, use high gain
signals. directional antenna or perform image processing at a
shorter range. These simple measures have been adopted
in the initial practical experiment carried out to acquire
data for image generation purposes. A high gain and
directional antenna for acquisition of reflected GPS
signals was employed and the range was kept short. This
resulted in strong reflections from the spherical target
covered with aluminum foil and the reflected signal was

Direct Signal

Correlation
power 24 dB

Reflected Signal
Tails of Auto Correlation
Function

0 ( time )

Fig. 6. Correlation of direct and reflected GPS signal

Journal of Global Positioning Systems (2009)
Vol.8, No 1: 88-100
comparable in magnitude with the direct signal. However, Define imaging parameters
in case of extended data acquisition ranges where weak (Integration time, sampling
frequency, IF and tow etc)
reflected signals are received, this issue has to be kept in Load ephemeris data
mind. Generate Candidate
position for the target
Generate PRN sequence
for reference satellite
Other researchers working in related field also came Define sub-segment size
(size of linear array)
Define Segment size, Skip
across this issue and have referred to it in terms of direct Time, Dwell Time
path interference (DPI) and cross correlations. M Initiate for loop
Generate local signal

Cherniakov et al. explained that taking into account only (corresponds to skip time) Initiate for loop
(Number of segments = Number
one satellite as the source of DPI, traditional adaptive side of samples / sub-Segment size)
Load data file
lobe cancellation methods can be applied for DPI
Mix signal to base band
suppression [20]. In a related study E P Glennon provided Update IF to resolve
a brief description of the GPS cross correlation problem Update phase and Doppler ambiguity
(software GPS receiver)
frequency after each
and sketched a comparison between the cross correlation iteration
problem as experienced by GPS and the CDMA cellular Calculate delay at Mid
communications. He examined several different types of point of the sub-segment

cross correlation mitigation techniques, ranging from

Determine delay offset
successive interference cancellation through to subspace
projection methods and showed that the mitigation
Reflected signal Calculate delay
techniques are also applicable to the removal of
autocorrelation side lobes [21]. Calculate delay at sub- Calculate delay at sub-
segment finish segment finish Direct signal
3.2. Accuracy of the ADC (Analogue to Digital Calculate delay at sub- Calculate delay at sub-
Converter) segment start segment start

In order to perform necessary signal processing the GPS Construct template for Construct template for
reflected signal direct signal
signal is to be converted to digital format after down
conversion from RF of 1.575 GHz to a manageable IF of Cross correlate template Determine phase
for reflected signal Cross correlate template
a few MHz. This process is to be performed by the with direct signal

Analogue to Digital Converter (ADC). The best possible

bit accuracy for the ADC was calculated by Matlab® Interpolate cross
correlated results (Linear
simulations. It is imperative to consider that a higher than Interpolatoin)

required accurate device will lead to unnecessary expensive Update tow, sub-segment
and complex hardware and a burden on the computation Increament tow and point count

power in the simulation process. On the other hand an to the requisite place in the
file
inferior ADC will result in signal degradation. It was Generate image

deduced that at least two bits are required as the minimum Fig. 7. Block diagram of the reconstruction code
accuracy. This was an important result and proved that a
commercial off-the-shelf GPS receiver can be used as RF
front end in the imaging hardware. The reason for results modular approach has been employed. The GPS signal
obtained was due to the fact that variance of quantisation was generated separately in C++ environment using the
noise (σ2) = q2 / 12, where q is the step size or (Least Borland C++ Builder. The program simulated all the
significant bit) LSB and hence Standard deviation (σ) = q / functions of the Maxim MAX 2741 IC used as the RF
√12 = q / 3.5. For two bits there are four level (2N-1) and for front end in the hardware. This approach speeded up at
one bit only two levels are present, hence quantization noise least the signal generation part of the code by a factor of
will become double and comparable to the thermal noise in about 9.
the signal.
The reconstruction code developed in Matlab® was also
modified in an effort to make it more suitable for a
3.3. Improvement in Computational Efficiency
practical scenario. In order to improve the computational
It was observed that the Matlab® simulations were very
efficiency the calculations were performed on linear
expensive in terms of computation time and memory and
arrays (or vectors) rather that multi-dimensional arrays.
will be feasible only for off-line applications. One of the
The block diagram of the reconstruction algorithm is
many methods of improving the computational efficiency
depicted in figure (7).
is by switching from interpretive Matlab ® to compiled
code. In order to improve the speed of the program a
Journal of Global Positioning Systems (2009)
Vol.8, No 1: 88-100
As evident from the diagram the two signals (locally
GPS Satellite
generated and received GPS signals) are multiplied and
added at each segment and a different value of correlation
is obtained for each instant in time. The value is stored in
the array by suitable Matlab® commands. However, the Receiver start Receiver finish
Position Position
values of the delay and offset PRN for the locally
generated or reference signal are to be calculated first.
These parameters are determined by measuring the
distance between reference satellite to candidate position
and receiver to candidate position. The delays are
calculated not for a particular target position, but for all

(0, 0, 6378000)
candidate positions in the area of interest. The simulation
calculated the delay from satellite to the receiver for each
box or segment as shown in the figure (8).

Initially, this represented a two-dimensional array

representing delays for each candidate position. The code Reference One pixel or element
also calculated the delay offset from the centre and
Linear Array of delays
considered it as a reference, utilizing it for comparison of
the offset values. First of all, the code calculated the Minimum Delay Maximum Delay
Offset
delays at the midpoint of the segment as a reference and Offset

stored the value for each candidate position. The delay- Fig. 8. Calculation of Delays and Interpolation for Linear Array
offset range is also determined as an offset to delay at
centre or midpoint of segment. w(n   )  (1   ) w( n)  ( ) w( n  1) (6)
A template was generated in the form of a linear array or
vector of delay offsets relative to the image centre with where n is the integer part of the original index value and
initial value as minimum value of delay offset and final η is the fractional part [22].
value as the maximum of delay offset. The template for
direct signal was represented as follows: Thus by calculating for each location, an array or vector
was formed which contained all the linear representations
Sdir  ddir (t )e jwLt (4) of the range of possible wrapped signals. Owing to the
one dimensional array, there are fewer calculations while
where ddir(t) is the C/A code for direct signal after offset the correlation is performed, but the simulation still took a
has been introduced and wL is GPS L1 frequency in long time as the new calculations are more complex.
radians / seconds. The template was cross-correlated with Never-the-less even with improved resolution as
the direct received signal to determine the phase mentioned in the next section there was a further
difference. improvement of computational efficiency by a factor of 2
to 4.
Similarly the template for reflected signal was represented
as follows: Another method to improve the computational efficiency
was signal
Sref  dref (t )e jwLt (5)

where dref(t) is the C/A code for reflected signal after

offset has been introduced.

The reflected signal was cross-correlated with the

template for reflected signal and results were interpolated
to form a linear array. The method used was linear
interpolation, which can be represented by the following
equation.

Fig. 9. Raw reconstructed image

Journal of Global Positioning Systems (2009)
Vol.8, No 1: 88-100
conspicuous conclusions, firstly the targets appeared to
have the potential of being resolved to the size of one bin,
and secondly there was significant blurring. It is pertinent
to note that despite the targets exhibiting a star-like PSF,
achieving a resolution comparable to the wavelength of
GPS L1 frequency is a novel achievement. According to
the author’s knowledge it has not been documented
anywhere in the research literature.

3.5. Image Restoration

The image-processing task at hand was to reconstruct the
unfocused target signature into a focused image. A
technique of de-convolution of the signal with its PSF
using Wiener Filter has been employed to correct the
Fig. 10. Point Spread Function blurred image obtained form two isolated point reflectors.
In our simulation we know the signal and also calculate
segmentation, which was incorporated in the code to the PSF of the signal in order to utilize this filter. The
simulate low velocity and static or stationary GPS point spread function is calculated during the course of
receiver. Only selected amount of data was processed, simulation and used as the degradation function and is
while the rest was omitted and the phase, frequency and shown in figure (10).
TOW (time of week) were updated accordingly.
However, the method can be used for any imaging
arrangement and the fact that only part of the signal
acquired is processed for image generation, the technique
resulted in reduced computational time but comparable
image quality.

There is a simple method to improve computational

efficiency, which is sub-sampling or under sampling. In
the reconstruction algorithm a built in Matlab® function
resample was used to perform sub-sampling (up to a
factor of 10), thus reducing amount of data to process and
in turn improve computational efficiency, without
compromising image quality. Fig. 11. Deconvolved and reconstructed image

3.4 Resolution During the simulation white noise is assumed, whose

In order to test the achievable resolution the simulation power spectrum is constant, which simplifies things
space was represented by an array of 400 x 400 bins with considerably [23]. The Wiener Filter in Matlab® is based
pixel dimensions of 25 x 25 centimeters. Four point on the following equation with value of K selected as
reflector targets ware placed on the corners of a 20 m 2 1016.
area, having coordinates (10, 10, 6378000), (-10, 10,
6378000), (-10, -10, 6378000) and (10, -10, 6378000).
 
Each conducting point reflector had an effective radius of  (u , v )  G (u , v )  H *(u , v) 
F
10 centimeters and reflection coefficient of 1. This 2
particular value for the radius was selected to demonstrate  H (u, v)  K  (7)
that it is possible to detect a target comparable to the 20
cm wavelength of GPS L1 frequency. It was assumed that Where F(u,v) is the input image, G(u,v) is the degraded
the antenna for reflected signals has a gain of 30 dB and image and H(u,v) is the spatial representation of the
that the direct signals can be suppressed by 70 dB. The degradation function in the frequency domain and H*(u,v)
coherent integration time was set to 10 seconds. is the complex conjugate of H(u,v). The restored image in
the spatial domain is given by the inverse Fourier
The generation of image with the modified reconstruction transform of the frequency domain estimate. The
algorithms (receiver velocity of 200 m/s and sampling deconvolved and reconstructed image is shown in figure
frequency of 2 MHz), as shown in figure (9), revealed two
Journal of Global Positioning Systems (2009)
Vol.8, No 1: 88-100
(11). Further details in this regard have been documented GPS Satellite

in [24]. Direct Path Hardware Reconstruction Algorithms

Direct Signal

Correlation
GPS Front End for Aquisition and
and Matched Image
RHCP Antenna Direct Signals Tracking of Direct
Filter Generation
Signals
Procesing

Reflected Signals
GPS Front End for measured as
LHCP Antenna Reflected Signals compared to
direct signal
Reflected Path

Reflected Signal

Fig. 13. The Imaging Scenario

4. The Imaging Hardware

The imaging scenario is shown in figure (13). It is
apprized that the reflected signal is not expected to have
enough SNR to permit successful signal acquisition. So
Fig. 12. Simulation results for 7200 seconds the direct signal from a specific satellite received by the
RHCP antenna has been selected, locked and used as a
reference for the reflected signal. Once successful lock is
3.6. Simulations for Static GPS receiver achieved, the same local oscillator can be employed for
As mentioned earlier for static receiver configuration the the direct as well as the reflected signal. The GPS signals
change in geometry to generate the synthetic aperture is
are received, amplified, down-converted and digitized
provided by the moving GPS satellites, which have an into near base band samples, which are than processed
average orbit altitude of 20,200 km above the surface of
using software routines to acquire and track the direct
Earth and complete one orbit in approximately 11 hours GPS signal and later on matched filter processing with
and 58 minutes [1]. For the orbiting GPS satellite to indirect signals culminate the task of image generation.
provide the requisite synthetic aperture, GPS data in the
order of hours is required. However, storage of such huge In order to perform the field trails a custom LHCP helical
amount of data is not practical. Fortunately, not all antenna and two-channel GPS front end and data
acquired data was utilized and it was necessary during capturing device for acquisition of GPS signals has been
image reconstruction to process selected amount of data prepared. The antenna furnished a compact size and high
while omitting other and updating the phase, frequency gain solution for signal acquisition. More importantly
and TOW (time of week) accordingly. This resulted in
being of reverse in polarization as compared to direct
dividing the GPS IF data into smaller parts or chunks
signals and pointed away from GPS satellites, it delivered
denoted as skip time (Tsk) and dwell time (Tdw). The
former represented the time interval during which data is excellent immunity against direct signal interference.
not processed or omitted in order for the change in
geometry to become effective or build up. The latter The design requirement for the data collection device of
denoted the time during which the reconstruction the imaging system was a low-noise, dual-input receiver
algorithms process the GPS data for correlation of direct operating at the frequency (L1) used by the GPS satellite
and reflected signals in order to generate an image. signal. Initial design reviews led to the adoption of two
MAX2741 IC's as the basis of the design. The IC is an
Detailed simulations have revealed that at least one hour L1-band dual-conversion GPS receiver which down
(3600 seconds) of integration time is required to converts the 1575.42MHz L1-GPS signal to a 37.38MHz
successfully resolve the targets with acceptable spatial first (Intermediate Frequency) IF, and then a more
resolution. However, two hours of data furnished a sub manageable second IF of a few MHz. The on-chip
meter level range resolution, as depicted in figure (12). Analogue to Digital Converter (ADC) samples the down-
The image quality is directly proportional to the physical converted GPS signal at the 2nd IF with the rate of 19.2
change in geometry of the GPS satellite, thus the MSPS. Sampled output is provided in either 2-bit (1-bit
resolution and quality of image depends upon the total magnitude, 1-bit sign) or 3-bit (2-bit magnitude, 1-bit
integration time. sign) formats, as determined by the ADC mode
configuration bit. As depicted in figure (14), one front end
of the device receives the direct GPS signals, while the
Journal of Global Positioning Systems (2009)
Vol.8, No 1: 88-100
other receives the reflected signals through the LHCP inline amplifier of gain 30 dB was used to amplify the
antenna. Further details regarding the imaging hardware weak reflected GPS signals.
have been documented in [25].
As expected a very strong signal was received with the
RHCP GPS antenna, yielding good acquisition even for
few ms of integration time. The comparatively weaker
reflected signals acquired by the LHCP antenna required
much longer integration times (200 ms) to achieve
comparable SNR. Longer acquisition time results in the
cancellation of uncorrelated noise, thus improving the
SNR. The Doppler frequency of both signals is the same
but code offset is slightly different corresponding to the
extra distance that the reflected signals have to travel.

The length of one GPS C/A code is 1023 chips and is

transmitted with a frequency of 1.023 MHz. Taking into
account the speed of light the length of one chip can be
(a)
calculated to be 300 m. As mentioned above the signal is
sampled at 19.2 MHz. Thus each code sample
corresponds to about 15 meters. It is possible to
distinguish between the direct and reflected signal if the
path length between direct and reflected signal is a
multiple of 15 meters.

During the experiment performed the difference in code

samples of direct and reflected signal was 3 or 4 which
came out to be about 45 to 60 meters and corresponds to
the round trip distance between antenna and the large
brick building or reflective surface. Figure (15a)
compares the direct and reflected signal correlation peaks
clearly depicting this path difference. The scale for the
(b) correlation value (y-axis) is different for both signals on
Fig. 15. Comparison of direct and reflected GPS signals for (a) 50 maccount of the varying acquisition time.
(b) 4 m (round trip distances)

5. GPS Data Acquisition and Imaging

In order to test both channels of the GPS front end and
data capturing device and to confirm that the signal
received by the LHCP antenna is in fact the reflected CP)
RH
signal an experiment was performed in front of a large te nna ( ignals
Target An irect s
d
brick building. The LHCP antenna was positioned so as to for
Transmitter
receive the signal bouncing off the building. A Wi-Sys Ant
e nn
refl a (LH
ecte C
d si P) for
gna
ls
GPS Satellite

Fig. 16. Location of both antennas and the target

Fig. 14. Two-channel GPS front end and data capturing device
Journal of Global Positioning Systems (2009)
Vol.8, No 1: 88-100
In order to verify the results, GPS IF data was collected at environment to have a clear view of the sky and thus
a position of only two meters away from the building or increasing the chances of receiving more GPS satellites
reflecting surface. A path difference of only one code and for extended duration. The target can be seen in the
samples was observed among the direct and reflected middle of the diagram, as the antenna's main lobe was
signal as shown in figure (15b). A near specular GPS aimed roughly towards the target center. Unfortunately,
reflected signal was received suggesting an optimum an increase in the integration time also resulted in more
geometry for reception of reflected signals. Thus the signal being received by the highly directional antenna
experimental results substantiated that the signal present and some of the signal was not received from the target,
at the LHCP antenna is in fact the reflected signal and but bounced off from nearby objects and building front.
construction of the hardware has been successful.

After successfully simulating the static receiver imaging

scenario, modifying the reconstruction algorithms to
process longer duration GPS signals and acquiring the
reflected GPS signals. The imaging hardware including
the two antennas and data acquisition device was
positioned in front of the University building, the location
of the direct and reflected antennas and target are shown
in figure (16). A 0.5m2 spherical target wrapped in
aluminum foil was placed in front of the university
building.

During acquisition of data for imaging purposes the

nearby buildings formed an urban canyon type
environment thus limiting a clear of the view of the sky.
The most suitable satellite (in terms of signal strength and
visibility) was GPS BIIRM-2 (PRN 31) and was therefore
selected as the reference satellite.

A comparison has been made between the acquisition

diagram for 5 ms of direct signal and 100 ms of reflected
signal in figure (17). The coarse frequency of both signals
is the same but code offset is slightly different
corresponding to the extra distance that the reflected
signals have to travel. During analysis of individual data
files, this difference fluctuated between one or two code
samples (may be due to the quantization effect of the
ADC) corresponding to a path difference of about 20
meters between direct and reflected GPS signals. (a)

In order for the GPS satellite to provide the requisite

change is geometry, 80 files were down loaded at an
interval of about 30 seconds each. The length of
individual file was 4 seconds, but only the first few
milliseconds of each file were used during the
reconstruction process. During simulations it was deduced
that a 30 seconds interval among the acquired files is just
enough to resolve the targets.

The GPS data provided about 2400 seconds for change in

geometry, just enough to identify the target, but
compromising the resolution. Figure (18) compares the
image by reconstructing 2400 seconds of actual data with a
simulated signal of 2400 seconds. In future it is
(b)
recommended to perform signal acquisition in an open
Fig. 17. Comparison of acquisition plot (PRN 31) for acquisition
time (a) of 5 ms for direct signal (RHCP antenna) (b) of 100 ms
for reflected signal (LHCP antenna)
Journal of Global Positioning Systems (2009)
Vol.8, No 1: 88-100
available during the experiments. The target has been
detected in a hostile environment with the help of
extremely weak reflected GPS signals that are
omnipresent, but exhibit an appalling SNR. It is further
apprized that the change in geometry to process the data
with the help of SAR technique was provided by the
orbiting GPS satellite. This particular method has so far
not been utilized in a practical environment for imaging
purposes.

Although steady progress was made during research

endeavor and an image of the area of interest, roughly
identifying a 0.5 m2 sphere target was generated.
Extensive field experimentation with the imaging device
needs to be performed. In order to confirm the detection
of target, more data is to be acquired at a suitable location
(a) that ensures visibility of the reference GPS satellite for at
least two to three hours for a high resolution image. The
resolution can be improved by correlating for longer
periods of time and correcting the relative phases of direct
and reflected GPS signals.

6. Conclusions and Recommendations

During the course of current research endeavour, the
entire GPS signal generation and image reconstruction
process has been simulated. The smeared image obtained
due to poor quality PSF was improved by means of a
Wiener Filter based deconvolution method. In this context
the deconvolution method is rather simplistic and this
issue will form the basis of further research. It will be
(b)
interesting to explore other deconvolution techniques like
blind de-convolution.
Fig. 18. Comparison of image obtained by reconstructing 2400
sec of (a) Actual signal (b) Simulated signal The results have demonstrated that this method can be
used in high resolution terrain and change monitoring
applications to obtain low cost and detailed images of a
fixed location. However, it is understood that in a
As evident front the results, although it is roughly
practical situation the resultant image might be severely
possible to distinguish the target, there is excessive
degraded, the objective of further research efforts is to
backscatter or clutter and noise exhibited in the images.
improve and bring the image as close to the ideal as
The situation was exacerbated due to the fact that the high
possible. Some of the factors having an adverse effect on
gain LHCP antenna was placed in front of a large brick
the signal include noise, strength of reflected signals,
building. The antenna was designed to acquire the weak
direct signal interference, atmospheric effects, radar
reflected GPS signal, but returns from objects around the
clutter or backscatter and hardware limitations
target were also received and displayed in the image. As
mentioned earlier that due to hardware limitations, the
data files were acquired after an interval of every 30
References
seconds each, such low frequency for acquiring temporal
samples induced an aliasing effect. In future it is
recommended to modify the imaging hardware to acquire 1. Kaplan, E.D., Understanding GPS Principles and application:
Artech House Publishers.
GPS data after every three seconds. 2. Martin-Neira, M., “A Passive Reflectometry and Interferometry
System (PARIS): Application to Ocean Altimetry,”. ESA J, 1993.
The exhibited images may seem very primitive, but it has vol. 17: p. 331-335.
to be kept in mind that no dedicated radar transmitter was
Journal of Global Positioning Systems (2009)
Vol.8, No 1: 88-100
3. Glennon, EP, Dempster, AG and Rizos, CR., Feasibility of air
target detection using GPS as a bistatic radar. Journal of Global
Positioning Systems, 5(1-2), 119-126.
4. Cherniakov, M.N., D. Kubik, K., Air target detection via bistatic
radar based on LEOS communication signals. IEE Proceedings-
Radar Sonar and Navigation, 2002: p. 33-38.
5. James L Garrison, A.K., Valery U Zavorotny, and Stephen J
Katzberg, Wind Speed Measurement using Forward Scattered
GPS Signals. IEEE Transactions on Geo-Science And Remote
Sensing, 2002. Vol. 40, (No. 1).
6. D Masters, S.K., P Axelrad. Airborne GPS Bistatic Radar Soil
Moisture Measurements during SMEX02. in Proceedings of IEEE
International IGARSS. 2003.
7. A Komjathy, J.A.M., V U Zavorotny, P Axelrad and S J Katzberg.
Towards GPS Surface Reflection Remote Sensing of Sea Ice
Conditions. in Sixth Int. Conf. on Remote Sensing for Marine &
Coastal Environments. 2000.
8. Stephen J. Katzberg and James L. Garrison, J.L., Utilizing GPS To
Determine Ionospheric Delay Over the Ocean, Research Center
Hampton: Virginia, USA.
9. Yonghong Li, C.R., Eugene Donskoi, John Homer, Bijan
Mojarrabi, “3D Multistatic SAR System for Terrain Imaging based
on Indirect GPS Signals” Journal of Global positioning system,
2002. Vol. 1(No. 1): p. 34-39.
10. Griffiths, H.D., Bistatic and Multistatic Radar, in IEE Military
Radar Seminar. 2004: Shrivenham, UK.
11. Cherniakov, M.N., D. Kubik, K., Air target detection via bistatic
radar based on LEOS communication signals. IEE Proceedings-
Radar Sonar and Navigation, 2002. 149(1): p. 33-38.
12. M.Cherniakov, X.H., T. Zeng, Signal detectability in SS-BSAR
with GNSS non-cooperative transmitter. Iee Proceedings-Radar
Sonar and Navigation, 2005. 152: p. 124-132.
13. Soumekh, M., Synthetic Aperture Radar signal processing with
Mat lab Algorithms. 1999: John Wiley & Sons, Inc.
14. Cardillo, G.P. On the use of Gradient to determine bistatic SAR
resolution. in IEEE Antennas and Propagation Society
International Symposium, 1990. AP-S. 'Merging Technologies for
the 90's'. Digest. 1990. Dallas, TX, USA.
15. Cherniakov, M., T. Zeng, and T. Long, Generalized approach to
resolution analysis in BSAR. IEEE Transactions on Aerospace and
Electronic Systems, 2005. 41(2): p. 461-474.
16. Cherniakov, M., Zeng, T. and Plakidis, E., Ambiguity Function for
Bistatic SAR and its Application in SS-BSAR Performance
Analysis, in Radar Conference. 2003: Adelaide, Australia.
17. Peyton Z. Peebles, J., Radar Principles: John Wiley & Sons, Inc.
18. Sklar, B., Digital Communication Fundamentals and Application:
Prentice Hall International Edition.
19. M. Usman, D W Armitage., A Remote Imaging System Based on
Reflected GPS Signals, in ICAST 2006. 2006: Islamabad, Pakistan.
20. He, X., M. Cherniakov, and T. Zeng, Signal detectability in SS-
BSAR with GNSS non-cooperative transmitter. Iee Proceedings-
Radar Sonar and Navigation, 2005. 152(3): p. 124-132.
21. E. P. Glennon, A.G.D., A Review of GPS Cross Correlation
Mitigation Techniques.
22. Smyth, T., Waveshaping Synthesis. February 18, 2008, School of
Computing Science, Simon Fraser University.
23. R. C. Gonzalez, R.E.W., “Digital Image Processing”. Second
Edition ed. 2002: Prentice-Hall Inc.
24. M Usman, D W Armitage, De-Convolution of images obtained by
Correlation of Direct and Reflected GPS Signals, in 6th
International Bhurban Conference on Applied Sciences and
Technology, 2008, Islamabad, Pakistan.
25. M Usman, D W Armitage Acquisition of Reflected GPS Signals
for Remote Sensing Applications. in Proceedings of ICAST 2008.
2008. Islamabad.