%test measurements
%for lab
%16/01/07

clear;%clear the workspace - place where variables are stored
clc; %clear the command window
%close all; %close all figures



%ANALYSIS PARAMETERS
N = 1024*8;         %DFT length
Fs = 44100;
f0 = 1000; 
% f0 = Fs/N*187;                % -see the difference when 'the puriest possible sinusoid is generated'  

%SINUSOID AMPLITUDE
a = 10^(-18/20);      % - 18dB
% a = 10^(-6/20);     % - 6dB        - you can also try different values
% a = 10^(-0/20);     % - 0dB


%generate the test signal
[x] = sinwave('test.wav',f0,Fs,a,1);    


%/////////////////BLACK BOX
%
%
%RUN THROUGH BLACK BOX
!./black_box test.wav 
%
%
%/////////////////////////////

%read the output file
[y,Fs] = wavread('test_out.wav');

%traspose input file ( it is a vector) - this is not important
x = x';


%///////////////////////////ANALYSIS

%take N samples for analysis
xn = x(1:N);
yn = y(1:N);


%APPLY WINDOWS
%window
% w = rectwin(N);     %rectangular window
w = hann(N);       %von hann window
% w = hamming(N);     %hamming window
% w = blackman(N);    %blackman window


%windowing and normalization
xwin = xn.*w/sum(w);        %input signal
ywin = yn.*w/sum(w);        %output signal    


%PERFROM DFT
Y = fft(ywin);          %size N
X = fft(xwin);          %size N    


Y = 2*Y(1:N/2);         %take non-begative freqs nad scaling
X = 2*X(1:N/2);

%generate frequency vector
f = Fs*(0:N/2-1)/N;


%%%PLOTS
%magnitude vectors against frequency vectors 

%output magnitude spectrum  - linear scale
figure(3);
plot(f,abs(Y),'k');
%output magnitude spectrum  - dB scale
figure(4);
plot(f,20*log10(abs(Y)),'r','LineWidth',2);


%##########################################################################






































% %remove dc
% load filt2;
% 
% K = size(C,1);
% 
% xfilt = G(1)*x;
% 
% for i=1:K
%     xfilt = filter(C(i,1:3),C(i,4:6),xfilt);
%     xfilt = xfilt*G(i+1);
% end;
% 
% 
% 
% 
% %remove 1kHz
% load filt;
% 
% K = size(C,1);
% 
% xfilt = G(1)*xfilt;
% 
% for i=1:K
%     xfilt = filter(C(i,1:3),C(i,4:6),xfilt);
%     xfilt = xfilt*G(i+1);
% end;
% 
% xfn = xfilt(1:N);
% xfwin = xfn.*w/sum(w);
% 
% 
% Xf = fft(xfwin);
% Xf = 2*Xf(1:N/2);
% 
% 
% figure(6);
% plot(f,abs(Y),'r--',f,abs(Xf),'k');









