__FIELD OF THE INVENTION__
The present invention relates generally to positioning
systems. More particularly, the present invention pertains to the application of
a filter, such as the extended Kalman filter, in such positioning systems, or in
systems that provide estimates of one or another aspect besides position of the
motion of an object, using information acquired from cellular base stations and
from satellites.

__BACKGROUND OF THE INVENTION__
A satellite positioning system (SATPS) receiver generally
determines its positions by triangulating its Line of Sight (LOS) range to several
satellites or space vehicles. A GPS receiver, for example, computes a 4-dimensional
solution involving latitude, longitude, altitude, and time using the LOS ranges
to as few as four satellites. The accuracy of the solution is a direct function
of the accuracy of the range measurements.

SATPS receivers are growing rapidly in popularity and application.
GPS receivers, for example, are now common in aviation, marine, and terrestrial
applications. An increasingly common terrestrial application for GPS receivers is
in automobiles. In the automotive context, the vehicle's location is typically displayed
on an electronic display of a street map. It is vital in this context, therefore,
to provide the driver with continuously updated position solutions, collectively
called a "ground track," that accurately track the vehicle's movement from one moment
to the next. Experience shows that consumers consider ground-track fidelity as one
of the most important criteria in choosing a receiver. It is extremely important,
therefore, that the ground-track displayed on the GPS receiver's electronic map
not have spurious jumps, stair steps, spikes, jigs, or jogs that are unrelated to
the vehicle's actual path.

There are a number of factors, however, that may cause
discontinuities in the position solutions used to determine the ground-track of
an automotive SATPS receiver. One source of position solution discontinuity is "Selective
Availability" (SA), which restricts the accuracy of civilian GPS receivers to roughly
100 meters. SA is intentionally used by the U.S. Government for purposes of national
security. The Department of Defense (DOD) implements SA by purposely injecting error
into the satellite range signals.

Another common source of position solution discontinuity
is due to the phenomenon known as multi-path, where the true LOS signal from a given
satellite reaches the GPS receiver's antenna, along with additional signals providing
supposedly the same information, the additional signals caused by reflection from
nearby objects, such as buildings or cliffs. The multi-path phenomenon is particularly
troublesome for automotive receivers because they are frequently used in cities
and surrounded by tall buildings. This environment is sometimes called an "urban
canyon" due to the canyon-like setting it resembles. Regardless of source, multi-path
can be a very vexing problem because the additional signals may be very strong,
but very wrong.

Yet another source of position solution discontinuity is
that the constellation of satellites used by a SATPS receiver can change; the SATPS
receiver may see a different constellation of satellites from one moment to the
next. If the GPS receiver is located in an urban canyon environment, for example,
individual satellites may become blocked and later unblocked as the receiver moves
past different buildings. The discontinuity arises in this situation because the
error in a position solution is based on the constellation of satellites used. (Two
satellites located in approximately the same direction will provide position information
with larger error than two satellites in very different directions, all other things
being equal.) If the position solution provided by a GPS receiver is suddenly based
on a different constellation, the different error may cause a jump or discontinuity
in position.

It is known in the art to use a Kalman filter to account
for the uncertainties in measurement data provided to a positioning system receiver.
Fig. 1 is a simplified flow diagram of a conventional GPS-type positioning system
10 including an RF antenna 11, a measurement engine 12 and a Kalman filter 14, providing
a position estimate *x̂*(*k*) for position at time instant k. The
measurement engine 12 receives RF signals from a plurality of orbiting satellites
via the antenna 11 and provides the Kalman filter 14 with measured position and
velocity, i.e. measured state information as opposed to the predicted state information
provided by the Kalman filter based on the measured values.

The construction of the measurement engine 12 varies from
application to application. Generally, the measurement engine 12 contains the analog
electronics (e.g. preamplifiers, amplifiers, frequency converters, etc.) needed
to pull in the RF signal, and further includes a code correlator for detecting a
particular GPS code corresponding to a particular satellite. The measurement engine
12 estimates the line of sight (LOS) range to a detected satellite using a local,
onboard GPS clock and data from the satellite indicating when the satellite code
was transmitted. The LOS ranges determined this way are called *pseudo-ranges*
because they are only estimates of the actual ranges, based on the local detection
time. In the positioning system 10 of Fig. 1, the measurement engine 12 converts
the pseudo-ranges it acquires over time to measurements z(k) of the state of the
process, i.e. to a position and velocity of the moving object whose position is
being determined.

In estimating the state x(k) of a process (such as the
motion of a vehicle), a (standard) Kalman filter relies on the assumption that the
process evolves over time according to a linear stochastic difference equation,
such as:
$$\mathrm{x}\left(\mathrm{k},+,1\right)=\mathrm{A}\left(\mathrm{k}\right)\mathrm{x}\left(\mathrm{k}\right)+\mathrm{B}\mathrm{u}\left(\mathrm{k}\right)+\mathrm{w}\left(\mathrm{k}\right)$$

where w(k) is the process noise, A(k) is an *n×n* matrix relating the
state at time step k to the state at time step k+1 in the absence of a driving function,
u(k) is a control input to the state x, and B is an nxl matrix that relates the
control input to the state x. A standard Kalman filter further relies on the assumption
that the measurements used by the Kalman filter in estimating the state of the process
are linearly related to the state of the process, i.e. that a measurement z(k) (i.e.
the measured position at time instant k) is corresponds to the state x(k) (i.e.
the actual position at time instant k) according to an equation of the form:
$$\mathrm{z}\left(\mathrm{k}\right)=\mathrm{H}\left(\mathrm{k}\right)\mathrm{x}\left(\mathrm{k}\right)+\mathrm{s}\left(\mathrm{k}\right)$$

