Using the Fast Fourier Transform (FFT)

to approximate the Fourier Series

Due: Wednesday, January 31st, 5pm, EECE 229

email: pattichis@eece.unm.edu

Consider the periodic triangular-wave function shown below.

Note that the triangular-wave function exhibits quarter-wave symmetry. Furthermore, the function is periodic with period .

Problem 1 (20 points). Compute the Fourier Series of . Make sure to express your answer in the exponential form of the Fourier Series.

Problem 2 (45 points total, CORRECTED). Verifying the Computation of the Fourier Series.

2 (a) (5 points) Use the triangle_wave function provided (look at the end of this assignment) to compute f at 16 points:

period = 4;

amplitude = 1;

no_of_points = 16;

[t, f] = triangle_wave (amplitude, no_of_points, period);

2(b) (10 points) Define a column vector c_column_vector that contains the harmonic coefficients:

% A vector of 16 elements, only c₈ is multiplied by 2.

c_column_vector = [c₀ c₁ c₂ ... c₇ 2c₈ c_-7 c_-6 ... c_-1]; % The C_i come from Problem 1.

c_column_vector = c_column_vector (:);

In Matlab, use 1i to represent j. For example, 3+5j is represented using 3+5*1i.

To generate the sum of harmonics efficiently, without using for loops, use:

% Generate: n*m

N = no_of_points;

new_time = 0:no_of_points-1;

harmonics_num = 0:no_of_points-1;

time_matrix = new_time(:)*harmonics_num;

% generate a matrix of all the harmonics along the columns.

% exp(1i*2*pi*(n*m)/N)

harmonics_matrix = exp((2*pi*1i/N)*time_matrix);

% Matrix multiplication adds up all the harmonic contributions.

f_true_harms = harmonics_matrix * c_column_vector;

% Remove the imaginary part. Imaginary part due to

% numerical approximation.

f_true_harms = real(f_true_harms);

2(c) (5 points) Verify that your answer "looks" correct by plotting the sum of your harmonics against the actual function using:

% PRINT this plot and submit with your assignment.

figure(1);

clf;

plot(t, f, '-', t, f, '*', t, f_true_harms, '-', t, f_true_harms, '+');

2(d) (5 points) Find the location where the maximum absolute error occurs:

% Fix dimensions:

f_true_harms = f_true_harms(:);

f = f(:);

% Find error:

[max_error, error_location] = max(abs(f_true_harms - f));

t_error = t(error_location);

disp(sprintf('Maximum error=%f occurs for t = %f', ...

max_error, t_error));

Submit (via cut and paste) the output after you run this code.

2(e) (5 points) Compute the mean-square error:

% Find mean-absolute error:

disp(sprintf('Mean-square error = %f for %d points.', ...

mean(abs(f_true_harms-f).^2), no_of_points));

Submit (via cut and paste) the output after you run this code.

2(f) (15 points) Redo parts 2(a)-2(d) for no_of_points=32. Your mean-square error should be less.

The only real change is in typing the 32 coefficients in 2(b).

Problem 3 (35 points total) (CORRECTED). Compute the absolute-difference between the true values and the ones obtained by using the fft.

3 (a) (5 points) Compute the fft of the time-domain sequence:

% Do only for no_of_points=16.

f_fft = fft (f);

3(b) (25 points) Carefully derive an approximate relationship for the exponential form of the Fourier series.

Apply this relationship to place the approximate coefficients in an array c_approx as given in 2(b).

% Compute the approximation

correction_factor = ????? % You must compute the correction factor.

c_approx = f_fft*correction_factor;

Hint: You need to take into account the summation approximation to the integral. Type help fft to see what Matlab

computes and compare what Matlab does with the summation that we did in class.

3(c) (5 points) Evaluate the error between the actual and estimated harmonic coefficients:

% Compute the error:

disp('FFT error is:');

disp(transpose(abs(c_approx(:) - c_column_vector(:))));

Submit (via cut and paste) the output after you run this code.

Appendix: triangle_wave.m

function [t, f] = triangle_wave (A, len, period)

%

% [t, f] = triangle_wave (A, len, period)

%

% Given:

% len: the number of points to use for sampling the triangle wave.

% (must be a multiple of 4).

% period: the period of the function.

%

% Computes:

% t: an array of the time values.

% f: an array of the function values.

%

% Example:

% A = 1;

% period = 4;

% len = 16;

% [t, f] = triangle_wave (A, len, period);

% plot(t,f); hold on; plot(t, f, '*');

%

t = linspace(0, period, len+1); % generate all the points.

t = t(1:len); % remove the last point so it does not repeat.

f = A*t(1:len/4+1) / (period/4); % Equation of a line from (0,0) to (T0/4, A).

f = [f f(len/4:-1:1)]; % Mirror the first line to get to (T0/2, 0).

f = [f -f(len/2:-1:2)]; % Mirror and negate the first half.