matlab实现基于QPSK调制分析4种不同的信道（高斯信道、瑞利信道、莱斯信道、Natagami信道）在不同信噪比下的误码率性...

%% % copyright yh % e-mail: [email protected] % If you have any question, you could contact me by wechat. % 现代通信技术project1仿真 2018年10月30号 % 这个程序仿真了QPSK在大尺度和小尺度衰落情况下的不同SNR下的BER % 程序运行matlab版本：R2014b 运行前请检查安装目录下toolbox是否存在comm通信工具箱 %% %--------------------------系统初始化参数设置 % clc; clear all; snr = -25:2:15; %信噪比 KFactor = 3; % Rician K-factor Ts=1e-4; % sampling time in second Fd=100; % doppler frequency in Hz Tau=[0 1.5e-4 2.5e-4]; % delay for the three paths PdB=[-3, -9, -10]; % power in each of the three paths chan = ricianchan(Ts, Fd, KFactor, Tau, PdB); %瑞利信道参数 chanraylei = rayleighchan(Ts, Fd, Tau, PdB); chan.ResetBeforeFiltering = 0; %----------------------------------- its = 1; %蒙特卡洛试验次数 ser1=[];ser2=[];ser3=[];ser4=[]; errors1=zeros(length(snr),1); errors2=zeros(length(snr),1); errors3=zeros(length(snr),1); errors4=zeros(length(snr),1); Z_in_constell=[];Z_in_constell1=[];Z_in_constell2=[];Z_in_constell3=[]; Z_qd_constell=[];Z_qd_constell1=[];Z_qd_constell2=[];Z_qd_constell3=[]; %星座图初始化 hScope1 = comm.ConstellationDiagram('ReferenceConstellation', [1+1i -1+1i -1-1i 1-1i],... 'XLimits', [-2 2], 'YLimits', [-2 2],... 'Position', figposition([1.5 72 17 20]), ... 'Name','发送信号星座图',... 'ShowTrajectory','true',... 'ColorFading','true');%'SamplesPerSymbol',length(yrician),... 'SymbolsToDisplaySource','Property',... hScope2 = comm.ConstellationDiagram('ReferenceConstellation', [1+1i -1+1i -1-1i 1-1i],...%[0.5+0.5i -0.5+0.5i -0.5-0.5i 0.5-0.5i],... 'XLimits', [-1.5 1.5], 'YLimits', [-1.5 1.5],... 'Position', figposition([1.5 50 17 20]), ... 'Name','高斯信道',... 'ShowTrajectory',false,... 'ColorFading',false);%'SamplesPerSymbol',length(yrician),... 'SymbolsToDisplaySource','Property',... hScope3 = comm.ConstellationDiagram('ReferenceConstellation', [1+1i -1+1i -1-1i 1-1i],...%[0.5+0.5i -0.5+0.5i -0.5-0.5i 0.5-0.5i],... 'XLimits', [-1.5 1.5], 'YLimits', [-1.5 1.5],... 'Position', figposition([1.5 28 17 20]), ... 'Name','莱斯信道',... 'ShowTrajectory',false,... 'ColorFading',false); hScope4 = comm.ConstellationDiagram('ReferenceConstellation', [1+1i -1+1i -1-1i 1-1i],...%[0.5+0.5i -0.5+0.5i -0.5-0.5i 0.5-0.5i],... 'XLimits', [-1.5 1.5], 'YLimits', [-1.5 1.5],... 'Position', figposition([1.5 6 17 20]), ... 'Name','瑞利信道',... 'ShowTrajectory',false,... 'ColorFading',false); hScope5 = comm.ConstellationDiagram('ReferenceConstellation', [1+1i -1+1i -1-1i 1-1i],...%[0.5+0.5i -0.5+0.5i -0.5-0.5i 0.5-0.5i],... 'XLimits', [-1.5 1.5], 'YLimits', [-1.5 1.5],... 'Position', figposition([19 72 17 20]), ... 'Name','nakagami-m信道',... 'ShowTrajectory',false,... 'ColorFading',false); for si=1:length(snr) for j=1:its % data=[0 1 0 1 1 1 0 0 1 1]; % information Number_of_bit=750; data=randi([0 1],Number_of_bit,1); figure(1) stem(data, 'linewidth',3), grid on; title('未调制原始数据流'); axis([ 0 100 0 1.5]); data_NZR=2*data-1; % Data Represented at NZR form for QPSK modulation s_p_data=reshape(data_NZR,2,length(data)/2); % S/P convertion of data br=10^9; %Let us transmission bit rate 1GHz fc=br; % minimum carrier frequency T=1/br; % bit duration t=T/99:T/99:T; % Time vector for one bit information %---------------------------------------- %QPSK调制 %生成随机数据符号流 y=[]; y_in=[]; y_qd=[]; for j=1:length(data)/2 y1=s_p_data(1,j)*cos(2*pi*fc*t); % inphase component y2=s_p_data(2,j)*sin(2*pi*fc*t) ;% Quadrature component y_in=[y_in y1]; % inphase signal vector y_qd=[y_qd y2]; %quadrature signal vector y=[y y1+y2]; % modulated signal vector end Tx_sig=y; % transmitting signal after modulation tt=T/99:T/99:(T*length(data))/2; %---------------------------------------- %画发送信号星座图 for k = 1:1:length(data)/2 y_constell(k) = complex(s_p_data(1,k),s_p_data(2,k)); end y_constell = y_constell.'; step(hScope1, y_constell) % Step scope y_constell = y_constell.'; y_constell = []; %---------------------------------------- figure(2) subplot(3,1,1); plot(tt,y_in,'linewidth',3), grid on; title('QPSK同相调制分量'); xlabel('time(sec)'); ylabel(' amplitude(volt)'); subplot(3,1,2); plot(tt,y_qd,'linewidth',3), grid on; title('QPSK正交调制分量'); xlabel('time(sec)'); ylabel(' amplitude(volt)'); subplot(3,1,3); plot(tt,Tx_sig,'r','linewidth',3), grid on; title('QPSK调制信号 （同相和正交分量的和）'); xlabel('time(sec)'); ylabel(' amplitude(volt)'); %% %------------------------------------------------------- %通过自由空间以及阴影衰落信道 d0=100; sigma=3; distance=1000; Gt=[1 1 0.5]; Gr=[1 0.5 0.5]; Exp=[2 3 6]; %其余的为备选参数 lamda=3e8/fc; PL= -20*log10(lamda/(4*pi*d0))+10*Exp(1)*log10(distance/d0); ydb = PL + sigma*randn(size(distance)); Tx_sig_dBm = 10*log10(abs(1000*Tx_sig)) - ydb; tx_sig_real = Tx_sig.*10.^(-ydb/10); % tx_sig_real = Tx_sig; figure(3) plot(tt,tx_sig_real,'r','linewidth',3), grid on; title('QPSK信号经过对数正态阴影衰落'); xlabel('time(sec)'); ylabel(' amplitude(volt)'); %------------------------------------------------------- %瑞利信道&莱斯信道&natagami-m信道 if (si==1) chan.AvgPathGaindB = chan.AvgPathGaindB + snr(si); chanraylei.AvgPathGaindB = chanraylei.AvgPathGaindB + snr(si); else chan.AvgPathGaindB = chan.AvgPathGaindB + snr(si) - snr(si-1); chanraylei.AvgPathGaindB = chanraylei.AvgPathGaindB + snr(si) - snr(si-1); end % 调整这里的信噪比的单位关系 yrician = awgn(filter(chan,tx_sig_real),snr(si),'measured'); yraylei = awgn(filter(chanraylei,tx_sig_real),snr(si),'measured'); ynakagami = (yrician + yraylei)/2; % yrician = filter(chan,Tx_sig_dBm); % yraylei = filter(chanraylei,Tx_sig_dBm); % yrician = filter(chan,tx_sig_real); % yraylei = filter(chanraylei,tx_sig_real); % ynakagami = (yrician + yraylei)/2; yawgn = awgn(tx_sig_real,snr(si),'measured'); %AWGN信道 figure(4) subplot(4,1,1); plot(tt,yrician,'r','linewidth',3), grid on; title('QPSK信号经过对数正态阴影衰落 + 莱斯信道'); xlabel('time(sec)'); ylabel('amplitude(volt)'); subplot(4,1,2); plot(tt,yraylei,'b','linewidth',3), grid on; title('QPSK信号经过对数正态阴影衰落 + 瑞利信道'); xlabel('time(sec)'); ylabel(' amplitude(volt)'); subplot(4,1,3); plot(tt,ynakagami,'g','linewidth',3), grid on; title('QPSK信号经过对数正态阴影衰落 + nakagami-m信道'); xlabel('time(sec)'); ylabel(' amplitude(volt)'); subplot(4,1,4); plot(tt,yawgn,'b','linewidth',1), grid on; hold on; title('QPSK信号经过对数正态阴影衰落 + AWGN信道'); xlabel('time(sec)'); % axis([0 0.5*10^-7 -6 6]); ylabel(' amplitude(volt)'); %-------------------------------------------------------- %接收端增益放大 yrician = yrician .* 10^(+ydb/10); yraylei = yraylei .* 10^(+ydb/10); ynakagami = ynakagami