 |
Description  |
|
|
BACKGROUND OF THE INVENTION
The present invention relates to a filtering system and more particularly
to an equalization system and method for utilizing adaptive digital
filters with non-linear frequency resolution.
Quality audio products are designed with the goal of reproducing as
accurately as possible at the listener's ears the acoustic signal
originally recorded or broadcast. Yet despite the many improvements that
have been made in audio technology in the past several years, there still
remains at least one major obstacle to achieving that goal.
Room boundaries can have a significant effect on the sound radiated by a
loudspeaker and eventually perceived by the listener. Reflections off
walls and furniture combine at the listener's ears in a complex manner
such that the various frequency components in the music are unbalanced,
influencing the sound to a greater extent than any other component in the
system. The problem is difficult to deal with because the extent of the
problem can only be assessed by making a measurement of the system at the
exact listening position.
A solution to the problem is to use adaptive digital filters to develop
inverse filters for the loudspeaker/room response. The invertability of
these system responses has been studied and numerous solutions have been
proposed. Several of these solutions involve frequency domain transform
techniques to design finite impulse response ("FIR") filters. Various
configurations of time domain adaptive FIR filters have also been
developed. Although adaptive infinite impulse response ("IIR") filters are
also applicable, problems of convergence and stability generally make them
less practical.
In most systems using inverse filters for equalization, one or more
loudspeakers are located in a small to medium sized room such as a studio
control room or domestic living room. During a separate calibration mode,
a test signal is output through a loudspeaker and received at a
microphone, located near, but by necessity not coincident to, the location
of the listener's ears. Inverse filter coefficients are then generated
from the measured transfer function. These coefficients are then
transferred to a fixed digital filter for use in the playback mode, at
which time the system processes the audio source material in real time.
Frequency response anomalies in a room are the results of the reinforcement
and cancellation which occur when sound waves from various sources (i.e.,
direct and reflected) add together in and out of phase. It has been found
that the average distance between pressure maxima in a room is about 0.9
times the wavelength. It follows that the level of high frequency sounds,
with short wavelengths, will vary significantly between nearby points,
while that of the longer wavelength low frequencies will be less position
sensitive.
This issue is generally ignored in applying equalization, but becomes
extremely important as the resolution and accuracy of the inverse
correction filter improves. Since a listener must use two ears, generally
separated by about 20 cm, it is not physically possible to provide perfect
equalization across the full audio band at both ears simultaneously. This
clearly suggests that while low frequencies may be corrected quite
accurately, higher frequencies must be treated differently.
To further complicate the issue, the auditory system has the ability to
discriminate direct sound from later reflections, as well as the ability
to detect the direction from which a sound is coming. It also perceives
tones on a logarithmic frequency scale, rather than the linear range in
which adaptive filters operate. While the ability to generate very
accurate equalization filters is one of the goals of the known adaptive
systems, i.e., minimize an error in the least mean squared sense, it is
not necessarily correct from psychoacoustic criteria.
From the above discussions, it is clear that improvements over the current
state of the art require a means for more effectively controlling filter
accuracy as a function of frequency and space. Existing techniques
approach the problem by either controlling the resolution of the filter
directly, as in U.S. Pat. No. 4,628,530 to Philips, or with a multi-band
approach, using high resolution digital filters for low frequencies, and
low resolution analog or digital filters at high frequencies. Current
means of implementing this multi-band approach require that the signal to
be equalized be split into two or more frequency bands and operated on by
parallel filters. This has several disadvantages. First, since there is no
interaction between the filters in the various bands, it makes the
adaptation of the filters difficult. Secondly, the additional processing
steps of band splitting and recombination distort the signal and introduce
noise into the system.
Referring to FIG. 1, the basic structure of a known adaptive FIR filter is
shown. The output of the filter y(n) is the linear combination of current
and delayed signal values x(n-i) scaled by the filter coefficients ai,
where 0.ltoreq.i.ltoreq.N-1, with N being the number of signal values,
i.e.
y(n)=a.sub.0 x(n)+a.sub.1 x(n-1)+. . . +a.sub.N-1 x(n-N+1)
The filter coefficients a0 to aN-1 are updated based on an error signal
e(n), which is the difference between the filter's output y(n) and a
reference signal d(n). Any known method may be used for performing this
update including those described in "Adaptive Signal Processing", edited
by L. H. Sibul, 1987 IEEE Press, New York. Such known methods typically
attempt to minimize some function of the error signal e(n). The
coefficient update equation for the LMS algorithm is (with K being the
convergence factor):
a.sub.i (n+1)=a.sub.i (n)+K e(n) x(n-i)
By choosing an appropriate input and reference signal, this technique can
be used to adaptively design digital filters with responses matched to the
given signals.
FIG. 2 illustrates a predictive filter structure used for equalization,
which is described in U.S. Pat. No. 4,458,362 issued to Berkovitz et al.
Here the filter A(n) is updated based on the error between the current
input signal sample and its predicted value y(n). If x(n) is the output of
a system driven by a "white" sequence, it can be shown that the resulting
filter is an inverse of the system response, with properties suitable for
use as an equalizer.
A digital filter's frequency resolution is directly proportional to its
length. If we define resolution to be the minimum 3 dB feature bandwidth,
and assume that an adapted FIR filter represents a rectangularly windowed
version of an optimal Infinite Impulse Response (IIR) filter, then
f.sub.res =0.89.times.f.sub.s /N
where fres is the resolution in Hertz, f.sub.s is the sampling frequency in
Hertz, and N is the total number of FIR filter coefficients.
For high quality audio applications, a common sampling frequency is 44,100
H.sub.z. For the filter to have better than 20 H.sub.z resolution, which
would be needed for satisfactory equalization at low frequencies, N must
be greater than 1960. However, in typical audio applications, this
resolution is required only at very low frequencies, which is a small
fraction of the total signal bandwidth f.sub.s /2. For most adaptive
algorithms such a long filter introduces computational difficulties. In
the case of algorithms such as the LMS, the problem is one of convergence
and "misadjustment" due to algorithm noise. This leads to a significant
disparity between the adaptive filter and the theoretical optimal filter.
It is therefore a principal object of the present invention to provide
inverse filters whose accuracy can be controlled as a function of
frequency.
An additional object of this invention is to provide a method for
accurately adapting an FIR filter when a large number of filter taps are
required to obtain adequate frequency resolution.
Yet another object of this invention is to provide an efficient means for
implementing long FIR filters which do not introduce amplitude or phase
distortions into the response by band splitting and recombination.
SUMMARY OF THE INVENTION
The system and method of the present invention provides a means for
designing a single fixed FIR filter adaptively from measured data, in a
manner whereby the filter's frequency and time resolution can be
controlled. The resulting filter exhibits properties which allow it to be
efficiently implemented in various multi-rate configurations.
The system and method of the present invention exploits several properties
of FIR filters to solve the problems associated with prior art systems.
The underlying principle is the duality which exists between the time and
frequency responses, i.e., narrow time domain events have broad frequency
domain components, while narrow frequency domain features are
correspondingly long time domain events. Specifically, the system and
method produce an FIR filter with high resolution at low frequencies by
having a large number of coefficients, but reduces resolution at higher
frequencies by allowing only a fraction of the coefficients to adapt to
the high frequency part of the signal. This is accomplished by using a
multi-rate, segmented adaptation procedure, such that resolution and
bandwidth are controlled independently at each stage. If desired, the
resulting filter can be made to approximate constant Q resolution. In
addition, by adapting only a short part of the filter at a time,
misadjustment is minimized.
These and other features and objects of the present invention will be more
fully understood from the following detailed description which should be
read in light of the accompanying drawings in which corresponding
reference numerals refer to corresponding parts throughout the several
views.
BRIEF DESCRIPTION OF THE DRAWINGS
FIG. 1 is a schematic view of an adaptive FIR filter.
FIG. 2 is a schematic view of a linear predictor equalizer.
FIG. 3(a)-3(d) are graphs of time domain and corresponding frequency domain
results from operation of the method of the present invention.
FIG. 4 is a schematic diagram showing the steps in the removal of the
sampling artifacts from data.
FIG. 5 is a schematic diagram of an adaptive FIR filter of the invention.
FIG. 6 is an impulse response of an equalization filter using system and
method of the present invention.
FIG. 7 is a parallel implementation of the filter of the present invention.
FIG. 8 is a two stage multi rate implementation of the filter of the
present invention.
FIG. 9 is a schematic diagram of a three stage multi rate implementation of
the filter of the present invention.
FIG. 10 is a graph of the band-limited loudspeaker/room response for which
equalization is to be adapted using the system and method of the present
invention.
FIG. 11 is a graph showing the impulse response of a 320 coefficient
Equalization filter adapted using a conventional LMS approach.
FIG. 12 is a graph showing frequency response of the filter whose time
response is shown in FIG. 11.
FIG. 13 is a graph of the impulse response of a 320 coefficient
Equalization filter using the system and method of the present invention.
FIG. 14 is a graph of the frequency response of the filter whose time
response is shown in FIG. 13.
FIG. 15 is a flow chart for the method of the present invention performing
half band processing.
FIG. 16 is a schematic diagram of the components of the present invention.
DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS
The basic principle set forth in the description that follows is the
segmentation of an adaptive filter in both the time and frequency domains
for the purpose of controlling the resolution of the resulting filter.
This is accomplished by sequentially adapting portions of the final filter
to data at different sampling rates, and using interpolation techniques to
make the transition between different rates. The method is meant to be
applied where adaptive filtering techniques are used to design
time-invariant filters from measured data, such as equalizers.
For purposes of illustration, the system and method of the present
invention will be described in an embodiment using the LMS adaptive
algorithm in a predictive equalization filter configuration. Those skilled
in the art will recognize that the technique can be applied to other
filter configurations and adaptive algorithms in a straightforward
extensions of the basic principles involved. The method as described
requires several operations to be performed sequentially, but except for
the sampling of signals, they need not be performed in real time. In the
preferred embodiment, half-band decimation and interpolation are used when
discussing sampling rate conversions, since these lend themselves well to
producing the quasi-constant-Q resolution desired. However, any other
integer or non-integer ratio can be equally applicable, and is not meant
to be excluded from other embodiments.
Referring now to FIG. 16, the embodiment of the hardware used to implement
the adaptive filter and method of the present invention will now be
described. A ROM or RAM 52 is used for storing test signals and is down
loaded from a host computer (not shown). A signal RAM 56 is used for
storing signals during processing and a suitable RAM would be the Motorola
MCM 6164. A clock 58 operates at a frequency which is a multiple of the
desired system sampling frequencies and a program divider 60 generates
various sampling frequencies from the clock 58. A digital signal
processing microcomputer 62 controls the operation of the system of the
present invention by executing various signal processing operations,
control and input output functions. Example of a suitable digital signal
processor is the Motorola DSP 56001. Secondary digital signal processors
64 operate in parallel to provide the necessary computational power and a
suitable processor would be the Motorola DSP 56200. A Sigma-Delta type
analog digital converter 66 operates at various sampling frequencies and
an antiliasing filter and sample and hold circuit are an inherent part of
the Sigma-Delta design and a suitable device can be obtained from Motorola
under the designation 56ADC16. A digital to analog conversion system 68
includes a digital interpolation filter such as the Philips NPC5803 and
the Analog Devices AD1864. A Bessel type low pass filter 70 is used to
remove the high frequency components from the output signal. An optional
interface 54 to the host computer 62 may also be used.
The system and method operate as follows. Given that the equalization
filter is to operate over a total bandwidth of .OMEGA. Hz, then the
sampling rate, f.sub.samp, is chosen such that
f.sub.samp .gtoreq.2.OMEGA.
in accordance with the well known Nyquist criterion for sampled signals.
The total filter length N is determined by the minimum desired low
frequency resolution, f.sub.0res of the filter according to the
relationship:
N=(0.89.times.f.sub.samp)/f.sub.0res
The filter is then constructed from a set of adaptively-derived segments,
such that the frequency resolution varies as a function of frequency. The
term "segment" is used here instead of "band" to differentiate the method
of the present invention from that of conventional frequency band
splitting, since the present filter is being limited in both the time and
frequency domains. The filter may then be implemented in a conventional
FIR filter structure, or in a preferred embodiment, as a parallel
structure which exploits the special time-frequency relationship inherent
in the coefficients, as discussed below and shown in FIG. 7.
The table below illustrates the relationships between the various parameter
of adaptation for the most general case.
______________________________________
Seg- Sampling Segment
Combined
Res-
ment Freq. Range
Frequency Length Length olution
______________________________________
0 0 to f.sub.1 Hz
2f.sub.1 Hz
p.sub.O
n.sub.0 f.sub.0res
1 f.sub.1 to f.sub.2 Hz
2f.sub.2 Hz
p.sub.1
n.sub.1 f.sub.1res
2 f.sub.2 to f.sub.3 Hz
2f.sub.3 Hz
p.sub.2
n.sub.2 f.sub.2res
: : : : : :
i f.sub.i to f.sub.i + 1
f.sub.isamp Hz
p.sub.i
n.sub.i f.sub.ires
: : : : : :
J f.sub.j to f.sub.samp /2
f.sub.samp
p.sub.j
n.sub.j f.sub.jres
______________________________________
The segment i is defined in terms of frequency range and resolution, and
those parameters determine minimum values for the segment adaptation
length p.sub.i, total length n.sub.i, and the sampling frequency
f.sub.isamp. The segment adaptation length p.sub.i is determined as
follows:
##EQU1##
Note that for ease of implementation, each intermediate sampling frequency
should be an integer submultiple of fsamp. Further computational
efficiency can be achieved by using half-band processing at each stage,
which results in the parameter relationships shown in the following table.
__________________________________________________________________________
Frequency Sampling
Segment
Combined
Seg.
Range Frequency
Length
Length Resolution
__________________________________________________________________________
0 0 to f.sub.samp /2.sup.J + 1
f.sub.samp /2.sup.J
p.sub.0
n.sub.0
f.sub.0res
1 f.sub.samp /2.sup.J + 1 to f.sub.samp /2.sup.J
f.sub.samp /2.sup.J - 1
p.sub.1
2n.sub.0
f.sub.1res
2 f.sub.samp /2.sup.J to f.sub.samp /2.sup.J - 1
f.sub.samp /2.sup.J - 2
p.sub.2
4n.sub.0
f.sub.2res
: : :
: : : :
J f.sub.samp /4 to f.sub.samp /2
f.sub.samp
p.sub.J
2.sup.J n.sub.0
f.sub.Jres
__________________________________________________________________________
Briefly, the following steps are involved in the adaptation process. Both
the time domain and corresponding frequency domain results for these steps
are shown graphically in FIGS. 3a-3d. In the first step, all n.sub.0
coefficients of the first segment are set to zero, as shown in FIG. 3a.
Typically, n.sub.0 is equal to p.sub.0 in the first segment adapted. The
coefficients are adapted according to any chosen adaptive algorithm,
operating on input data sampled at the appropriate sample rate
f.sub.0samp. The resulting n.sub.0 coefficients, shown in FIG. 3b,
comprise the first segment. The filter is then interpolated by a factor L
to length n.sub.1 and scaled, forming the n.sub.1 total filter
coefficients at that stage, comprising the p.sub.1 coefficients to be
readapted, plus the interpolated tail representing the lower frequency
response of resolution f.sub.0res, as shown in FIG. 3c. The second segment
consists of the first p.sub.1 values readapted according to the chosen
algorithm, this time with input data obtained at L times the previous
sample rate. The results of this are shown in FIG. 3d. Subsequent segments
are obtained by repeating this sequence of basic operations; the entire
adaptive process is described below in more detail for the half-band case
following the discussion of several related issues.
For all segments but the last, it is necessary that the sampled data be
free of artifacts such as DC bias introduced by the analog to digital
converter, and other sampling related artifacts which would not appear if
the signal is sampled at the f.sub.samp rate. One such important artifact
is caused by the antialiasing filters in the Analog to Digital converters.
This introduces an attenuation of the high frequency signal content of the
band and a corresponding phase shift. It must not be present in the data
used by the adaptive filter, because the filter will adapt to those
characteristics and consequently be in error in subsequent segments, where
such artifacts are not present. If the adaptation is performed in real
time, then it is necessary that the Analog to Digital converter oversample
the data by a factor of two, the resulting digitized samples processed
through a half band filter and decimated by a factor of 2, with the output
used by the adaptive filter. This is shown schematically in FIG. 4. If the
adaptation procedure is being performed in non-real time using previously
sampled and stored data, then the data may be sampled at the higher rate
and post processed in a similar manner to obtain equivalent results.
The length of the sampled data sequences which are processed by the filter
during adaptation must be sufficient to allow the filter to converge for
the given filter length chosen adaptation algorithm. The parameters which
affect this are algorithm specific, but well known to those skilled in the
art. For the example described here using the LMS algorithm, the
convergence factor K is initially chosen to be near its maximum allowable
value for stable behavior, a value which can be calculated from the total
energy in the signal and the number of adapting coefficients in the filter
using well-known relationships. Then, in order to minimize misadjustment,
K is scaled by a constant less than but very close to 1.0 after each
sampling interval. When the value of K becomes small enough that
adaptation of the coefficients essentially ceases, then the process may be
terminated, and the adapted filter coefficients used as the results of the
current segment. If the adaptation is performed in non-real time using
stored samples, a significant savings in memory may be achieved if a
minimum data set is initially sampled and stored, then cycled through
repeatedly while the convergence constant is reduced as described above.
In this case the data sample size should be at least twice the total
filter length n.sub.i of the segment being adapted.
When the filter is used in the embodiment of an equalizer, the test signal
must be chosen carefully so that a correct inverse filter can be
generated. Specifically, the test signal must have a white spectrum and be
uncorrelated, and in particular, this property should be ensured in a
finite number of samples. This can be achieved by using a maximum length
sequence, or preferably obtained by creating a sequence with the desired
properties using an inverse Fourier Transform method. In this latter
method, the magnitude of the frequency transform will be specified, and
the phase component of the transform will be generated from a sequence of
random numbers between +.pi. and -.pi.. In order to generate a real
function, the real part of the transform is made symmetric about 0, and
the imaginary part anti-symmetric. By applying an inverse transform to
this data, such as is available using well known Fast Fourier Transform
techniques, a time domain signal which is a power of 2 in length and has
the desired spectral properties will be generated. A further advantage of
using such a method is its ability to generate any arbitrary spectral
shape, which can be used to superimpose a desired spectral shape on the
equalization filter. This is accomplished by designing the test signal
noise to have the inverse characteristic other than white, specifically
that of the desired post equalization spectrum.
Since sample-by-sample adaptive algorithms such as the LMS never actually
achieve the optimal solution for the coefficients at any given sampling
interval, but rather do so only in the mean, performance may be improved
by obtaining the mean using the following processing steps. Once algorithm
convergence has been achieved, R additional samples are processed with a
value of K several times greater than its final previous value. After each
of the samples are processed and coefficients updated in the chosen
manner, the filter coefficients are copied from the adaptive filter and
added into a suitable external memory array for the purposes of averaging.
After all R samples have been processed, the coefficient sums in the
memory array are each divided by R to obtain their average value, and
these are then used as the coefficient results of the current segment.
Once filter coefficients are obtained from the adaptation of the current
segment, they must be further processed in order to be used as initial
conditions for the next segment of the filter. First the coefficients are
interpolated to the new sampling rate f.sub.2 =L.times.f.sub.1 using well
known techniques employing linear phase FIR low pass filters with
normalized cutoff frequency of .pi./L and unity passband gain. For this
example L=2 (half-band case) as noted earlier. The lowpass interpolation
filter may be designed using any of the well-known techniques for
designing FIR filters. A particularly straightforward method is that of
the window-function technique, described in most texts on digital filter
design. For this application, the length of the interpolating filter M
must be odd and the quantity (M-1)/(2L) an integer in order for the
interpolation filter delay to be removed. Examples of typical values for M
are 41, 61 and 253 for L=2.
Interpolation of the segment filter coefficients is accomplished by
inserting (L-1) zero value coefficients between each of the original
coefficients, and appending [(M-1)/2]-1 zeroes to th | | |
