Cs 312 Project 2: Signal Processing

CS 312: SIGNAL PROCESSING LAB CS 312 TA

1. Objectives • To become familiar with the convolution operator and how it applies to signal processing. • To apply the fast fourier transform to a real world problem. • To understand how a divide and conquer algorithm can increase efficiency.

2. Background Signal processing deals with the analysis and manipulation of signals. A signal is any time varying quantity, and can be thought of as a function over time. An example of a signal is a codified message being sent accross a communication channel where at each time t, a certain symbol is received. A sound recording is an example of a signal which has a specific amplitude (or volumn) at each time t. The pitch heard is determined by the frequency of the changing signal. Signals can be analyzed and manipulated just like functions. For example two signals a and b can be added together to produce a + b = (a0 + b0 , a1 + b1 , . . . , an−1 + bn−1 ). One way of combining two signals is the convolution. Suppose you have two vectors or signals of length m and n respecively: a = (a0 , a1 , . . . , am−1 ), b = (b0 , b1 , . . . , bn−1 ). The convolution of a and b, written a ∗ b, will be a vector of length m + n − 1, where coordinate k is equal to the sum of the products a and b who’s indecies sum to k. X

ai bj .

(i,j):i+j=k

For example, if c = a ∗ b, then c4 = a0 b4 + a1 b3 + a2 b2 + a3 b1 + a4 b0 because each index of a and b sum to 4. So, a∗b = (a0 b0 , a0 b1 +a1 b0 , a0 b2 +a1 b1 +a2 b0 , . . . , am−2 bn−1 +am−1 bn−2 , am−1 bn−1 ). 1

2

CS 312 TA

One way to get the convolution is to imagine the matrix a0 b0 a1 b0 a2 b0 .. .

a0 b1 a1 b1 a2 b1 .. .

... ... ... .. .

a0 bn−2 a1 bn−2 a2 bn−2 .. .

a0 bn−1 a1 bn−1 a2 bn−1 .. .

am−1 b0 am−1 b1 . . . am−1 bn−2 am−1 bn−1 then each coordinate of the convolution vector is found by summing along the diagonals, starting on the left entry of column one and summing diagonnaly up and to the right. The convolution operator is used widely in signal processing for many different applications. When two signals are convolved they produce a new signal that is a combination of the other two. With sound signals this can produce many different effects like echo, noise cancelation, and various types of noise filters. In order to compute the convolution of two vectors a and b the following sum needs to be calculated X ai bj (i,j):i+j=k

which is just the definition of convolution. Here we will describe the “reverseand-slide” method for computing this sum. Let c = a ∗ b. Start by listing the vector a in reverse order and then list the vector b underneath it so that the first entry of b is right under the first entry of a as shown below. am−1 am−2 . . . a2 a1 a0 b0 b1 b2 . . . bn−2 bn−1 Now the first index of the convolution, c0 , is found by multiplying columns together and adding the results together, treating the empty spots as 0. So c0 = a0 b0 . Now to get the next index shift the bottom row to the left by one and repeat the sum of the column multiplication. am−1 am−2 . . . a2 a1 a0 b0 b1 b2 . . . bn−2 bn−1 This results in c1 = a1 b0 + a0 b1 . After the next iteration am−1 am−2 . . . a2 a1 a0 b0 b1 b2 . . . bn−2 bn−1 You get c2 = a2 b0 + a1 b1 + a0 b2 . Continue in this manner until the last two iterations which look like this b0 b1 b2 . . .

am−1 am−2 . . . a2 a1 a0 bn−2 bn−1

b0 b1 b2 . . . bn−2

am−1 am−2 . . . a2 a1 a0 bn−1

CS 312: SIGNAL PROCESSING LAB

3

