

% Demo program that shows DF detection using FFT method V1 July 29 2004
% developed by Michael Gaspa and Peng Wie Hsia

% ****************** Part I  **************************
% Compose 5 sine functions with know frequencies in f1, f2, f3, f4 and f5 with
% amplitude of +-A1, A2, A3, A4 and A5 respectively
% The corresponding phases are phi1, phi2, phi3, phi4 and phi5 respectively
% Display the composite wave in Figure 1 (Comp_func)

% pause

% ****************** Part II **************************
% Use FFT to decompose Comp_func so individual frequency is extracted
% Display the frequency spectrum in Figure 2
% where is the peak frequency this is the dominant frequency

% pause

% Program to alter frequencies for different sine waves, and make composite wave.


%---------------------Wave Frequencies-------------------------

f1 = 3;         % Frequency 1
f2 = 6;         % Frequency 2
f3 = 9;         % Frequency 3
f4 = 12;        % Frequency 4
f5 = 15;        % Frequency 5

%-----------------------Phase angles----------------------------
% All phase angles initialized to zero (must be between -pi and pi)

phi1 = 0;           % Phase angle 1 
phi2 = 0;           % Phase angle 2
phi3 = 0;           % Phase angle 3 
phi4 = 0;           % Phase angle 4 
phi5 = 0;           % Phase angle 5

%-----------------------Amplitudes------------------------------
% All phase angles initialized to one
  
A1 = 1;     % Amplitude 1
A2 = 1;     % Amplitude 2    
A3 = 1;     % Amplitude 3
A4 = 1;     % Amplitude 4    
A5 = 1;     % Amplitude 5

%---------------------------------------------------------------

Samp_rate = 1000; % Sampling Rate
N = 0:999; % Vector of number of samples
n = length(N); % number of samples


%---------------------------------------------------------------
% Calculate discrete frequency from sampling rate

d_freq1 = f1 / Samp_rate;
d_freq2 = f2 / Samp_rate;
d_freq3 = f3 / Samp_rate;
d_freq4 = f4 / Samp_rate;
d_freq5 = f5 / Samp_rate;

%---------------------------------------------------------------   
%Calculate and graph function for each frequency and label axes
function1 = A1 * sin(2 * pi * N * d_freq1 + phi1);
subplot(6, 1, 1);
plot(N, function1);
%axis([1 1000 -3 3]);
xlabel(''); 
ylabel([num2str(f1), 'hz' ]);
title('Composite Function Generator', 'FontSize', 13);
function2 = A2 * sin(2 * pi * N * d_freq2 + phi2);
subplot(6, 1, 2);
plot(N, function2);
xlabel('');
ylabel([num2str(f2), 'hz' ]);
function3 = A3 * sin(2 * pi * N * d_freq3 + phi3);
subplot(6, 1, 3);
plot(N, function3);
xlabel('');
ylabel([num2str(f3), 'hz' ]);
function4 = A4 * sin(2 * pi * N * d_freq4 + phi4);
subplot(6, 1, 4);
plot(N, function4);
xlabel('');
ylabel([num2str(f4), 'hz' ]);
function5 = A5 * sin(2 * pi * N * d_freq5 + phi5);
subplot(6, 1, 5);
plot(N, function5);
xlabel('');
ylabel([num2str(f5), 'hz' ]);


%---------------------------------------------------------------
% Sum all functions to make composite function (and plot)

Comp_func = function1 + function2 + function3 + function4 + function5;
subplot(6, 1, 6);
plot(N, Comp_func);
xlabel('Sample #', 'FontSize', 11);
ylabel('Sum');
pause;

%--------------------------------------------------------------
% Perform Fast Fourier Transform on composite function
Fast_four_trans = (abs(fft(Comp_func)))/n;
figure;
subplot(3, 1, 1);                           % Plot FFT on 3 different scaled x-axes
plot(N, Fast_four_trans);
axis([1 1000 0 1]);
title('Frequency Spectrum', 'FontSize', 14);

subplot(3, 1, 2);
plot(N, Fast_four_trans);
axis([1 500 0 1]);

subplot(3, 1, 3);
plot(N, Fast_four_trans,'r:o');
axis([1 50 0 1]);


