msk的调制解调MATLAB源代码function out = delay(data,n,sample_number)%data:延迟的数据%n:延迟码元个数%sample_number:码元采样个数out = zeros(1,length(data));out(n*sample_number+1:length(data)) =data(1:length(data)-n*sample_number);function [data_diff] = difference(data)%差分编码%************************************************************************* *%data 输入信号%data_diff 差分编码后信号%************************************************************************* *%--------------------------------------------------------------------------data_diff = zeros(1,length(data));data_diff(1) = 1 * data(1); %1为差分编码的初始参考值for i = 2:length(data)data_diff(i) = data_diff(i-1) * data(i);end%************************************************************************* *function [signal_out,I_out,Q_out] =mod_msk(data,data_len,sample_number,Rb)%MSK基带调制%************************************************************************* *% data 调制信号% data_len 码元个数% sample_number 每个码元采样点数% Rb 码元速率% signal_out 基带调制输出% I_out I路输出% Q_out Q路输出%************************************************************************* *% data_len = 10; %码元个数% sample_number = 8; %采样点数% Rb = 16000; %码元速率% data1 = randint(1,data_len);% data = 2*data1-1; %传输的序列Tb = 1/Rb; %码元时间fs = Rb*sample_number; %采样速率%--------------------------------------------------------------------------%差分编码[data_diff] = difference(data);%************************************************************************* *%--------------------------------------------------------------------------%并串转换，延时I(1) = 1; %fai0 = 0,cos(fai0) = 1for i = 1:2:data_lenQ(i) = data_diff(i);Q(i+1) = data_diff(i);endfor i = 2:2:data_lenI(i+1) = data_diff(i);I(i) = data_diff(i);endfor i = 1:sample_numberI1(i:sample_number:data_len*sample_number) = I(1:data_len);Q1(i:sample_number:data_len*sample_number) =Q(1:data_len);end%************************************************************************* *%--------------------------------------------------------------------------%乘加权函数t=1/fs:1/fs:data_len*Tb;I_out = I1 .* cos(pi*t/2/Tb);Q_out = Q1 .* sin(pi*t/2/Tb);%************************************************************************* *%--------------------------------------------------------------------------%调制信号产生signal_out = I_out + j*Q_out;%************************************************************************* *% %--------------------------------------------------------------------------% %画图% subplot(221)% plot(data,'.-');title('MSK传输的数据');xlabel('时间');ylabel('幅度') % subplot(222)% plot(data_diff,'.-');title('差分后的数据');xlabel('时间');ylabel('幅度') % subplot(223)% plot(I1,'.-');title('加权前I路');xlabel('时间');ylabel('幅度');% subplot(224)% plot(Q1,'.-');title('加权前Q路');xlabel('时间');ylabel('幅度');%% figure(2)% subplot(221)% plot(cos(pi*t/2/Tb),'.-');title('加权函数cos(πt/(2Tb))');xlabel('时间');ylabel('幅度')% subplot(222)% plot(sin(pi*t/2/Tb),'.-');title('加权函数sin(πt/(2Tb))');xlabel('时间');ylabel('幅度')% subplot(223)% plot(I_out,'.-');title('加权后I路');xlabel('时间');ylabel('幅度');% subplot(224)% plot(Q_out,'.-');title('加权后Q路');xlabel('时间');ylabel('幅度');% %********************************************************************** ****function [signal_out,I_out,Q_out,phase] =mod_msk2(data,data_len,sample_number,Rb)%MSK基带调制%************************************************************************* *% data 调制信号% data_len 码元个数% sample_number 每个码元采样点数% Rb 码元速率% signal_out 基带调制输出% I_out I路输出% Q_out Q路输出%************************************************************************* *% data_len = 10; %码元个数% sample_number = 8; %采样点数% Rb = 16000; %码元速率% data1 = randint(1,data_len);% data = 2*data1-1; %传输的序列Tb = 1/Rb; %码元时间fs = Rb*sample_number; %采样速率%--------------------------------------------------------------------------%采样for i = 1:sample_numberdata_sample(i:sample_number:data_len*sample_number) = data;end%************************************************************************* *%--------------------------------------------------------------------------%计算相位phase = zeros(1,data_len*sample_number);phase(1) = data_sample(1) * pi/2/sample_number;for i = 2:data_len*sample_numberphase(i) = phase(i-1) + data_sample(i-1) * pi/2/sample_number; end%************************************************************************* *%--------------------------------------------------------------------------%I/QI_out = cos(phase);Q_out = sin(phase);*%--------------------------------------------------------------------------%调制信号产生signal_out = I_out + j*Q_out;%************************************************************************* *%MSK调制，差分解调方法一clear allclose all%--------------------------------------------------------------------------%参数设置data_len = 10000; %码元个数sample_number = 8; %采样个数Rb = 24000; %码元速率fc = 96000; %载波频率*%--------------------------------------------------------------------------%随机产生传输信号data=rand_binary(data_len);%************************************************************************* *%--------------------------------------------------------------------------%MSK基带调制[signal_out,I_out,Q_out] =mod_msk(data,data_len,sample_number,Rb);%************************************************************************* *%--------------------------------------------------------------------------%中频搬移multi = fc/Rb;I_temp=interp(I_out,multi);Q_temp=interp(Q_out,multi);Fs=fc*sample_number;t=1/Fs:1/Fs:length(I_temp)*1/Fs;signal_i=I_temp.*cos(2*pi*fc*t);signal_q=Q_temp.*sin(2*pi*fc*t);signal_mod=I_temp.*cos(2*pi*fc*t)-Q_temp.*sin(2*pi*fc*t);%************************************************************************* *%--------------------------------------------------------------------------%加噪声for SNR = 0:8signal_mod1 = awgn(signal_mod,SNR);%--------------------------------------------------------------------------%去载波N=300; % 滤波器的阶数为(N+1) F=[0,fc-1000,fc+1000,Fs/2]*2/Fs;A=[1,1,0,0];lpf=firls(N,F,A);[amp_lpf,w]=freqz(lpf);I_dem=signal_mod1.*cos(2*pi*fc*t)*2;I_dem=conv(I_dem,lpf);I_dem=I_dem(N/2+1:N/2+length(I_temp));Q_dem=signal_mod1.*sin(2*pi*fc*t)*2;Q_dem=conv(Q_dem,lpf);Q_dem=-Q_dem(N/2+1:N/2+length(I_temp));I_dem_out=zeros(1,length(I_dem)/multi); % 抽取Q_dem_out=zeros(1,length(Q_dem)/multi);for i=1:length(I_dem_out)I_dem_out(i)=I_dem(multi*(i-1)+1);Q_dem_out(i)=Q_dem(multi*(i-1)+1);end;%********************************************************************** ****%--------------------------------------------------------------------------%差分解调demod_data = zeros(1,data_len);demod_data(1) = Q_dem_out(sample_number);for i = 2:data_lendemod_data(i) =Q_dem_out(i*sample_number)*I_dem_out((i-1)*sample_number) -I_dem_out(i*sample_number)*Q_dem_out((i-1)*sample_number);end%********************************************************************** ****%--------------------------------------------------------------------------%判决demod_data = demod_data>0;demod_data = 2*demod_data-1;%********************************************************************** ****%--------------------------------------------------------------------------%计算误码率[num,ber(SNR+1)]=symerr(demod_data,data);%********************************************************************** ****end%************************************************************************* *%--------------------------------------------------------------------------%误码率曲线semilogy([0:8],ber,'r*-');%************************************************************************* *%--------------------------------------------------------------------------%误码率理论值snr = 0:0.1:8;for i = 1:length(snr)snr1(1,i) = 10^(snr(1,i)/10);ps(1,i) = 1/2 * erfc(sqrt(snr1(1,i)));pe(1,i) = 2 * ps(1,i);endhold onsemilogy([0:.1:8],pe);%************************************************************************* *%MSK调制，差分解调方法二clear allclose all%--------------------------------------------------------------------------%参数设置data_len = 10000; %码元个数sample_number = 8; %采样个数Rb = 24000; %码元速率fc = 96000; %载波频率%************************************************************************* *%--------------------------------------------------------------------------%随机产生传输信号data=rand_binary(data_len);%************************************************************************* *%--------------------------------------------------------------------------%MSK基带调制[signal_out,I_out,Q_out] =mod_msk(data,data_len,sample_number,Rb);%************************************************************************* *%--------------------------------------------------------------------------%中频搬移multi = fc/Rb;I_temp=interp(I_out,multi);Q_temp=interp(Q_out,multi);Fs=fc*sample_number;t=1/Fs:1/Fs:length(I_temp)*1/Fs;signal_i=I_temp.*cos(2*pi*fc*t);signal_q=Q_temp.*sin(2*pi*fc*t);signal_mod=I_temp.*cos(2*pi*fc*t)-Q_temp.*sin(2*pi*fc*t);%************************************************************************* *%--------------------------------------------------------------------------%加噪声for SNR = 0:8signal_mod1 = awgn(signal_mod,SNR);%--------------------------------------------------------------------------%去载波N=300; % 滤波器的阶数为(N+1) F=[0,fc-1000,fc+1000,Fs/2]*2/Fs;A=[1,1,0,0];lpf=firls(N,F,A);[amp_lpf,w]=freqz(lpf);I_dem=signal_mod1.*cos(2*pi*fc*t)*2;I_dem=conv(I_dem,lpf);I_dem=I_dem(N/2+1:N/2+length(I_temp));Q_dem=signal_mod1.*sin(2*pi*fc*t)*2;Q_dem=conv(Q_dem,lpf);Q_dem=-Q_dem(N/2+1:N/2+length(I_temp));I_dem_out=zeros(1,length(I_dem)/multi); % 抽取Q_dem_out=zeros(1,length(Q_dem)/multi);for i=1:length(I_dem_out)I_dem_out(i)=I_dem(multi*(i-1)+1);Q_dem_out(i)=Q_dem(multi*(i-1)+1);end;%********************************************************************** ****%--------------------------------------------------------------------------%差分解调demod_data = zeros(1,data_len);demod_data(1) = Q_dem_out(sample_number);for i = 2:2:data_lendemod_data(i) =-I_dem_out(i*sample_number)*Q_dem_out((i-1)*sample_number);endfor i = 3:2:data_lendemod_data(i) =Q_dem_out(i*sample_number)*I_dem_out((i-1)*sample_number);end%********************************************************************** ****%--------------------------------------------------------------------------%判决demod_data = demod_data>0;demod_data = 2*demod_data-1;%********************************************************************** ****%--------------------------------------------------------------------------%计算误码率[num,ber(SNR+1)]=symerr(demod_data,data);%********************************************************************** ****end%************************************************************************* *%--------------------------------------------------------------------------%误码率曲线semilogy([0:8],ber,'r*');%************************************************************************* *%--------------------------------------------------------------------------%误码率理论值snr = 0:0.1:8;for i = 1:length(snr)snr1(1,i) = 10^(snr(1,i)/10);ps(1,i) = 1/2 * erfc(sqrt(snr1(1,i)));pe(1,i) = 2 * ps(1,i);endhold onsemilogy([0:.1:8],pe);%************************************************************************* *%MSK调制，解调clear allclose all%--------------------------------------------------------------------------%参数设置data_len = 30000; %码元个数sample_number = 8; %采样个数Rb = 24000; %码元速率fc = 96000; %载波频率%************************************************************************* *%--------------------------------------------------------------------------%随机产生传输信号data = rand_binary(data_len);%************************************************************************* *%--------------------------------------------------------------------------%MSK基带调制[signal_out,I_out,Q_out] =mod_msk(data,data_len,sample_number,Rb);%************************************************************************* *%--------------------------------------------------------------------------%中频搬移multi = fc/Rb;I_temp=interp(I_out,multi);Q_temp=interp(Q_out,multi);Fs=fc*sample_number;t=1/Fs:1/Fs:length(I_temp)*1/Fs;signal_i=I_temp.*cos(2*pi*fc*t);signal_q=Q_temp.*sin(2*pi*fc*t);signal_mod=I_temp.*cos(2*pi*fc*t)-Q_temp.*sin(2*pi*fc*t);% signal_mod1=real((I_temp+j*Q_temp).*exp(j*2*pi*fc*t));%************************************************************************* *%--------------------------------------------------------------------------%加噪声for SNR = 0:8signal_mod1 = 0.01 * awgn(signal_mod,SNR);%--------------------------------------------------------------------------%去载波N=300; % 滤波器的阶数为(N+1) F=[0,fc-1000,fc+1000,Fs/2]*2/Fs;A=[1,1,0,0];lpf=firls(N,F,A);[amp_lpf,w]=freqz(lpf);I_dem=signal_mod1.*cos(2*pi*fc*t)*2 .* cos(pi*t*Rb/2);I_dem=conv(I_dem,lpf);I_dem=I_dem(N/2+1:N/2+length(I_temp));Q_dem=signal_mod1.*sin(2*pi*fc*t)*2 .* sin(pi*t*Rb/2);Q_dem=conv(Q_dem,lpf);Q_dem=-Q_dem(N/2+1:N/2+length(I_temp));I_dem_out=zeros(1,length(I_dem)/multi); % 抽取Q_dem_out=zeros(1,length(Q_dem)/multi);for i=1:length(I_dem_out)I_dem_out(i)=I_dem(multi*(i-1)+1);Q_dem_out(i)=Q_dem(multi*(i-1)+1);end;%********************************************************************** ****%--------------------------------------------------------------------------%解调，判决demod_data = zeros(1,data_len);for i = 1:data_lenI_dem1(i) = I_dem_out(i*sample_number)>0;Q_dem1(i) = Q_dem_out(i*sample_number)>0;endfor i = 1:2:data_lenm1(i) = Q_dem1(i);endfor i = 2:2:data_lenm1(i) = I_dem1(i);endm1 = 2*m1-1;%********************************************************************** ****%--------------------------------------------------------------------------%差分解码demod_data = zeros(1,data_len);demod_data(1) = 1 * m1(1);for i = 2:data_lendemod_data(i) = m1(i-1) * m1(i);end%********************************************************************** ****%--------------------------------------------------------------------------%计算误码率[num,ber(SNR+1)]=symerr(demod_data,data);%********************************************************************** ****end % /for snr%************************************************************************* *%--------------------------------------------------------------------------%误码率曲线semilogy([0:8],ber,'r*-');%************************************************************************* *%--------------------------------------------------------------------------%误码率理论值snr = 0:0.1:8;for i = 1:length(snr)snr1(1,i) = 10^(snr(1,i)/10);ps(1,i) = 1/2 * erfc(sqrt(snr1(1,i)));pe(1,i) = 2 * ps(1,i);endhold onsemilogy([0:.1:8],pe);%************************************************************************* *function [data_binary,data_binary1]=rand_binary(data_len);%随机产生一个二进制序列作为仿真用的消息序列%************************************************************************* *%data 序列长度%data_binary 产生的二进制序列%************************************************************************* *%--------------------------------------------------------------------------data1=randn(1,data_len);data_binary=zeros(1,data_len);data_binary1=zeros(1,data_len);for i=1:data_lenif data1(i)<0data_binary(i)=-1;data_binary1(i)=0;elsedata_binary(i)=1;data_binary1(i)=1;end end。