Loren on the Art of MATLAB

Turn ideas into MATLAB

Note

Loren on the Art of MATLAB has been archived and will not be updated.

Vuvuzela Denoising with Parametric Equalizers

Today I’d like to introduce a guest blogger, Charu, who is a Technical Marketing Manager at The MathWorks. Charu regularly presents webinars on signal processing topics, and most recently presented a webinar highlighting R2010a capabilities in the area of MATLAB signal processing.

Contents

Introduction

Hi there - soccer fever is in full swing, and the incessant sounds of the vuvuzela are a big part of the soccer viewing experience. This South African instrument has become so popular that it now has a dedicated Facebook page! Well, 2 weeks ago, Loren featured a very interesting MATLAB Central submission by Vincent Choqueuse on vuvuzela denoising using the spectral subtraction method. This cool demo inspired me to try out other filtering techniques for the vuvuzela denoising problem.

A quick search on “vuvuzela filtering” suggested that notch filters and adaptive filters were the most common methods, so I decided to try a method that uses notch parametric equalizers. Parametric equalizers allow the user to control the bandwidth, gain and center frequency of the filter. Since the game commentary is a streaming broadcast, I also wanted to see if I could apply stream processing techniques using new System objects in MATLAB.

This example demonstrates how you can design parametric equalizer filters in MATLAB and use them to filter out vuvuzela sounds from an audio stream. This demo uses Signal Processing Toolbox, Filter Design Toolbox and Signal Processing Blockset.

Set up input and output streams

First, let us look at how we can create streaming inputs and outputs in MATLAB. To get the data, download 'vuvuzela.wav' from this MATLAB Central entry.

We need an input audio stream that can read data from a WAV file. To do this, we can use the MultimediaFileReader System object.

Create a System object to read an audio stream from the WAV file

hSound = signalblks.MultimediaFileReader(...
    'Vuvuzela.wav',...
    'SamplesPerAudioFrame', 1024, ...
    'AudioOutputDataType', 'double');

Get the input sampling frequency

hInfo = info(hSound);
Fs = hInfo.AudioSampleRate;

Create a System object to play back audio to the sound card

hPlayer = signalblks.AudioPlayer('SampleRate', Fs);

Create a System object to write the demo output to a WAV file.

hOutputFile = signalblks.MultimediaFileWriter(...
    'with_then_without_vuvuzela.wav',...
    'SampleRate', Fs);

Create System objects for logging pure noise and filtered output signals.

hLogNoise = signalblks.SignalLogger;
hLogOutput = signalblks.SignalLogger;

Extracting the noise spectrum

The first 5 seconds of the input do not contain any significant data, ie, no noise even, so let’s advance the input by 5 seconds, i.e., about 100 frames (@Fs = 22050Hz, 1024 samples per frame).

for frame_index=1:100
    step(hSound);
end

To obtain our “pure noise” signal, we use a 1-second portion of the signal that has only the vuvuzela sound, but no speech.

for noise_frames = 1:20
    noise = step(hSound);
    step(hLogNoise, noise);
end

Plot the noise spectrum.

pwelch(hLogNoise.Buffer, [], [], 8192, Fs);

Identifying vuvuzela frequencies

Let’s zoom into the plot to identify noise peaks.

The fundamental and the first 4 harmonic frequencies of the vuvuzela are at: 235Hz, 476Hz, 735Hz, 940Hz, 1180Hz.

Designing the parametric equalizers

Our next task is to design notch parametric equalizers that will selectively filter out these frequencies. Since the notch filters may also remove components of the speech signal at the above frequencies, we will apply a peak filter to boost the speech signal at the output of the notch filters.

To design parametric equalizers, we need to specify filter order, center frequency, and the bandwidth or Q-factor of the filter. To determine the bandwidth or Q-factor of the filter, let’s take a closer look at the noise spectrum. At the fundamental and 1st harmonic frequeuncies (235Hz, 465Hz), we see that the peaks are sharp. Beyond the 1st harmonic, the noise peaks are spread across a wider band. So, for the notch filters at the lower frequencies, we specify a narrow stopband, i.e., a high Q, and for higher harmonics, we specify a wider stopband, i.e., a low Q. After some trial and error, I arrived at the following parameters:

Here are the specifications for notch parametric equalizers.

F01 = 235; F02 = 476; F03 = 735; F04 = 940; F05 = 1180;