where s(k) is the measurement noise, and the mxn matrix H relates the state x(k)
to the measurement z(k). If either of these assumptions do not apply, a standard
Kalman filter cannot be used, or at least not used without first deriving from information
provided (by e.g. satellites) about the process measurement data that is, at least
approximately, linearly related to the state of the process, and at least not without
taking measurements at close enough intervals that the effects of any non-linearity
in the evolution of the process does not become important from one measurement to
the next.

In the case of a positioning system in a vehicle where
the positioning system uses information from satellites as the basis for position
measurements, the dependence between the state (position, but in general also velocity)
of the process (motion of the vehicle hosting the positioning system), and the position
measurements is non-linear, i.e. instead of equation (2), the satellite provides
a pseudo-range:
$$\mathrm{\&rgr;}\left(k\right)=\sqrt{{\displaystyle \sum _{i=1}^{3}{\left[{x}_{i},\left(k\right),-,{x}_{i}^{s},\left(k\right)\right]}^{2}}}+c{t}_{o}$$

in which x_{i}(k) is, e.g. the i^{th} component of the three-dimensional
position of the vehicle, *x*^{5}, is the same for the satellite, t_{o}
is a clock offset between the positioning system clock in the vehicle and the satellite
clock, and c is used to designate the speed of light.

It is known in the art, generally, that in case either
of the assumptions on which a standard Kalman filter relies is violated, a so-called
extended Kalman filter can be used. (Sometimes what is called a converted-measurement
Kalman filter (CMKF) is used, instead of an extended Kalman filter, in case of non-linear
measurements.) In such a filter, there is a linearizing about the current mean and
covariance of the estimated state of a process and for the measurements used in
estimating the state of the process. This linearizing makes possible the use of
a Kalman-type filter in estimating the state of a process, an approach that is recursive
(the solution at a preceding instant, along with a current measurement, leading
to the solution at the next instant), as compared with other approaches, such as
iterative or direct solution, that rely on the current measurement (alone) for solution.

Although extended Kalman filters are known in the art,
the prior art does not teach using such a filter, or using other types of filters,
in a positioning system based on information provided by cellular base stations,
or, ideally by both cellular base stations and satellites, so that at any one instant
of time, information either from only cellular base stations, or from both cellular
base stations and satellites, are used as one measurement of the state of the system.

What is further needed is a way for a positioning system
to combine multiple such measurements, at any instant of time, so as to further
increase the reliability of the positioning system. Since cellular base station
data does not suffer from the same errors as affect satellite data, a positioning
system can be made more reliable if it uses both kinds of data separately and simultaneously,
providing therefore at least two measurements at the same instant of time, one purely
based on satellite data and one purely based on cellular data. (The main sources
of noise and inaccuracy in cellular data are not the same as for satellites; instead
of atmosphere, multi-path, and selective availability, error in a cellular system
is caused by obstacles in the signal path, range measurement quantization, and interfering
frequencies.)

Even though the use of two different types of sources of
positioning measurements can provide an inherently more reliable positioning system,
it is still true that at any instant of time a particular measurement or a particular
piece of information on which a measurement is based can be greatly impaired by
one or another source of error (not usually the selective diversity, which is intended
only to "dither" the correct measurement). Thus, it remains advantageous, even in
case of a hybrid positioning system (hybrid in the sense of different kinds of sources
of measurement), to systematically eliminate some measurements or some pieces of
information (e.g. a pseudo-range value) from which a measurement is determined (by
triangulation of several pseudo-range values), the eliminating based on
*statistical* measures of the error of the measurements or pieces of information.

__SUMMARY OF THE INVENTION__
Accordingly, the present invention provides a filter for
a generalized positioning system, and also a corresponding generalized positioning
system and method. All of these are based on using information from cellular base
stations, at least part of the time, in providing a succession of estimates of one
or another aspect of the state of motion of an object, such as estimates of position
or velocity or both. The generalized positioning system includes: a filter, responsive
to an initial estimate of a state of motion of the object, and an initial value
of state covariance, and further responsive to a succession of measurements of the
state, each corresponding to a different instant of time, for statistically determining
the succession of state estimates; and a measurement engine, responsive to information
provided by cellular base stations, from which a measurement of the state can be
determined, for providing the succession of measurements of the state. In some aspects
of the invention, the filter is a non-linear filter, such as an extended Kalman
filter, and in other aspects of the invention the filter is a linear filter.

In a further aspect of the invention, the generalized positioning
system also includes: means for determining association probabilities for a plurality
of state measurements all at a same instant of time, each association probability
corresponding to a particular state measurement at the instant of time, and also
means for determining an association probability providing a value for the probability
that none of the plurality of measurements is correct; means for combining into
a single measurement innovation each individual innovation using the association
probabilities as weightings; and means for determining the covariance of the updated
state based on the individual association probabilities, the association probability
for providing a value for the probability that none of the measurements is correct,
and the combined measurement innovation.

In another aspect of the invention, the filter of the generalized
positioning system is responsive to measurements provided by a measurement engine
receiving information from both satellites and cellular base stations.

In yet another aspect of the invention, the generalized
positioning system also includes an output display, wherein both the measurement
engine and the output display are co-located with the object for which the state
information is to be measured, and further wherein the filter, the means for determining
association probabilities, the means for combining into a single measurement innovation
each individual innovation, and the means for determining the covariance of the
updated state are all hosted by a non-local facility, one that is separate and remote
from the object for which the state information is to be measured, and wherein the
non-local facility and the measurement engine and the output display are in radio
communication.

In yet still another aspect of the invention, the generalized
positioning system filter is responsive to a finite gate width &ggr; of a validation
region, and rejects a measurement if it falls outside the validation region.

