Contents
Home Work 1 - Low Pass Filtering a Signal 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:07:45
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-1)/2)+1), len_x-((M-1)/2)]; 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-((M-1)/2)); % preassign space for y so it is not growing in a loop. for i = (((M-1)/2)+1):len_x-(M-1)/2 y(i) = 0; for j = (-(M-1)/2):((M-1)/2); 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(((M-1)/2)+1:len_x-(M-1)/2); t_sampled = t((M-1)/2+1:len_x-(M-1)/2); % 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 4 and goes to sample number 98
len_t =
95
len_y =
95
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')
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')
