Final Project of Digital Signal Processing

 

Implementation and Analysis of Artificial Reverberation

 

Wen-Hung Lo

Jen-Chieh Hsiao

Abstract

 

  Artificial reverberation can be implemented in many ways. However, some manners could not generate natural reverberations. In the following cases, we try to compare each filter and simulate a better filter, which can make a natural reverberation almost like a real reverberant environment.


Motivation


Acoustic energy is reflected from various boundary surfaces, so that the listener can hear the multitude of reflection from many directions. If symphony orchestra is recorded in an anechoic room, it will sound terrible because of no reverberation. Through DSP manner, we can convert a natural sound into one with reverberation effect. For instance, if we have a Karaoke in a small room and want echo sound effect, in this DSP manner, we can create artificial reverberation effect for music or singing. During the simulation, we can find many ways, which can simulate reverberation. So, we try to find what’s the difference among those filters and analyze the pros and cons of each spectrum generated by those filters.

 

Features Overview of Reverberation

 

What is reverberation?

 

  1. Reflected waves of the sound

  2. Same waveform

  3. Lower level and delayed

  4. Reverberation Time: A typical listening room measuring about 2.5 X 6.5 X 3.75 metres (HLW) will have a Reverberation Time of about 0.4 seconds

 

Major concert halls have a far longer RT in the region of 1.5 to 2.5 seconds.

 


Filter Design Method

 
In the below examples, we try to simulate the typical listening room. So we take R=5000, which can have around 0.4 seconds of reverberation time.

 

1.      Single Echo Filter:


In this manner, we try to make delay units by using FIR transfer function. The transfer function is  y[n]=x[n]+
αx[n-R],  |α|<1,  The design block diagram is as Fig.1-1. In this case, we take α=0.8 and R=5000 (Refer to the Matlab code as Appendix part a.)  As a result, if we try to give an impulse as an input and get the impulse response as in Fig 1-2. We can find magnitude response as Fig 1- 3 like a comb filter and the decaying amplitude is NOT natural.

 

 

  • Fig 1-1. Single Echo FIR Filter Structure

 

We use the Matlab to get the original voice signal as Fig 1-4 and its spectrum as Fig 1-5.
  Then, input the wave file,
mssp1.wav, into this filter and get the signal diagram as Fig 1-6 and spectrum diagram as Fig 1-7. In Fig1-7, we can find the spectrum is cut in much frequency because of the comb filter effect. By the way, the decaying magnitude is not a natural mode.  Therefore, it’s not a very good way to implement artificial reverberator like this.


 

 

  • Fig1-2. (a) Typical impulse response for R=5000 and a=8

  • Fig1-3. (b) Signals of origial samples

  • Fig1-4. (c) Signals after single echo FIR filter processing
  • Fig1-5. (d) Magnitude response for R=5000 and a=8
  • Fig1-6. (e) Spectrogram of original sample
  • Fig1-7. (f) Spectrogram of single echo filter

 

2. Multiple echo FIR filter:

Fig 2-1 shows the design block diagram of Multiple Echo FIR Filter. (Refer to Matlab code in Appendix part c.)  In this case, we take R as 1000 due to too much in calculation in R=5000

 

 

  • Fig 2-1. Multiple FIR Echo Filter

 

The Transfer function is

 

,|α|<1

In Fig 2-2, we can find the impulse response like a natural decay.
In Fig 2-3, the frequency response is like a mixed comb filter.

  • Fig2-2.Impulse Response (left)

  • Fig2-3. Magnitude Response with α=0.8 for N=6 and R=4 (right)

 

Original
R=1000
R=2000
R=3000
R=4000
R=5000

 

 

3.      IIR Filter Structure:

Fig 3-1 shows the design block diagram of Multiple Echo IIR Filter. (Refer to Matlab Code in Appendix part c.) In this case, we also take R as 1000.



Fig 3-1 Multiple Echo IIR Filter

 

The Transfer Function is

, |α|<1

Fig 3-2 shows more natural decays in impulse response.
But Fig 3-3 result in many peaks in frequency response. If we don't take a big value in R, it may cause output frequency different from the original one.

 

 

  • Fig3-2.Impulse Response (left)

  • Fig3-3. Magnitude Response with α=0.8 for R=8 (right)

 

4.    Allpass Reverberator Structure
Fig 4-1 shows the design block diagram of Allpass Reverberator.(Matlab code in Appendix part d.)

 

  • Fig 4-1. Allpass Reverberator

    

 

The transfer Function is

, |α| <1
  