In yet even a further aspect of the invention, the generalized
positioning system also includes: an output display, wherein both the measurement
engine and the output display are co-located with the object for which the state
information is to be measured, and further wherein the filter is hosted by a non-local
facility, one that is separate and remote from the object for which the state information
is to be measured, and wherein the non-local facility and the measurement engine
and the output display are in radio communication.

In one particular application, the generalized positioning
system filter is an extended Kalman filter (EKF).

__BRIEF DESCRIPTION OF THE DRAWINGS__
The above and other objects, features and advantages of
the invention will become apparent from a consideration of the subsequent detailed
description presented in connection with accompanying drawings, in which:

- Fig. 1 is a flow diagram/ block diagram of a positioning system based on a Kalman
filter;
- Fig. 2 is a flow diagram/ block diagram of a positioning system based on an
extended Kalman filter, according to the present invention, with information acquired
from both satellites and cellular base stations;
- Fig. 3 is a flow diagram/ block diagram of a positioning system based on an
extended Kalman filter, according to the present invention, modified to incorporate
probability data association (PDA), with statistical pre-processing of inputs, and
with information acquired from both satellites and cellular base stations;
- Fig. 4 is a flow diagram/ block diagram of a distributed positioning system
similar to the system of Fig. 3 but with some components remote from the object
making use of the positioning system; and
- Fig. 5 is a flowchart of a method of filtering measurements of state information,
according to the present invention.

__BEST MODE FOR CARRYING OUT THE INVENTION__
In the preferred embodiment, as is described below, the
present invention involves the use of an extended Kalman filter (EKF) as part of
a positioning (navigation) system. It is to be understood, however, that the present
invention is not restricted to using only an EKF or any other non-linear filter
in such a system, but instead can use any kind of filter in a position-measuring
system, or in a system that measures another aspect of the motion of an object,
such as velocity, and so is here called a generalized positioning system (because
it does not necessarily measure position).

*Background on Standard Kalman Filter and Notation Used*
Referring again to Fig. 1, in a positioning system 10 according
to the prior art, a Kalman filter 24 is sometimes used to estimate the state x(k)
of a process, such as the position and velocity of a vehicle in motion, based on
measurements related to the process. The (unknown and never perfectly knowable)
state x(k) of a discrete-time controlled process is assumed governed by equation
(1), and measurement z(k) is assumed to correspond to the state x(k) according to
equation (2). The measurement z(k) is provided by a measurement engine 22 based
on information received by an antenna 21.

The Kalman filter uses the measurements z(k), along with
knowledge of the dynamical equation (eq. 1), to determine an *a posteriori*
estimate *x̂*(*k*|*k*-1) that takes into account the measurements
z(k) more or less, depending on a Kalman gain, K(k), calculated as described below.
The estimate *x̂*(*k*|*k*-1) is said to be a *posteriori*
in that it is made after the measurements z(k), taking those measurements into account.

A central assumption in using a Kalman filter is that the
process noise and the measurement noise are both white, and with normal probability
distributions:
$$\mathrm{p}\left(\mathrm{w},\left(\mathrm{k}\right)\right)\sim \mathrm{N}\left(0,,,\mathrm{Q},\left(\mathrm{k}\right)\right)$$
$$\mathrm{p}\left(\mathrm{s},\left(\mathrm{k}\right)\right)\sim \mathrm{N}\left(0,,,\mathrm{R},\left(\mathrm{k}\right)\right)$$

where Q(k) is the process noise covariance (matrix), and R(k) is the measurement
error covariance (matrix). (Here, N(a,b) is the normal distribution with mean a
and covariance b.) In the implementation of a Kalman filter, the R(k) and Q(k) for
the filter are sometimes measured, or values are assumed, once and for all before
using the filter. The selection of values to be used for R(k) and Q(k) is said to
constitute "tuning" the filter, and is sometimes performed off-line, using another
Kalman filter. In some applications, R(k) and Q(k) are taken to be constant.

Referring again to Fig. 1, the operation of a Kalman filter,
in general, is in two parts: a time update part and a measurement update part. In
the time update part, the output of the Kalman filter (i.e. the a *posteriori*
estimated state *x̂*(*k*)) is used in the assumed dynamical equation
(1) to determine an a *priori* estimate *x̂*(*k+*1|*k*)
of the next state of the process:
$$\widehat{x}\left(k,+,1,|k\right)=\mathrm{A}\left(\mathrm{k}\right)\widehat{x}\left(k\right)+\mathrm{B}\mathrm{u}\left(\mathrm{k}\right).$$
Here, the notation *x̂*(*k*+1|*k*) is used to show expressly
that the value obtained is based on assuming the a *posteriori* estimate
*x̂*(*k*).

The time update part also determines the so-called a
*priori* estimate error covariance using the equation:
$$\mathrm{P}\left(\mathrm{k}|,\mathrm{k},-,1\right)=\mathrm{A}\left(\mathrm{k}\right)\mathrm{P}\left(\mathrm{k}\right){\mathrm{A}}^{\mathrm{T}}\left(\mathrm{k}\right)+\mathrm{Q}\left(\mathrm{k}\right).$$

where the a *priori* estimate error covariance P(k|k-1) is based on the a
*priori* estimate error x(k) *- x̂*(*k*|*k -* 1), and
is defined as:
$$\mathrm{P}\left(\mathrm{k}|,\mathrm{k},-,1\right)=\mathrm{E}\left[\left(\mathrm{x},\left(\mathrm{k}\right),-,\widehat{x},,\left(k|,k,-,1\right)\right),,{\left(\mathrm{x},\left(\mathrm{k}\right),-,\widehat{x},,\left(k|,k,-,1\right)\right)}^{\mathrm{T}}\right]$$

in which E(a) is the expectation of the random variable *a*, and where P(k)
is the corresponding a *posteriori* estimate error covariance, based on the
a *posteriori* estimate error x(k) - *x̂*(*k*).

