[data,Fs,Nbits] = WAVREAD('sample.wav');		%get raw sample from file
%plot(data);
%length(range)
%data = data(1:256);		%this extracts specific window of sample file
DATA1 = fft(data);		%performs fft
range1 = linspace(-Fs/2, Fs/2, length(DATA1));		%sets x-axis range
DATA2 = fft(data, 512);
range2 = linspace(-Fs/2, Fs/2, length(DATA2));
DATA3 = fft(data, 1024);
range3 = linspace(-Fs/2, Fs/2, length(DATA3));

%Magnitude Phase plot
figure(1);
subplot(2,1,1);
plot(range1, fftshift(20*log10(abs(DATA1))));		%plot spectrum in dB and shift result
title(['Frequency Spectrum of sample.wav: ',num2str(length(DATA1)), '-point FFT']);
xlabel('Frequency Content (Hz)');
ylabel('Amplitude (dB)');
subplot(2,1,2);
prange = linspace(-pi,pi,length(DATA1));		%plot phase
plot(prange, angle(DATA1));
title('Phase of sample.wav'); 

%Spectral/Density resolution comparison plot
%Spectral --
figure(2);
subplot(3,1,1);
plot(range1, fftshift(20*log10(abs(DATA1))));
title(['Frequency Spectrum of sample.wav w/ ',num2str(length(DATA1)), '-point FFT performed']);
xlabel('Frequency Content (Hz)');
ylabel('Amplitude (dB)');
grid;
subplot(3,1,2);
plot(range2, fftshift(20*log10(abs(DATA2))));
title(['Frequency Spectrum of sample.wav w/ ',num2str(length(DATA2)), '-point FFT performed']);
xlabel('Frequency Content (Hz)');
ylabel('Amplitude (dB)');
grid;
subplot(3,1,3);
plot(range3, fftshift(20*log10(abs(DATA3))));
title(['Frequency Spectrum of sample.wav w/ ',num2str(length(DATA3)), '-point FFT performed']);
xlabel('Frequency Content (Hz)');
ylabel('Amplitude (dB)');
grid;

%Density -- 
wl1 = 0.25;
wh1 = 0.80;
FDATA1 = fir1(10,[wl1 wh1]);
h1 = freqz(FDATA1,1,512);
f = [0 0.37 0.49 0.6 0.79 1]; a = [0 0 0.2 1 0 0];
b = remez(10,f,a);
[h2,w2] = freqz(b,1,512);
c = butter(10,[wl1 wh1]);
h3 = freqz(c,1,512);
figure(3);
subplot(2,2,1);
plot(linspace(0,1,512),abs(h1));
title('FIR Filter w/ Hamming window');
subplot(2,2,2);
plot(linspace(0,1,512),abs(h2));
title('Chebyshev Filter');
subplot(2,2,3);
plot(linspace(0,1,512),abs(h3));
title('Butterworth Filter');
subplot(2,2,4);



%Comments
%As the number of zeros padded to the signal file are increased, the spectral resolution 
%increases. This can be seen in figure 2 with the fft take on different sample sizes. The 
%plots become more jagged as the sample size increases. As with the spectral resolution, the 
%Density resolution increases in proportion to the number of samples taken.

% Wasnt sure how to filter the samples with the filter funtion. All the funtions are looking 
% for a numberator and denominator. So, I do not have any plots of the Density resolution

%Summary Report
disp(['Total number of samples used in raw wav file: ', num2str(length(data))]);
disp(['Sampling frequency of wav file: ', num2str(Fs),' Hz']);
disp(['DFT performed: ', num2str(length(DATA1)), '-point FFT']);

% MP3 technology has become very popular as the popularity of the internet and digial music have
% exploded over the last few years. The MP3 encoder employs a number of steps to reduce the size
% of a file by an order of magnitude or more. The encoder doesnt just apply a filter blindly to
% the file as was done in this project. Rather, it takes advantage of certain shortcoming of the
% human auditory system. For instance, if there is a loud sound at some frequency in time, humans cannot
% hear any other frequencies during that same time. Also, the human ear's frequency responce is 
% quite a bit smaller then what is normally recorded on a CD. So, these frequencies are filtered
% just as in the project. There are also a number of other shortcomings of the human auditory 
% system that are take advantage of to reduce to size of the file. But, the encoder is more then
% just a filter, it is a data compressor. MP3 files can be expanded to retreive most of the 
% content of the original file.