In this case , we take α as 0.8 and R as 1000.

  

Fig 4-2 shows an impulse response like a more natural decay than Multiple Echo IIR filter because in natural environment the original may be absorbed.

Fig 4-3 sows frequency response like a more complex comb filter which can make all frequency pass through the filter.

   

 

 

  • Fig4-2.Impulse Response (left)

  • Fig4-3. Magnitude Response with α=0.8 for R=4 (right)

 


5.          A Natural Sounding Reverberator Schem:


In this case we try to simulate with a combination of 4 IIR filters and 2 Allpass filters. We can have a general model of standard reverberation filter. With a suitable choice of decay factor,
α, and delay units, we can have a more natural reverberator.

Fig 5-1 shows the design block diagram of Allpass Reverberator. (Refer to the Matlab code as the Appendix part. e )


 

  • Fig 5-1. Natural Sounding Reverberator

1950's Style Echo

Echo Hall

   The 4 IIR filters are to simulate the reverberation from many directions, like reflection , scattering and refraction. After going into the following 2 allpass filters, we can generate more natural reverberation because of the characteristics of allpass and natural decaying impulse response from allpass filters.

 


Conclusion



 1. We compare difference of FIR and IIR filters. FIR is like a comb filter and IIR like a filter with many peaks. Therefore, FIR can retain more original frequency than IIR but IIR can have morenatural decay characteristics. It’s a trade-off. Depends on what kind of requirement that we want to implement.
 

2. In the case 5, a natural sounding reverberator, we try to simulate a reverberation paradigm like a network. By combining IIR and Allpass filters with suitable decay factors and delay units, we can simulate reverberant environment in any possible location, such as a small room or concert hall.

Appendix


(a) Single Echo FIR Filter

%Single echo FIR filter: Echos are simply generated by delay units
%Read in the sound data
[d,r]=wavread('mssp1.wav');
%R sampling periods delay generated by the FIR filter
R = 5000;
%Difference equation: y[n]=x[n]+ax[n-R]
%Equivalently, by the transfer function H(z)=1+az^(-R)
num=[1,zeros(1,R-1),0.8];
den=[1];
%The output of the FIR filter is computed using the function 'filter'
d1 = filter(num,den,d);
soundsc(d1,r);
%Typical impulse response for R=5000 and a=8
subplot(3,2,1);
i =[1,zeros(1,R)];
d2= filter(num,den,i);
stem(d2);grid
xlabel('Sample index');
ylabel('Amplitude');
title(['(a)Typical impulse response for R=5000 and a=8']);
%Signals of original samples
subplot(3,2,3);
plot(abs(d));
axis([0,110250,0,0.8]);
xlabel('Sample index');
ylabel('Amplitude');
title(['(b)Signals of original samples']);
%Signals after single echo FIR filter processing
subplot(3,2,5);
plot(abs(d1));
axis([0,110250,0,0.8]);
xlabel('Sample index');
ylabel('Amplitude');
title(['(c)Signals after single echo FIR filter processing']);
 
%Magnitude response for R=5000 and a=8
subplot(3,2,2);
tfe(d,d1,1024,[],[],512);
[h1,w] = freqz(num,den,512);
plot(w/pi,20*log10(abs(h1)));grid
xlabel('Normalized Frequency');
ylabel('Magnitude');
title(['(d)Magnitude response for R=5000 and a=8']);
%Spectrogram of original sample
subplot(3,2,4);
specgram(d,1024,r);
xlabel('Time index');
ylabel('Frequency');
title(['(e)Spectrogram of original sample']);
colorbar;
%Spectrogram of single echo FIR filter
subplot(3,2,6);
specgram(d1,2048,r);
xlabel('Time index');
ylabel('Frequency');
title(['(f)Spectrogram of single echo FIR filter']);
colorbar;

(b) Multiple Echo FIR Filter

%Using FIR filter to generate multiples echoes spaced R sampling periods 
%apart with exponentially decaying amplitudes.
%Multiple echo FIR filter with a=0.8 for N=6 and R=1000
[d,r]=wavread('mssp1.wav');
N=6;
a=0.8;
b=a^N;
R=1000;
num=[1,zeros(1,R*N-1),-b];
den=[1,zeros(1,R-1),-a];
d1=filter(num,den,d);
soundsc(d1,r);
%Typical impulse response with a=0.8 for N=6 and R=4
subplot(1,2,1);
RI=4;
I =[1,zeros(1,N*RI-1)];
numi=[1,zeros(1,RI*N-1),-b];
deni=[1,zeros(1,RI-1),-a];
d2= filter(numi,deni,I);
stem(d2);grid;
axis([0,N*RI,0,1]);
xlabel('Sample index');
ylabel('Amplitude');
title(['Typical impulse response with a=0.8 for N=6 and R=4']);
%Magnitude response with a=0.8 for N=6 and R=4
subplot(1,2,2);
tfe(I,d2,1024,[],[],512);
[h1,w] = freqz(numi,deni,512);
plot(w/pi,20*log10(abs(h1)));grid
xlabel('Normalized Frequency');
ylabel('Magnitude');
title(['Magnitude response with a=0.8 for N=6 and R=4']);