In the second, measurement update part of the operation
of a Kalman filter, the Kalman gain K(k+1) is calculated, using:
$$\mathrm{K}\left(\mathrm{k},+,1\right)=\mathrm{P}\left(\mathrm{k}+1|,\mathrm{k}\right){\mathrm{H}}^{\mathrm{T}}\left(\mathrm{k},+,1\right){\left[\mathrm{H},,\left(\mathrm{k},+,1\right),,\mathrm{P},,\left(\mathrm{k}+1|,\mathrm{k}\right),,{\mathrm{H}}^{\mathrm{T}},,\left(\mathrm{k},+,1\right),+,\mathrm{R},,\left(\mathrm{k},+,1\right)\right]}^{-1}$$
and is then used to determine the a posteriori state estimate *x̂*(*k*+1),
according to:
$$\widehat{x}\left(k,+,1\right)=\widehat{x}\left(k+1|,k\right)+\mathrm{K}\left(\mathrm{k},+,1\right)\mathrm{\&ugr;}\left(\mathrm{k},+,1\right),$$

where u (k+1) is the so-called measurement innovation (or measurement residual),
defined as:
$$\mathrm{\&ugr;}\left(\mathrm{k},+,1\right)=\mathrm{z}\left(\mathrm{k},+,1\right)-\mathrm{H}\left(\mathrm{k},+,1\right)\widehat{x}\left(k+1|,k\right).$$
The estimate *x̂*(*k*+1) is then output to an output display 16
as the Kalman filter's estimate of the state of the process as the time corresponding
to k+1, and is also used in the next first, time update part of the operation of
the Kalman filter. Finally, the first, update part of the operation also requires
a new a *posteriori* estimate error covariance, P(k+1), which is calculated
using:
$$\mathrm{P}\left(\mathrm{k},+,1\right)=\left[\mathrm{I},-,\mathrm{K},,\left(\mathrm{k},+,1\right),,\mathrm{H},,\left(\mathrm{k},+,1\right)\right]\mathrm{P}\left(\mathrm{k}+1|,\mathrm{k}\right).$$
The operation of the Kalman filter must of course be started with initial estimates
for *x̂*(*k*+1|*k*) and for P(k+1|k), and the prior art includes
suggestions for how to select initial values.

*Extended Kalman Filter with Cellular Measurement Information*
Referring now to Fig. 2, in the present invention, a positioning
system 20 is based on an extended Kalman filter 24 (EKF), so as to enable using
as a (single) measurement at an instant of time, information acquired from cellular
base stations, and preferably, also from satellites. From such information, as in
a standard Kalman filter, a time update of the state of the system (i.e. the position
and, optionally, velocity of the object hosting the positioning system) and then
a measurement update is performed. Measurement z(k) is provided by a heterogeneous
measurement engine 22 at an instant k of time, based on information received by
different antennae 21a and 21b, which are different in kind in that one is intended
to receive information from one kind of source of information, such as satellites,
and the other is intended to receive information from another kind of source of
information, such as cellular base stations.

At some times in the operation of a positioning system
using the present invention, information from cellular bases stations may be so
less reliable than information from satellites that it is advantageous to use only
satellite data. In the best mode, a positioning system according to the present
invention, besides having the ability to automatically choose to disregard the cellular
data (or the satellite data), would allow a user to command that the cellular data
(or the satellite data) not be used.

In performing the state and time updates, and thus ultimately
providing an estimate *x̂*(*k*) of the state of the system based
on measurements z(k) and previous estimate *x̂*(*k*-1), the EKF
does not assume the linear stochastic difference equations (1) and (2), but instead
assumes that the process is governed by a non-linear stochastic difference equation,
indicated as:
$$x\left(k,+,1\right)=f\left[x,\left(k\right),,,u,\left(k\right),,,w,\left(k\right)\right]$$
(where again u(k) is a driving function, and w(k) is a zero mean process noise)
and, correspondingly, that the measurement z(k) and state x(k) have a non-linear
relationship, indicated as:
$$z\left(k\right)=h\left[x,\left(k\right),,,s,\left(k\right)\right]$$
(where s(k) is again the zero-mean measurement noise).

With these assumptions and knowledge of the functions f
and h, operation of the EKF is based on approximating the state and measurement
vectors (the approximation of the measurement vector for use in determining the
measurement innovation) using equations (12) and (13) but using a value of zero
for both the process and measurement noise, w and s, i.e.
$$\tilde{x}\left(k,+,1\right)\approx f\left[\widehat{x},\left(k\right),,,u,\left(k\right),,,0\right],\hspace{0.17em}\mathrm{and}$$
$$\tilde{z}\left(k\right)=h\left[\tilde{x},\left(k\right),,,0\right]$$
and asserting that the approximated state vector is the a *priori* (predicted)
state vector, i.e.
$$\widehat{x}\left(k+1|,k\right)=\tilde{x}\left(k,+,1\right).$$

The non-linear equations (12) and (13) are then linearized
by performing essentially a Taylor series expansion about *x̃* and
*z̃*, as given by equations (14) and (15). This linearizing process produces
the equations of the EKF (see below) in terms of four Jacobians, which are each
matrices, having components:
$${A}_{i,j}=\frac{\partial {f}_{\left[i\right]}}{\partial {x}_{\left[j\right]}}\left[\widehat{x},\left(k\right),,,u,\left(k\right),,,0\right],$$
$${W}_{i,j}=\frac{\partial {f}_{\left[i\right]}}{\partial {w}_{\left[j\right]}}\left[\widehat{x},\left(k\right),,,u,\left(k\right),,,0\right],$$
$${H}_{i,j}=\frac{\partial {h}_{\left[i\right]}}{\partial {x}_{\left[j\right]}}\left[f,\left(\widehat{x},\left(k\right),,,u,\left(k\right),,,0\right),,,0\right],\hspace{0.17em}\mathrm{and}$$
$${V}_{i,j}=\frac{\partial {h}_{\left[i\right]}}{\partial {s}_{\left[j\right]}}\left[f,\left(\widehat{x},\left(k\right),,,u,\left(k\right),,,0\right),,,0\right]$$

