clc
clear all
close all
M=4; % M=4 for QPSK & M=16 for 16-QAM
bits=log2(M);
%if M==4
E=2;
%elseif M==16
E=10;
%end
FFTsize = 1024;
numSymbols =FFTsize;
CPsize = 128; % Cyclic prefix length
SNR = 0:50;
no_of_frames = 100;
numRuns = 10^7;
Ng=CPsize;
Nt=128;
Ternary=randi([0 2],1,Nt)-1;
pedAchannel = [1 10^(-9.7/20) 10^(-22.8/20)];
pedAchannel = pedAchannel/sqrt(sum(pedAchannel.^2));
vehAchannel = [1 0 10^(-1/20) 0 10^(-9/20) 10^(-10/20) 0 0 0 10^(-15/20) 0 0 0 10^(-20/20)];
vehAchannel = vehAchannel/sqrt(sum(vehAchannel.^2));
idenChannel = 1;
%channel = idenChannel;
channel = pedAchannel;
%channel = vehAchannel;
%equalizerType ='ZERO';
equalizerType ='MMSE';
H_channel = fft(channel,FFTsize);
%**************************************************************************
fs=4069; %samples frequency
T = 1/fs;
t = (0:numSymbols-1)/fs;
F = fs * (0 : numSymbols/2) / numSymbols; %Since the fft result is symmetrical, only the positive half is sufficient for spectral representation
ofdmSymbol=zeros(1,no_of_frames*(numSymbols+Ng));
x=randi([0 M-1],1,no_of_frames*numSymbols);
offset=zeros(1,no_of_frames*(2*numSymbols+Ng));
offset_index=zeros(1,no_of_frames*(2*numSymbols+Ng));
shift_index=zeros(1,no_of_frames*(2*numSymbols+Ng));
for n=1:no_of_frames %sequence of frames
x1=x(1+(n-1)*numSymbols:n*numSymbols); %1 1025 2049...(start of each frame)
inputSymbols=modulate(modem.qammod('M',M,'SymbolOrder','Gray'),x1); %modulated data
TxSamples = sqrt(numSymbols)*ifft(inputSymbols);
TxSamples(end-Ng+1:end-Ng+Nt)=TxSamples(end-Ng+1:end-Ng+Nt)+Ternary; %cyclic prefix containg ternary code
ofdmSymbol(1+(n-1)*(numSymbols+Ng):n*(numSymbols+Ng))=[TxSamples(end-Ng+1:end) TxSamples];
end
%**************************************************************************
%channel
%**************************************************************************
RxSamples = filter(channel,1,ofdmSymbol);
complexNoise = (randn(1,length(ofdmSymbol))+1i*randn(1,length(ofdmSymbol)))/sqrt(2);
for k=1:length(SNR)
% noisePower = 10.^(-SNR(k)/10);
% RxSamples_n = RxSamples + (sqrt(noisePower)*complexNoise);
RxSamples_n = awgn(RxSamples,SNR(k),'measured');
%**************************************************************************
for n=1:no_of_frames
EstSymbols = RxSamples_n(1+(n-1)*(numSymbols+Ng):n*(numSymbols+Ng));
for s=1:Ng
Shift(s)=sum(EstSymbols(s:s+Nt-1).*Ternary) ;
end
w = 2/numSymbols* abs(fft(Shift,numSymbols));
[value shift_index(n,k)]=max(abs(Shift));
recSymbols=EstSymbols(shift_index(n,k)+Ng:shift_index(n,k)+Ng+numSymbols-1); %%%sometimes gives error%%%
Y = fft(recSymbols, FFTsize);
if equalizerType == 'ZERO'
Y = Y./H_channel;
elseif equalizerType == 'MMSE'
C = conj(H_channel)./(conj(H_channel).*H_channel + 10^(-SNR(k)/10));
Y = Y.*C;
end
Y1=sqrt(numSymbols)*ifft(Y,FFTsize);
Y1(end-Ng+1:end-Ng+Nt)=Y1(end-Ng+1:end-Ng+Nt)-Ternary;
Y2=fft(Y1,FFTsize);
f=demodulate(modem.qamdemod('M',M,'SymbolOrder','Gray','OutputType','integer'),Y2);
xr(1+(n-1)*numSymbols:n*numSymbols)=f;
end
d(k)=biterr(x,xr)/(no_of_frames*numSymbols*bits);
end
for n = 1: no_of_frames-2
for theta = ((1+(n-1))*(numSymbols+Ng)):(n+2)*(numSymbols+Ng)-1 %1152:n*3200=2048
p=(conj(RxSamples_n((theta-Ng+1:theta)-numSymbols)));
o=(RxSamples_n(theta-Ng+1:theta));
offset(theta)=(1/2*pi)*atan(sum(imag(o.*p))./sum(real(o.*p)));
[value offset_index]=max(offset);
end
end
figure(1)
plot(F(1:256),w(1:256))
save fftdata SNR d
figure (2)
semilogy(SNR,d)
xlabel('SNR in dB')
ylabel('Bit Error Rate')
grid
评论0