Transcript Room Impulse Response Measurement

Digital sound processing
Convolution
Digital Filters
FFT
Fast Convolution
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 1
Sampling
 Sampling an electrical signal means capturing its value
repeatedly at a constant, very fast rate.
 Sampling frequency (fs) is defined as the number of samples
captured per second
 The sampled value is known with finite precision, given by the
“number of bits” of the analog-to-digital converter, which is limited
(typically ranging between 16 and 24 bits)
 Of consequence, in a time-amplitude chart, the
analog waveform is approximated by a sequence
of points, which lye in the knots of a lattice, as
both time and amplitude are integer multiplies of
small “sampling units” of time and amplitude
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 2
Time/frequency discretization
DV
Dt
Analog signal (true)
Digital signal (sampled)
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 3
Fidelity of sampled signals
 Can a sampled digital signal represent faithfully
the original analog one?
 YES, but only if the following “Shannon theorem”
is true:
“Sampling frequency must be at least twice of
the largest frequency in the signal being
sampled”
A frequency equal to half the sampling frequency is named
the “Nyquist frequency”– for avoiding the presence of
signals at frequencies higher than the Nyquist’s one, an
analog low-pass filter is inserted before the sampler.
It is called an “anti Aliasing” filter.
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 4
Common cases
 CD audio – fs = 44.1 kHz – discretization = 16 bit
Nyquist frequency is 22.05 kHz, the anti-aliasing starts at 20 kHz, so that at
22.05 kHz the signal is already attenuated by at least 80 dB. Hence the filter
is very steep, causing a lot of artifacts in time domain (ringing, etc.)
 DAT recorder – fs = 48 kHz – discretization = 16 bit
Nyquist frequency is 24 kHz, the anti-aliasing starts at 20 kHz, so that at 24
kHz the signal is already attenuated by at least 80 dB. Now the filter is less
steep, and the time-domain artifacts are almost gone.
 DVD Audio – fs = 96 kHz – discretization = 24 bit
Nyquist frequency is 48 kHz, but the anti-aliasing starts around 24 kHz, with
a very gentle slope, so that at 48 kHz the signal is attenuated by more than
120 dB.
 Such a gentle filter is very “short” in time domain, hence there are virtually no
time-domain artifacts.
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 5
Impulse Response
Time of flight
Direct sound
Early reflections
Reverberant tail
System
under test
Unit pulse d
09.11.2012 Angelo Farina
System’s impulse response
UNIPR | All Rights Reserved | Confidential
| Page 6
A simple linear system
Real-world system (one input, one output)
CD player
Amplifier
Loudspeaker
Microphone
Analyzer
“SYSTEM”
Block diagram
x(t)
h(t)
y(t)
Input signal
System’s Impulse
Response
(Transfer function)
Output signal
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 7
FIR Filtering (Finite Impulse Response)
x(t)  x(i  Dt)
h(t)  h(i  Dt)
y(t)  y(i  Dt)
The effect of the linear system h on the signal x passing through it is
described by the mathematical operation called “convolution”, defined by:
N 1
y(i)   x i  j  h j
j0
This “sum of products” is also called FIR filtering, and models
accurately any kind of linear systems.
This is usually written, in compact notation, as:
y(i)  xi   h j
“convolution” operator
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 8
IIR Filtering (Infinite Impulse Response)
x(t)  x(i  Dt)
a ( j) 


b( j) 
y(t)  y(i  Dt)
Alternatively, the filtering caused by a linear system can also be described
by the following recursive formula:
y( i ) 
N 1
N 1
j0
j1
 xi  j  a  j   yi  j  b j
•
In practice, the filter is computed not only from the input samples x, but also as a
function of the output samples y, obtained at the previous time steps.
•
In many cases, this method allows for representinging faithfully the behaviour of the
system with a much smaller number of coefficients than when employing FIR
filtering.
•
However, modern algorithms on fast computers make
FIR filtering preferable and even faster
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 9
The FFT Algorithm
 The Fast Fourier Transform (FFT) is often employed in Acoustics, with
two goals:
► Performing spectral analysis with constant bandwidth
► Fast FIR filtering
 FFT transforms a segment of time-domain data in the corresponding
spectrum, with constant frequency resolution, starting at 0 Hz (DC) up to
Nyquist frequency (which is half of the sampling frequency)
 The longer the time segment, the narrower will be the frequency