where x and z are actual state and measurement vectors respectively, *x̂*(*k*|*k*-1)
is an a *priori* state estimate, *x̂*(*k*) is an a
*posteriori* state estimate, and w and s again represent the process and measurement
noise. Thus, e.g. the i-j^{th} element of the Jacobian matrix A is the partial
derivative of the i^{th} component of the vector function f with respect
to the j^{tn} component of the state vector x(k), evaluated at w(k)=0 and
at x(k)=*x̂*(k). (The state x can have more than three, spatial coordinates;
it can e.g. have six components, three for spatial position and three for velocity.)

Based on linearizing the governing equations (12) and (13),
and using the Jacobians defined by equations (17)-(20), an EKF operates as follows.
Referring again to Fig. 2, first an initial estimate is provided for *x̂*(*k*|*k*-1)
(i.e. a value for each component of the state vector is provided for the predicted
state at the first instant of time), and also for the covariance P(k|k-1) (i.e.
a value is provided for the covariance of each component of the state vector at
the first instant of time). The prior art includes suggestions for how to select
initial values.

Next, a measurement update is performed in which the estimate
of the state vector output by the positioning system is made. This is accomplished
first by computing the so-called Kalman gain, K(k), according to:
$$\mathrm{K}\left(\mathrm{k}\right)=\mathrm{P}\left(\mathrm{k}|,\mathrm{k},-,1\right){\mathrm{H}}^{\mathrm{T}}\left(\mathrm{k}\right){\left[\mathrm{H},\left(\mathrm{k}\right),,\mathrm{P},,\left(\mathrm{k}|,\mathrm{k},-,1\right),,{\mathrm{H}}^{\mathrm{T}},\left(\mathrm{k}\right),+,\mathrm{V},\left(\mathrm{K}\right),,\mathrm{R},\left(\mathrm{k}\right),,{\mathrm{V}}^{\mathrm{T}},\left(\mathrm{k}\right)\right]}^{-1}$$
and is then used to determine the a *posteriori* state estimate
*x̂*(*k*+1), according to:
$$\widehat{x}\left(k\right)=\widehat{x}\left(k|,k,-,1\right)+\mathrm{K}\left(\mathrm{k}\right)\mathrm{\&ugr;}\left(\mathrm{k}\right),$$

where u (k) is the measurement innovation (or measurement residual), which for the
EKF is calculated using:
$$\mathrm{\&ugr;}\left(\mathrm{k}\right)=\mathrm{z}\left(\mathrm{k}\right)-\mathrm{h}\left[\widehat{x},,\left(k|,k,-,1\right),,,0\right]$$
based on defining the measurement innovation as z(k)-*z̃*(*k*),
using equation (15) for *z̃*(*k*) and asserting equation (16). The
estimate *x̂*(*k*) is then output to an output display 16 as the
Kalman filter's estimate of the state of the process as the time corresponding to
k+1, and is also used in the next first, time update part of the operation of the
Kalman filter. (The matrix V is the identity matrix in the case that the noise is
assumed to be white.) Finally, as in operation of the standard Kalman filter, the
first, measurement update part of the operation also requires a new *a posteriori*
estimate error covariance, P(k), which is again (as in the standard Kalman filter)
calculated using:
$$\mathrm{P}\left(\mathrm{k}\right)=\left[\mathrm{I},-,\mathrm{K},\left(\mathrm{k}\right),,\mathrm{H},\left(\mathrm{k}\right)\right]\mathrm{P}\left(\mathrm{k}|,\mathrm{k},-,1\right).$$

In the next, time update part, the output of the Kalman
filter (i.e. the a *posteriori* estimated state *x̂*(*k*))
is used in equation (14) and using the asserted equation (16), to determine an a
*priori* estimate *x̂*(*k*+1|*k*) of the next state of
the process:
$$\widehat{x}\left(k+1|,k\right)\approx f\left[\widehat{x},\left(k\right),,,u,\left(k\right),,,0\right].$$
The time update part also again (as in the standard Kalman filter) determines the
so-called a *priori* estimate error covariance, here using the equation:
$$\mathrm{P}\left(\mathrm{k}|,\mathrm{k},-,1\right)=\mathrm{A}\left(\mathrm{k}\right)\mathrm{P}\left(\mathrm{k}\right){\mathrm{A}}^{\mathrm{T}}\left(\mathrm{k}\right)+\mathrm{W}\left(\mathrm{k}\right)\mathrm{Q}\left(\mathrm{k}\right){\mathrm{W}}^{\mathrm{T}}\left(\mathrm{k}\right)$$

where the *a priori* estimate error covariance P(k|k-1) is based on the a
*priori* estimate error x(k)-*x̂*(*k*|*k*-1), and is defined
as:
$$\mathrm{P}\left(\mathrm{k}|,\mathrm{k},-,1\right)=\mathrm{E}\left[\left(\mathrm{x},\left(\mathrm{k}\right),-,\widehat{x},,\left(k|,k,-,1\right)\right),,{\left(\mathrm{x},\left(\mathrm{k}\right),-,\widehat{x},,\left(k|,k,-,1\right)\right)}^{\mathrm{T}}\right]$$

in which E(a) is the expectation of the random variable *a*, and where P(k)
is the corresponding a *posteriori* estimate error covariance, based on the
a *posteriori* estimate error x(k)*-x̂*(*k*). (The matrix W
is the identity matrix in the case that the noise is assumed to be white.)

