function[pb,ps]=cm_sm32(snr_in_dB)
%[pb,ps]=cm_sm32(snr_in_dB)
%        CM_SM32求出QPSK以dB為單位的給定信噪比的
%        比特誤碼率和符號誤碼率
N=10000;
E=1;                      %每符號能量
snr=10^(snr_in_dB/10);    %信噪比
sgma=sqrt(E/(3*2*snr));       %噪聲標準方差
%信號映射
s000=[1,0];
s001=[cos(pi/4),sin(pi/4)];
s010=[cos(3*pi/4),sin(3*pi/4)];
s011=[0,1];
s100=[cos(7*pi/4),sin(7*pi/4)];
s101=[0,-1];
s110=[-1,0];
s111=[cos(5*pi/4),sin(5*pi/4)];

%數據源的產生
for i=1:N,
	temp=rand;         %在區間（0，1 )間的一個均勻隨機變量
	if  (temp<0.125),        %信源輸出為“000”的概率為1/8
	    dsource1(i)=0;    
	    dsource2(i)=0;
            dsource3(i)=0;
	elseif (temp<0.25),      %信源輸出為“001”的概率為1/8
	    dsource1(i)=0;
	    dsource2(i)=0;
            dsource3(i)=1;
	elseif (temp<0.375),     %信源輸出為“010”的概率為1/8
	    dsource1(i)=0;
            dsource2(i)=1;
	    dsource3(i)=0;
        elseif (temp<0.5),     %信源輸出為“011”的概率為1/8
	    dsource1(i)=0;
            dsource2(i)=1;
	    dsource3(i)=1;
        elseif (temp<0.625),     %信源輸出為“100”的概率為1/8
	    dsource1(i)=1;
            dsource2(i)=0;
	    dsource3(i)=0;
        elseif (temp<0.75),     %信源輸出為“101”的概率為1/8
	    dsource1(i)=1;
            dsource2(i)=0;
	    dsource3(i)=1;
        elseif (temp<0.875),     %信源輸出為“110”的概率為1/8
	    dsource1(i)=1;
            dsource2(i)=1;
	    dsource3(i)=0;
	else                    %信源輸出為“111”的概率為1/8
    	    dsource1(i)=1; 
	    dsource2(i)=1;
            dsource3(i)=1;
	end;
end;
%判決、誤碼率的計算
numofsymbolerror=0;
numofbiterror=0;
for i=1:N,
    %在判決器處的接收信號，對于第i個符號為：
	n(1)=gngauss(sgma);
	n(2)=gngauss(sgma);
	if ((dsource1(i)==0) & (dsource2(i)==0) & (dsource3(i)==0)),
		r=s000+n;
	elseif ((dsource1(i)==0) & (dsource2(i)==0) & (dsource3(i)==1)),
		r=s001+n;
    elseif ((dsource1(i)==0) & (dsource2(i)==1) & (dsource3(i)==0)),
		r=s010+n;
    elseif ((dsource1(i)==0) & (dsource2(i)==1) & (dsource3(i)==1)),
		r=s011+n;
    elseif ((dsource1(i)==1) & (dsource2(i)==0) & (dsource3(i)==0)),
		r=s100+n;
    elseif ((dsource1(i)==1) & (dsource2(i)==0) & (dsource3(i)==1)),
		r=s101+n;
    elseif ((dsource1(i)==1) & (dsource2(i)==1) & (dsource3(i)==0)),
		r=s110+n;
	else 
		r=s111+n;
	end;
    %以下為計算互相關量度
	c000=dot(r , s000);
	c001=dot(r , s001);	
	c010=dot(r , s010);
	c011=dot(r , s011);
    c100=dot(r , s100);
	c101=dot(r , s101);	
	c110=dot(r , s110);
	c111=dot(r , s111);
    %第i個符號的判決按如下方式進行
	c_max=max([c000 c001 c010 c011 c100 c101 c110 c111]);
	if (c000==c_max)
		decis1=0;decis2=0;decis3=0;
	elseif (c001==c_max)
		decis1=0;decis2=0;decis3=1;
    elseif (c010==c_max)
		decis1=0;decis2=1;decis3=0;
    elseif (c011==c_max)
		decis1=0;decis2=1;decis3=1;
    elseif (c100==c_max)
		decis1=1;decis2=0;decis3=0;
    elseif (c101==c_max)
		decis1=1;decis2=0;decis3=1;
	elseif (c110==c_max)
		decis1=1;decis2=1;decis3=0;
	else
		decis1=1;decis2=1;decis3=1;
	end;
	%如果判決不對，計錯器加1
	symbolerror=0;
	if (decis1~=dsource1(i)),
		numofbiterror=numofbiterror+1;
		symbolerror=1;
	end;
	if (decis2~=dsource2(i)),
		numofbiterror=numofbiterror+1;
		symbolerror=1;
	end;
    if (decis3~=dsource3(i)),
		numofbiterror=numofbiterror+1;
		symbolerror=1;
	end;
	if (symbolerror==1),
		numofsymbolerror=numofsymbolerror+1;
	end;
end;
ps=numofsymbolerror/N;     %因為總共有N個符號
pb=numofbiterror/(3*N);    %因為發送了3N個比特