N = 2; %Filter order
Gref = 0; %Reference gain
G0 = -20; %Attenuate by 20dB
Qa = 10;  %Higher Q-factor for fundamental and 1st harmonic
Qb = 5;   %Lower Q-factor for higher harmonics

f1 = fdesign.parameq('N,F0,Qa,Gref,G0',N,F01,Qa,Gref,G0,Fs);
f2 = fdesign.parameq('N,F0,Qa,Gref,G0',N,F02,Qa,Gref,G0,Fs);
f3 = fdesign.parameq('N,F0,Qa,Gref,G0',N,F03,Qb,Gref,G0,Fs);
f4 = fdesign.parameq('N,F0,Qa,Gref,G0',N,F04,Qb,Gref,G0,Fs);
f5 = fdesign.parameq('N,F0,Qa,Gref,G0',N,F05,Qb,Gref,G0,Fs);

Here are specifications for peak parametric equalizer.

F06 = 480;  %Center frequency to boost

N_peak = 8; %Use higher filter order for peak filter
Qc = 0.5;   %Use very low Q-factor for boost equalizer filter
G1 = 9;     %Boost gain by 9dB

f6 = fdesign.parameq('N,F0,Qa,Gref,G0',N_peak,F06,Qc,Gref,G1,Fs);

Visualizing magnitude response of the filters

Hp1 = design(f1, 'butter');
Hp2 = design(f2, 'butter');
Hp3 = design(f3, 'butter');
Hp4 = design(f4, 'butter');
Hp5 = design(f5, 'butter');
Hp6 = design(f6, 'butter');
  hFV = fvtool([Hp1 Hp2 Hp3 Hp4 Hp5 Hp6]);
  legend(hFV,'NotchEQ #1', 'NotchEQ #2', 'NotchEQ #3',...
        'NotchEQ #4','NotchEQ #5','Peak EQ');

Implementing the parametric equalizers

To implement each filter, we use a System object that realizes an SOS (second-order-sections) structure. Using a System object to implement the filter enables us to apply the filter in a stream processing loop, without having to worry about updating filter states between each iteration.

hEQ1 = signalblks.BiquadFilter('SOSMatrix',Hp1.sosMatrix,'ScaleValues',Hp1.ScaleValues);
hEQ2 = signalblks.BiquadFilter('SOSMatrix',Hp2.sosMatrix,'ScaleValues',Hp2.ScaleValues);
hEQ3 = signalblks.BiquadFilter('SOSMatrix',Hp3.sosMatrix,'ScaleValues',Hp3.ScaleValues);
hEQ4 = signalblks.BiquadFilter('SOSMatrix',Hp4.sosMatrix,'ScaleValues',Hp4.ScaleValues);
hEQ5 = signalblks.BiquadFilter('SOSMatrix',Hp5.sosMatrix,'ScaleValues',Hp5.ScaleValues);
hEQ6 = signalblks.BiquadFilter('SOSMatrix',Hp6.sosMatrix,'ScaleValues',Hp6.ScaleValues);

Filtering the input in a stream processing loop

To filter out the vuvuzela noise, we use a stream processing loop that reads in one frame of input data from the WAV file, applies the series of 6 parametric equalizers, and then plays back the filtered output through the sound card.

for frames = 1:400
    %Step through input WAV file one frame at a time
    input = step(hSound);

    %Apply notch EQ filters first, then peak EQ filter
    out1 = step(hEQ1, input);
    out2 = step(hEQ2, out1);
    out3 = step(hEQ3, out2);
    out4 = step(hEQ4, out3);
    out5 = step(hEQ5, out4);
    filt_output = step(hEQ6, out5);

    %Play 200 frames of input signal, then 200 frames of filtered output
    %to compare original and filtered signals
    if frames < 200
      step(hPlayer, input);
    else
       step(hPlayer, filt_output);
    end

   %Log filtered output to buffer, WAV file
   step(hLogOutput, filt_output);
   step(hOutputFile, filt_output);

end

Visualizing and listening to output results

Listen to the noisy and filtered signals here.



Visualize the spectrum of the filtered signal. If you zoom into the plot, you can see that the vuvuzela frequencies have been attenuated by our notch parametric equalizers.

    pwelch(hLogOutput.Buffer, [], [], 8192, Fs);

Cleanup

Close input and output System objects to ensure that we leave all objects in their proper states.

close(hSound);
close(hPlayer);
close(hOutputFile);

Things to Try

This demo was a lot of fun to work on and I hope you’ll find it interesting. Feel free to download the code and try your own equalizer parameters. Have fun!




Published with MATLAB® 7.10


  • print