*Extended* Kalman *Filter with Probabilistic* Data Association
The positioning system based on an EKF, as described above,
with inputs from cellular base stations, or with inputs from cellular base stations
and also satellites, proceeds on the basis of a single measurement z(k) at each
instant of time. (The single measurement z(k) uses information from several sources,
at least some of them cellular base stations.) The positioning system can be made
more reliable by combining several measurements z_{i}(k) at each instant
of time, at least some of which are based on information provided by cellular base
stations. Such a positioning system is here called a multi-source *positioning*
system.

Referring now to Fig. 3, in another aspect of the present
invention, a positioning system 30 uses an EKF 34 modified to use probabilistic
data association (PDA), so as to enable using a combination of measurements at each
instant of time, at least one measurement being based at least partially on information
received from a cellular base station. In the preferred embodiment of this aspect
of the present invention, up to m different measurements z_{1}(k), z_{2}(k),
..., z_{m}(k) are provided by a modified heterogeneous measurement engine
32 at each instant *k* of time, based on information received by different
antennae 21a and 21b, which are different in kind in that one is intended to receive
information from one kind of source of information, such as satellites, and the
other is intended to receive information from another kind of source of information,
such as cellular base stations. The modified heterogeneous measurement engine 32
here, in this multi-source positioning system that is based on what is here called
a *multi-source EKF,* is fundamentally the same as the heterogeneous measurement
engine 22 of the single-measurement-at-a-time EKF 24, except that it provides the
m measurements at each instant k of time.