(c) Multiple Echo FIR Filter
%IIR filter generating an infinite number of echoes with a=0.8 for R=1000
[d,r]=wavread('mssp1.wav');
num=[0,zeros(1,1000),1]
den=[1,zeros(1,1000),-0.8];
d1=filter(num,den,d);
soundsc(d1,r);
%Typical impulse response with a=0.8 for R=8
subplot(1,2,1);
I =[1,zeros(1,120)];
numi=[0,zeros(1,7),1]
deni=[1,zeros(1,7),-0.8];
d2= filter(numi,deni,I);
stem(d2);grid;
axis([0,120,0,1]);
xlabel('Sample index');
ylabel('Amplitude');
title(['Typical impulse response with a=0.8 for R=7']);
%Magnitude response with a=0.8 for R=8
subplot(1,2,2);
tfe(I,d2,1024,[],[],512);
[h1,w] = freqz(numi,deni,512);
plot(w/pi,20*log10(abs(h1)));grid
xlabel('Normalized Frequency');
ylabel('Magnitude');
title(['Magnitude response with a=0.8 for R=7']);

(d) Allpass Reverberator
%Allpass reveberator with a=0.8 for R=8
[d,r]=wavread('mssp1.wav');
num=[0.8,zeros(1,999),1];
den=[1,zeros(1,999),0.5];
d1=filter(num,den,d);
soundsc(d1,r);
%Typical impulse response with a=0.8 for R=4
subplot(1,2,1);
I =[1,zeros(1,60)];
numi=[0.8,zeros(1,3),1];
deni=[1,zeros(1,3),0.8];
d2= filter(numi,deni,I);
stem(d2);grid;
axis([0,60,-1,1]);
xlabel('Sample index');
ylabel('Amplitude');
title(['Typical impulse response with a=0.8 for R=4']);
%Magnitude response with a=0.8 for R=4
subplot(1,2,2);
tfe(I,d2,1024,[],[],512);
[h1,w] = freqz(numi,deni,512);
plot(w/pi,20*log10(abs(h1)));grid
xlabel('Normalized Frequency');
ylabel('Magnitude');
title(['Magnitude response with a=0.8 for R=4']);

(e) Natural Sound Reverberator
%A proposed natural sounding reveberator
a1=0.6;
a2=0.4;
a3=0.2;
a4=0.1;
a5=0.8;
a6=0.8;
a7=0.8;
R1=2000;
R2=4000;
R3=1000;
R4=500;
R5=2000;
R6=2000;
[d,r]=wavread('mssp1.wav');
num1=[0,zeros(1,R1-1),1];
den1=[1,zeros(1,R1-1),-a1];
d1=filter(num1,den1,d);
num2=[0,zeros(1,R2-1),1];
den2=[1,zeros(1,R2-1),-a2];
d2=filter(num2,den2,d);
num3=[0,zeros(1,R3-1),1];
den3=[1,zeros(1,R3-1),-a3];
d3=filter(num3,den3,d);
num4=[0,zeros(1,R4-1),1];
den4=[1,zeros(1,R4-1),-a4];
d4=filter(num4,den4,d);
dIIR=d1+d2+d3+d4;
num5=[a5,zeros(1,R5-1),1];
den5=[1,zeros(1,R5-1),a5];
dALL1=filter(num5,den5,dIIR);
num5=[a6,zeros(1,R6-1),1];
den5=[1,zeros(1,R6-1),a6];
dALL2=filter(num5,den5,dALL1);
dTOTAL=d+a7*dALL2;
soundsc(dTOTAL,r);




Reference:

1.      Sanjit K. Mitra, Digital Signal Processing Second Edition, 2000
2.      Applications of Digital Signal Processing, Alan V. Oppenheim, 1978
3.      Speech And Audio Signal Processing, Ben Gold and Nelson Morgan, 2000
4.      An Analysis/Synthesis Approach To Real-Time Artificial Reverberation, 1992 IEEE