 |
Description  |
|
|
BACKGROUND OF THE INVENTION
The invention described in this patent application relates to the problem
of estimation of constant parameters of multiple signals received by an
array of sensors in the presence of additive noise. There are many
physical problems of this type including direction finding (DF) wherein
the signal parameters of interest are the directions-of-arrival (DOA's) of
wavefronts impinging on an antenna array (cf. FIG. 1), and harmonic
analysis in which the parameters of interest are the temporal frequencies
of sinusoids contained in a signal (waveform) which is known to be
composed of a sum of multiple sinusoids and possibly additive measurement
noise. In most situations, the signals are characterized by several
unknown parameters all of which need to be estimated simultaneously (e.g.,
azimuthal angle, elevation angle and temporal frequency) and this leads to
a multidimensional parameter estimation problem.
High resolution parameter estimation is important in many applications
including electromagnetic and acoustic sensor systems (e.g., radar, sonar,
electronic surveillance systems, and radio astronomy), vibration analysis,
medical imaging, geophysics, well-logging, etc. In such applications,
accurate estimates of the parameters of interest are required with a
minimum of computation and storage requirements. The value of any
technique for obtaining parameter estimates is strongly dependent upon the
accuracy of the estimates. The invention described herein yields accurate
estimates while overcoming the practical difficulties encountered by
present methods such as the need for detailed a priori knowledge of the
sensor array geometry and element characteristics. The technique also
yields a dramatic decrease in the computational and storage requirements.
The history of estimation of signal parameters can be traced back at least
two centuries to Gaspard Riche, Baron de Prony, (R. Prony, Essai
experimental et analytic, etc.L'Ecole Polytechnique, 1:24-76, 1795) who
was interested in fitting multiple sinusoids (exponentials) to data.
Interest in the problem increased rapidly after World War II due to its
applications to the fast emerging technologies of radar, sonar and
seismology. Over the years, numerous papers and books addressing this
subject have been published, especially in the context of direction
finding in passive sensor arrays.
One of the earliest approaches to the problem of direction finding is now
commonly referred to as the conventional beamforming technique. It uses a
type of matched filtering to generate spectral plots whose peaks provide
the parameter estimates. In the presence of multiple sources, conventional
beamforming can lead to signal suppression, poor resolution, and biased
parameter (DOA) estimates.
The first high resolution method to improve upon conventional beamforming
was presented by Burg (J. P. Burg, Maximum entropy spectral analysis, In
Proceedings of the 37th Annual International SEG Meeting, Oklahoma City,
Okla., 1967). He proposed to extrapolate the array covariance function
beyond the few measured lags, selecting that extrapolation for which the
entropy of the signal is maximized. The Burg technique gives good
resolution but suffers from parameter bias and the phenomenon referred to
as line splitting wherein a single source manifests itself as a pair of
closely spaced peaks in the spectrum. These problems are attributable to
the mismodeling inherent in this method.
A different approach aimed at providing increased parameter resolution was
introduced by Capon (J. Capon, High resolution frequency wave number
spectrum analysis, Proc. IEEE, 57:1408-1418, 1969). His approach was to
find a weight vector for combining the outputs of all the sensor elements
that minimizes output power for each look direction while maintaining a
unit response to signals arriving from this direction. Capon's method has
difficulty in multipath environments and offers only limited improvements
in resolution.
A new genre of methods were introduced by Pisarenko (V. F. Pisarenko, The
retrieval of harmonics from a covariance function, Geophys. J. Royal
Astronomical Soc., 33:347-366, 1973) for a somewhat restricted formulation
of the problem. These methods exploit the eigenstructure of the array
covariance matrix. Schmidt made important generalizations of Pisarenko's
ideas to arbitrary array/wavefront geometries and source correlations in
his Ph.D. thesis titled A Signal Subspace Approach to Multiple Emitter
Location and Spectral Estimation, Standford University, 1981. Schmidt's
MUltiple SIgnal Classification (MUSIC) algorithm correctly modeled the
underlying problem and therefore generated superior estimates. In the
ideal situation where measurement noise is absent (or equivalently when an
infinite amount of measurements are available), MUSIC yields exact
estimates of the parameters and also offers infinite resolution in that
multiple signals can be resolved regardless of the proximity of the signal
parameters. In the presence of noise and where only a finite number of
measurements are available, MUSIC estimates are very nearly unbiased and
efficient, and can resolve closely spaced signal parameters.
The MUSIC algorithm, often referred to as the eigenstructure approach, is
currently the most promising high resolution parameter estimation method.
However, MUSIC and the earlier methods of Burg and Capon which are
applicable to arbitrary sensor array configurations suffer from certain
shortcomings that have restricted their applicability in several problems.
Some of these are:
Array Geometry and Calibration--A complete characterization of the array in
terms of the sensor geometry and element characteristics is required. In
practice, for complex arrays, this characterization is obtained by a
series of experiments known as array calibration to determine the so
called array manifold. The cost of array calibration can be quite high and
the procedure is sometimes impractical. Also, the associated storage
required for the array manifold is 2 ml.sup.g words (m is the number of
sensors, l is the number of search (grid) points in each dimension, and g
is the number of dimensions) and can become large even for simple
applications. For example, a sensor array containing 20 elements,
searching over a hemisphere with a 1 millirad resolution in azimuth and
elevation and using 16 bit words (2 bytes each) requires approximately 100
megabytes of storage! This number increases exponentially as another
search dimension such as temporal frequency is included. Furthermore, in
certain applications the array geometry may be slowly changing such as in
light weight spaceborne antenna structures, sonobuoy and towed arrays used
in sonar etc., and a complete characterization of the array is never
available.
Computational Load--In the prior methods of Burg, Capon, Schmidt and
others, the main computational burden lies in generating a spectral plot
whose peaks correspond to the parameter estimates. For example, the number
of operations required in the MUSIC algorithm in order to compute the
entire spectrum, is approximately 4 m.sup.2 l.sup.g. An operation is
herein considered to be a floating point multiplication and an addition.
In the example above, the number of operations needed is approximately
4.times.10.sup.9 which is prohibitive for most applications. A powerful 10
MIP (10 million floating point instructions per second) machine requires
about 7 minutes to perform these computations! Moreover, the computation
requirement grows exponentially with dimension of the parameter vector.
Augmenting the dimension of the parameter vector further would make such
problems completely intractable.
The technique described herein is hereafter referred to as Estimation of
Signal Parameters using Rotational Invariance Techniques (ESPRIT). ESPRIT
obviates the need for array calibration and dramatically reduces the
computational requirements of previous approaches. Furthermore, since the
array manifold is not required, the storage requirements are eliminated
altogether.
SUMMARY OF THE INVENTION
ESPRIT is an alternative method for signal reception and source parameter
estimation which possesses most of the desirable features of prior high
resolution techniques while realizing substantial reduction in computation
and elimination of storage requirements. The basic properties of the
invention may be summarized as follows:
1. ESPRIT details a new method of signal reception for source parameter
estimation for planar wavefronts.
2. The method yields signal parameter estimates without requiring knowledge
of the array geometry and sensor element characteristics, thus eliminating
the need for sensor array calibration.
3. ESPRIT provides substantial reduction in computation and elimination of
storage requirements over prior techniques. Referring to the previous
example, ESPRIT requires only 4.times.10.sup.4 computations compared to
4.times.10.sup.9 computations required by prior methods, and reduces the
time required from 7 minutes to under 4 milliseconds. Furthermore, the 100
megabytes of storage required is completely eliminated.
4. A feature of the invention is the use of an array of sensor pairs or
doublets (used synonymously herein) where the sensors in each pair are
identical and each group of pairs has a common displacement vector.
Briefly, in accordance with the invention, an array of signal sensor pairs
is provided in which groups of sensor pairs have a uniform relative
displacement vector within each group, but the displacement vector for
each group is unique. The sensors in each pair must be matched, however
they can differ from other sensor pairs. Moreover, the characteristics of
each sensor and the array geometry can be arbitrary and need not be known.
Within each group, the sensor pairs can be arranged into two subarrays, X
and Y, which are identical except for a fixed displacement (cf. FIG. 2).
For example, in order to simultaneously perform temporal frequency and
spatial angle estimation, one group of sensor pairs would share a common
spatial displacement vector while the second group would share a common
temporal displacement. In general, for each additional type of parameter
to be estimated, a sensor group sharing a common displacement is provided.
Furthermore, the number of sensor pairs in each group must be more than
the number of sources whose parameters are to be estimated.
Having provided an array of sensors which meets the specifications outlined
above, signals from this array of sensor pairs are then processed in order
to obtain the parameter estimates of interest. The procedure for obtaining
the parameter estimates in accordance with one embodiment employing
standard least-squares estimation techniques may be outlined as follows:
1. Using the array measurements from a group of sensor pairs, determine the
auto-covariance matrix R.sub.xx of the X subarray in the group and the
cross-covariance matrix R.sub.xy between the X and Y subarrays in the
group.
2. Determine the smallest eigenvalue of the covariance matrix R.sub.xx and
then subtract it out from each of the elements on the principal diagonal
of R.sub.xx. The results of the subtraction are referred to hereinafter as
C.sub.xx.
3. Next, the generalized eigenvalues (GE's) .gamma.i of the matrix pair
C.sub.zz, R.sub.xy are determined. A number d of the GE's will lie on or
near the unit circle and the remaining m-d noise GE's will lie at or near
the origin. The number of GE's on or near the unit circle determine the
number of sources, and their angles are the phase differences sensed by
the sensor doublets in the group for each of the wavefronts impinging on
the array. These phase differences are directly related to the directions
of arrival.
4. The procedure is then repeated for each of the groups, thereby obtaining
the estimates for all the parameters of interest (e.g. azimuth, elevation,
temporal frequency).
Thus, the number of sourses and the parameters of each source are the
primary quantities determined.
In another embodiment of the invention, the processing of signal
measurements from the two subarrays to identify the number of sources and
estimate parameters thereof utilizes a total least-squares estimation
technique. The total least-squares algorithm represents an improvement and
simplification of the least squares algorithm.
ESPIRIT can be further extended to the problem of determining the array
geometry a posteriori, i.e., obtaining estimates of the sensor locations
given the measurements. Source powers and optimum weight vectors for
solving the signal copy problem, a problem involving estimation of the
signals received from the sources one at a time eliminating all others,
can also be estimated in a straightforward manner as follows:
1. The optimum weight vector for signal copy for the i.sup.th signal is the
generalized eigenvector (GV) e.sub.i corresponding to the i.sup.th GE
.gamma.i.
2. For the case when the sources are uncorrelated, the direction vector
a.sub.i for the i.sup.th wavefront is given by R.sub.xy e.sub.i. With
these direction vectors in hand, the array geometry can be estimated by
solving a set of linear equations.
3. Using the direction vectors a.sub.i, the signal powers can also be
estimated by solving a set of linear equations.
The invention and objects and features thereof will be more readily
apparent from the following example and appended claims.
BRIEF DESCRIPTION OF THE DRAWINGS
FIG. 1 is a graphic representation of a problem of direction-of-arrival
estimation in which two sources are present and being monitored by a
three-element array of sensors.
FIG. 2 is a graphic representation of a similar problem in which the two
signals are now impinging on an array of sensors pairs in accordance with
the invention.
FIG. 3 is a graphic illustration of the parameter estimates from a
simulation performed in accordance with the invention in which three
signals were impinging on an array of eight sensor doublets and
directions-of-arrival were being estimated.
DETAILED DESCRIPTION OF THE DRAWINGS
As indicated above, the invention is directed at the estimation of constant
parameters of signals received by an array of sensor pairs in the presence
of noise. The problem can be visualized with reference to FIG. 1 in which
two signals (s.sub.1 and s.sub.2) are impinging on an array of three
sensors (r.sub.1, r.sub.2, r.sub.3). It is assumed in this illustrated
example that the sources and sensors lie in a plane; thus only two
parameters need be identified, the azimuth angle of the two signals.
Heretofore, techniques such as MUSIC have been able to accurately estimate
the DOA's of the two signals; however the characteristics of each sensor
must be known as well as the overall array geometry. This leads to
exceedingly large storage requirements when the array must be calibrated,
and a correspondingly large computation time in the execution of the
algorithms.
In accordance with the present invention, array (manifold) calibration is
not required in ESPRIT as long as the array is comprised of (groups of)
matched sensor pairs sharing a common displacement vector. This is
illustrated in FIG. 2 in which the two signals (s.sub.1 and s.sub.2) are
sensed by receiver pairs (r.sub.1, r'.sub.1 ; r.sub.2, r'.sub.2 ; and
r.sub.3, r'.sub.3). The only requirements of the array are that the
sensors in each pair are offset by the same vector as indicated, and that
the number of sensor pairs exceeds the number of sources as is the case in
this example.
This figure illustrates only a single group; the extension to several
groups requires adding sensor pairs with a displacement vector different
from the displacement vectors of the single group.
The performance of the invention is graphically illustrated in FIG. 3 which
presents the results of a simulation performed according to the
specifications of ESPRIT. The simulation consisted of an array with 8
doublets. The elements in each of the doublets were spaced a quarter of a
wavelength apart. The array geometry was generated by randomly scattering
the doublets on a line 10 wavelengths in length such that the doublet axes
were all parallel to the line. Three planar and weakly correlated signal
wavefronts impinged on the array at angles 20.degree., 22.degree., and
60.degree., with SNRs of 10, 13 and 16 db relative to the additive
uncorrelated noise present at the sensors. The covariance estimates were
computed from 100 snapshots of data and several simulation runs were made
using independent data sets.
FIG. 3 shows a plot of the GE's obtained from 10 independent trials. The
three small circles on the unit circle indicate the locations of the true
parameters and the pluses are the estimates obtained using ESPRIT. The
GE's on the unit circle are closely clustered and the two sources
2.degree. apart are easily resolved.
As illustrated, accurate estimates of the DOA's are obtained. Furthermore,
ESPRIT has several additional features which are enumerated below.
1. ESPRIT appears to be very robust to errors in estimating the minimum
eigenvalue of the covariance R.sub.xx. It is also robust to the numerical
properties of the algorithm used to estimate the generalized eigenvalues.
2. ESPRIT does not require the estimation of the number of sources prior to
source parameter estimation as in the MUSIC algorithm, where an error in
the estimate of the number of sources can invalidate the parameter
estimates. In accordance with the invention, ESPRIT simultaneously
estimates the signal parameters and the number of sources.
APPLICATIONS
There are a number of applications that exploit one or more of the
important features of ESPRIT, i.e., its insensitivity to array geometry,
low computational load and no storage requirements. Some of these are
described below.
1. Direction-of-Arrival Estimation
(a) Space Antennas--Space structures are necessarily light weight, very
large and therefore fairly flexible. Small disturbances can cause the
structure to oscillate for long periods of time resulting in a sensor
array geometry which is time-varying. Furthermore, it is nearly impossible
to completely calibrate such an array as the setting up of a suitable
facility is not practical. On the other hand, the use of matched pairs of
sensor doublets whose directions are constantly aligned by a low-cost
star-tracking servo results in total insensitivity to the global geometry
of the array. Note that signal copy can still be performed, a function
which is often a main objective of such large spaceborne antenna arrays.
In fact, a connected structure for the array is not required! Rather, only
a collection of relatively small antenna doublets is needed, each
possessing a star-tracker or earth-based beacon tracker for alignment.
Ease of deployment, maintenance, and repair of such disconnected arrays
can have significant cost and operational benefits (for example, a
defective unit can be merely transported to a space station or back to the
earth for repair).
(b) Sonobuoys--Sonobuoys are air-dropped and scatter somewhat randomly on
the ocean surface. The current methods of source location require complete
knowledge of the three dimensional geometry of the deployed array. The
determination of the array geometry is both expensive and undesirable
(since it involves active transmission thus alerting unfriendly
elements!). Using ESPRIT, vertical alignment of doublets can be achieved
using gravity as a reference. Horizontal alignment can be obtained via a
small servo and a miniature magnetic sensor (or even use an acoustic
spectral line radiated from a beacon or the target itself). Within a few
minutes after the sonobuoys are dropped, alignment can be completed and
accurate estimates of DOA's become available. As before, signal copy
processing is also feasible. Furthermore, the sonobuoy array geometry can
itself be determined should this be of interest.
(c) Towed Arrays--These consist of a set of hydrophones placed inside an
acoustically transparent tube that is towed well behind a ship or
submarine. The common problem with towed arrays is that the tube often
distorts from the assumed straight line geometry due to ocean and tow-ship
induced disturbances. Therefore, prior array calibration becomes invalid.
In the new approach, any translational disturbance in the doublets is of
no consequence. Therefore by selective use of doublets (whose orientation
can be easily sensed) that are acceptably co-directional, reliable source
DOA estimates can still be obtained.
(d) Mobile DF and Signal Copy Applications--Often, mobile (aircraft, van
mounted) direction finding (DF) systems cannot meet the vast storage and
computational requirements of the prior methods. ESPRIT can drastically
reduce such requirements and still provide good performance. This has
particular applicability in the field of cellular mobile communications
where the number of simultaneous users is limited due to finite bandwidth
constraints and cross-talk (interchannel interference). Current techniques
for increasing the number of simultaneous users exploit methods of signal
separation such as frequency, time and code division multiplexing apart
from the area multiplexing inherent to the cellular concept. Using
directional discrimination (angle division multiplexing), the number of
simultaneous users could be increased significantly. ESPRIT provides a
simple and relatively low cost technique for performing the signal copy
operation through angular signal separation. The estimation (possibly
recursively) of the appropriate generalized eigenvector is all that is
needed in contrast to substantially more complex procedures required by
prior methods.
2. Temporal Frequency Estimation--There are many applications in radio
astronomy, modal identification of linear systems including structural
analysis, geophysics sonar, electronic surveillance systems, analytical
chemistry etc., where a composite signal containing multiple harmonics is
present in additive noise. ESPRIT provides frequency estimates from
suitably sampled time series at a substantially reduced level of
computation over the previous methods.
3. Joint DOA-Frequency Estimation--Applications such as radio astonomy may
require the estimation of declination and right ascension of radio sources
along with the frequency of the molecular spectral lines emitted by them.
Such problems also arise in passive sonar and electronic surveillance
applications. As previously noted, ESPRIT has particularly important
advantages in such multi-dimensional estimation problems.
Having concluded the summary of the invention and applications, a detailed
mathematical description of the invention is presented.
PROBLEM FORMULATION
The basic problem under consideration is that of estimation of parameters
of finite dimensional signal processes given measurements from an array of
sensors. This general problem appears in many different fields including
radio astronomy, geophysics, sonar signal processing, electronic
surveillance, structural (vibration) analysis, temporal frequency
estimation, etc. In order to simplify the description of the basic ideas
behind ESPRIT, the ensuing discussion is couched in terms of the problem
of multiple source direction-of-arrival (DOA) estimation from data
collected by an array of sensors. Though easily generalized to higher
dimensional parameter spaces, the discussion and results presented deal
only with single dimensional parameter spaces, i.e., azimuth only
direction finding (DF) of far-field point sources. Furthermore, narrowband
signals of known center frequency will be assumed. A DOA/DF problem is
classified as narrowband if the sensor array width is small compared to
the inverse of the transit time of a wavefront across the array. The
generality of the fundamental concepts on which ESPRIT is based makes the
extension to signals containing multiple frequencies straightforward as
discussed later. Note that wideband signals can also be handled by
decomposing them into narrowband signal sets using comb filters.
Consider a planar array of arbitrary geometry composed of m matched sensor
doublets whose elements are translationally separated by a known constant
displacement vector as shown in FIG. 2. The element characteristics such
as element gain and phase pattern, polarization sensitivity, etc., may be
arbitrary for each doublet as long as the elements are pairwise identical.
Assume there are d<m narrowband stationary zero-mean sources centered at
frequency .omega..sub.0, and located sufficiently far from the array such
that in homogenous isotropic transmission media, the wavefronts impinging
on the array are planar. Additive noise is present at all the 2 m sensors
and is assumed to be a stationary zero-mean random process that is
uncorrelated from sensor to sensor.
In order to exploit the translational invariance property of the sensor
array, it is convenient to describe the array as being comprised of two
subarrays, X and Y, identical in every respect although physically
displaced (not rotated) from each other by a known displacement vector.
The signals received at the i.sup.th doublet can then be expressed as:
##EQU1##
where s.sub.k (.) is the k.sup.th signal (wavefront) as received at sensor
1 (the reference sensor) of the X subarray, .theta..sub.k is the direction
of arrival of the k.sup.th source relative to the direction of the
translational displacement vector, .alpha..sub.i (.theta..sub.k) is the
response of the i.sup.th sensor of either subarray relative to its
response at sensor 1 of the same subarray when a single wavefront impinges
at an angle .theta..sub.k, .DELTA. is the magnitude of the displacement
vector between the two arrays, c is the speed of propagation in the
transmission medium, n.sub.x.sbsb.i (.) and n.sub.y.sbsb.i (.) are the
additive noises at the elements in the i.sup.th doublet for subarrays X
and Y respectively.
Combining the outputs of each of the sensors in the two subarrays, the
received data vectors can be written as follows:
##EQU2##
The vector s(t) is a d.times.1 vector of impinging signals (wavefronts) as
observed at the reference sensor of subarray X. The matrix .PHI. is a
diagonal d.times.d matrix of the phase delays between the doublet sensors
for the d wavefronts, and can be written as:
.PHI.=diag[e.sup.j.omega..sbsp.0.sup..DELTA. sin .theta..sbsp.1.sup./e, . .
. , e.sup.j.omega..sbsp.0.sup..DELTA. sin .theta..sbsp.d.sup./e ].(4)
Note that .PHI. is a unitary matrix (operator) that relates the
measurements from subarray X to those from subarray Y. In the complex
field, .PHI. is a simple scaling operator. However, it is isomorphic to
the real two-dimensional rotation operator and is herein referred to as a
rotation operator. The m.times.d matrix A is the direction matrix whose
columns {a(.theta..sub.k), k=1, . . . , d} are the signal direction
vectors for the d wavefronts.
a.sup.T (.theta..sub.k)=[.alpha..sub.1 (.theta..sub.k), . . . ,
.alpha..sub.m (.theta..sub.k)]. (5).
The auto-covariance of the data received by subarray X is given by:
R.sub.xx =E[x(t)x*(t)]=ASA*+.sigma..sup.2 I, (6)
where S is the d.times.d covariance matrix of the signals s(t), i.e.,
S=E[s(t)s(t)*], (7)
and .sigma..sup.2 is the covariance of the additive uncorrelated white
noise that is present at all sensors. Note that (.)* is used herein to
denote the Hermitean conjugate, or complex conjugate transpose operation.
Similarly, the cross-covariance between measurements from subarrays X and
Y is given by:
R.sub.xy =E[x(t)y(t)*]=AS.PHI.*A*. (8).
This completes the definition of the signal and noise model, and the
problem can now be stated as follows:
Given measurements x(t) and y(t), and making no assumptions about the
array geometry, element characteristics, DOA's, noise powers, or the
signal (wavefront) correlation, estimate the signal DOA's.
ROTATIONALLY INVARIANT SUBSPACE APPROACH
The basic idea behind the new technique is to exploit the rotational
invariance of the underlying signal subspaces induced by the translational
invariance of the sensor array. The following theorem provides the
foundation for the results presented herein.
Theorem: Define T as the generalized eigenvalue matrix associated with the
matrix pencil {(R.sub.xx -.lambda..sub.min I), R.sub.xy } where
.lambda..sub.min is the minimum (repeated) eigenvalue of R.sub.xx. Then,
if S is nonsingular, the matrices .PHI. and T are related by
##EQU3##
to within a permutation of the elements of .PHI..
Proof: First it is shown that ASA* is rank d and R.sub.xx has a
multiplicity (m-d) of eigenvalues all equal to .sigma..sup.2. From linear
algebra,
.rho.(ASA*)=min(.rho.(A),.rho.(S)) (10)
where .rho.(.) denotes the rank of the matrix argument. Assuming that the
array geometry is such that there are no ambiguities (at least over the
angular interval where signals are expected), the columns of the m.times.d
matrix A are linearly independent and hence .rho.(A)=d. Also, since S is a
d.times.d matrix and is nonsingular, .rho.(S)=d. Therefore, .rho.(ASA*)=d,
and consequently ASA* will have m-d zero eigenvalues. Equivalently
ASA*+.sigma..sup.2 I will have m-d minimum eigenvalues all equal to
.sigma..sup.2. If {.lambda..sub.1 >.lambda..sub.2 >. . . >.lambda..sub.m }
are the ordered eigenvalues of R.sub.xx, then
.lambda..sub.d+1 =. . . =.lambda..sub.m =.sigma..sup.2. (11)
Hence,
R.sub.xx -.lambda..sub.min I=R.sub.xx -.sigma..sup.2 I=ASA*.(12)
Now consider the matrix pencil
C.sub.xx -.gamma.R.sub.xy
=ASA*-.gamma.AS.PHI.*A*=AS(I-.gamma..PHI.*)A*;(13)
where C.sub.xx .apprxeq.R.sub.xx -.lambda..sub.min.sup.xx I. By inspection,
the column space of both ASA* and AS.PHI.*A* are identical. Therefore,
.rho.(ASA*-.gamma.AS.PHI.*A*) will in general be equal to d. However, if
.gamma.=e.sup.j.omega..sbsp.0.sup..DELTA. sin .theta..sbsp.i.sup./e,(14)
the i.sup.th row of (I-e.sup.j.omega..sbsp.0.sup..DELTA. sin
.theta..sbsp.i.sup./e .PHI.) will become zero. Thus,
.rho.(I-e.sup.j.omega..sbsp.0.sup..DELTA. sin .theta..sbsp.i.sup./e
.PHI.)=d-1. (15)
Consequently, the pencil (C.sub.xx -.gamma.R.sub.xy) will also decrease in
rank to d-1 whenever .gamma. assumes values given by (14). However, by
definition these are exactly the generalized eigenvalues (GEV's) of the
matrix pair {C.sub.xx,R.sub.xy }. Also, since both matrices in the pair
span the same subspace, the GEV's corresponding to the common null space
of the two matrices will be zero, i.e., d GEV's lie on the unit circle and
are equal to the diagonal elements of the rotation matrix .PHI., and the
remaining m-d (equal to the dimension of the common null space) GEV's are
at the origin. This completes the proof of the theorem.
Once .PHI. is known, the DOA's can be calculated from:
.theta..sub.k =arc sin {c.PHI..sub.k k/.omega..sub.0 .DELTA.}.(16)
Due to errors in estimating R.sub.xx and R.sub.xy from finite data as well
as errors introduced during the subsequent finite precision computations,
the relations in (9) and (11) will not be exactly satisfied. At this
point, a procedure is proposed which is not globally optimal, but utilizes
some well established, stepwise-optimal techniques to deal with such
issues.
SUBSPACE ROTATION ALGORITHM (ESPRIT)
The key steps of the algorithm are:
1. Find the auto- and cross-covariance matrix estimates R.sub.xx and
R.sub.xy from the data.
2. Compute the eigen-decomposition of R.sub.xx and R.sub.xy and then
estimate the number of sources d and the noise variance .sigma..sup.2.
3. Compute rank d approximations to ASA* and AS.PHI.*A* given
.sigma..sup.2.
4. The d GEV's of the estimates of ASA* and AS.PHI.*A* that lie close to
the unit circle determine the subspace rotation operator .PHI. and hence,
the DOA's.
Details of the algorithm are now discussed.
Covariance Estimation
In order to estimate the required covariances, observations x(t.sub.j) and
y(t.sub.j) at time instants t.sub.j are required. Note that the subarrays
must be sampled simultaneously. The maximum likelihood estimates (assuming
no underlying data model) of the auto- and cross-covariance matrices are
then given by
##EQU4##
The number of snapshots, N, needed for an adequate estimate of the
covariance matrices depends upon the signal-to-noise ratio at the array
input and the desired accuracy of the DOA estimates. In the absense of
noise, N>d i | | |