According to the present invention, it is assumed that
the probability density function of the object's state estimate *x̂*(*k*)
(i.e. the object's position and velocity), as provided by the EKF of the present
invention, is normally distributed when conditioned on the set Z(k-1) of all previous
measurements up to the (k-1)th instant of time, i.e. that
$$\mathrm{p}\left[\widehat{x}\left(k\right)|,\mathrm{Z},,\left(\mathrm{k},-,1\right)\right]=\mathrm{N}\left[\widehat{x},\left(k\right),;,\widehat{x},,\left(k|,k,-,1\right),,,\mathrm{P},,\left(\mathrm{k}|,\mathrm{k},-,1\right)\right],$$
which corresponds to the assumptions indicated as equations (3) and (4). It is
further assumed that all measurement innovations (residuals) &ugr;_{j}(k)
are also normally distributed when conditioned on all previous measurements Z(k-1),
i.e. that
$$\mathrm{p}\left[{\mathrm{\&ugr;}}_{\mathrm{j}}\left(\mathrm{k}\right)|,\mathrm{Z},,\left(\mathrm{k},-,1\right)\right]=\mathrm{N}\left[{\mathrm{\&ugr;}}_{\mathrm{j}},\left(\mathrm{k}\right),;,0,,,\mathrm{P},,\left(\mathrm{k}|,\mathrm{k},-,1\right)\right],$$
the mean taken to be zero since the measurement noise is assumed to be white. Under
these assumptions, association probabilities (i.e. associated with measurements)
are determined according to:
$${\mathrm{\&bgr;}}_{\mathrm{j}}\left(\mathrm{k}\right)={e}_{j}\left(k\right){\left[b,\left(k\right),+,{\displaystyle \sum _{j=1}^{m\left(k\right)}{e}_{j}\left(k\right)}\right]}^{-1}$$
for j=1 ... m(k), where m(k) is the number of measurements; where
$${\mathrm{e}}_{\mathrm{j}}\left(\mathrm{k}\right)={{\mathrm{P}}_{\mathrm{G}}}^{-1}\mathrm{N}\left[{\mathrm{\&ugr;}}_{\mathrm{j}},\left(\mathrm{k}\right),;,0,,,\mathrm{S},\left(\mathrm{k}\right)\right],$$

in which &ugr;_{j}(k) is the measurement innovation for the j^{th}
measurement and S(k) is the innovation covariance at instant k, i.e. *S*(*k*)*=H*(*k*)*P*(*k*|*k-*1)*H*^{T}
(*k*)*+R*(*k*); where
$$\mathrm{b}\left(\mathrm{k}\right)=\mathrm{m}\left(\mathrm{k}\right)\left(1,-,{\mathrm{P}}_{\mathrm{D}},,{\mathrm{P}}_{\mathrm{G}}\right){\left[{\mathrm{P}}_{\mathrm{D}},,{\mathrm{P}}_{\mathrm{G}},,\mathrm{V},\left(\mathrm{k}\right)\right]}^{-1}$$

in which V(k) is the volume of a validation region, described below (equation 40);
and where a special association probability is defined corresponding to the probability
that none of the measurements is correct (or, equivalently, the probability that
at least one measurement is correct, subtracted from one), given by
$${\mathrm{\&bgr;}}_{0}\left(\mathrm{k}\right)=b\left(k\right){\left[b,\left(k\right),+,{\displaystyle \sum _{j=1}^{m\left(k\right)}{e}_{j}\left(k\right)}\right]}^{-1}.$$

Note that the Kalman gain K(k), given by equation 8, can
also be expressed in terms of the innovation covariance *S*(*k*)*=H*(*k*)*P*(*k*|*k-*1)*H*^{T}
(*k*)*+R*(*k*)*, as K*(*k*)*=P*(*k*|*k-*1)*H*^{T}
(*k*)*S*
^{-1}(*k*).

In equation (31) the quantity P_{D} is the so-called
probability of detection, and is pre-determined; it is the probability that the
correct measurement is obtained. It is specified a *priori* by a user of the
positioning system. There are many ways to evaluate P_{D}; it can be evaluated
statistically, analytically, or heuristically (trial and error). Since a measurement
will always have some error associated with a condition impeding the measurement,
such as an inability to see all of the satellites or base stations all of the time
because of blockage or equipment failures including receiver malfunctions, in practice
P_{D} will never be equal to unity.

The quantity P_{G} in equations (30) and (31) is
a so-called gate probability. It is a factor used to restrict the normal distribution
of the measurement innovation based on the width of a validation gate. If there
is no pre-validation of measurements, so that any measurement is provided to the
EKF, then P_{G} has a value of one, corresponding to a gate of infinite
width. P_{G} has a value of less than one when measurements too far from
a predicted measurement are ejected. Using such a value for P_{G} amounts
to eliminating the edges of the normal distribution of the measurements. Values
of P_{G} of less than one are found from a chi-square tail probability tables
as a function of the measurement dimension n_{z} (i.e. the number of degrees
of freedom of the measurement, which for only a range measurement is one, and for
a combined range and Doppler measurement is two) and the width *&ggr;* of
the validation gate (taking *&ggr;* as the value of the tail probability).
Using g to denote √*&ggr;*, (g therefore equal to a number of standard
deviations), g can be pre-determined to have a value limiting the accepted measurement
to only 2/3 of any set of measurements, on average. In the present invention, as
indicated in Fig. 2, it is the value of g that is pre-determined, and from which
the value of P_{G} is inferred.

The volume V(k) of the validation region is defined as:
$$\mathrm{V}\left(\mathrm{k}\right)={c}_{{n}_{z}}{g}^{{n}_{z}}{\left|S,\left(k\right)\right|}^{1/2}$$
(using g=√&ggr;) in which
$${c}_{{n}_{z}}=\frac{{\mathrm{\&pgr;}}^{{n}_{z}/2}}{\mathrm{\&Ggr;}\left({n}_{z},/,2,+,1\right)},$$

where &Ggr;(*&agr;*) is the gamma function (so that, e.g. &Ggr;(n+1)=n!
where n is an integer). *c*_{nz}
can also be interpreted to as the volume of the n_{z}-dimensional unit
hypersphere. Thus, e.g. c_{1}=2; c_{2}=&pgr; ; c_{3}=4&pgr;/3,
etc. Equation (34) is the general expression for *c*_{nz}.

With the association probabilities &bgr;_{i}(k)
determined according to equations (29) and (32), a combined innovation u (k) is
determined, according to:
$$\mathrm{\&ugr;}\left(\mathrm{k}\right)={\displaystyle \sum _{i=1}^{m\left(k\right)}{\mathrm{\&bgr;}}_{i}\left(k\right){\mathrm{\&ugr;}}_{i}\left(k\right)}$$
and used to update the state of the process, according to equation (9), i.e. as
in a standard Kalman filter.

Although the update to the state of the process is calculated
as usual, the covariance associated with the updated state, i.e. the *a posteriori*
estimate error covariance, P(k+1), is calculated differently than in a standard
Kalman filter, i.e. not according to equation (11). Instead, it is obtained using:
$$\mathrm{P}\left(\mathrm{k}\right)={\mathrm{\&bgr;}}_{0}\left(\mathrm{k}\right)\mathrm{P}\left(\mathrm{k}|,\mathrm{k},-,1\right)+\left[1,-,{\mathrm{\&bgr;}}_{0},\left(\mathrm{k}\right)\right]{\mathrm{P}}^{\mathrm{C}}\left(\mathrm{k}\right)+\tilde{P}\left(k\right)$$

where P^{c}(k) is the covariance of the state updated with the assumed to
be correct measurement (i.e. it is the covariance of the state evaluated with &bgr;_{0}(*k*)
is set to zero), and is given by
$${\mathrm{P}}^{\mathrm{C}}\left(\mathrm{k}\right)=\mathrm{P}\left(\mathrm{k}|,\mathrm{k},-,1\right)+\mathrm{K}\left(\mathrm{k}\right)\mathrm{S}\left(\mathrm{k}\right){\mathrm{K}}^{\mathrm{T}}\left(\mathrm{k}\right),$$
and where *P̃*(*k*) is the spread of innovations, given by
$$\tilde{P}\left(k\right)=K\left(k\right)\left[{\displaystyle \sum _{i=1}^{m\left(k\right)}{\mathrm{\&bgr;}}_{i}\left(k\right){\mathrm{\&ugr;}}_{i}\left(k\right){\mathrm{\&ugr;}}_{i}^{T}\left(k\right)},-,\mathrm{\&ugr;},\left(k\right),,{\mathrm{\&ugr;}}^{T},\left(k\right)\right]{K}^{T}\left(k\right).$$

Thus, the weighting of measurements in a Kalman filter,
given by equation (36) and intended to correspond to the reliability of the measurements,
is different in the PDA-modified EKF of the present invention, given by equation
(11). In fact, the weighting differs from the standard Kalman filter weighting in
that it is more responsive to changes in the reliability of the measurements. This
greater responsiveness by a PDA-modified EKF is prudent because the PDA-modified
EKF will reject measurements that are too far from what is expected, so that giving
greater credence to the remaining measurements is sensible. The modified weighting
is thus more appropriate for measurements provided by different measurement sources,
such as both satellites and cellular base stations. To significantly alter the weighting
in a standard Kalman filter, it is usually necessary to use some means external
to the EKF, or else to reinitialize the EKF. In the preferred embodiment in which
the PDA-modified EKF is implemented as software, the cost of the greater responsiveness
is a slightly greater computational load, compared to a standard EKF.

*Extended* Kalman *Filter-Based Positioning System with Pre-Validation*
Still referring to Fig. 3, in another aspect of a positioning
system according to the present invention, a pre-validation of measurements is performed
by a pre-validation module using a pre-determined (finite) value of g corresponding
to a finite gate width &ggr;, using the same notation as above. In performing
the pre-validation, the measurements are assumed to be normally distributed, having
a probability density function,
$$\mathrm{p}\left(\mathrm{z}\left(\mathrm{k},+,1\right)|,\mathrm{Z},\left(\mathrm{k}\right)\right)=\mathrm{N}\left[\mathrm{z},,\left(\mathrm{k},+,1\right),;,\widehat{z},,\left(\mathrm{k}+1|,\mathrm{k}\right),,,\mathrm{S},,\left(\mathrm{k},+,1\right)\right],$$

where again Z(k) is the set of all previous measurements, of any dimension, up to
the k^{th} instant, and where *ẑ*(*k*+1|*k*), called
the measurement prediction, is merely shorthand for the term *H*(*k+*1)*x̂*(*k+*1|*k*),
as in. e.g. equation (10). Next, a validation region is defined, based on the above-described
gate parameter gamma, that changes with each instant k of time; it is given by,
$$V\left(k,+,1,,,\mathrm{\&ggr;}\right)=\left\{z,,\left(k,+,1\right),:,{\left[z,,\left(k,+,1\right),-,\widehat{z},,\left(k+1|,k\right)\right]}^{T},,{S}^{-1},,\left(k,+,1\right),,\left[z,,\left(k,+,1\right),-,\widehat{z},,\left(k+1|,k\right)\right],\le ,\mathrm{\&ggr;}\right\}$$

where S^{-1}(k+1) is the inverse of the innovation covariance, and where
&ggr; is again, as above, a pre-defined threshold obtained from the one-sided
chi square distribution. Again denoting √&ggr; by g, equation (40) is also
sometimes written as,
$${\mathrm{\&ugr;}}^{T}\left(k,+,1\right){S}^{-1}\mathrm{\&ugr;}\left(k,+,1\right)\le {g}^{2}$$

where, as above, &ugr;(*k*+1) is the measurement innovation, i.e. z (k+1)
-*ẑ* (k+1|k).

There is a certain probability that a measurement will
fall in the validation region. The value of g is pre-determined to affect such a
probability: the greater the value of g, the greater the probability that a measurement
will fall in the validation region. If, for a particular measurement, the left-hand
side of equation (41) evaluates to an amount greater than g^{2}, the measurement
is declared an outlier, and is not used by the EKF in making its next estimate of
position (state of the process).

Thus, in the present invention, one assumes a Gaussian
probability density function for measurements, with a width determined by the value
of g, a value that can be set based on one or more of various considerations, as
discussed above.

It is important to note that some of the terms in the above
PDA and g-sigma gate equations are directly obtained from the extended EKF used
in estimating position. It is also important to note that the use of PDA weighting
does not exclude using signal to nose ratios (SNRs), bit error rates (BERs), or
other inputs for outlier rejection or measurement weighting.

*Distributed Filtering*
Although in the preferred embodiment the filtering is performed
in the positioning system itself, it is also important and in fact in some implementations
desirable to have one or more components of the filter implemented in a computing
facility that is separate and apart from the positioning system components moving
with the object being tracked, and may be in only radio communication with those
positioning system components.

Thus, referring now to Fig. 4, a distributed positioning
system according to the present invention is shown as including local components
41 moving with the object being tracked and also including non-local components
34, the filtering hardware modules, hosted by a remote facility. The non-local components
34 perform the filtering computations of the present invention, and are in only
radio communication with the local modules 41 of the positioning system, the heterogeneous
measurement engine 22 and the output display 16. The radio communication is indicated
by the dashed lines 42 representing radio communication signals by which the outputs
of the heterogeneous measurement system are provided to the filter modules 34 (i.e.
the non-local components), and by the dashed lines 43 representing radio communication
signals by which the outputs of the filter modules 34 are provided to the output
display 16. (The antennae and associated receiver and transmitter hardware needed
to accomplish the indicated radio communications are implied.)

In such an embodiment, the non-local components 34 perform
the filtering computations for not only a particular positioning system with local
components 41, but also for other positioning systems, not shown, also having only
some components on-board the object making use of the positioning system. In some
applications, the non-local components 34 would be part of a telecommunications
network that would communicate with each different positioning system using a protocol
similar to what is used in cellular telephone communications.

Referring now to Fig. 5, a flowchart of a method, according
to the present invention, of filtering measurements in a generalized positioning
system is shown.

It is to be understood that the above-described arrangements
are only illustrative of the application of the principles of the present invention.
In particular, as indicated above, the present invention comprehends not only filters
that might be referred to as extended Kalman filters, but any kind of filter in
a position-measuring system, e.g. any kind of filter that accepts measurements of
state information for a process that is governed by a non-linear stochastic difference
equation (equation 12). Further, a positioning system according to the present invention
will not, of course, always yield an estimate of position, since in some situations
there will not be adequate measurement data due, for example, to environmental affects
that temporarily prevent receiving data from the number of sources needed to estimate
position. In such a case, the positioning system might temporarily estimate only
a heading and a clock bias.

Thus, it should also be understood that the present invention
applies equally well to systems that are not, strictly, positioning systems, but
instead are systems that in general estimate one or another aspect of the motion
of an object, position being only one aspect. In case of estimating position, the
system may aptly be called a "positioning system," even though it may estimate values
for other aspects of the motion. But in other cases the aspect of the motion that
is most important may be other than position; it may, for example, be angular velocity,
or linear or angular acceleration, or simply linear velocity. In general, therefore,
the present invention is for use as part of a system that is most aptly described
as a generalized positioning system, and is to be understood to encompass estimating,
in some applications, position, but in other applications, other aspects of motion,
not necessarily including position.

In addition, numerous modifications and alternative arrangements
of the embodiments disclosed here may be devised by those skilled in the art without
departing from the spirit and scope of the present invention, and the appended claims
are intended to cover such modifications and arrangements.