山东大学qpsk monte carlo 仿真实验报告.doc

QPSK通信系统的 Monte Carlo仿真实验报告 2020年2月11日一、实验目的1二、实验原理及框图1三、实验内容2四、程序及其注释21、主函数（不包括Monte Carlo仿真部分）22、Monte Carlo仿真主函数33、测试（7，4）汉明码性能的主函数34、二进制的格雷码信源子函数45、二进制信号变为四进制符号子函数46、星座映射子函数57、产生两路正交高斯噪声信号子函数58、最小距离判决子函数59、最大投影点判决子函数610、四进制信号变为二进制子函数611、 误符号率计算子函数712、通过符号计算误比特率子函数713、星座图子函数714、汉明编码子函数815、汉明解码和纠错子函数816、 从汉明码中提取原二进制码不纠错子函数9五、 仿真结果和图形9六、 结果分析和结论15七、 遇到的问题、解决方法和收获151一、实验目的1、提高独立学习的能力；2、培养发现问题、解决问题和分析问题的能力；3、学习Matlab的使用；4、掌握4PSK通信系统的MonteCarlo仿真方法；5、掌握4PSK通信系统的的组成原理；6、比较编码信号和为编码信号在随机信道中的传输，加深对纠错编码的理解；二、实验原理及框图1.未加信道纠错码的QPSK调制通信系统 （见图（一）高斯随机数发生器均匀随机数发生器检测器+4PSK映射+ 2比特符号高斯随机数发生器比较符号差错计数器比特差错计数器图一 未加信道纠错码的QPSK调制通信系统框图2.信道纠错编码（7,4）汉明码+QPSK调制的通信系统（见图（二）输入数据QPSK解调信道解码（汉明码）输出数据QPSK调制信道编码（汉明码）AWGN信道图二 信道纠错编码（7,4）汉明码+QPSK调制的通信系统框图三、实验内容（1） 未加信道纠错编码的QPSK调制通信系统1） 最大投影点准则进行判决a，计算噪声方差分别为0、0.1、0.5、1.0时的符号差错概率和比特差错概率；b，画出在每种方差时，在检测器输入端1000个接收到的信号加噪声样本（星座图）；c，分别画出数据点为1000,5000,10000,100000,时的Monte Carlo仿真误比特率曲线和理论误比特率曲线，比较区别，分析数据点的数量对仿真结果的影响。2） 将检测器的判决准则改为最小距离法（星座图上符号间的距离），比较与上面结果的区别（2） 信道纠错编码（7,4）汉明码+QPSK调制的通信系统1）比较经过信道纠错编码的QPSK调制系统与未加信道调制编码的系统的传输性能。四、程序及其注释1、主函数（不包括Monte Carlo仿真部分）clear;clc;N=5000;sgm2=0.1; gray_in = source( N); %信源，用格雷码表示s_in=gray2symbol(gray_in);%变为四进制符号 c = mapping( s_in );%星座图映射n = orthogonal_gaussian( N,sgm2 );%产生两路正交的高斯噪声信号cn = awgn_out( c,n);%经过信道s_out_distance = minimum_distance_judge( cn );%最小距离法判决s_out_projection = maximum_projection_judge( cn );%最大投影判决gray_out = s2gray( s_out_projection );%把四进制码译为2进制格雷码ser_distance=SER(s_in,s_out_distance);%最小距离法误码率ser_projection=SER(s_in,s_out_projection);%最大投影法误码率s_difference=SER(s_out_projection,s_out_distance);%两种判决方法结果不同率ber_distance=BER(s_in,s_out_distance);%最小距离法误比特率，由四进制码计算二进制误比特率ber_projection=SER(gray_in,gray_out);%最大投影法误比特率，由二进制码套用误码率的函数计算b_difference=BER(s_out_projection,s_out_distance);%两种方法的不同比特率 constellaion_qpsk(s_in,cn);%星座图2、Monte Carlo仿真主函数clear;clc;N=10000; for i=1:106 SNR(i)=(i-1)/5-10;%信噪比范围 r=10(SNR(i)/10); sgm2(i)=1/(2*r);%通过信噪比计算出噪声的方差 gray_in= source( N); s_in=gray2symbol(gray_in);%变为四进制符号 c = mapping( s_in ); n = orthogonal_gaussian( N,sgm2(i); cn = awgn_out( c,n); s_out_projection = maximum_projection_judge( cn ); ber_projection(i)=BER(s_in,s_out_projection); %仿真的误码率 theory(i)=0.5*erfc(sqrt(r);%理论误码率 end semilogy(SNR,ber_projection,o);%画出log的图像hold on;semilogy(SNR,theory)grid on;xlabel(Eb/N0(dB);ylabel(BER);hold off;3、测试（7，4）汉明码性能的主函数clear;clc;N=10000;HN=7*N/4; for i=1:106 SNR(i)=(i-1)/5-10;%信噪比范围 r=10(SNR(i)/10); sgm2(i)=1/(2*r);%通过信噪比计算出噪声的方差 gray_in = source( N);%信源的格雷码hn = hamming_encode( gray_in );%汉明编码s_in=gray2symbol(hn);%变为四进制符号 c = mapping( s_in );%星座映射n = orthogonal_gaussian( HN,sgm2(i) );%产生两路正交的高斯信号cn = awgn_out( c,n);%经过信道 s_out_distance = minimum_distance_judge( cn );%最小距离法判决 gray_out_74 = s2gray( s_out_distance );%把四进制码译为2进制格雷码gray_out1=hamming_decode( gray_out_74 );%汉明解码和纠错gray_out2=none_hamming_decode( gray_out_74 );%同一个信号没有汉明纠错 hamming_ber(i)=SER(gray_in,gray_out1);%计算误码率none_hamming_ber(i)=SER(gray_in,gray_out2);end semilogy(SNR,hamming_ber,*);%画出两种方法的误码率图像hold on;semilogy(SNR,none_hamming_ber,ro);legend(采用纠错码的BER,未采用纠错码的BER);grid on;xlabel(Eb/N0(dB);ylabel(BER);hold off;4、二进制的格雷码信源子函数function a= source( N) %N为四进制符号长度，a为输出的二进制格雷码a=zeros(1,2*N);x=rand(1,N);for i=1:N if (x(i)=0.25)%根据四种符号产生格雷码 a(2*i-1)=0; a(2*i)=0; elseif (x(i)0.25) a(2*i-1)=0; a(2*i)=1; elseif (x(i)0.5) a(2*i-1)=1; a(2*i)=1; else a(2*i-1)=1; a(2*i)=0; endendend5、二进制信号变为四进制符号子函数function an =gray2symbol(bn)%bn为格雷码，输出an为四进制0到3N=length(bn);an=zeros(1,N/2);for i=1:N/2 if(bn(2*i-1)=0&bn(2*i)=0)%两位二进制码换成一个四进制符号 an(i)=0; else if(bn(2*i-1)=0&bn(2*i)=1) an(i)=1 ; else if(bn(2*i-1)=1&bn(2*i)=1) an(i)=2; else an(i)=3; end end endend6、星座映射子函数function a = mapping( s )%s为四进制信号，a为2*n的矩阵表示星座映射坐标for i=1:length(s) if s(i)=0 a(1,i)=1; a(2,i)=0; elseif s(i)=1 a(1,i)=0; a(2,i)=-1; elseif s(i)=2 a(1,i)=-1; a(2,i)=0; else a(1,i)=0; a(2,i)=1; endEnd end 7、产生两路正交高斯噪声信号子函数function n = orthogonal_gaussian( N,sgm2 ) %N为长度，sgm2为方差 for i = 1:N a=rand; r=sqrt(sgm2)*sqrt(2*log(1/(1-a); b=rand; nc(i)=r*cos(2*pi*b); ns(i)=r*sin(2*pi*b); endn=zeros(2,N);n(1,:)=nc;n(2,:)=ns;end8、最小距离判决子函数function s_out = minimum_distance_judge( an )%an为输入要判决的2*n矩阵，s_out为判决输出的四进制信号N=length(an);s_out=zeros(1,N);d=zeros(4,N);md=zeros(1,N);for i=1:N d(1,i)=(an(1,i)-1)2+an(2,i)2;%计算距离 d(2,i)=an(1,i)2+(an(2,i)+1)2; d(3,i)=(an(1,i)+1)2+an(2,i)2; d(4,i)=an(1,i)2+(an(2,i)-1)2; md(i)=min(d(:,i);%取最小点判决 if (md(i)=d(1,i) s_out(i)=0; end if (md(i)=d(2,i) s_out(i)=1; end if (md(i)=d(3,i) s_out(i)=2; end if (md(i)=d(4,i) s_out(i)=3; endend 9、最大投影点判决子函数function s_out = maximum_projection_judge( an )%an为需要判决的序列，s_out为判决输出四进制信号 N=length(an);s_out=zeros(1,N);d=zeros(4,N);md=zeros(1,N);for i=1:N d(1,i)=an(1,i); d(2,i)=-1*an(2,i); d(3,i)=-1*an(1,i); d(4,i)=an(2,i); md(i)=max(d(1,i),d(2,i),d(3,i),d(4,i); if (md(i)=d(1,i) s_out(i)=0; end if (md(i)=d(2,i) s_out(i)=1; end if (md(i)=d(3,i) s_out(i)=2; end if (md(i)=d(4,i) s_out(i)=3; end end10、四进制信号变为二进制子函数function bn = s2gray( an )%an为四进制信号，bn为译出的格雷码N=length(an);bn=zeros(1,2*N);for i=1:N if(an(i)=0) bn(2*i-1)=0; bn(2*i)=0; else if(an(i)=1) bn(2*i-1)=0; bn(2*i)=1; else if(an(i)=2) bn(2*i-1)=1; bn(2*i)=1; else bn(2*i-1)=1; bn(2*i)=0; end end end End end 11、 误符号率计算子函数（也可输入二进制计算误比特率）function c = SER( an,bn )%an为原符号序列，bn为判决输出的符号序列，c为误符号率N=length(an);cn=(an=bn);cn=double(cn);numwrong=cn*cn;c=numwrong/N;end 12、通过符号计算误比特率子函数function c = BER( an,bn )%an，bn为原符号序列和判决出的序列，c为误比特率N=length(an);cn=abs(an-bn);num=0;for i=1:N if (cn(i)=1)|(cn(i)=3) num=num+1; else if cn(i)=2 num=num+2; end endendc=num/(2*N);end 13、星座图子函数function constellaion_qpsk(an,bn)%an为原信号，bn为接收到的信号figurefor i=1:length(an) if (an(i)=0) plot(bn(1,i),bn(2,i),*); hold on; else if (an(i)=1) plot(bn(1,i),bn(2,i),R*); hold on; else if (an(i)=2) plot(bn(1,i),bn(2,i),Y*); hold on; else if (an(i)=3) plot(bn(1,i),bn(2,i),G*); hold on; end end end end endhold off;end 14、汉明编码子函数function bn = hamming_encode( an )%an为二进制格雷码，bn为74汉明码N=length(an);n=N/4;bn=zeros(1,7*n);for i=1:n bn(7*i-6)=xor(an(4*i),xor(an(4*i-2),an(4*i-3); bn(7*i-5)=xor(an(4*i),xor(an(4*i-1),an(4*i-3); bn(7*i-4)=xor(an(4*i),xor(an(4*i-1),an(4*i-2); bn(7*i-3)=an(4*i-3); bn(7*i-2)=an(4*i-2); bn(7*i-1)=an(4*i-1); bn(7*i)=an(4*i);End end15、汉明解码和纠错子函数function bn =hamming_decode( an )%an为收到的74汉明码，bn为译码输出的二进制码N=length(an)/7;bn=zeros(1,4*N);for i=1:N s1=xor(xor(an(7*i),an(7*i-1),xor(an(7*i-2),an(7*i-4); s2=xor(xor(an(7*i),an(7*i-1),xor(an(7*i-3),an(7*i-5); s3=xor(xor(an(7*i),an(7*i-2),xor(an(7*i-3),an(7*i-6); bn(4*i-3)=an(7*i-3); bn(4*i-2)=an(7*i-2); bn(4*i-1)=an(7*i-1); bn(4*i)=an(7*i); if (s1=0&s2=1&s3=1) bn(4*i-3)=1-an(7*i-3); else if (s1=1&s2=0&s3=1) bn(4*i-2)=1-an(7*i-2); else if (s1=1&s2=1&s3=0) bn(4*i-1)=1-an(7*i-1); else if (s1=1&s2=1&s3=1) bn(4*i)=1-an(7*i); end end end end end End16、 从汉明码中提取原二进制码不纠错子函数为了在同一次仿真过程中比较纠错和不纠错的差别，编写了这一子函数，直接取（7,4）汉明码中的思维二进制码出来。从而可以和纠错的进行误码率比较。function bn =none_hamming_decode( an )%输入an为汉明码，输出为后四位N=length(an)/7;bn=zeros(1,4*N);for i=1:N bn(4*i-3)=an(7*i-3); bn(4*i-2)=an(7*i-2); bn(4*i-1)=an(7*i-1); bn(4*i)=an(7*i);End end5、 仿真结果和图形噪声方差00.10.51.0误率比特BER00.01850.15850.2470误符号率SER00.03700.29200.4220表1 不同噪声方差时的误符号率和误码率，信号长度N=1000噪声方差200.10.51.0最大投影点判决误率比特BER00.01850.15850.2470误符号率SER00.03700.29200.4220最小距离判决误率比特BER00.01850.15850.2470误符号率SER00.03700.29200.4220表2 最大投影点和最小距离法的误码率误比特率比较图1.1 2=0的星座图，信号长度N=1000图1.2 2=0.1的星座图，信号长度N=1000图1.3 2=0.5的星座图，信号长度N=1000图1.4 2=1的星座图，信号长度N=1000图2.1 数据点为1000时的Monte Carlo仿真误比特率曲线和理论曲线图2.2 数据点为5000时的Monte Carlo仿真误比特率曲线和理论曲线图2.3 数据点为10000时的Monte Carlo仿真误比特率曲线和理论曲线图2.4 数据点为100000时的Monte Carlo仿真误比特率曲线和理论曲线图3.1 不同信噪比条件下使用纠错编码和未使用的误比特率比较，信号数目N=1000时图3.2 不同信噪比条件下使用纠错编码和未使用的误比特率比较，信号数目N=10000时6、 结果分析和结论1、 最大投影点准则在不同方差下的误码率和无比特率，观察表1可见 2=0时没有误码，因为没有噪声。观察图1.1可见星座点集中在四个点上，没有误码。 2=0.1时，有了误码，而且误比特率严格为误符号率的一半，说明每个符号都只误到了相邻符号，没有出现一个符号两位都错。星座图如图1.2，点已经开始分散，范围也变大，但还能大概看出边界，有的点越过了边界，产生误码。 2=0.5和1时，误比特率和误符号