N=10000; 
R=3;  
M=16; 
v=sqrt(-1); 
x=randint(1,N,M)+1;  
s_16apsk=[R*exp(v*2*pi*3/24); 
       R*exp(v*2*pi*21/24); 
       R*exp(v*2*pi*9/24); 
       R*exp(v*2*pi*15/24); 
       R*exp(v*2*pi*1/24); 
       R*exp(v*2*pi*23/24); 
       R*exp(v*2*pi*11/24); 
       R*exp(v*2*pi*13/24); 
       R*exp(v*2*pi*5/24); 
       R*exp(v*2*pi*19/24); 
       R*exp(v*2*pi*7/24); 
       R*exp(v*2*pi*17/24); 
       1/3*2.^0.5*R*exp(v*2*pi*3/24); 
       1/3*2.^0.5*R*exp(v*2*pi*21/24); 
       1/3*2.^0.5*R*exp(v*2*pi*9/24); 
       1/3*2.^0.5*R*exp(v*2*pi*15/24)].'; 
   
for n=1:N ;
    x_16apsk(n)=s_16apsk(x(n)); 
end 
figure(1); 
plot(real(x_16apsk),imag(x_16apsk),'.');title('16apsk信號星座圖');grid; 
for i=1:30
    x_16apsknoisy=awgn(x_16apsk,i,'measured'); 
    r_16apsk=abs(x_16apsknoisy);
    q_16apsk=angle(x_16apsknoisy); 
    for n=1:N
        if r_16apsk(n)>=(1/3*2.^0.5+1)/2*R
            if q_16apsk(n)<=4*pi/12&q_16apsk(n)>2*pi/12 
                x_de16apsk(n)=s_16apsk(1); 
            elseif q_16apsk(n)<=-2*pi/12&q_16apsk(n)>-4*pi/12 
                x_de16apsk(n)=s_16apsk(2); 
            elseif q_16apsk(n)<=10*pi/12&q_16apsk(n)>8*pi/12 
                x_de16apsk(n)=s_16apsk(3); 
            elseif q_16apsk(n)<=-8*pi/12&q_16apsk(n)>-10*pi/12 
                x_de16apsk(n)=s_16apsk(4); 
            elseif q_16apsk(n)<=2*pi/12&q_16apsk(n)>0*pi/12 
                x_de16apsk(n)=s_16apsk(5); 
            elseif q_16apsk(n)<=0*pi/12&q_16apsk(n)>-2*pi/12 
                x_de16apsk(n)=s_16apsk(6); 
            elseif q_16apsk(n)<=12*pi/12&q_16apsk(n)>10*pi/12 
                x_de16apsk(n)=s_16apsk(7); 
            elseif q_16apsk(n)<=-10*pi/12&q_16apsk(n)>-12*pi/12 
                x_de16apsk(n)=s_16apsk(8); 
            elseif q_16apsk(n)<=6*pi/12&q_16apsk(n)>4*pi/12 
                x_de16apsk(n)=s_16apsk(9); 
            elseif q_16apsk(n)<=-4*pi/12&q_16apsk(n)>-6*pi/12 
                x_de16apsk(n)=s_16apsk(10); 
            elseif q_16apsk(n)<=8*pi/12&q_16apsk(n)>6*pi/12 
                x_de16apsk(n)=s_16apsk(11); 
            elseif q_16apsk(n)<=-6*pi/12&q_16apsk(n)>-8*pi/12 ;
                x_de16apsk(n)=s_16apsk(12); 
            end
        elseif r_16apsk(n)<(1/3*2.^0.5+1)/2*R
            if q_16apsk(n)<=pi/2&q_16apsk(n)>0*pi/4 
                x_de16apsk(n)=s_16apsk(13); 
            elseif q_16apsk(n)<=0&q_16apsk(n)>-pi/2 
                x_de16apsk(n)=s_16apsk(14); 
            elseif q_16apsk(n)<=pi&q_16apsk(n)>pi/2 
                x_de16apsk(n)=s_16apsk(15); 
            else 
                x_de16apsk(n)=s_16apsk(16); 
            end 
        end 
    end 
    [num,rt]=symerr(x_16apsk,x_de16apsk); 
    ber_16apsk(i)=rt; 
end
n=1:30;
figure(2); 
semilogy(n,ber_16apsk,'*');  
hold on
semilogy(n,ber_16apsk,'-'); 
title('16apsk信噪比與誤碼率的關系曲線'); 
grid;