resolution:
[N sampled points in time] => [N/2+1 frequency bands]
(the +1 represents the band at frequency 0 Hz, that is the DC component
– but in acoustics, this is always with zero energy…)
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 10
The FFT Algorithm
The number of points in the time block must be a power of 2 – for example:
4096, 8192, 16384, etc.
Time signal (64 points)
IFFT
FFT
The inverse transform is also possible
(from frequency to time)
09.11.2012 Angelo Farina
Frequency spectrum
(32 bands + DC)
UNIPR | All Rights Reserved | Confidential
| Page 11
Complex spectrum, autospectrum
 FFT yields a complex spectrum, at every frequency we get a value made of
a real and an imaginary parts (Pr, Pi), or, equivalently, by modulus and
phase
 In many cases the phase is considered meaningless,and only the
magnitude of the spectrum is plotted in dB:
{p1, p2 , p3, p4 ,...,p N }  [ FFT]  {P0 , P1, P2 ,...PN / 2}
 Pr f 2  Pi f 2 
 Pf   P' f 
Lp f   10  log10 
  10  log10 

2
2
po
 po



•
The second version of the formula contains the definition of the
Autospectrum, that is the product, at every frequency, of the spectral
complex number P(f) with its complex conjugate P’(f)
09.11.2012 Angelo Farina
Elaborazione numerica
suono
UNIPRdel
| All
Rights Reserved | Confidential
| Page 12
Complex spectrum, autospectrum
In other cases, also
the phase
information is
relevant, and is
charted separately
(mainly when the
FFT is applied to an
impulse response).
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 13
“leakage” and “windows”
 One of the assumptions of Fourier analysis is that the time-segment
analysed represents a complete period of a periodic waveform
 This is generally UNTRUE: the imperfect connection of the end of a
segment with the beginning of the next one (identical, as the signal is
assumed to be periodic), causes a “click”, which produces a wide-band
“white noise”, contaminating the whole spectrum (“leakage”):
Leakage
Theoretical
spectrum
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 14
“leakage” and “windows”
 If we want to analyze a generic, aperiodic signal, we need to “window” the
signal inside the block being analyzed, bringing it to zero at both ends
 To this purpose, many differnet types of “windows” are used, named
“Hanning”, “Hamming”, “Blackmann”, “Kaizer”, “Bartlett”, “Parzen”, etc.
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 15
Window overlapping
 The problem is that events occurring near the ends of two adjacent blocks
are substantially not analyzed
 To avoid this loss of information, instead of shifting the analysis window by
one whole block, we need to analyze partially-overlapped blocks, with at
least 50% overlapping, usually overalpped at 75% or even more
09.11.2012 Angelo Farina
Block 1
Window
FFT
Block 2
Window
FFT
Block 3
Window
FFT
UNIPR | All Rights Reserved | Confidential
| Page 16
Averaging, waterfall, spectrogram
 Once a sequence of FFT spectra is obtained, we can average them either
exponentially (Fast, Slow) o linearly (Leq), emulating a SLM
 Alternatively, we can visualize how the spectrum changes over time, by two
graphical representations called “waterfall” and “spectrogram” (or “sonogram”)
17
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 17
Fast FIR filtering with FFT
Convolution is signifcantly faster if performed in
frequency domain:
x(n)
FFT
x(n)  h(n)
y(n)
y(n)
Problems
Solution
X(k)
X(k)  H(k)
IFFT
Y(k)
• The whole lenght of signal must be recorded before being processed
• if N is large, a lot of memory is required.
• “Overlap & Save”Algorithm
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 18
Overlap & Save
Convoluzione veloce FFT con Overlap & Save (Oppenheim &
Shafer, 1975):
xm(n)
FFT
N-point
x
h(n)
IFFT
Xm(k)H(k)
Select last
N–Q+1
samples
FFT
N-point
Append to
y(n)
Problems
Solution
• Excessive processing latency between input and output
• If N is large, a lot of memory is still required
• “uniformly-partitioned Overlap & Save”
09.11.2012 Angelo Farina
UNIPR | All Rights Reserved | Confidential
| Page 19
Uniformly Partitioned Overlap & Save
Filter’s impulse
response h(n) is also
partitioned in a
number of blocks
1st block
2nd block
09.11.2012 Angelo Farina
3rd block
4th block
Now latency and
memory occupation
reduce to the length
of one single block
UNIPR | All Rights Reserved | Confidential
| Page 20