Giving you cm+n−2 = am−1 bn−2 + am−2 bn−1 and cm+n−1 = am−1 bn−1 for the last two terms. Analyzing this algorithm yields a run time of O(m ∗ n) or O(n2 ) if the signals are the same length. Go ahead and prove it to yourself. There is another way to compute the convolution of two vectors or signals that will run in O(n log n). This requires using the Fast Fourier Transform. We will treat our vectors/signals as polynomials. The vector a = (a0 , a1 , . . . , am−1 ) becomes the polynomial A(x) = a0 + a1 x + a2 x2 + . . . + am−1 xm−1 and vector b becomes the polynomial B(x) = b0 + b1 x + b2 x2 + . . .+bn−1 xn−1 . Now C(x) = A(x)B(x) is a polynomial where the coefficients form a vector c which is equal to the convolution of a and b. The algorithm proceeds as follows: (1) First we choose m + n values x1 , x2 , . . . , xm+n and evaluate A(xj ) and B(xj ) for each of j = 1, 2, . . . , m + n. (2) We can now compute C(xj ) for each j very easily: C(xj ) = A(xj )B(xj ). (3) Finally, we recover C from its values on x1 , x2 , . . . , xm+n . By using polynomial interpolation any polynomial of degree d can be reconstructed from its values on any set of d + 1 or more points. Since A and B each have degree at most m − 1 and n − 1 respectively, their product C has degree at most m + n − 2 so it can be reconstructed from m + n − 1 or more points. 3. Signal Input The files that you will be given as signals to convolve are .raw files which are essentially .wav files that had the headers removed. They are all standardized to the following standard, Sample rate: 22050, 16 bit word sizes per sample, with linear 2’s compliment encoding, and 1 channel (mono as opposed to stereo). Which means that each sample of the sound file is a 16 bit 2’s compliment number, so you can just read in 16 bits and interpret them as an integer (or short). Since it is a binary file there are no spaces or parsing, just a stream of 16 bit words. A helpful class you may want to use for reading in binary data is a BinaryReader. To output the convolved signal the program should just output 16 bit integers (or shorts) to a stream one after another. Make sure that you don’t just output standard 32 bit integers or there will be extra data in the sound file and it will sound garbled. Again, a helpful class you may want to use is the BinaryWriter. Once you output the raw sound file you can listen to it if you have a player that can play raw sound files. One that I know of is GoldWave. (Trial available at http://www.goldwave.com/release.php) When opening a raw file it asks for the file format (sample rate, word size, encoding, and channel) and when you enter the format (as specified above: 22050 sample rate, 16 bit word size, linear 2’s comp. encoding, 1 channel) it will then play the raw file. Alternatively you can convert the file to a wav file by using

4

CS 312 TA

Sox (downloadable at http://sox.sourceforge.net/). To convert a file from .raw to .wav using our format you just enter this command in the command prompt. sox -r 22050 -s -w rawfilename.raw -r 22050 -s -w newfilename.wav You can also convert any .wav file to our .raw file format by entering the command sox wavfilename.wav -r 22050 -s -w newfilename.raw In these commands -r is the option to set the sample rate, -s is for linear 2’s complement encoding, and -w is for 16 bit word size. These options must always be specified for .raw files, whereas in .wav files this information is contained in the header. 4. To Do Your assignment for this project is to implement the Fast Fourier Transform algorithm in C#. Your program will read in two different signals from two different files. It will then produce a new signal that is the convolution of the two input signals and output this new signal to a new file. The filenames for both input signals and output signal will need to be entered on the command line as parameters. Running your program from the command line must look like this C:\>convolution.cs inputfile1.raw inputfile2.raw outputfile.raw Here are some guidelines/helps for implementing Fast Fourier Transform (1) Follow the algorithm as it shows in the slides. (2) Use the complex roots of unity to decide which points to evaluate A(xj ) and B(xj ) at. (3) Your program should use divide and conquer to evaluate the functions A and B in O(n log n) time. (4) For polynomial interpolation define a new function D(x) as shown in class and use your divide and conquer algorithm to evaluate D(x) at the 2nth roots of unity. (5) Once you have D(x) find all the coefficients of C(x) using cs = 1 2πji/k . 2n D(ω2n−s,2n ) where ωj,k is the complex number e (6) One easy way to do the inverse FFT is as follows: (a) take the conjugate of each point of the signal. (b) take the forward FFT. (c) take the pointwise conjugate again. (d) divide each point by n. 5. paper Once you have your program completed, you need to convolve the five pairs of signals found in table 1 which show five different applications of convolution. Complete a lab write-up containing the following.

CS 312: SIGNAL PROCESSING LAB

5

Table 1. Signal and Impulse Response pairs that should be convolved together to pass off the project. Signals Impulse Responses Desired Result BlakesPiano.raw StNicolaesChurch.raw Blake’s piano in St. Nicolaes Church. Classical.raw lowpass.raw Filter out the high pitched sounds. Classical.raw allpass.raw Filter out high and low sounds. (telephone) Classical.raw highpass.raw Filter out the low sounds. NoisySong.raw notch60hz.raw Should remove the buzzing sound. • An estimate of the time it took to complete the lab along with a brief description of anything cool/extra that you implemented in the project. • We described the “reverse-and-slide” method of computing the convolution in the background section. Compare and contrast this algorithm with the FFT algorithm implemented in this project. Email your report, your code, an executable of your program that can be run using the command line description above, and the 5 files containing the convolution of the 5 pairs of convolved signals in .wav format to the TA Paul using blackboard.