Contents

Home Work 1 - Low Pass Filtering in Real Time Using Difference Equations

% Programmer: Aaron Klapheck
% Date: 24-Feb-09
clear, clc, home
fprintf('The date and time: %s \n \n', datestr(now))

The date and time: 03-Mar-2009 11:08:02 
 

The Original Signals

t = 0:.01:1; % Sample time
x = 5*sin(10*t) + sin(60*t) + .5*sin(200*t) + .5*sin(300*t); % Signal with noise
f_s = 100; % Hz. Sample rate. Given by (c-a)/b, where t = a:b:c. .
f_nyquist = f_s/2;  % By definition

% give a top signal frequency:
% sample frequency => 2*(signal frequency)

plot(t, x)

The theory of difference equations.

% The average arround a point (eqn 1):
% y[i] = 1/M*(summation from j=-(M-1)/2 to (M-1)/2 of (x[i+j])).

% Note: because future values are being used this method will only work for
% post-processed data, not data being collected in real-time. For the
% formula above M must be an odd number; ergo (M-1). i must be greater than
% (M-1)/2


% The average difference equation (eqn 2):
% y[i] = 1/M*(summation from j=-(M-1) to 0 of (x[i+j])).

% Note: i must be equal to or greater than M. The equation above applies if
% the sample x[i] is given. If the data is sampled in such a way that
% x[i-1] is the larges value given then simply use j=-M and i then must be
% greater than M. This equation can be used on data that is being aquired
% in real time.


% Meaning of i's, j's, and M's:
%   All must be integers.
%   i - The input coordinate that coorospondes to the output y[i].
%   M - The number of points that will be averaged to filter the signal.
%   j - Is responsible for taking the M number of points from the signal
%       that will be averaged.

Using difference equations to filter the signal using equation 1.

% Use equation 1 (post processing):
% y[i] = 1/M*(summation from j=-(M-1)/2 to (M-1)/2 of (x[i+j]));

M = 7; % Numbr of points to average.
len_x = length(x);

A = [M, len_x];
fprintf('\ni, or signal values sampled, starts at sample number %1.0f and goes to sample number %3.0f\n', A)

y = zeros(1,len_x); % preassign space for y so it is not growing in a loop.


for i = 1:len_x-M
    y(i) = 0;
    for j = 0:(M-1)
        y(i) = y(i) + x(i+j);
    end
    y(i) = y(i)/M;
end

% get rid of the fist few data points which are zero.
y = y(1:len_x-M);

t_sampled = t(1:len_x-M);

% Make sure the length of y matches the lenght of t.
len_t = length(t_sampled)
len_y = length(y)

plot(t_sampled, y)

i, or signal values sampled, starts at sample number 7 and goes to sample number 101

len_t =

    94


len_y =

    94

Compairson of Sampled Signal to Filtered Signal

plot(t_sampled, y, t, x), xlabel('Time'), ylabel('Amplitude'), ...
    legend('Filtered Signal','Sampled Signal','Location','North'), ...
    title('Compairson of Sampled Signal to Filtered Signal')

% Notice that using the difference equation based on past data values
% causes the filtered signal to be phase shifted behind the orignial
% signal. This magnitude of this phase shift is related to the number of
% points used to take the average.

Perform an FFT of both the original and the filtered signal

% Decide which exponant of 2 to use to make the fft go faster.
% This involves changing the sample frequecy to a larger number.
% This change of the sample frequency will change the look of the fft plot
% but the peaks which tell us where the natural frequencies are present will
% not change depending on fft_freq. To test this out set fft_freq equal to
% 100 or any other similar value and view the graph before and after.
exp = log(f_s)./log(2);
fft_freq = 2^(ceil(exp));


x_fft = fft(x, fft_freq);
x_fft = abs(x_fft);     % Eliminate immaginary componant
x_fft = x_fft(1:(fft_freq/2)); % Plot from 0 to .5*f_s

frequency_vector = f_nyquist*linspace(0,1,(fft_freq/2));


y_fft = fft(y, fft_freq);
y_fft = abs(y_fft);     % Eliminate immaginary componant
y_fft = y_fft(1:(fft_freq/2));

plot(frequency_vector, x_fft, frequency_vector, y_fft), ...
    xlabel('Frequency'), ylabel('Frequency Content'), ...
    title('Frequency Content in the Original and Filtered Signal'), ...
    legend('Original Signal', 'Filtered Signal')

Published with MATLAB® 7.5