% function [Ps,Pb] = pb_qam_awgn (EbN0db,M,method,map);
% ------------------------------------------------------------------------------
% EINGABE:
%   EbN0db: Signal/Noise pro Infobit [db] 
%           (Skalar, Spalten- oder Zeilenvektor)
%        M: Stufigkeit der Modulation (default=4)
%           (Skalar)
%   method: Loesungsverfahren, das fuer M-QAM angewendet werden soll
%           'approx' (default)   
%             Approximation ueber die Symbolfehlerwahrscheinlichkeit
%           'int'
%             Zurueckfuehrung auf die Berechnung der PAM
%           (string)
%      map: Art des Mappings
%           'gray' (default)   
%             Gray-Codierung
%           'nat'
%             natuerliches Mapping
%           (string)
%
% AUSGABE:
%   Ps: Symbolfehlerwahrscheinlichkeit (Spaltenvektor)
%   Pb: Bitfehlerwahrscheinlichkeit (Spaltenvektor)
%
% QUELLE: 
%   [Gla83] 
%
% AUTOR: Juergen Rinas,  08.12.1998
% Erweitert um Symbolfehlerraten von Volker Kuehn, 04.01.01
% ------------------------------------------------------------------------------

function [Ps,Pb] = pb_qam_awgn (EbN0db,M,method,map)

if (sscanf(version,'%d')<5) 
  disp(['# Achtung: ',mfilename,' wurde unter Matlab 5.1 entwickelt.']); 
end;
if (nargin<2) M=4; end;             % default: M=4
if (nargin<3) method='approx'; end; % default: method='approx'
if (nargin<4) map='gray'; end;      % default: map='gray'
K=log2(M);
if (K~=fix(K))
  disp(['# Achtung(',mfilename,') M ist keine Potenz von 2.']); 
end;
sqrtM=sqrt(M);
if (sqrtM~=fix(sqrtM))
  disp(['# Achtung(',mfilename,') sqrt(M) ist keine ganze Zahl.']); 
end;

EbN0=10.^(EbN0db(:)./10);
EsN0=K*EbN0;

% exakte Loesung fuer Ps
Ps=2*(sqrt(M)-1)/sqrt(M)*erfc(sqrt(EsN0*3/(2*(M-1)))) ...
   .*(1- (sqrt(M)-1)/(2*sqrt(M))*erfc(sqrt(EsN0*3/(2*(M-1)))) );

if ((upper(method(1))=='A') & (upper(map(1))=='G'))
  % Naeherung fuer grosse EbNo
  Pb=1/K*Ps;
else
  [tmp,Pb]=pb_pam_awgn(EbN0db,sqrt(M),'int',map);
end;

% EOF --------------------------------------------------------------------------