%% define the signal
dt = 1/500; % time step
t = 0:dt:100; % make the time signal
f_sig=5; % define the freq of the signal
sig = sin(2*pi*f_sig*t)+(0.1+0.1)*rand(1,length(t)); % make the signal and add random noise from [a b]; (b-a)*rand(1,length(t))
% sig = zeros(1,length(t)) + 0.00001;
%% Windowing
M = length(t); 
w = hanning(M); 
sig = w'.*sig; 
%% Take the fft
 % pads the signal with 0s till the length of the signal is a power of 2
 % I think this is done for speed
n = 2^nextpow2(length(sig));
FT = fft(sig,n);
f = (0:(n)/2)/((n)*dt); % create the freqs in freq domain
F = 2*abs(FT(1:(n)/2+1)/(n)); % normalizes and takes the magnitude of the FT
F(1)=0; % the first element corosponds to the signal offset. I don't care about that
% Obtaining the magnitude and frequency
[maxpos,m] = max(F);
peak = f(m);
%% plot stuff
subplot(2,1,1)
plot(t,sig,'.b')
hold on
plot(t,sig,'-b')
xlabel('Time')
ylabel('Windowed Signal')

subplot(2,1,2)
semilogy(f,F,'.b')
hold on
semilogy(f,F,'-b')
xlabel('freq')
ylabel('Amplitude fft')
title(['Peak freq. is at ',num2str(peak),' Actual freq. is ',num2str(f_sig)])

disp(' ')
disp(['Peak freq. is at ',num2str(peak)])
disp(['Actual freq. is ',num2str(f_sig)])
disp(['The error is ',num2str(100*(1-peak/f_sig)),' %'])
disp(' ')