- Email: [email protected]
A measurement method for the linear and nonlinear properties of electro-acoustic transmission systems
SIGNAL
PROCESSING ELSEVIER
Signal Processing 64 (1998) 49-60
A measurement method for the linear and nonlinear properties of electro-acoustic transmission systems F. Heinle”, R. Rabensteinb>*, A. Stengerb aPhilips Semiconductors, Technology Centre for Mobile Communications, Stromerstrasse 5-7, D-90443 Niirnberg, Germany bLehrstuhl ftir Nachn’chtentechnik, Universitiit Erlangen-Niirnberg, Cauerstrasse 7, D-91058 Erlangen, Germany Received 9 October 1997
Abstract Competitive audio consumer products require not only cheap signal processing hardware but also low-cost analog equipment and sound transducers. The nonlinear distortions produced by these electro-acoustic transmission systems cannot be described and analyzed by standard methods based on linear systems theory alone. In order to take the nonlinear properties into account, we present a measurement method for the linear and nonlinear transmission characteristics of almost arbitrary systems and show its application to the analysis of electro-acoustic systems. Examples demonstrate the measurement of the impulse response of a loudspeaker-enclosure-microphone-system with cheap analog equipment in the presence of background noise. 0 1998 Elsevier Science B.V. All rights reserved. Zusammenfassung Wettbewerbsfahige Audio-Produkte fiir den Konsumentenbereich erfordern nicht nur billige digitale Hardware zur Signalverarbeitung sondern such kostengiinstige analoge Bauteile und Schallwandler. Die nichtlinearen Verzerrungen, die von diesen elektroakustischen Ubertragungssystemen erzeugt werden, kiinnen nicht allein durch die Standardmethoden der linearen Systemtheorie beschrieben und analysiert werden. Zur Beriicksichtigung der Nichtlinearitaten stellen wir eine MeBmethode fur die linearen und nichtlinearen Ubertragungseigenschaften sehr allgemeiner Systeme vor und zeigen ihre Anwendung auf die Analyse elektroakustischer Systeme. Beispiele demonstrieren die Messung der Impulsantwort eines Lautsprecher-Raum-Mikrofon-Systems mit billigen analogen Komponenten und mit Hintergrund-Gerluschen. 0 1998 Elsevier Science B.V. All rights reserved. R&sum& Les produits audio de grande consommation compttitifs r&lament non seulement un materiel de traitement du signal peu cotiteux mais igalement un equipement analogique et des transducteurs de son de prix peu Bevt. Les distortions non-lintaires produites par ces systemes de transmission Clectro-acoustique ne peuvent Ctre d&c&es et analys6e.s par des methodes standards bastes sur la theorie des systtmes lineaires seule. Pour tenir compte des proprietes non-lineaires, nous presentons une methode de mesure des caracteristiques lintaires et non-line&es de transmission de systeme presqu’arbitraires et montrons son application a l’analyse des systemes tlectro-acoustiques. Des exemples sont don& sur la mesure de la reponse impulsionnelle d’un systeme haut-parleur plus microphone incorporant un Bquipement analogique de faible co%. 0 1998 Elsevier Science B.V. All rights reserved. Keywor&: Measurement;
Nonlinear systems; Electra-acoustic
systems; Frequency response; Power spectral density
*Corresponding author. Tel.: + 49 9131 85 8717/7101; fax: + 49 9131 85 8849, e-mail: [email protected]. 0165-1684/98/$19.00 0 1998 Elsevier Science B.V. All rights reserved. PZISO165-1684(97)00175-8
50
F. Heinle et al. / Signal Processing
1. Introduction Real-time signal processing used to be the most expensive part of acoustic echo and noise control applications. The additional cost of complementing the digital hardware with good-quality audio equipment was tolerable and, consequently, nonlinear effects of sound transducers could be neglected. The situation is changing with the availability of cheap computing power for real-time applications. Competitive audio consumer products require not only cheap signal processing hardware but also low-cost sound transducers. Software-only solutions of speech communication features for desktop computers must rely on built-in microphones and speakers, whatever their quality may be. Therefore, nonlinear distortions have to be taken into account in the design of low-cost electro-acoustic systems. Furthermore, the digital transmission of speech or audio signals are subject to nonlinear effects. Subband-coding with the least possible number of bits assigned to each band is a standard technique. Also low-cost analog-digital converters produce distortions, which cannot be described by linear effects alone. Many acoustic echo and noise control applications require measuring or estimating the properties of the loudspeaker-enclosure-microphone-system (LEMS) including any digital pre- and post-processing. The nonlinear effects described above may show up in many places of this electro-acoustic transmission chain, so that it is usually impossible to consider them by a proper theoretical analysis. Also measurements of the distortion factor of isolated components do not give a complete picture of the nonlinear behaviour of the overall system. Moreover, common measurement methods for the room-impulse response record only the linear transmission characterics and are blind to nonlinearities. This contribution presents a measurement method for the linear and nonlinear transmission characteristics of almost arbitrary systems and shows its application to the analysis of electroacoustic systems. The measurement principle described here is based on an analysis method for finite wordlength implementations of time-invariant digital filters [l,ll,lZ] and periodically time-vary-
64 (1998) 49-60
ing multirate systems [5,6,8]. It has been successfully applied to the simultaneous determination of the frequency response, the alias components, and the quantization noise of implemented filter banks. Section 2 describes the electro-acoustic transmission system under consideration. The theoretical foundations of the measurement method are covered by Section 3 along with two alternate interpretations given in Section 4. Finally, some examples for the measurement of room-impulse responses with cheap sound transducers and background noise are presented in Section 5.
2. Problem description 2.1. Overview Fig. 1 shows an electro-acoustic transmission system for the measurement of an LEMS, consisting of amplifiers, digital-to-analog and analog-todigital converters. A discrete measurement system MS excites the transmission system with a discrete input signal u(k) and records its response y(k). The design of the measurement system is the topic of this contribution.
2.2. Linear approximation The transmission of acoustic signals from the loudspeaker to the microphone is described by the linear wave equation [ 171. However, if other sound sources are present in the LEMS, the measurement signal may be superimposed by additive background noise. Amplifiers and sound transducers can be modelled as linear systems if the distortions produced by them are negligible. This assumption
Fig. 1. Electra-acoustic transmission system for the measurement of a loudspeaker-enclosure-microphone-system (LEMS).
F. Heinle et al. / Signal Processing 64 (1998) 49-60
v(k)0
hiraW
by(k)
+y(k)
v(k)0 Fig. 3. Nonlinear model.
Fig. 2. Linear model.
is only valid for good-quality (expensive) equipment. The converters are nonlinear systems by their very nature. They can be approximated by linear systems only for a high wordlength of the digital signals. If background noise, distortions of the analog equipment, and quantization effects are neglected, the electro-acoustic system from Fig. 1 can be replaced by a linear model with an overall impulse response hIi” between input u(k) and output y(k) (see Fig. 2). An estimation G*,,(k) of this impulse response can be obtained with the Discrete Fourier Transformation (DFT), e.g. in the most simple case as
51
y(k)
v(k)
Fig. 4. Weakly nonlinear model.
both systems have to be obtained experimentally from suitably chosen input signals v(k) and the corresponding responses y(k). A method for the determination of S,_ and SN from measurements of u(k) and y(k) is described in the following.
3. Description of the measurement method
A number of techniques exist, which differ in the choice of the input signal u(k), the DFT length M, and in more elaborate details of obtaining the estimate &,,,(k) [16]. However, all these techniques rely on the assumption of linearity and are blind to any of the nonlinear effects described above.
2.3. Weakly nonlinear approximation A detailed account of all possible nonlinear contributions in Fig. 1 leads to a general nonlinear model according to Fig. 3. Due to the variety and different nature of the nonlinear components, an exact nonlinear model S of the transmission system from Fig. 1 would be too complex to be handled efficiently. In order to circumvent this difhculty, we approximate the nonlinear system S by a weakly nonlinear model [ll] as shown in Fig. 4. It consists of a parallel arrangement of a linear system S,_ and a nonlinear system SN. Again, it is not feasible to derive S,_and SN analytically from the components of the transmission system. Instead, the characteristics of
3.1. Frequency response of the linear subsystem SL
The frequency response of the linear subsystem St is determined such that the contribution of the nonlinear subsystem SN to the total output signal y(k) becomes a minimum. This is achieved by expressing the output y,(k) of S, in Fig. 4 by a convolution with the impulse response h(k) of SL: yL(k) = h(k) * u(k). The minimization of the quadratic mean value of the output of the nonlinear system E(k) = B{ln(k)l’) is equivalent to determining the impulse response h(k) such that the expected value b{ .} s(k) = Qly(k) - h(k) * W21
(2)
becomes a minimum. By expanding the squared magnitude in Eq. (2) we obtain s(k) = QWy*(W
-
“c” h*(lc)&{y(k)u*(k - rc)>
K=-m
F. Heinle et al. / Signal Processing
52 +m
c
+
rcl=-m
x
+m
c
KS=-cc
‘c”kwJ,(
-
H(e@) at M discrete frequencies Sz,,= ,u27r/M, /l = 0, . ..) M - 1. The values H(ej”) of Eq. (6) can now be expressed by V(p) = DFT,{u(k)} and y(,4 = DFT,{y(k)} VI>
4m*(%)
b(u(k - q)u*(k
- Q))
- 4
K=-cc
where q,,(k), q,,(k) and q,,(k) are the auto- and cross-correlation sequences between input and output. The unknown values h(k) of the impulse response follow by differentiation of z(k) with respect to h*(p) from the minimum condition - W)
= 0,
64 (1998) 49-60
vp.
ah*(P)
The differentiation of a function f(x) with respect to its conjugate complex variable x* can be carried out in a formally correct way by the so-called Wirtinger calculus [2,13]. The result is (5) Finally, we obtain by Fourier transformation 9 an estimation H(ejR) = F{h(k)} of the frequency response of the linear subsystem
The expected values c?{.} are calculated from the measured input and output sequences by applying L sample sequences u,(k), I = 0, . . . , L - 1, of the same stochastic process to the transmission system and by averaging the results. Thus, the estimated value fi(ej’p) for the frequency response of the linear subsystem is obtained
The procedure is simplified, if the DFT-spectrum of the sample sequences has a constant magnitude IV,(p)1= IV(& This can be achieved by deriving the sample sequences u,(k) from an arbitrary periodic, deterministic signal u(k) by adding a stochastic phase term VA(k)= uA(k + M) = u(k)ejvA.
(9)
The phase qpl is a stochastic variable, equally distributed in the interval [ - rr, + n). For input signals of this kind, the computation of the estimated value simplifies to fi(ejDu)
=
Fig. 5 shows the structure of the resulting measurement arrangement. In our application, a complex chirp signal
@&j’) = ~"~(cP~~(P)~and @,A’“,= ~{cP,~P>> are
u(k) = ej(VM)k’
the power spectral density (PSD) and the crosspower density of u(k) and y(k), respectively. For simplicity, we furthermore assume that @JejT # 0 for all frequencies of interest. Since these power densities are not known in advance, they have to be determined from the signals u(k) and y(k). In general, a wide range of possible input signals and a variety of corresponding spectral estimation methods could be utilized. Here, we will focus on periodic signals u(k) and y(k) with period M. Then it is sufficient to determine
Fig. 5. Measurement the linear subsystem
(11)
arrangement for the frequency &. (-)* denotes the conjugate
response value.
of
F. Heinle et al. J Signal Processing 64 (1998) 49-60
has been used as the deterministic component u(k) of the input. This signal achieves the minimum possible crest factor.
3.2. Power spectral density at the output of the nonlinear subsystem SN Once the frequency response of the linear subsystem SL is estimated, we can use this result to obtain the PSD @,,,(ej’r) of the nonlinear subsystem SN. For the estimated output of the nonlinear system follows from Fig. 4 that f?(k)= y(k) - 6(k) * v(k),
(12)
where h(k) is the estimated impulse response of SL. The estimation of Qnn(ejn*)is performed in the frequency domain by evaluation of
L-l
53
Thus, the PSD estimate is smeared since the spectrum R(ejT of r(k) has high sidelobes. Many spectral estimation procedures such as the WelchBartlett method or the Blackman-Tukey method try to circumvent this problem by using different windows, e.g. the Bartlett, Hamming or Hanning window. Nevertheless, all these procedures only use possibly overlapping signal blocks of length M for the estimation. As the window length M remains unchanged the sidelobes can only be reduced on the expense of the width of the mainlobe. Recently, in [9] a new windowing technique was proposed which overcomes most of these problems. Here, we will only give a brief survey of some major points of this new method. Obviously, the spectral properties of a window sequence can be improved by increasing the window length. Therefore, the method uses iV subsequent blocks of nl(k) and one common windowf(k) of length NM for estimating PSD at M frequencies the discrete 52, = 2x,u/M, p = 0, . . . ,M - 1. Thus, the PSD estimate (14) has to be modified as follows:
x 1 CIW412- lfi@~~)121W~121. (13)
1
A=0
L-l
IN.+-1
12
4n(P)= M(LL_ l,Azo) kTof(hW-‘2~‘k’M( 3.3. A new windowing technique The proposed method for estimating the PSD was originally developed for the investigation of fixed-point digital filters under the assumption of periodic quantization noise. However, it is wellknown that averaged periodograms &nn(ejop) =
L-l x
M-IN-1
2
c
c
1
A=0
k=O
v=O
f(k+vM)n,(k+vM)e-j2”“k’M
v
1
M(L1_ l)L$ IWP12 a
0
are biased in case of nonperiodic noise sequences nn(k), i.e. the estimate does not equal the true PSD even for an infinite number L of single measurements [7]. Due to the implicit windowing of nl(k) with a rectangular window r(k) =
1 = M(L - 1)
With F(e@) = Fv(k)}, yields
the PSD estimate
then
~,n(ej(Pu-r))lF(e’“)(2d[.
(18)
1, k=O ,..., M-l, 0,
else,
the expected value of Eq. (14) yields
(15) Now, it is possible to deduce the following desirable properties of the window sequence f(k) from the above equations: The PSD estimate at Sz, = 27cp/M should comprise all spectral components of the continous l
54
F. Heinle et al. / Signal Processing
PSD in the frequency bin 27c(p-4)/M ,< 0 < 27c@+ 4)/M. The optimum window therefore has the spectrum
(19) Consequently, the PSD estimate is the integral over the respective frequency bin at sl,. Obviously, F&ejo) cannot be realized with a finite window length. Therefore, we have to design a window f(k) which approximates the optimum window under the given constraints. The proposed method is capable of measuring the frequency response and the PSD simultaneously. For efficiency reasons, it is therefore necessary that the windowing does not cause any interference of periodic output signal components which result from different periodic input signal components. This can be easily achieved by prescribing the following zeroes of the window spectrum: F(ej”s)=O,
p=l,...,
M-l.
(20)
Finally, the windowing should preserve the total power, i.e.
From this constraint, the last condition immediately results in M-l
“z.
I~(&(R-2”“‘M’)12 ; M2,
(22)
Evidently, the window sequence has to be a socalled Mth-band filter which is well-known from filter bank theory [15]. In [9] a highly efficient method for designing such windows was introduced. Our PSD estimation procedure which is based on these window sequences shows an excellent performance even in case of a modest N. This will also be illustrated by the examples at the end of this paper.
64 (1998) 49-60
4. Alternate interpretations The derivation of the frequency response of the linear subsystem in Section 3.1 is quite straightforward and self-contained. In order to show interrelations with other results, we give here two alternate interpretations of the obtained estimate from the point of view of optimal filtering theory and estimation theory.
4. I. Wiener Jilter The minimization of a quadratic error term according to Eq. (2) is also used for the derivation of the optimal reconstruction filter (Wiener filter). However, the scope of Wiener filtering is different from the measurement method considered here. Nevertheless, the measurement problem and the optimal filter problem can be discussed in the same framework. Both lead to the same minimization problem and consequently the solutions are formally identical. At first we consider the general optimal reconstruction problem as shown in Fig. 6. For later reference to the measurement problem, the signals in Fig. 6 have been denoted with the same letters as in the weakly nonlinear model from Fig. 4. However, the meaning of the signals in the optimal filter problem is different. y(k) is an unknown input signal to the system S. Only the output v(k) is available by observation. The problem is to find a linear system H such that its output yL(k) is a reconstruction of the unknown signal y(k). The frequency response H(ejy of the system H is determined by least-squares minimization of the difference
n(k)= y(k)- ydk) = Y(k)- h(k)* 40,
Fig. 6. Optimal filter problem.
(23)
F. Heinle et al. / Signal Processing 64 (1998) 49-60
with respect to the coefficients h(k) of the impulse response of H. The result H(ej”) is expressed by the power spectral density Q&j”) of the observed signal v(k) and the cross-power density @JejT, which reflects the available knowledge of the system S as
In order to establish a relation to the measurement problem considered previously, we redraw Fig. 6 and arrive at the modified arrangement of Fig. 7. Note that the problem description is identical to Fig. 6. The shaded area indicates that the available knowledge of the reconstruction problem consists of the correlation properties of the signals u(k) and y(k) represented by Qi”“(ejn)and @Jej@). Now we return to the measurement problem as described in Section 2. Fig. 8 shows the arrangement for the detection of the linear subsystem of the weakly nonlinear model according to Fig. 4. The signal designations in Figs. 8 and 4 correspond to each other. The derivation of the frequency response H(ejO) of the linear subsystem by
Fig. 7. Modified arrangement of the optimal filter problem.
Fig. 8. Measurement problem.
55
minimization of n(k) = y(k) - y,(k) = y(k) - h(k) * a(k)
(25)
has been shown in Section 3.1. Here, we focus on the similarities and differences of the measurement problem according to Fig. 8 and the modified arrangement of the optimal filter problem according to Fig. 7. An obvious difference is that the role of v(k) and y(k) as input and output signals of S are interchanged. However, both problems share the property that the system S is described by the cross-power spectrum @&jr?). A second difference is the purpose of the linear system H. It serves to reconstruct an unknown signal in the optimal filter problem and it describes the linear part of an unknown system in the measurement problem. However, the minimization problems (23) and (25) are identical for both cases, irrespective of the different meaning of the signals u(k) and y(k) in Figs. 7 and 8. Since both the way of describing the system S and the minimization problems (23) and (25) formally agree, we can expect that also the solutions for the frequency responses H(eja) of the unknown systems H coincide. This is indeed the case, as is seen by comparison of the frequency response of the optimal filter according to Eq. (24) and the frequency response of the linear subsystem of Eq. (6). This shows that the foundation of the presented measurement method can also be interpreted by formal agreement with optimal filter theory.
4.2. Maximum likelihood method The frequency response estimate can even be derived in a quite different way which ensures optimality in terms of estimation theory. This approach is called maximum likelihood (ML) estimate. We will derive this estimate for a single measurement and subsequently generalize the results. We start from a frequency domain description of the weakly nonlinear system according to Fig. 9 (cf. Figs. 4 and 8). It can be formulated in terms of the DFT values of blocks of measured data with V@) = DFT{u(k)} and Y(p), N(p) alike. H(p) is the DFT of
F. Heinle et al. / Signal Processing 64 (1998) 49-60
56
Now, the PDF in the frequency domain can be determined easily using the DFT matrix W, and n = WffN/M. We obtain
Fig. 9. Additive noise description system.
of the weakly nonlinear
(26)
Under the assumption that the probability density function (PDF) of the disturbing signal component N(p) is known, we try to determine a frequency response estimate ET,,(p) which most probably caused the observed output Y(p) of the system under test. For convenience, we introduce vectors n, N, V, Y, H and & which consist of the corresponding time- or frequency-domain samples. The ML estimate of the frequency response then can be found by maximizing the conditional PDF p&YlH) with respect to the frequency response .H, i.e. we have to determine
Instead of maximizing p&V) it is convenient introduce the so-called log-likelihood function L(H) = ln(n”ldet @,,I *p&V))
Due to Eq. (26) the output Y(p) is distributed around the deterministic value H(p)V(p) with the PDF of the disturbing signal N(p) [14], i.e. (28)
Thus, we have to maximize the noise PDF p,(N) with respect to the frequency response estimate. The required frequency domain PDFp,(N) can be determined from its time-domain equivalent in a straightforward manner. For this purpose, the joint PDF of M subsequent samples of the underlying noise-like process n(k) = n,(k) + j,,(k) is needed. Under quite general conditions, this process is normally distributed and has an autocorrelation cp..(l) = (Pi.&) + V,,,(A) with am,,. = cp&~). With the autocorrelation matrix en,, = b{nnH) the required joint PDF then can be written as P&l) =
1 exp - nH@,‘n. r?(det @,,,I
(29)
to
(31)
and maximize it with respect to H(p). With N(p) = Y(p) - H(p)&) the unknown noise sequence can be eliminated and the log-likelihood function L(H) = - h”il
Mi’ c,,(Y*(v) v-o jJ=o
- H*(v) v*(v))( Y(p) - H(cL)V(p))
(32)
follows. From the maximum condition
-
ptiH(YIH) = P&’ - diag@(p)) V = p&V.
c (30)
the impulse response of the unknown linear subsystem: Y(cl) = H(p))(p) + N(p).
=:
&LW@))L0,
(33)
the ML estimate immediately results: (34)
Y(p) = &.(1L)Q).
If we furthermore restrict ourselves to signals with a deterministic magnitude, i.e.
Iw~12= ql%412h
(35)
this estimate can be generalized measurements and we obtain
&Id4 = GLnh)~ =
&
for multiple
g{ Y(p)J’*(,+
(36)
Obviously, this is equivalent to the aforementioned estimate (7) if Eq. (35) is taken into account. Furthermore, it can be shown that the estimation error d(p) = @VCL)IQ)1 is given by (37)
F. Heinle et al. / Signal Processing 64 (1998) 49-60
with constant factors upV.This is a generalization of the Cramer-Rao bound [10,14] for the joint estimation of multiple parameters. Due to the equality in Eq. (37), the ML estimate (34) meets this bound. Thus, it is an efficient estimate in terms of estimation theory. Nevertheless, it has to be emphasized that this result only holds for multiple measurements if input signals with a deterministic magnitude according to Eq. (35) are used.
5. Results This section presents measurements of an LEMS both with high- and low-quality equipment. The LEMS was an anechoic chamber containing the measurement gear which caused some reflections. The high-quality equipment was a studio microphone, a preamplifier (Bryston Mod. PB4) and an active loudspeaker (GENELEC Mod. 1013A). For the low-quality measurements, preamplifier and loudspeaker were replaced by a one-chip amplifier (TBA 820M) with a cardboard mounted noname speaker (6 cm diameter). Fig. 10 shows the frequency response of the LEMS, measured with high quality audio equipment. Here, we used a standard measurement procedure according to Eq. (1) which assumes a strictly linear system behaviour. Fig. 11 shows the same measurement using lowcost equipment and again a standard measurement method. Due to this equipment a different frequency response has to be expected even though the room impulse response was not affected by the change. Nevertheless, the measured frequency response shows an unexpectedly strong noise-like behaviour. This is due to the fact that conventional
Fig. 10. Frequency response using high-quality audio equipment and standard measurement method.
51
measurement procedures are not able to distinguish between linear and nonlinear system components, the latter being caused by distortions of the amplifier and speaker. Thus, the result is a superposition of both components. The frequency response measurement shown in Fig. 12 was performed with an additional broadband background noise 4 dB below the signal level. It was produced by a noise generator fed into a guitar amplifier. The detrimental distortion of the measurement result by the superimposed noise source is obvious. In contrast, Fig. 13 shows the measurement results with our proposed method for the low-cost equipment with background noise as in Fig. 12. This method is capable of recording linear and nonlinear effects separately. The method yields the frequency response of the linear subsystem as well as the PSD of nonlinear distortions including noise. Obviously, the true frequency response deviates from the original one (Fig. 10) due to additional linear distortions of the low-cost equipment. The remaining set of results demonstrates the suitability of the windowing technique described in Section 3.3 for the discrimination between LEMS response and narrowband background noise.
b E-10 0-20. a -300 ’
1
3
4
Fig. 11. Frequency response using low-cost audio equipment (standard measurement method).
Fig. 12. Frequency response using low-cost audio equipment with background noise (standard measurement method).
58
F. Heinle et al. / Signal Processing 64 (1998) 49-60
Fig. 13. Frequency response and PSD using low-cost audio equipment with background noise (proposed measurement method).
To this end, we repeat the measurement with high-quality equipment similar to Fig. 10. During the measurement, an additional noise source was present in the LEMS which emitted a narrowband random signal with a signal-to-noise-ratio (SNR) of 3.5 dB. The noise signal was obtained from a broadband noise source by filtering with a digital FIR filter with a passband between 1 and 1.5 kHz. Fig. 14 shows the frequency response and the PSD obtained with the proposed measurement method and a rectangular data window. As a reference, the squared magnitude of the FIR form filter which was used for the generation of the narrowband noise is shown. Although averaging with L = 100 sample sequences was used (see Eqs. (8) and (lo)), the measured PSD still differs significantly from the true one. This deviation can be attributed to the nonperiodicity of the noise signal. The same system was measured with the Mthband window described in Section 3.3 and with L = 10 and sample sequences, respectively. Two effects can be observed from the results in Figs. 15 and 16. First, the corruption of the frequency response is evidently restricted to the actual frequency range of the bandpass noise. Second, the PSD of the nonlinear disturbances closely resembles the true PSD of the noise signal. The latter is given again by the squared magnitude of the FIR form filter. The obvious improvement in performance is due to the use of an appropriately designed data window.
Fig. 14. Frequency response and PSD using high-quality audio equipment with bandpass noise (proposed measurement method, rectangular window, L = 100 sample sequences).
N
4
.
_-1;;. m 20. x30.
R-40a -501 I
-@O
I”?
form filter
I .
I
1
1
3
4
f&W Fig. 15. Frequency response and PSD using high-quality audio equipment with bandpass noise (proposed measurement method, M&band window, L = 10 sample sequences), squared magnitude of form filter.
Finally, we repeat the same experiment with the bandpass noise replaced by a sine wave of 1 kHz and an SNR of 4 dB. Fig. 17 shows the results for L = 100 sample sequences. The distortion of the measured frequency response of the linear subsystem by the sine wave is hardly visible. On the other hand, we see that the PSD gives a perfect picture of the sinusoidal disturbance of the measurement. The absence of any spectral smearing or leakage clearly demonstrates the suitability of the Mth-band window for spectral estimation.
F. Heinle et al. / Signal Processing 64 (1998) 4%60
Fig. 16. Frequency response and PSD using high-quality audio equipment with bandpass noise (proposed measurement method, Mth-band window, L = 100 sample sequences), squared magnitude of form filter.
59
linear and nonlinear properties of electro-acoustic transmission systems. The fundamentals of the measurement method can be derived from a weakly nonlinear system model by minimization of a quadratic error criterion. Nonlinear effects are characterized by the PSD of an additive noise signal. A specially adapted windowing technique provides extremely sharp estimates also for narrowband noise sources. The interrelations with optimal filtering theory and estimation theory have been shown by alternate interpretations of the measurement principle using the Wiener filter and the maximum likelihood approach. The results presented show the ability to clearly distinguish between the linear and nonlinear properties of an electro-acoustic transmission system. The frequency response of the linear subsystem and the PSD which describes nonlinear distortions and background noise can be obtained separately.
E g-10. $20
*
N -300”
1
3
\J 4
f&z] 01
q
I
I
1 ,
.
1
1
3
4
f&I Fig. 17. Frequency response and PSD using high-quality audio equipment with sinusoidal noise (proposed measurement method, Mth-band window, L = 100 sample sequences).
6. Conclusions We have presented a method for measuring the properties of the loudspeaker-enclosure-microphone-system in the presence of nonlinear distortions and background noise. These disturbances have to be considered in the evaluation of the performance of low-cost audio products in office or home environments. The measurement method is based on the analysis of implemented digital filters and multirate systems. Here, it has been extended to the investigation of the
Acknowledgements The authors wish to thank Eberhard Hlnsler for valuable hints, Thomas Klinger for preparing the results and Ivan Selesnick for his help with the manuscript.
References [l] Y. Dong, Ein neues Verfahren zur Messung der Eigenschaften schwach nichtlinearer Systeme, Ph.D. Thesis, Universitat Erlangen-Ntimberg, Lehrstuhl fiir Nachrichtentechnik, Erlangen, 1990. 121 W. Fischer, I. Lieb, Ausgewiihlte Kapitel aus der Funktionentheorie, Vieweg, Braunschweig, 1988. c31B. Girod, R. Rabenstein, A. Stenger, Einfuhrung in die Systemtheorie, Teubner, Stuttgart, 1997. c41E. Hansler, Statistische Signale, Springer, Berlin, 1997. PI F. Heinle, H.W. SchiiBler, An enhanced method for measuring the performance of multirate systems, in: Proc. Internat. Conf. on DSP, Limassol, Cyprus, 1995. pp. 182-187. VI F. Heinle, H.W. SchiiBler, Measuring the performance of implemented multirate systems, in: Proc. Internat. Conf. Acoust. Speech Signal Process. (ICASSP 96), Atlanta, USA, 1996, pp. 2754-2757. Fl SM. Kay, Modem Spectral Estimation - Theory & Application, Prentice-Hall, Englewood Cliffs, NJ, USA, 1988.
60
F. Heinle et al. / Signal Processing 64 (1998) 49-60
[8] R. Reng, H.W. SchiiSler, Measurement of aliasing distortions and quantization noise in multirate system, in: Proc. Internat. Symp. on Circuits and Systems, 1992, pp. 2328-2331. [9] H. Repp, Fensterung beim Rauschklirnnegverfahren, Ph.D. Thesis, Universitlt Erlangen-Niimberg, Lehrstuhl fur Nachrichtentechnik, Erlangen, 1998. [lo] L.L. Scharf, Statistical Signal Processing: Detection, Estimation, and Time Series Analysis, Addison-Wesley, Reading, USA, 1991. [ll] H.W. SchiiSler, Y. Dong, A new method for measuring the performance of weakly nonlinear and noisy systems, FREQUENZ 44 (1990) 82-87. [12] H.W. SchiiSler, F. Heinle, Measuring the properties of implemented digital systems, FREQUENZ 48 (1994) 3-7.
[13] H. Tietz, Funktionentheorie, in: R. Sauer, I. Szabb (Eds.), Mathematische Hilfsmittel des Ingenieurs, Vol. I, Springer, Berlin, 1970. [14] H.L. Van Trees, Detection, Estimation, and Modulation Theory, Part I, Wiley, New York, 1968. [lS] P.P. Vaidyanathan, Multirate Systems and Filter Banks, Prentice-Hall, Englewood Cliffs, NJ, USA, 1993. [16] N. Xiang, Using M-sequences for determining the impulse responses of LTI-systems, Signal Processing 28 (2) (1992) 139-159. [I73 L.J. Ziomek, Fundamentals of Acoustic Field Theory and Space-Time Signal Processing, CRC Press, Boca Raton, 1995.
PROCESSING ELSEVIER
Signal Processing 64 (1998) 49-60
A measurement method for the linear and nonlinear properties of electro-acoustic transmission systems F. Heinle”, R. Rabensteinb>*, A. Stengerb aPhilips Semiconductors, Technology Centre for Mobile Communications, Stromerstrasse 5-7, D-90443 Niirnberg, Germany bLehrstuhl ftir Nachn’chtentechnik, Universitiit Erlangen-Niirnberg, Cauerstrasse 7, D-91058 Erlangen, Germany Received 9 October 1997
Abstract Competitive audio consumer products require not only cheap signal processing hardware but also low-cost analog equipment and sound transducers. The nonlinear distortions produced by these electro-acoustic transmission systems cannot be described and analyzed by standard methods based on linear systems theory alone. In order to take the nonlinear properties into account, we present a measurement method for the linear and nonlinear transmission characteristics of almost arbitrary systems and show its application to the analysis of electro-acoustic systems. Examples demonstrate the measurement of the impulse response of a loudspeaker-enclosure-microphone-system with cheap analog equipment in the presence of background noise. 0 1998 Elsevier Science B.V. All rights reserved. Zusammenfassung Wettbewerbsfahige Audio-Produkte fiir den Konsumentenbereich erfordern nicht nur billige digitale Hardware zur Signalverarbeitung sondern such kostengiinstige analoge Bauteile und Schallwandler. Die nichtlinearen Verzerrungen, die von diesen elektroakustischen Ubertragungssystemen erzeugt werden, kiinnen nicht allein durch die Standardmethoden der linearen Systemtheorie beschrieben und analysiert werden. Zur Beriicksichtigung der Nichtlinearitaten stellen wir eine MeBmethode fur die linearen und nichtlinearen Ubertragungseigenschaften sehr allgemeiner Systeme vor und zeigen ihre Anwendung auf die Analyse elektroakustischer Systeme. Beispiele demonstrieren die Messung der Impulsantwort eines Lautsprecher-Raum-Mikrofon-Systems mit billigen analogen Komponenten und mit Hintergrund-Gerluschen. 0 1998 Elsevier Science B.V. All rights reserved. R&sum& Les produits audio de grande consommation compttitifs r&lament non seulement un materiel de traitement du signal peu cotiteux mais igalement un equipement analogique et des transducteurs de son de prix peu Bevt. Les distortions non-lintaires produites par ces systemes de transmission Clectro-acoustique ne peuvent Ctre d&c&es et analys6e.s par des methodes standards bastes sur la theorie des systtmes lineaires seule. Pour tenir compte des proprietes non-lineaires, nous presentons une methode de mesure des caracteristiques lintaires et non-line&es de transmission de systeme presqu’arbitraires et montrons son application a l’analyse des systemes tlectro-acoustiques. Des exemples sont don& sur la mesure de la reponse impulsionnelle d’un systeme haut-parleur plus microphone incorporant un Bquipement analogique de faible co%. 0 1998 Elsevier Science B.V. All rights reserved. Keywor&: Measurement;
Nonlinear systems; Electra-acoustic
systems; Frequency response; Power spectral density
*Corresponding author. Tel.: + 49 9131 85 8717/7101; fax: + 49 9131 85 8849, e-mail: [email protected]. 0165-1684/98/$19.00 0 1998 Elsevier Science B.V. All rights reserved. PZISO165-1684(97)00175-8
50
F. Heinle et al. / Signal Processing
1. Introduction Real-time signal processing used to be the most expensive part of acoustic echo and noise control applications. The additional cost of complementing the digital hardware with good-quality audio equipment was tolerable and, consequently, nonlinear effects of sound transducers could be neglected. The situation is changing with the availability of cheap computing power for real-time applications. Competitive audio consumer products require not only cheap signal processing hardware but also low-cost sound transducers. Software-only solutions of speech communication features for desktop computers must rely on built-in microphones and speakers, whatever their quality may be. Therefore, nonlinear distortions have to be taken into account in the design of low-cost electro-acoustic systems. Furthermore, the digital transmission of speech or audio signals are subject to nonlinear effects. Subband-coding with the least possible number of bits assigned to each band is a standard technique. Also low-cost analog-digital converters produce distortions, which cannot be described by linear effects alone. Many acoustic echo and noise control applications require measuring or estimating the properties of the loudspeaker-enclosure-microphone-system (LEMS) including any digital pre- and post-processing. The nonlinear effects described above may show up in many places of this electro-acoustic transmission chain, so that it is usually impossible to consider them by a proper theoretical analysis. Also measurements of the distortion factor of isolated components do not give a complete picture of the nonlinear behaviour of the overall system. Moreover, common measurement methods for the room-impulse response record only the linear transmission characterics and are blind to nonlinearities. This contribution presents a measurement method for the linear and nonlinear transmission characteristics of almost arbitrary systems and shows its application to the analysis of electroacoustic systems. The measurement principle described here is based on an analysis method for finite wordlength implementations of time-invariant digital filters [l,ll,lZ] and periodically time-vary-
64 (1998) 49-60
ing multirate systems [5,6,8]. It has been successfully applied to the simultaneous determination of the frequency response, the alias components, and the quantization noise of implemented filter banks. Section 2 describes the electro-acoustic transmission system under consideration. The theoretical foundations of the measurement method are covered by Section 3 along with two alternate interpretations given in Section 4. Finally, some examples for the measurement of room-impulse responses with cheap sound transducers and background noise are presented in Section 5.
2. Problem description 2.1. Overview Fig. 1 shows an electro-acoustic transmission system for the measurement of an LEMS, consisting of amplifiers, digital-to-analog and analog-todigital converters. A discrete measurement system MS excites the transmission system with a discrete input signal u(k) and records its response y(k). The design of the measurement system is the topic of this contribution.
2.2. Linear approximation The transmission of acoustic signals from the loudspeaker to the microphone is described by the linear wave equation [ 171. However, if other sound sources are present in the LEMS, the measurement signal may be superimposed by additive background noise. Amplifiers and sound transducers can be modelled as linear systems if the distortions produced by them are negligible. This assumption
Fig. 1. Electra-acoustic transmission system for the measurement of a loudspeaker-enclosure-microphone-system (LEMS).
F. Heinle et al. / Signal Processing 64 (1998) 49-60
v(k)0
hiraW
by(k)
+y(k)
v(k)0 Fig. 3. Nonlinear model.
Fig. 2. Linear model.
is only valid for good-quality (expensive) equipment. The converters are nonlinear systems by their very nature. They can be approximated by linear systems only for a high wordlength of the digital signals. If background noise, distortions of the analog equipment, and quantization effects are neglected, the electro-acoustic system from Fig. 1 can be replaced by a linear model with an overall impulse response hIi” between input u(k) and output y(k) (see Fig. 2). An estimation G*,,(k) of this impulse response can be obtained with the Discrete Fourier Transformation (DFT), e.g. in the most simple case as
51
y(k)
v(k)
Fig. 4. Weakly nonlinear model.
both systems have to be obtained experimentally from suitably chosen input signals v(k) and the corresponding responses y(k). A method for the determination of S,_ and SN from measurements of u(k) and y(k) is described in the following.
3. Description of the measurement method
A number of techniques exist, which differ in the choice of the input signal u(k), the DFT length M, and in more elaborate details of obtaining the estimate &,,,(k) [16]. However, all these techniques rely on the assumption of linearity and are blind to any of the nonlinear effects described above.
2.3. Weakly nonlinear approximation A detailed account of all possible nonlinear contributions in Fig. 1 leads to a general nonlinear model according to Fig. 3. Due to the variety and different nature of the nonlinear components, an exact nonlinear model S of the transmission system from Fig. 1 would be too complex to be handled efficiently. In order to circumvent this difhculty, we approximate the nonlinear system S by a weakly nonlinear model [ll] as shown in Fig. 4. It consists of a parallel arrangement of a linear system S,_ and a nonlinear system SN. Again, it is not feasible to derive S,_and SN analytically from the components of the transmission system. Instead, the characteristics of
3.1. Frequency response of the linear subsystem SL
The frequency response of the linear subsystem St is determined such that the contribution of the nonlinear subsystem SN to the total output signal y(k) becomes a minimum. This is achieved by expressing the output y,(k) of S, in Fig. 4 by a convolution with the impulse response h(k) of SL: yL(k) = h(k) * u(k). The minimization of the quadratic mean value of the output of the nonlinear system E(k) = B{ln(k)l’) is equivalent to determining the impulse response h(k) such that the expected value b{ .} s(k) = Qly(k) - h(k) * W21
(2)
becomes a minimum. By expanding the squared magnitude in Eq. (2) we obtain s(k) = QWy*(W
-
“c” h*(lc)&{y(k)u*(k - rc)>
K=-m
F. Heinle et al. / Signal Processing
52 +m
c
+
rcl=-m
x
+m
c
KS=-cc
‘c”kwJ,(
-
H(e@) at M discrete frequencies Sz,,= ,u27r/M, /l = 0, . ..) M - 1. The values H(ej”) of Eq. (6) can now be expressed by V(p) = DFT,{u(k)} and y(,4 = DFT,{y(k)} VI>
4m*(%)
b(u(k - q)u*(k
- Q))
- 4
K=-cc
where q,,(k), q,,(k) and q,,(k) are the auto- and cross-correlation sequences between input and output. The unknown values h(k) of the impulse response follow by differentiation of z(k) with respect to h*(p) from the minimum condition - W)
= 0,
64 (1998) 49-60
vp.
ah*(P)
The differentiation of a function f(x) with respect to its conjugate complex variable x* can be carried out in a formally correct way by the so-called Wirtinger calculus [2,13]. The result is (5) Finally, we obtain by Fourier transformation 9 an estimation H(ejR) = F{h(k)} of the frequency response of the linear subsystem
The expected values c?{.} are calculated from the measured input and output sequences by applying L sample sequences u,(k), I = 0, . . . , L - 1, of the same stochastic process to the transmission system and by averaging the results. Thus, the estimated value fi(ej’p) for the frequency response of the linear subsystem is obtained
The procedure is simplified, if the DFT-spectrum of the sample sequences has a constant magnitude IV,(p)1= IV(& This can be achieved by deriving the sample sequences u,(k) from an arbitrary periodic, deterministic signal u(k) by adding a stochastic phase term VA(k)= uA(k + M) = u(k)ejvA.
(9)
The phase qpl is a stochastic variable, equally distributed in the interval [ - rr, + n). For input signals of this kind, the computation of the estimated value simplifies to fi(ejDu)
=
Fig. 5 shows the structure of the resulting measurement arrangement. In our application, a complex chirp signal
@&j’) = ~"~(cP~~(P)~and @,A’“,= ~{cP,~P>> are
u(k) = ej(VM)k’
the power spectral density (PSD) and the crosspower density of u(k) and y(k), respectively. For simplicity, we furthermore assume that @JejT # 0 for all frequencies of interest. Since these power densities are not known in advance, they have to be determined from the signals u(k) and y(k). In general, a wide range of possible input signals and a variety of corresponding spectral estimation methods could be utilized. Here, we will focus on periodic signals u(k) and y(k) with period M. Then it is sufficient to determine
Fig. 5. Measurement the linear subsystem
(11)
arrangement for the frequency &. (-)* denotes the conjugate
response value.
of
F. Heinle et al. J Signal Processing 64 (1998) 49-60
has been used as the deterministic component u(k) of the input. This signal achieves the minimum possible crest factor.
3.2. Power spectral density at the output of the nonlinear subsystem SN Once the frequency response of the linear subsystem SL is estimated, we can use this result to obtain the PSD @,,,(ej’r) of the nonlinear subsystem SN. For the estimated output of the nonlinear system follows from Fig. 4 that f?(k)= y(k) - 6(k) * v(k),
(12)
where h(k) is the estimated impulse response of SL. The estimation of Qnn(ejn*)is performed in the frequency domain by evaluation of
L-l
53
Thus, the PSD estimate is smeared since the spectrum R(ejT of r(k) has high sidelobes. Many spectral estimation procedures such as the WelchBartlett method or the Blackman-Tukey method try to circumvent this problem by using different windows, e.g. the Bartlett, Hamming or Hanning window. Nevertheless, all these procedures only use possibly overlapping signal blocks of length M for the estimation. As the window length M remains unchanged the sidelobes can only be reduced on the expense of the width of the mainlobe. Recently, in [9] a new windowing technique was proposed which overcomes most of these problems. Here, we will only give a brief survey of some major points of this new method. Obviously, the spectral properties of a window sequence can be improved by increasing the window length. Therefore, the method uses iV subsequent blocks of nl(k) and one common windowf(k) of length NM for estimating PSD at M frequencies the discrete 52, = 2x,u/M, p = 0, . . . ,M - 1. Thus, the PSD estimate (14) has to be modified as follows:
x 1 CIW412- lfi@~~)121W~121. (13)
1
A=0
L-l
IN.+-1
12
4n(P)= M(LL_ l,Azo) kTof(hW-‘2~‘k’M( 3.3. A new windowing technique The proposed method for estimating the PSD was originally developed for the investigation of fixed-point digital filters under the assumption of periodic quantization noise. However, it is wellknown that averaged periodograms &nn(ejop) =
L-l x
M-IN-1
2
c
c
1
A=0
k=O
v=O
f(k+vM)n,(k+vM)e-j2”“k’M
v
1
M(L1_ l)L$ IWP12 a
0
are biased in case of nonperiodic noise sequences nn(k), i.e. the estimate does not equal the true PSD even for an infinite number L of single measurements [7]. Due to the implicit windowing of nl(k) with a rectangular window r(k) =
1 = M(L - 1)
With F(e@) = Fv(k)}, yields
the PSD estimate
then
~,n(ej(Pu-r))lF(e’“)(2d[.
(18)
1, k=O ,..., M-l, 0,
else,
the expected value of Eq. (14) yields
(15) Now, it is possible to deduce the following desirable properties of the window sequence f(k) from the above equations: The PSD estimate at Sz, = 27cp/M should comprise all spectral components of the continous l
54
F. Heinle et al. / Signal Processing
PSD in the frequency bin 27c(p-4)/M ,< 0 < 27c@+ 4)/M. The optimum window therefore has the spectrum
(19) Consequently, the PSD estimate is the integral over the respective frequency bin at sl,. Obviously, F&ejo) cannot be realized with a finite window length. Therefore, we have to design a window f(k) which approximates the optimum window under the given constraints. The proposed method is capable of measuring the frequency response and the PSD simultaneously. For efficiency reasons, it is therefore necessary that the windowing does not cause any interference of periodic output signal components which result from different periodic input signal components. This can be easily achieved by prescribing the following zeroes of the window spectrum: F(ej”s)=O,
p=l,...,
M-l.
(20)
Finally, the windowing should preserve the total power, i.e.
From this constraint, the last condition immediately results in M-l
“z.
I~(&(R-2”“‘M’)12 ; M2,
(22)
Evidently, the window sequence has to be a socalled Mth-band filter which is well-known from filter bank theory [15]. In [9] a highly efficient method for designing such windows was introduced. Our PSD estimation procedure which is based on these window sequences shows an excellent performance even in case of a modest N. This will also be illustrated by the examples at the end of this paper.
64 (1998) 49-60
4. Alternate interpretations The derivation of the frequency response of the linear subsystem in Section 3.1 is quite straightforward and self-contained. In order to show interrelations with other results, we give here two alternate interpretations of the obtained estimate from the point of view of optimal filtering theory and estimation theory.
4. I. Wiener Jilter The minimization of a quadratic error term according to Eq. (2) is also used for the derivation of the optimal reconstruction filter (Wiener filter). However, the scope of Wiener filtering is different from the measurement method considered here. Nevertheless, the measurement problem and the optimal filter problem can be discussed in the same framework. Both lead to the same minimization problem and consequently the solutions are formally identical. At first we consider the general optimal reconstruction problem as shown in Fig. 6. For later reference to the measurement problem, the signals in Fig. 6 have been denoted with the same letters as in the weakly nonlinear model from Fig. 4. However, the meaning of the signals in the optimal filter problem is different. y(k) is an unknown input signal to the system S. Only the output v(k) is available by observation. The problem is to find a linear system H such that its output yL(k) is a reconstruction of the unknown signal y(k). The frequency response H(ejy of the system H is determined by least-squares minimization of the difference
n(k)= y(k)- ydk) = Y(k)- h(k)* 40,
Fig. 6. Optimal filter problem.
(23)
F. Heinle et al. / Signal Processing 64 (1998) 49-60
with respect to the coefficients h(k) of the impulse response of H. The result H(ej”) is expressed by the power spectral density Q&j”) of the observed signal v(k) and the cross-power density @JejT, which reflects the available knowledge of the system S as
In order to establish a relation to the measurement problem considered previously, we redraw Fig. 6 and arrive at the modified arrangement of Fig. 7. Note that the problem description is identical to Fig. 6. The shaded area indicates that the available knowledge of the reconstruction problem consists of the correlation properties of the signals u(k) and y(k) represented by Qi”“(ejn)and @Jej@). Now we return to the measurement problem as described in Section 2. Fig. 8 shows the arrangement for the detection of the linear subsystem of the weakly nonlinear model according to Fig. 4. The signal designations in Figs. 8 and 4 correspond to each other. The derivation of the frequency response H(ejO) of the linear subsystem by
Fig. 7. Modified arrangement of the optimal filter problem.
Fig. 8. Measurement problem.
55
minimization of n(k) = y(k) - y,(k) = y(k) - h(k) * a(k)
(25)
has been shown in Section 3.1. Here, we focus on the similarities and differences of the measurement problem according to Fig. 8 and the modified arrangement of the optimal filter problem according to Fig. 7. An obvious difference is that the role of v(k) and y(k) as input and output signals of S are interchanged. However, both problems share the property that the system S is described by the cross-power spectrum @&jr?). A second difference is the purpose of the linear system H. It serves to reconstruct an unknown signal in the optimal filter problem and it describes the linear part of an unknown system in the measurement problem. However, the minimization problems (23) and (25) are identical for both cases, irrespective of the different meaning of the signals u(k) and y(k) in Figs. 7 and 8. Since both the way of describing the system S and the minimization problems (23) and (25) formally agree, we can expect that also the solutions for the frequency responses H(eja) of the unknown systems H coincide. This is indeed the case, as is seen by comparison of the frequency response of the optimal filter according to Eq. (24) and the frequency response of the linear subsystem of Eq. (6). This shows that the foundation of the presented measurement method can also be interpreted by formal agreement with optimal filter theory.
4.2. Maximum likelihood method The frequency response estimate can even be derived in a quite different way which ensures optimality in terms of estimation theory. This approach is called maximum likelihood (ML) estimate. We will derive this estimate for a single measurement and subsequently generalize the results. We start from a frequency domain description of the weakly nonlinear system according to Fig. 9 (cf. Figs. 4 and 8). It can be formulated in terms of the DFT values of blocks of measured data with V@) = DFT{u(k)} and Y(p), N(p) alike. H(p) is the DFT of
F. Heinle et al. / Signal Processing 64 (1998) 49-60
56
Now, the PDF in the frequency domain can be determined easily using the DFT matrix W, and n = WffN/M. We obtain
Fig. 9. Additive noise description system.
of the weakly nonlinear
(26)
Under the assumption that the probability density function (PDF) of the disturbing signal component N(p) is known, we try to determine a frequency response estimate ET,,(p) which most probably caused the observed output Y(p) of the system under test. For convenience, we introduce vectors n, N, V, Y, H and & which consist of the corresponding time- or frequency-domain samples. The ML estimate of the frequency response then can be found by maximizing the conditional PDF p&YlH) with respect to the frequency response .H, i.e. we have to determine
Instead of maximizing p&V) it is convenient introduce the so-called log-likelihood function L(H) = ln(n”ldet @,,I *p&V))
Due to Eq. (26) the output Y(p) is distributed around the deterministic value H(p)V(p) with the PDF of the disturbing signal N(p) [14], i.e. (28)
Thus, we have to maximize the noise PDF p,(N) with respect to the frequency response estimate. The required frequency domain PDFp,(N) can be determined from its time-domain equivalent in a straightforward manner. For this purpose, the joint PDF of M subsequent samples of the underlying noise-like process n(k) = n,(k) + j,,(k) is needed. Under quite general conditions, this process is normally distributed and has an autocorrelation cp..(l) = (Pi.&) + V,,,(A) with am,,. = cp&~). With the autocorrelation matrix en,, = b{nnH) the required joint PDF then can be written as P&l) =
1 exp - nH@,‘n. r?(det @,,,I
(29)
to
(31)
and maximize it with respect to H(p). With N(p) = Y(p) - H(p)&) the unknown noise sequence can be eliminated and the log-likelihood function L(H) = - h”il
Mi’ c,,(Y*(v) v-o jJ=o
- H*(v) v*(v))( Y(p) - H(cL)V(p))
(32)
follows. From the maximum condition
-
ptiH(YIH) = P&’ - diag@(p)) V = p&V.
c (30)
the impulse response of the unknown linear subsystem: Y(cl) = H(p))(p) + N(p).
=:
&LW@))L0,
(33)
the ML estimate immediately results: (34)
Y(p) = &.(1L)Q).
If we furthermore restrict ourselves to signals with a deterministic magnitude, i.e.
Iw~12= ql%412h
(35)
this estimate can be generalized measurements and we obtain
&Id4 = GLnh)~ =
&
for multiple
g{ Y(p)J’*(,+
(36)
Obviously, this is equivalent to the aforementioned estimate (7) if Eq. (35) is taken into account. Furthermore, it can be shown that the estimation error d(p) = @VCL)IQ)1 is given by (37)
F. Heinle et al. / Signal Processing 64 (1998) 49-60
with constant factors upV.This is a generalization of the Cramer-Rao bound [10,14] for the joint estimation of multiple parameters. Due to the equality in Eq. (37), the ML estimate (34) meets this bound. Thus, it is an efficient estimate in terms of estimation theory. Nevertheless, it has to be emphasized that this result only holds for multiple measurements if input signals with a deterministic magnitude according to Eq. (35) are used.
5. Results This section presents measurements of an LEMS both with high- and low-quality equipment. The LEMS was an anechoic chamber containing the measurement gear which caused some reflections. The high-quality equipment was a studio microphone, a preamplifier (Bryston Mod. PB4) and an active loudspeaker (GENELEC Mod. 1013A). For the low-quality measurements, preamplifier and loudspeaker were replaced by a one-chip amplifier (TBA 820M) with a cardboard mounted noname speaker (6 cm diameter). Fig. 10 shows the frequency response of the LEMS, measured with high quality audio equipment. Here, we used a standard measurement procedure according to Eq. (1) which assumes a strictly linear system behaviour. Fig. 11 shows the same measurement using lowcost equipment and again a standard measurement method. Due to this equipment a different frequency response has to be expected even though the room impulse response was not affected by the change. Nevertheless, the measured frequency response shows an unexpectedly strong noise-like behaviour. This is due to the fact that conventional
Fig. 10. Frequency response using high-quality audio equipment and standard measurement method.
51
measurement procedures are not able to distinguish between linear and nonlinear system components, the latter being caused by distortions of the amplifier and speaker. Thus, the result is a superposition of both components. The frequency response measurement shown in Fig. 12 was performed with an additional broadband background noise 4 dB below the signal level. It was produced by a noise generator fed into a guitar amplifier. The detrimental distortion of the measurement result by the superimposed noise source is obvious. In contrast, Fig. 13 shows the measurement results with our proposed method for the low-cost equipment with background noise as in Fig. 12. This method is capable of recording linear and nonlinear effects separately. The method yields the frequency response of the linear subsystem as well as the PSD of nonlinear distortions including noise. Obviously, the true frequency response deviates from the original one (Fig. 10) due to additional linear distortions of the low-cost equipment. The remaining set of results demonstrates the suitability of the windowing technique described in Section 3.3 for the discrimination between LEMS response and narrowband background noise.
b E-10 0-20. a -300 ’
1
3
4
Fig. 11. Frequency response using low-cost audio equipment (standard measurement method).
Fig. 12. Frequency response using low-cost audio equipment with background noise (standard measurement method).
58
F. Heinle et al. / Signal Processing 64 (1998) 49-60
Fig. 13. Frequency response and PSD using low-cost audio equipment with background noise (proposed measurement method).
To this end, we repeat the measurement with high-quality equipment similar to Fig. 10. During the measurement, an additional noise source was present in the LEMS which emitted a narrowband random signal with a signal-to-noise-ratio (SNR) of 3.5 dB. The noise signal was obtained from a broadband noise source by filtering with a digital FIR filter with a passband between 1 and 1.5 kHz. Fig. 14 shows the frequency response and the PSD obtained with the proposed measurement method and a rectangular data window. As a reference, the squared magnitude of the FIR form filter which was used for the generation of the narrowband noise is shown. Although averaging with L = 100 sample sequences was used (see Eqs. (8) and (lo)), the measured PSD still differs significantly from the true one. This deviation can be attributed to the nonperiodicity of the noise signal. The same system was measured with the Mthband window described in Section 3.3 and with L = 10 and sample sequences, respectively. Two effects can be observed from the results in Figs. 15 and 16. First, the corruption of the frequency response is evidently restricted to the actual frequency range of the bandpass noise. Second, the PSD of the nonlinear disturbances closely resembles the true PSD of the noise signal. The latter is given again by the squared magnitude of the FIR form filter. The obvious improvement in performance is due to the use of an appropriately designed data window.
Fig. 14. Frequency response and PSD using high-quality audio equipment with bandpass noise (proposed measurement method, rectangular window, L = 100 sample sequences).
N
4
.
_-1;;. m 20. x30.
R-40a -501 I
-@O
I”?
form filter
I .
I
1
1
3
4
f&W Fig. 15. Frequency response and PSD using high-quality audio equipment with bandpass noise (proposed measurement method, M&band window, L = 10 sample sequences), squared magnitude of form filter.
Finally, we repeat the same experiment with the bandpass noise replaced by a sine wave of 1 kHz and an SNR of 4 dB. Fig. 17 shows the results for L = 100 sample sequences. The distortion of the measured frequency response of the linear subsystem by the sine wave is hardly visible. On the other hand, we see that the PSD gives a perfect picture of the sinusoidal disturbance of the measurement. The absence of any spectral smearing or leakage clearly demonstrates the suitability of the Mth-band window for spectral estimation.
F. Heinle et al. / Signal Processing 64 (1998) 4%60
Fig. 16. Frequency response and PSD using high-quality audio equipment with bandpass noise (proposed measurement method, Mth-band window, L = 100 sample sequences), squared magnitude of form filter.
59
linear and nonlinear properties of electro-acoustic transmission systems. The fundamentals of the measurement method can be derived from a weakly nonlinear system model by minimization of a quadratic error criterion. Nonlinear effects are characterized by the PSD of an additive noise signal. A specially adapted windowing technique provides extremely sharp estimates also for narrowband noise sources. The interrelations with optimal filtering theory and estimation theory have been shown by alternate interpretations of the measurement principle using the Wiener filter and the maximum likelihood approach. The results presented show the ability to clearly distinguish between the linear and nonlinear properties of an electro-acoustic transmission system. The frequency response of the linear subsystem and the PSD which describes nonlinear distortions and background noise can be obtained separately.
E g-10. $20
*
N -300”
1
3
\J 4
f&z] 01
q
I
I
1 ,
.
1
1
3
4
f&I Fig. 17. Frequency response and PSD using high-quality audio equipment with sinusoidal noise (proposed measurement method, Mth-band window, L = 100 sample sequences).
6. Conclusions We have presented a method for measuring the properties of the loudspeaker-enclosure-microphone-system in the presence of nonlinear distortions and background noise. These disturbances have to be considered in the evaluation of the performance of low-cost audio products in office or home environments. The measurement method is based on the analysis of implemented digital filters and multirate systems. Here, it has been extended to the investigation of the
Acknowledgements The authors wish to thank Eberhard Hlnsler for valuable hints, Thomas Klinger for preparing the results and Ivan Selesnick for his help with the manuscript.
References [l] Y. Dong, Ein neues Verfahren zur Messung der Eigenschaften schwach nichtlinearer Systeme, Ph.D. Thesis, Universitat Erlangen-Ntimberg, Lehrstuhl fiir Nachrichtentechnik, Erlangen, 1990. 121 W. Fischer, I. Lieb, Ausgewiihlte Kapitel aus der Funktionentheorie, Vieweg, Braunschweig, 1988. c31B. Girod, R. Rabenstein, A. Stenger, Einfuhrung in die Systemtheorie, Teubner, Stuttgart, 1997. c41E. Hansler, Statistische Signale, Springer, Berlin, 1997. PI F. Heinle, H.W. SchiiBler, An enhanced method for measuring the performance of multirate systems, in: Proc. Internat. Conf. on DSP, Limassol, Cyprus, 1995. pp. 182-187. VI F. Heinle, H.W. SchiiBler, Measuring the performance of implemented multirate systems, in: Proc. Internat. Conf. Acoust. Speech Signal Process. (ICASSP 96), Atlanta, USA, 1996, pp. 2754-2757. Fl SM. Kay, Modem Spectral Estimation - Theory & Application, Prentice-Hall, Englewood Cliffs, NJ, USA, 1988.
60
F. Heinle et al. / Signal Processing 64 (1998) 49-60
[8] R. Reng, H.W. SchiiSler, Measurement of aliasing distortions and quantization noise in multirate system, in: Proc. Internat. Symp. on Circuits and Systems, 1992, pp. 2328-2331. [9] H. Repp, Fensterung beim Rauschklirnnegverfahren, Ph.D. Thesis, Universitlt Erlangen-Niimberg, Lehrstuhl fur Nachrichtentechnik, Erlangen, 1998. [lo] L.L. Scharf, Statistical Signal Processing: Detection, Estimation, and Time Series Analysis, Addison-Wesley, Reading, USA, 1991. [ll] H.W. SchiiSler, Y. Dong, A new method for measuring the performance of weakly nonlinear and noisy systems, FREQUENZ 44 (1990) 82-87. [12] H.W. SchiiSler, F. Heinle, Measuring the properties of implemented digital systems, FREQUENZ 48 (1994) 3-7.
[13] H. Tietz, Funktionentheorie, in: R. Sauer, I. Szabb (Eds.), Mathematische Hilfsmittel des Ingenieurs, Vol. I, Springer, Berlin, 1970. [14] H.L. Van Trees, Detection, Estimation, and Modulation Theory, Part I, Wiley, New York, 1968. [lS] P.P. Vaidyanathan, Multirate Systems and Filter Banks, Prentice-Hall, Englewood Cliffs, NJ, USA, 1993. [16] N. Xiang, Using M-sequences for determining the impulse responses of LTI-systems, Signal Processing 28 (2) (1992) 139-159. [I73 L.J. Ziomek, Fundamentals of Acoustic Field Theory and Space-Time Signal Processing, CRC Press, Boca Raton, 1995.