Journal of Global Positioning Systems (2005)
Vol. 4, No. 1-2: 95-105
Real-time Doppler/Doppler Rate Derivation for Dynamic Applications
Jason Zhang
School of Mathematical and Geospatial Sciences, RMIT University
e-mail: jzhang@eos-aus.com Tel: + 61-02-6222 7949
Kefei Zhang, Ron Grenfell, Yong Li and Rod Deakin
School of Mathematical and Geospatial Sciences RMIT University
Received: 6 November 2004 / Accepted: 15 October 2005
Abstract. Precise GPS velocity and acceleration
determination relies on Doppler and/or Doppler rate
observations. There are no direct Doppler rate
measurements in GPS. Although every GPS receiver
measures Doppler shifts, some receivers output only
“raw” Doppler shift measurements and some don’t output
any at all. In the absence of raw Doppler and Doppler rate
measurements, a differentiator is necessary to derive
them from other GPS measurements such as the carrier
phase observations. For real-time dynamic applications,
an ideal differentiator should have a wideband frequency
response to cover all the dynamics. It should also have a
group delay as short as possible. In addition, a low-order
differentiator is more favourable for easy implementation.
This paper provides an overview of methods in
differentiator design for applications of GPS velocity and
acceleration determination. Low-order Finite Impulse
Response (FIR) differentiators proposed by Kavanagh are
introduced. A class of first-order Infinite Impulse
Response (IIR) differentiators are developed on the basis
of Al-Alaoui’s novel differentiator. For noise attenuation,
it is proposed to selectively use Kavangagh’s FIR
differentiators, and the first-order IIR filters derived for
adaptation to different dynamics.
Key words: GPS velocity determination, GPS
acceleration determination, differentiator design, FIR
filter, IIR filter, Doppler.
1. INTRODUCTION
Previously proposed methods for GPS velocity and
acceleration determination fall in two categories. One is
to derive velocity and acceleration directly from GPS
determined positions, another is based on the Doppler
shift method. The latter has several advantages: it doesn’t
rely on the precision of the positions from GPS, nor will
the accuracy dramatically degrade with an increase in
sampling rate (say 10Hz or more). Since there is no direct
Doppler rate observation in GPS measurements, as a
“virtual” observable, it must be derived in order for the
formulae presented by Jekeli and Garcia (1997) to be
applied directly in the Doppler shift method.
Every GPS receiver measures Doppler shifts. However,
this is primarily an intermediate process to obtain
accurate carrier phase measurements. Thus the quality of
Doppler shift output varies from receiver to receiver
depending on manufacturer. The Trimble 5700™
geodetic receiver, for instance, has a measurement
precision of ±1mm/s. The observed Doppler is from a
tracking loop that is updated at a very high rate. This also
enables the receiver to sense phase accelerations (Harvey
2004). Unfortunately the sensed phase acc-elerations and
the Doppler shift on L2 are discarded. Some other GPS
receivers, for example the Superstar II™ from NovAtel,
have only code and L1 phase outputs (SuperstarII 2004)
and the Doppler shifts are masked out of the
measurements. For our purposes to obtain accurate
velocity and acceleration using these types of receivers, it
is necessary to derive the Doppler shifts, i.e. the change
rates of the carrier phase from the measured carrier phase
measurements.
Differentiators are required to get the Doppler rate
“observable” for any type of receiver, or to get the
Doppler shift from the carrier phase. In real-time and
dynamic applications it is also desirable that the designed
96 Journal of Global Positioning Systems
differentiator should have a wideband frequency response
to cover the system dynamics. It should also have a group
delay as short as possible so as to get the Doppler shift or
Doppler rate instantaneously. For those receivers that
output only “raw” Doppler shifts, the derivation of
precise Doppler from the carrier phase plays a key role in
precise velocity and acceleration determination. This is
because the precision of carrier phase observables can be
fully exploited. The objective of this paper is to explore
the techniques to derive Doppler rate from GPS
measurements, or to derive precise Doppler shift from the
carrier phase in real time and in dynamic situations.
Several investigations have been conducted for this
purpose in the GPS measurement domain, and the
proposed methods can be categorised into:
(1) Curve fitting (Fenton and Townsend, 1994);
(2) Kalman smoother/filtering (Hebert, Keith et al.
1997);
(3) Taylor series approximation (Hebert, Keith et al.
1997; Cannon, Lachapelle et al. 1998; Bruton,
Glennie et al. 1999);
(4) Finite Impulse Filter (FIR) by using Fourier series
with window techniques (Bruton, Glennie et al.
1999); and
(5) FIR optimal design using the Remez exchange
algorithm (ibid).
The FIR filtering technique based on Taylor series
approximations was recently adopted to derive phase
accelerations by Kennedy (2003).
This paper briefly describes the digital differentiator
theory and states the design problems in real-time
dynamic GPS applications. It is followed by a
comprehensive literature review on each method referred
to in the above section. By comparing the various
differentiator designs, a series of first-order Infinite
Impulse Filters (IIR) are presented which are capable of
delivering the derivatives from input signals in real-time
dynamic situations. An adaptive scheme is also proposed
for noise attenuation.
2. Digital f il ter and digital differentiator Design
2.1. Digi tal Filtering
Suppose there is a discrete signal sequence of xn (n is an
integer) with a sampling period of T. A digital filter can
be regarded as a linear combination of the discrete
samples xn-k, together with the previous output yn-k., which
can be defined by the following formula (Hamming 1977,
p2):
∑∑ =
−∞=
−∞=
+=
M
1m
mnk
k
k
knkn ydxcy (1)
where yn is the output of the filter, ck and dk are filter
coefficients which are referred to as the impulse response
of the filter, which is the filter response for a unit input
signal pulse (ibid).
The coefficients of a filter completely define the property
of the filter and selectively suppress or enhance particular
parts of signals. When the coefficients of the second term
on the right hand side of Eq. (1) are nonzero, the filter is
referred to as a recursive filter since the output of yn-k has
been used recursively. The filter coefficients ck and dk are
usually time-invariant in classical filter designs. Their
values are carefully chosen to achieve the desired
filtering result. However, their values can be assigned on-
line to respond to the change of situations in the so-called
adaptive filter design. For practical applications, the
length of a realisable digital filter is always finite.
2.1.1. Transfer function
The transfer function of a discrete filter is defined as the
Z-transform of the filter output signal over the Z-
transform of the input signal, i.e.
=
∞=
−∞=
=
∞=
−∞=
⋅−
=
⋅−
=≡ M
1n
n
k
k
k
k
k
M
1m
mm
k
k
kk
zd1
zc
)y(Zd1
)x(Zc
)z(X
)z(Y
)z(H (2)
where z is a complex variable, and the Z-transform is a
linear transform whereby a discrete-time signal value of
xn is defined as
n
k
knkn xz)z(X)x(Z
−− == (3)
and where z-1 serves as a unit delay operator. The transfer
function is most important in filter design and analysis.
With the transfer function having been determined, one
can directly write out the impulse response of the filter
(filter coefficients), and further analyse the performance
of the filter either in the time domain or in the frequency
domain.
2.1.2. Frequ ency and a mp l itude respons e
The frequency response of a filter is defined as the
discrete Fourier transform of output signals over the
discrete Fourier transform of input signals
Zhang et al.: Real-Time Doppler/Doppler Rate Derivation for Dynamic Applications 97
=
ω−
∞=
−∞=
ω−
⋅−
=
ω
ω
M
1m
mj
m
k
k
kj
k
ed1
ec
)(X
)(Y
)(H (4)
The above frequency response function is obtained by
simply replacing the variable z in the transfer function by
the Fourier transform variable ejw. The frequency
response function allows us to evaluate the frequency
response of a filter on the unit cycle.
Factoring the magnitude of the frequency response into
the following form
)(j
e)(G)(H ωΘ
⋅ω=ω (5)
gives the amplitude response G(ω) which is the gain of
the filter. The phase response Θ(ω) shows the radian
phase shift experienced by each sinusoidal component of
the input signal. The phase and group delays of a filter
give the time delay in seconds experienced by each
sinusoidal component of the input signals:
)(
d
d
delaygroup
)(
delayphase ωΘ
ω
ω
ωΘ
(6)
In the case of a filter that has a linear phase response, the
group delay and the phase delay are identical, for
example whenω⋅
π
=ωΘ 2
)( .
2.1.3. Noise a mpl i fi c ation
A digital filter is a linear combination of input signals
that are usually contaminated by noise. For simplicity we
assume that the noise is Gaussian white, and thus the
error propagation law applies. This allows us to estimate
the noise amplification of the filter. Assume that the noise
of a series of L1 carrier phase measurements
nnn xx
ε
+= 0 is Gaussian white, where 0
n
x stands for
the true value of xn, and then the outcome of the finite
non-recursive filter is
∑∑=
−=
=
−=
−−
=
−=
ε+==
Kk
Kk
Kk
Kk
knk
0
knk
Kk
Kk
kn
kn cxcxcy (7)
and the variance of the filter can be evaluated by
(Hamming 1977, p14)
[]
∑∑∑ =
−=
=
−=
=
−=
σ=ε=
ε
Kk
Kk
2
k
2
x
Kk
Kk
2
kn
2
k
2
Kk
Kk
knk cEccE (8)
This shows that the sum of the squares of each coefficient
of a filter determines the noise amplification of the
filtering process.
Supposing that the variance of a recursive filter is 2
yn
σ,
and applying the preceding procedures, we have
2
y
M
1m
2
kn
k
k
2
k
2
x
2
ymnn dc
σ+σ=σ ∑∑ =
∞=
−∞=
(9)
Let us further assume that 2
y
2
y
2
y
2
ynMn2n1n σ=σ=σ=σ −−−",
and then the variance of the filter can be estimated by
2
x
M
1m
2
mn
k
k
2
k
2
y
d1
c
nσ⋅
≈σ
=
∞=
−∞= (10)
This indicates that we can either roughly estimate the
variance of the recursive filter or “precisely” calculate the
filter variance by computing the initial variance of the
recursive filter using Eq. (10), and then estimating the
variance of the filtered signals using Eq. (9).
2.2. Statement of Problem of Differentiator Design
Differentiator design has been the subject of extensive
investigation in digital signal processing. A main issue is
that a differentiator amplifies noise at high-frequencies
(Carlsson, Ahlen et al. 1991). As GPS signals are of low
98 Journal of Global Positioning Systems
Fig. 1: Power spectral densities for the 1Hz and 10Hz carrier phase signals
frequency character, (see Fig. 1), it is suggested that a
low pass filter would be suitable for the design of
differentiators. However, the change of dynamics in a
system is normally of high frequency. Hence we have to
deal with the complicated high frequencies with a
broad/full band differentiator. Another complication
arises from the signal correlation. It is shown that the
GPS carrier signals can be regarded as Gaussian white
only when the sampling rate is lower than 1Hz; when the
sampling rate goes higher, time correlations must be
considered (Bona 2000; Borre and Tiberius 2000).
Thirdly, the differentiation may be affected due to lack of
information on future signals since the application is real
time oriented. Finally there might be aliasing problems
due to sampling.
So the problem is to get the derivative from GPS
observations where both the signals and noise have
random characteristics. In the case of corrupting noise
being wideband white and the signal being a Gauss-
Markov process (mostlikely for GPS applications), it is
apparent that no differentiator is going to be perfect in
passing the desired derivative whilst suppressing the
noise (Brown and Hwang 1992, p172 ). This is a typical
Wiener filter problem (ibid). The solution is a
compromise between good differentiation and low noise
sensitivity to achieve a small total error.
The Kalman filter is a space-state solution of the Wiener
filter problem (ibid), which is formulated by using the
minimum mean-square-error estimation criterion in a
two-step recursive procedure. By assuming that both the
process driving noise and the measurement noise are
Gaussian white and there is no correlation between them,
it first predicts the signal state using the system dynamic
equation, and then updates the prediction with
measurements to get estimates. A successful Kalman
filter is subject to proper modelling of system dynamics
and the associated stochastic random process. It is
suggested that the less than satisfactory performance of
the Kalman filter in the case of Heber et al. (1997) is not
due to the Kalman filter approach itself, but due to the
improper modelling of the system state when it is highly
dynamic.
When the sampling rate is high, the theoretical difficulties
in Kalman filtering are mainly in the determination of the
random process of system driving noise, and the handling
of correlations of measurement noise and the cross-
correlation between the measurement and signal noises.
Another associated practical problem is the heavy
computational load in real-time data processing. Finally
the outcome of a Kalman filter is a smooth, band-limited
solution (Bruton, Glennie et al. 1999). Therefore, it is
reasonable to find solutions in the frequency domain
rather than in the state space using Kalman filters.
The digital differentiator design oriented in the frequency
domain should still consider the variance of the output.
Thus the criteria of the differentiator may be summarised
as follows:
the magnitude of frequency response is accurate in
low frequencies and is as close to the ideal
differentiator
ω
=
ω
j)(H (Stearns 2003, p127) as
possible in a broad band sense depending on the
system dynamics;
the phase response is linear or approximately linear;
the group delay is acceptably small;
the sum of the squares of filter coefficients can be
minimized; and
easy to be implemented in real time , i.e. to be causal
and low order since there are cycle slips and loss-of-
lock of signals.
Zhang et al.: Real-Time Doppler/Doppler Rate Derivation for Dynamic Applications 99
3. Taylor series approximations
Taylor series approximations have been widely used to
derive differentiators. The differentiators used by
(Cannon, Lachapelle et al. 1998), Hebert (1997) and
Kennedy (2002; 2003) are of low order Taylor series.
They are all in the form of central difference
approximations such as
−=
⋅=
N
Nk
kkn xcy (11)
where N is the order of Taylor series approximation. Fig.
2 depicts the frequency responses of some low order
central difference Taylor series approximations.
00.10.2 0.3 0.40.5 0.6 0.7 0.80.91
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1Normalized M agni tude Res ponse: Low order Taylor series approx i m ations
1st
3rd
5th
ideal
Fig. 2: Frequency response of low order central difference Taylor series
approximation
It is apparent that the higher the order, the closer that a
Taylor series approximation is to the ideal differentiator.
This suggests that broad band differentiators can be
designed based on Taylor series, and this can be observed
in Khan and Ohba (1999), who gave the explicit
coefficients ck by
)!kN()!kN(k
!N)1(
cc
0c
21k
kk
0
+−
==
=
+
(12)
00.10.2 0.3 0.40.5 0.6 0.7 0.8 0.91
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1Normaliz ed M agnitude Responses
5th
10th
20th
30
40
50
100
150
Fig.3: Frequency responses for arbitrary order Taylor series
approximations
As can be seen from Fig. 3, this type of differentiator is
characterised as having zero amplitude response in both
ω=0 and ω =1 (Nyquist frequency). Actually this is the
property of type III FIR filters (Chen 2001,p299 ) which
will be discussed later.
4. Curve fittin g with wi ndow
Jekeli and Garcia (1997) used fifth-order B-splines to
derive phase accelerations, and Fenton and Townsend
(1994) adopted parabolic functions to obtain the precise
Doppler. The referenced curve fitting techniques use
sliding windows wherein the data are fitted into
polynomials using the least squares approach. The
derivative of the central point of a window is obtained by
differentiating the polynomials with respect to time
accordingly.
Bruton (1999) gave an in-depth review of the curve
fitting differentiators. It is concluded that whether a curve
fitting uses a polynomial, a parabola, or a cubic spline,
the resultant differentiator approaches the ideal only at
lower frequencies. Since it is band-limited and lowpass, it
is suitable only for low dynamic or static applications.
Furthermore, performing the least squares estimation
involves intensive computation. Moreover, to obtain the
current derivative at t0, the curve fitting with window
requires the input at tk, which is a signal in the future.
Therefore we may conclude that the windowed curve
fitting approach is inappropriate for real-time dynamic
applications.
5. FIR filters
A Finite Impulse Response (FIR) filter consists of a
series of multiplications followed by a summation. The
FIR filter operation can be represented by the following
equation (Hamming 1977)
−=
⋅=
K
Ki
inin xcy (13)
A filter in this form is named FIR because the response to
an impulse dies away in a finite number of samples. Note
that this form is non-causal and unrealisable. In order to
present a causal FIR differentiator, changing the form is
required. This leads to
=
⋅=
N
0i
Nninxcy (14)
The Fourier series with window are classical in the design
of FIR filters where the impulse response is calculated by
the inverse discrete Fourier transform of the transfer
function, i.e. (Chen 2001)
100 Journal of Global Positioning Systems
ωω
π
πω
πω
ω
dejncMnj
d⋅⋅= =
−=
−− )(
2
1
][
2
)(
])sin[(
)(
])cos[(
Mn
Mn
Mn
Mn
=
π
ππ
(15)
where M=N/2 and the infinite length of Fourier terms is
truncated into finite terms. The truncation may cause a
discontinuity at the edges of the window and leads to
residual oscillations named Gibbs oscillations (ripples in
the amplitude response against frequency). Different
window methods can be used to smooth the glitches,
truncate the filter coefficients, and sharpen the frequency
response. Fig.4 gives the comparison between direct
truncation and applying the Kaiser window technique.
5.1. Type III FIR Differentiator Design
A FIR filter of type III has an odd length and anti-
symmetric impulse response. In this case, the
differentiator’s coefficients are
=
π−
=
Mnfor0
Mnfor
)Mn(
])Mncos[(
)n(cd (16)
where the sine term in Eq. (15) vanishes. To eliminate the
Gibbs phenomenon due to the finite truncation, a window
function is required. Among many windows that are
available, the Kaiser window is most popular. It can be
evaluated to any desired degree of accuracy using the
rapidly converging series of the zero-order Bessel
function of the first kind (Farlex 2004). The ripple of the
stopband can also be controlled by an adjustable variable
α to meet the optimal criteria given by Kumar and Roy
(1988) and Selesnick (2002). With the above procedures,
one can also design FIR differentiators with different cut-
off frequencies.
00.10.2 0.3 0.40.5 0.6 0.7 0.8 0.91
0
0.2
0.4
0.6
0.8
1
1.2 Normal i zed M agni t ude Respons es w/wt K ai ser window
ideal
direc t t runcat i on
Kai ser window
Fig. 4: Magnitude responses of FIR differentiators based on the window
technique
Theoretically, FIR filters of type III can be designed to
meet requirements at nearly all frequencies, as long as we
increase the filter order. However, since the frequency
response to the Nyquist frequency is zero, it is impossible
to design a full band type III differentiator. Although
such filters are causal and are linear in phase, the actual
derivative obtained is with respect to time t-(N/2)T. This
means that the more taps in a FIR filter, the longer the
group delay will be. This property of the FIR filter is
detrimental to the real-time requirements. However it can
be alleviated if the sampling period T is small. The
difficulty is that increasing the sampling frequency will
result in more noisy derivatives. Therefore trade-off and
compromise must be made to introduce this type of FIR
in real-time applications.
5.2. Type IV FIR Differentiators
Since a FIR filter of type III has the limitation that the
amplitude response must go to zero at the Nyquist
frequency, it is impossible to get a full band differentiator
using a finite number of coefficients. This can be shown
in Fig. 3 where transition frequency range of 0.85~1.0 is
associated with the 150th order (length of 301) central
difference Taylor series differentiator.
A FIR filter of type IV has an even length and anti-
symmetric impulse response. The type IV FIR is
preferable to a type III as a differentiator in terms of the
frequency response. This can be evidenced by the
simplest FIR differentiator of yn=xn-xn-1, which has a
frequency response of
2
j
j
1
e
2
sin2je1)(H
z1
)z(X
)z(Y
)z(H
ω
ω−
ω
⋅=−=ω⇒
−==
(17)
The corresponding amplitude response against low order
Taylor series approximations is shown in Fig. 5
00.10.2 0.3 0.4 0.5 0.6 0.7 0.8 0.91
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1Normaliz ed M agnitude Respons e
2-pnt differentiator
1st order
3rd order
5th order
ideal
Fig. 5: Frequency response of the simplest IV differentiator against low
order Taylor Series FIR filters
It can be seen that even though the differentiator is the
simplest form, it is closest to the ideal at low frequencies
(<0.2). It has a better amplitude response for the rest of
frequency band than its type III counterpart of first-order.
It also has a linear frequency response and therefore has a
constant group delay at half the sampling period. The
Zhang et al.: Real-Time Doppler/Doppler Rate Derivation for Dynamic Applications 101
type IV FIR differentiators are superior to the type III
FIR differentiators in terms of the frequency response,
since they have no disadvantageous characteristic of
being zero at ω=1.
Details of type IV differentiator design are referred to
Chen (2001,p332 ). An example differentiator of length 8
(7th order) is given with the transfer function of
7654
321
z0260.0z0509.0z1415.0z2732.1
z2732.1z1415.0z0509.00260.0)z(H
−−−−
−−−
+−+−
+−+−=
(18)
which minimizes
ωω−= π
π−
ωω dej)e(HE
2
2/jNj (19)
Therefore it is an optimal differentiator in the sense of
least squares with an excellent frequency response at high
frequency band. The noise amplification can be
calculated from Eq. (10) as σ2yn=3.2887, which is acc-
eptable so far.
It may be expected that a type IV FIR obtained from the
Remez exchange algorithm (Parks and McClellan 1972)
would be able to deliver a better performance. This is
because the Remez exchange algorithm is a minimax
optimal, i.e. minimize {maximum [Hideal(ω)-Hdisigned(ω)]}
for all frequencies, and is more difficult to
mathematically compute, but guarantees that the worst
case error has been reduced to a quantifiable value. To
verify this, the frequency responses have been depicted in
Fig. 6 for the 7th and 25th-order filters respectively by the
Remez algorithm
Fig. 6: Frequency response of type IV FIR filters by the Remez
exchange algorithm
The FIR filter design by the Remez algorithm is referred
to as the equal ripple design. This is because the method
can suppress the ripples from the Gibbs phenomenon
(Antoniou 1993) to a certain level and turn them into
equal ripples in both the passband and stopband.
It seems that type IV FIR differentiators using the Remez
exchange algorithm will give us a closing solution.
However, the resultant filters provide the first derivative
with system biases and higher level of noise.
Type IV FIR differentiators based on Taylor series
(Khan, Ohba et al. 2000) have also been tested in this
research. It has been found that wideband type IV
differentiators are associated with heavy noise
amplifications and big biases. Our investigation of type
IV FIR filters for differentiator design is still at an early
stage and continuing.
5.3. Other FIR Differe n t iators
In a series of publications, Kumar and Dutta (1988; 1988;
1989; 1989) presented optimal and maximally linear FIR
differentiators for low-frequency, mid-frequency, and
around specific frequency respectively. They gave the
explicit formulae and efficient recursive algorithms to
calculate the impulse response of filters. Their
contributions are highly appreciated, for example, as the
state of art differentiators by Al-Alaoui (1993). In the
case of signals that have low frequency components
contaminated by wideband noise, FIR differentiators of
optimum white-noise attenuation are desired. Kavanagh
(2001) investigated the impact of quantization noise on
signal from systems with low-frequency rates of change.
It is showed that the differentiator proposed by Vainio et
al. (1997)
1Nn0
)1N(N
)n21N(6
h2
n−≤≤
−−
= (20)
has an optimum white-noise attenuation and a constant
group delay. Kavanagh also proposed a better
differentiator for the rate experiencing slow changes
−=
−≤≤
=
=
1Nn
1N
1
2Nn00
0n
1N
1
hn (21)
This differentiator has the characteristic of minimising
the worst-case error. Clearly when N=2 (type IV), this
becomes the simplest two-point differentiator and when
N=3 (type III), this turns into the three-point first-order
differentiator of a Taylor series approximation.
00.2 0.4 0.6 0.8 1
0
0.2
0.4
0.6
0.8
1Magnit u de Res po nse of t he 25t h O rder Rem ez type IV Different i at or
00.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1Normali zed M agni tude Res ponse: 7th order Rem ez
102 Journal of Global Positioning Systems
6. IIR filters
There is another category of filters known as the Infinite
Impulse Filter (IIR). A causal IIR filter is represented by
∑∑ =
=
+=
M
1m
mnm
N
0k
knkn ydxcy (22)
where the output signal at a given instant is obtained as
the weighted sum of the signal xn-k, and the past outputs
of yn-m. As suggested by its name, an impulse input has a
response that lasts forever since the output will be
recursively used. It is the recursive characteristic that
allows IIR filters to be implemented with a lesser order
and better performance when compared with FIR filters.
Thus IIR filters are attractive for real-time applications.
An IIR filter is unstable if its response to a transient input
increases without bound. Poles and zeros are used to
analyse the stability of an IIR filter. The poles are the
roots of the denominator and the zeros are the roots of the
numerator in the transfer function. An IIR filter is stable
if and only if, all poles of H(z) are inside the unit circle
on the z-plane (Stearns, 2003, p83).
The IIR filter cannot be designed by calculating the
impulse response from the known frequency response as
is the case in FIR designs. Many IIR filters can be
derived from the analogue filter designs and then
transformed into the sampled z-plane. Another popular
method is the bilinear transform. The IIR differentiator
design has been of considerable interest (Rabiner and
Steiglitz, 1970). Among various recursive differentiator
designs, Al-Alaoui’s second order IIR family (1992;
1993; 1994) has been highly acknowledged and widely
used, for example (Chen and Lee, 1995).The novel
approach of designing digital differentiators by Al-Alaoui
is an extension of the method in designing analogue
differentiators by using integrators. That is, in the
analogue signal processing, differentiators are often
obtained by inverting the transfer functions of analogue
integrators.
The general procedures to derive the Al-Alaoui family
are as follows
design an integrator that has the same range and
accuracy as the desired differentiator;
invert the obtained transfer function of the
integrator;
reflect the poles that lie outside the unit circle to
inside, in order to stabilise the resultant transfer
function; and
compensate the magnitude using the reciprocals
of the poles that lie outside the circle.
6.1. Al-Alaoui’s First-Order Differentiator
A first-order IIR differentiator was developed by Al-
Alaoui (1993) with an effective range 0.78 of the Nyquist
frequency based on a non-minimum phase digital
integrator. The integrator is a synthesis of the rectangular
integrator and the trapezoidal integrator. By assigning
weighting factors of ¾ and ¼ to the transfer functions of
the integrators respectively, the ideal integrator, which
has the following transfer function, is approximated
1z
7z
8
T
)1z(2
)1z(T
4
1
1z
T
4
3
)z(H
4
1
)z(H
4
3
)z(HTRI
+
⋅=
+
+
⋅=
+=
(23)
Reflecting the zero z=-7 with its reciprocal -1/7, and
compensating the magnitude by multiplying r=7, results
in a minimum phase digital integrator with the transfer
function
1z
7
1
z
8
T7
)z(HI
+
= (24)
Inverting the above transfer function yields the Al-
Alaoui’s stabilized IIR differentiator of the first order
7
1
z
1z
T
7
8
)z(HD+
⋅= (25)
The characteristics of this differentiator is shown in Fig. 7.
This differentiator is able to approximate the ideal
differentiator up to 0.78 of the full band, and has an
outstanding “linear phase” response. Al-Alaoui reported
that within the effective frequency range, it has a less
than 2.0% magnitude error. Since it is of first-order, the
delay of the filter is just half of the sample thus it meets
every requirement to be used in real-time.
Zhang et al.: Real-Time Doppler/Doppler Rate Derivation for Dynamic Applications 103
00.10.2 0.3 0.40.5 0.6 0.7 0.80.91
0
0.5
1Normali zed M agni tude Response
Ideal
A l-Alaoui
00.10.2 0.3 0.40.5 0.6 0.7 0.80.91
0
50
100 P hase response of Al-A l aoui 's Di fferentiatorDif(Al -Alaoui )
Ideal
-5 -4-3 -2 -101 2 345
-1
0
1
Real Part
Im aginary Part
Zero/P ol e response of A l -A l aoui 's different i ator
Zeros at 1
Poles at -1/7
Fig. 7: Characteristics of the first order IIR differentiator
6.2. First-order IIR Differentiator Family
Al-Alaoui contributes the above differentiator as an
individual. However, a family of such first order
differentiators can be derived following his methodology.
That is, while Al-Alaoui designates the weighting factors
of ¾ and ¼ empirically, we may get the optimal weights
experimentally. To achieve this, a variable
α
is
introduced to adjust the weighting factor in the way of
)1z(2
1
1
z)1(T
)1z(2
)1z(T
)1(
1z
T
)z(H)1()z(H)z(H TRI
α−
α+
+α−
=
+
α−+
⋅α=
α−+α=
(26)
where 0<
α
<1 serves as a tuner to adjust the integrator so
that it better closes to the ideal.
α
=¾ can be used as a
good reference to refine the integrator in the desired
range of frequencies. Obviously it has a zero outside the
unit circle. Applying Al-Alaoui’s procedure to reflect the
zero with its reciprocal and to compensate the magnitude,
a variable integrator is obtained as
)1z(2
1
1
z)1(T
)z(HI
α+
α−
+α+
= (27)
Inverting the transfer function gives a new set of
differentiators with transfer functions as
α+
α−
+α+
=
1
1
z)1(T
)1z(2
)z(HD (28)
Since 1-
α
<1+
α
, the pole is well located inside the unit
circle and the resultant differentiators are, therefore,
stable. Setting
α
=¾ gives the transfer function proposed
by Al-Alaoui, and slightly changing
α
around ¾ results in
differentiators which outperform in target bandwidth. The
noise amplification of this kind of differentiator can be
evaluated using Eq. (10), which is only slightly noisier
than the simplest two-point differentiator.
7. Conclusions
The general theory on digital filter design has been
introduced. The aim of this research is to find appropriate
differentiators that can be used to derive Doppler
shifts/Doppler rates from GPS observables in real-time,
dynamic applications.
It is concluded that the differentiators obtained from both
curve fitting and Kalman filtering require intensive
computation and are lowpass. Thus they are not suitable
for real-time dynamic applications. Type III FIR
differentiators have the inherent nature in frequency
response of approaching zero at Nyquist frequency. To
extend the performance of type III FIR filters in the
higher frequency bands, one has to increase the filter taps.
This causes difficulties in managing the data since there
are cycle-slips and loss-of-lock signals. It also results in a
longer group delay that is detrimental for real-time
applications where instant response is desired.
Type IV FIR differentiators using Fourier series have
been found to have outstanding frequency response,
however, they are noisy and biased. It is found that only
the Kavanagh’s differentiators of type IV deliver good
first derivatives. However, they approximate the ideal
differentiator only at low frequencies (lower than 0.2 of
Nyquist frequency). Type III FIR filters can be used to
derive Doppler/Doppler rate “observables” in post
processing mode. Higher order central difference
approximations using Taylor series might outperform
windowed Fourier series since there is no truncation and
the associated Gibbs phenomenon.
It is demonstrated that IIR filters are more favourable for
real-time application. Since the outputs of the filter are
recursively used, they have much lower orders than the
FIR filters. The first-order IIR differentiator from Al-
Alaoui is ideal in terms of the frequency response, phase
linearity and half sample group delay. The proposed class
of first-order IIR differentiators allows us to choose an
optimum in the desired frequency range.
104 Journal of Global Positioning Systems
It is suggested that Kavanagh’s differentiators can be
used in static or in constant velocity modes. The proposed
first-order IIR differentiators can be adaptively used
when systems experience high dynamics.
ACKNOWLEDGEMENTS
The research described in this paper was partially
supported by the Cooperative Research Centre (CRC) for
microTechnology, Australia endorsed to A/Prof. Kefei
Zhang. Partial financial support from the Australia
Research Council Linkage Project (LP0455170) endorsed
to a research consortium led by A/Prof Kefei Zhang is
highly appreciated. Mr. Ravindra Babu at the SNAP,
UNSW is highly acknowledged for consultations and
discussions.
REFERENCES
Al-Alaoui M.A. (1992): Novel Approach to Designing Digital
Differentiators. Electronics Letters 28(15): 1376-1378.
Al-Alaoui M.A. (1993): Novel Digital Integrator and
Differentiator. Electronics Letters 29(4): 1376-1378.
Al-Alaoui M.A. (1994): Novel IIR Differentiator from the
Simpson Integration Rule. IEEE Transactions on Circuits
and System-I Fundamental , Theory and Applications
41(2).
Antoniou A. (1993): Digital Filters: Analysis, Design and
Appications. New York, McGraw-Hill 0-07-002121X.
Bona P. (2000): Precision, Cross Correlation, and Time
Correlation of GPS Phase and Code Observations. GPS
solutions 4(2): 3-13.
Borre K. and Tiberius C. (2000): Time Series Analysis of GPS
Observables. proceedings of the ION GPS 2000, Salt Lake
City, UT, 1885-1894.
Brown R.G. and Hwang P.Y.C. (1992): Introduction to
Random Signals and Applied Kalman Filtering. New
York, John Wiley &Sons, INC 0471-52573-1.
Bruton A.M., Glennie C.L. and Schwarz K.P. (1999):
Differentiation for High-Precision GPS Velocity and
Acceleration Determination. GPS solutions 2(4): 7-21.
Cannon M.E., Lachapelle G., Szarmes M.C., Hebert J.M., Keith
J. and Jokerst S. (1998): DGPS Kinematic Carrier Phase
Signal Simulation Analysis for Precise Velocity and
Position Determination. Journal of the Institute of
Navigation 44(2): 231-245.
Carlsson B., Ahlen A. and Sternad M. (1991): Optimal
Differentiation Based on Stochastic Signal Models. IEEE
Transactions on Signal Processing 39(2): 341-353.
Chen C.-K. and Lee J.-H. (1995): Design of High-Order
Digital Differentiators using L1 Error Criteria. IEEE
Transactions on Circuits and Systems II: Analogy and
Digital Signal Processing 42(4): 287-291.
Chen C.-T. (2001): Digital Signal Processing: Spectral
Computation and Filter Design. New York, Oxford,
Oxford University Press 0-19-513638-1.
Farlex (2004): TheFreeDictionaryCom: Kaiser Window. http://
encyclopedia.thefreedictionary.com/Kaiser%20window.
Fenton P. and Townsend B. (1994): NowAtel Communications
Ltd,-What's New? Proceedings of KIS'94, Banff, AB, CA,
25-29.
Hamming R.W. (1977): Digital Filters. Englewood Cliffs,
Prentice-Hall, Inc. 0-13-212571-4.
Harvey N. (2004): Doppler and Carrier Phase Measures in
Trimble 5700 GPS Receiver. Christchurch, New Zealand
(Personal communication).
Hebert C.J. (1997): High Accuracy GPS Velocity Using the
Carrier Phase Observable. Proceedings of KIS'97, Banff,
AB, Canada, 265-269.
Hebert C. J., Keith J., Ryan S., Szarmes M., Lachapelle G. and
Canon M. E. (1997): DGPS Kinematic Carrier Phase
Signal Simulation Analysis for Precise Aircraft Velocity
Determination. Proceedings of the ION GPS-97,
Albuquerque, NM,
Jekeli C. and Garcia R.. (1997): GPS Phase Accelerations for
Moving-base Vector Gravimetry. Journal of Geodesy
71(10): 630-639.
Kavanagh R.C. (2001): FIR Differentiators for Quantized
Signals. IEEE Transactions on Signal Processing 49(11):
2713-2720.
Kennedy S. (2002): Precise Acceleration Determination from
Carrier Phase Measurements. Proceeding of the ION
GPS-2002, Portalnd Oregon, USA, 962-972.
Kennedy S.L. (2003): Precise Acceleration Determination
from Carrier-Phase Measurements. Navigation, Journal
of the Institute of Navigation 50(1): 9-19.
Khan I.R. and Ohba R. (1999): New Design of Full Band
Differentiators Based on Taylor Series. IEEE
Proceedings-Vision, Image and Signal Processing 146(4):
185-189.
Khan I. R., Ohba R. and Hozumi N. (2000): Mathematical
Proof of Explicit Formulas for Tap-coefficients of Full-
band FIR Digital Differentiators. IEE Proceedings-
Vision, Image and Signal Processing 147(6): 553-555.
Kumar B. and Dutta S. C. R. (1988): Coefficients of Maximally
Linear, FIR Digital Differentiators for Low Frequencies.
Electronics Letters 24(9): 563-565.
Kumar B. and Dutta S.C.R. (1988): Design of Digital
Differentiators for Low Frequencies. Proceedings of the
IEEE 76(3): 287-289.
Kumar B. and Dutta S.C.R. (1989): Design of Efficient FIR
Digital Differentiators and Hilbert Transforms for
Midband Frequency Ranges. Internal 36(6).
Kumar B. and Dutta S.C.R. (1989): Maximally Linear FIR
Digital Differentiators for Midband Frequencies. Int. J.
Circuit Theory Appl. 17: 483-488.
Zhang et al.: Real-Time Doppler/Doppler Rate Derivation for Dynamic Applications 105
Parks T.W. and McClellan J.J. (1972): Chebyshev
Approximatin for Nonrecursive Digital Filters with
Linear Phase. IEEE Trans. Circuit Theory 19: 189-194.
Rabiner L. and Steiglitz K. (1970): The Design of Wide-band
Recursive and Nonrecursive Digital Differentiators. IEEE
Transactions on Audio and Electroacoustics 18(2): 204-
209.
Selesnick I.W. (2002): Maximally Flat Low-pass Digital
Differentiator. IEEE Transactions on Circuits and Systems
II: Analogue and Digital Signal Processing 49(3): 219.
Stearns S.D. (2003): Digital Signal Processing with Examples
in MATLAB. Boca Raton, London, New York,
Washington,D.C., CRC press 0-8493-1091-1.
SuperstarII (2004): Superstar II Product Brochure. 2004
http://www.seabed.nl/pdf/superstar.pdf.
Vainio O., Renfors M. and SaramT. (1997): Recursive
Implementation of FIR Differentiators with Optimum
Noise Attenuation. IEEE Transactions on Instrum. Meas.
46: 1202-1207.
Wolfe J.D., Williamson W.R. and Speyer J.L. (2003):
Hypothesis Testing for Resolving Integer Ambiguity in
GPS. Navigation, Journal of the Institute of Navigation
50(1): 45-56.