%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Creative Commons
% Attribution-Noncommercial 2.5 India
% You are free:
% to Share — to copy, distribute and transmit the work
% to Remix — to adapt the work
% Under the following conditions:
% Attribution. You must attribute the work in the manner 
% specified by the author or licensor (but not in any way 
% that suggests that they endorse you or your use of the work). 
% Noncommercial. You may not use this work for commercial purposes. 
% For any reuse or distribution, you must make clear to others the 
% license terms of this work. The best way to do this is with a 
% link to this web page.
% Any of the above conditions can be waived if you get permission 
% from the copyright holder.
% Nothing in this license impairs or restricts the author's moral rights.
% http://creativecommons.org/licenses/by-nc/2.5/in/

% Author	: Krishna
% Email		: krishna@dsplog.com
% Version	: 1.0
% Date		: 24 February 2008
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Script for computing the frequency offset from an OFDM short preamble
% constructed per symbol IEEE 802.11A  specifications. 

clear
fsMHz = 20; % sampling frequency

nFFTSize = 64;
% for each symbol bits a1 to a52 are assigned to subcarrier 
% index [-26 to -1 1 to 26] 
subcarrierIndex = [-26:-1 1:26];
inputFFTShortPreamble = [zeros(1,8) 1+j 0 0 0  -1-j 0 0 0 ... % [-32:-17]
 1+j 0 0 0  -1-j 0 0 0 -1-j 0 0 0   1+j 0 0 0 ...             % [-16:-1]
 0   0 0 0  -1-j 0 0 0 -1-j 0 0 0   1+j 0 0 0 ...             % [0:15]
 1+j 0 0 0   1+j 0 0 0  1+j 0 0 0   0   0 0 0 ];              % [16:31]

inputiFFT = sqrt(13/6)*fftshift(inputFFTShortPreamble); 

% taking ifft
outputiFFT = ifft(inputiFFT,nFFTSize); % generate 64 sample sequence

% concatenating multiple symbols to form 10short preamble
outputShortPreamble = [outputiFFT outputiFFT outputiFFT(1:32)]; 


% introducing frequency offset
fdeltakHz = 200;
outputWithFreqOffset = outputShortPreamble.*exp(j*2*pi*fdeltakHz*1e3*[0:length(outputShortPreamble)-1]/(fsMHz*1e6));

% estimating frequency offset
yt = outputWithFreqOffset;
ytDelayBuffer = zeros(1,0.8*fsMHz);
op = zeros(size(yt));

for ii = 1:length(yt)
   
   op(ii) = conj(yt(ii))*ytDelayBuffer(end);
   
   % shifting samples in the delay buffer
   ytDelayBuffer(2:end) = ytDelayBuffer(1:end-1);
   ytDelayBuffer(1) = yt(ii);

end

fdeltaEstimatekHz = -1*angle(op)/(2*pi*0.8*1e-6)/1000;

close all
plot(fdeltaEstimatekHz,'b','LineWidth',4)
hold on
plot(fdeltakHz*ones(size(yt)),'g','LineWidth',2);
legend('estimated','actual');
grid on
xlabel('sample number,n');
ylabel('frequency Estimate, kHz');
title('frequency offset estimation from 802.11a short preamble');
axis([0 160 0 250])
print('freqEstimate.png','-dpng','-S448,336');