%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 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/
% Script for simulating 16-PSK transmission and reception and compare the
% simulated and theoretical symbol error probability
% Checked for proper operation with Octave Version 3.0.0
% Author : Krishna
% Email : krishna@dsplog.com
% Version : 1.0
% Date : 16 February 2007
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% symbol error rate for 16-PSK modulation
clear
N = 2*10^5; % number of symbols
M = 16;
thetaMpsk = [0:M-1]*2*pi/M; % reference phase values
Es_N0_dB = [0:25]; % multiple Es/N0 values
ipPhaseHat = zeros(1,N);
for ii = 1:length(Es_N0_dB)
% symbol generation
% ------------------
ipPhase = randsrc(1,N,thetaMpsk);
ip = exp(j*ipPhase);
s = ip; % normalization of energy to 1
% noise
% -----
n = 1/sqrt(2)*[randn(1,N) + j*randn(1,N)]; % white guassian noise, 0dB variance
y = s + 10^(-Es_N0_dB(ii)/20)*n; % additive white gaussian noise
% demodulation
% ------------
% finding the phase from [-pi to +pi]
opPhase = angle(y);
% unwrapping the phase i.e. phase less than 0 are
% added 2pi
opPhase(find(opPhase<0)) = opPhase(find(opPhase<0)) + 2*pi;
% rounding the received phase to the closest
% constellation
ipPhaseHat = 2*pi/M*round(opPhase/(2*pi/M)) ;
% as there is phase ambiguity for phase = 0 and 2*pi,
% changing all phases reported as 2*pi to 0.
% this is to enable comparison with the transmitted phase
ipPhaseHat(find(ipPhaseHat==2*pi)) = 0;
% counting errors
nErr(ii) = size(find([ipPhase- ipPhaseHat]),2); % couting the number of errors
end
simBer = nErr/N;
theoryBer = erfc(sqrt(10.^(Es_N0_dB/10))*sin(pi/M));
close all
figure
semilogy(Es_N0_dB,theoryBer,'bs-','LineWidth',2);
hold on
semilogy(Es_N0_dB,simBer,'mx-','LineWidth',2);
axis([0 25 10^-5 1])
grid on
legend('theory', 'simulation');
xlabel('Es/No, dB')
ylabel('Symbol Error Rate')
title('Symbol error probability curve for 16-PSK modulation')
