pi/4dqpsk的matlab及FPGA仿真

#pi/4dqpsk的matlab及FPGA仿真| 来源: 网络整理| 查看: 265

基于八相偏移调制的原理，本着资料稀少的原则，分享一下自己的经验。

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % PIQPSKMODEM使用清单 % 卷积编码/编码的Ser/Ber~Snr曲线 % 仿真星座图，pi/4_QPSK信号的产生方法基带及调制； % 仿真差分解调 QPSK信号的过程及结果； % 符号速率＝1Mbps； % 基带成形滤波器滚降系数a=0.5； % 采样速率＝8Mbps； % 绘制Hilbert滤波器及普通带通滤波器的频率响应图 % 绘制 pi/4_QPSK信号的频谱及时域波形，绘制解调后的同相正交支路眼图 % 暂缺同步 close all; clear all; clc; ps=1*10^6; %码速率为1MHz a=0.8; %成形滤波器系数 B=(1+a)*ps; %中频信号处理带宽 Fs=8*10^6; %采样速率 fc=2*10^6; %载波频率 N=2000; %仿真数据的长度 t=0:1/Fs:(N*Fs/ps-1)/Fs;%产生长度为N,频率为fs的时间序列 s=randi(4,1,N); %产生随机四进制数据作为原始数据 A=[1i,2^0.5/2+2^0.5/2*1i,1,2^0.5/2-2^0.5/2*1i,-1i,-2^0.5/2-2^0.5/2*1i,-1,-2^0.5/2+2^0.5/2*1i]; Maptable=[0 0;0 1;1 1;1 0]; figure ; scatter(real(A),imag(A),'filled');title('pi/4QPSK星座图'); %将绝对码转化为二进制绝对码 bin_s=zeros(1,2*N); for i=1:1:N if(s(i)==1)bin_s(2*i-1)=0,bin_s(2*i)=0; else if(s(i)==2)bin_s(2*i-1)=0,bin_s(2*i)=1; else if(s(i)==3)bin_s(2*i-1)=1,bin_s(2*i)=1; else if(s(i)==4)bin_s(2*i-1)=1,bin_s(2*i)=0; end end end end end %对二进制绝对码进行卷积 Tre = poly2trellis(7,[133 171]); code_bin_s= convenc(bin_s,Tre); %将卷积码转换为四进制 code_s=zeros(1,2*N); for i=1:2*N if(code_bin_s(2*i-1)==0&&code_bin_s(2*i)==0)code_s(i)=1;end if(code_bin_s(2*i-1)==0&&code_bin_s(2*i)==1)code_s(i)=2;end if(code_bin_s(2*i-1)==1&&code_bin_s(2*i)==1)code_s(i)=3;end if(code_bin_s(2*i-1)==1&&code_bin_s(2*i)==0)code_s(i)=4;end end %将绝对码变换为相对码 xk=ones(1,N); yk=ones(1,N); for i=2:N if s(i)==1 xk(i)=xk(i-1)*cos(pi/4)-yk(i-1)*sin(pi/4); yk(i)=yk(i-1)*cos(pi/4)+xk(i-1)*sin(pi/4); elseif s(i)==2 xk(i)=xk(i-1)*cos(-pi/4)-yk(i-1)*sin(-pi/4); yk(i)=yk(i-1)*cos(-pi/4)+xk(i-1)*sin(-pi/4); elseif s(i)==3 xk(i)=xk(i-1)*cos(3*pi/4)-yk(i-1)*sin(3*pi/4); yk(i)=yk(i-1)*cos(3*pi/4)+xk(i-1)*sin(3*pi/4); elseif s(i)==4 xk(i)=xk(i-1)*cos(-3*pi/4)-yk(i-1)*sin(-3*pi/4); yk(i)=yk(i-1)*cos(-3*pi/4)+xk(i-1)*sin(-3*pi/4); end end %将卷积绝对码转化为卷积相对码 xkk=ones(1,2*N); ykk=ones(1,2*N); for i=2:2*N if code_s(i)==1 xkk(i)=xkk(i-1)*cos(pi/4)-ykk(i-1)*sin(pi/4); ykk(i)=ykk(i-1)*cos(pi/4)+xkk(i-1)*sin(pi/4); elseif code_s(i)==2 xkk(i)=xkk(i-1)*cos(-pi/4)-ykk(i-1)*sin(-pi/4); ykk(i)=ykk(i-1)*cos(-pi/4)+xkk(i-1)*sin(-pi/4); elseif code_s(i)==3 xkk(i)=xkk(i-1)*cos(3*pi/4)-ykk(i-1)*sin(3*pi/4); ykk(i)=ykk(i-1)*cos(3*pi/4)+xkk(i-1)*sin(3*pi/4); elseif code_s(i)==4 xkk(i)=xkk(i-1)*cos(-3*pi/4)-ykk(i-1)*sin(-3*pi/4); ykk(i)=ykk(i-1)*cos(-3*pi/4)+xkk(i-1)*sin(-3*pi/4); end end %对相对码数据以Fs频率采样 Ads_i=upsample(xk,Fs/ps); Ads_q=upsample(yk,Fs/ps); %加噪声 SNR=20; Ads_i=awgn(Ads_i,SNR); Ads_q=awgn(Ads_q,SNR); %设计平方根升余弦滤波器 n_T=[-2 2]; rate=Fs/ps; T=1; Shape_b = rcosfir(a,n_T,rate,T,'sqrt'); %对采样后的数据进行升余弦滤波; rcos_Ads_i=filter(Shape_b,1,Ads_i); rcos_Ads_q=filter(Shape_b,1,Ads_q); %产生同相正交两路载频信号 f0_i=cos(2*pi*fc*t); f0_q=sin(2*pi*fc*t); %产生PI/4_QPSK已调信号 piqpsk=rcos_Ads_i.*f0_i-rcos_Ads_q.*f0_q; dd=(xk+yk*1i)*(2^0.5)/2; Ps=1; Pn=Ps/SNR; z=sqrt(Pn/2)*(randn(size(dd))+1i*randn(size(dd))); y=dd+z; figure;scatter(real(y(1:N/2)),imag(y(1:N/2)),'filled');title('π/4QPSK接收星座图'); %设计Hilbert滤波器及相同阶数的普通带通滤波器 fpm=[0 0.25 1 3 3.75 4]*10^6*2/Fs; %firpm函数的频段向量 magpm=[0 0 1 1 0 0]; %firpm函数的幅值向量 n=30; %滤波器阶数 h_bpf=firpm(n,fpm,magpm,'hilbert滤波器') ;%Hilbert带通滤波器 bpf=firpm(n,fpm,magpm); %普通带通滤波器 figure; %绘制Hilbert滤波器及普通带通滤波器频率响应 freqz(h_bpf); freqz(bpf); %完成对PI/4_QPSK信号的Hilbert滤波及普通滤波 piqpsk_i=filter(bpf,1,piqpsk); piqpsk_q=filter(h_bpf,1,piqpsk); %对普通带通滤波后的数据进行一个符号周期延时处理 piqpsk_di=[zeros(1,Fs/ps),piqpsk_i(1:length(piqpsk_i)-Fs/ps)]; %实现差分解调 demod_mult_i=piqpsk_i.*piqpsk_di; demod_mult_q=piqpsk_q.*piqpsk_di; %对乘法运算后的同相正交支路滤波 demod_i=filter(Shape_b,1,demod_mult_i); demod_q=filter(Shape_b,1,demod_mult_q); %绘制解调后的同相正交支路眼图 eyediagram(demod_i,4*Fs/ps) eyediagram(demod_q,4*Fs/ps) %得到的解码信号demod_i=sing,g为相位，demod_q为cosg %绘制pi4_QPSK信号频谱、pi4_QPSK信号时域波形 figure; m_piqpsk=20*log10(abs(fft(piqpsk,1024)));m_piqpsk=m_piqpsk-max(m_piqpsk); %设置幅频响应的横坐标单位为MHz x_f=[0:(Fs/length(m_piqpsk)):Fs/2];x_f=x_f/10^6; %只显示正频率部分的幅频响应 mpiqpsk=m_piqpsk(1:length(x_f)); %设置时域波表的横坐标单位为us Len=100;%设置时域波形显示的点数 x_t=1:Len;%产生长度为Len的时间序列 x_t=x_t/Fs*10^6; %显示所需的频谱及时域波形 subplot(211); plot(x_f,mpiqpsk); legend('PI/4 QPSK信号频谱'); xlabel('频率(MHz)');ylabel('幅度(dB)');grid on; subplot(212);plot(x_t,piqpsk(101:Len+100)); legend('PI/4 QPSK时域信号波形'); xlabel('时间(us)');ylabel('幅度(V)');grid on; % figure; % subplot(211); % x_t=1:8:1000; % plot(x_t,demod_i(x_t));title("解调后同相基带波形"); % subplot(212); % plot(x_t,demod_q(x_t));title("解调后正交基带波形"); % % %进行信号解调，统计不同信噪比下的误码率与误符号率，并于理论结果对比验证 sum=2; H=cos(pi/4); for i=2:length(dd) B=[ real(dd(i-1))*H-imag(dd(i-1))*H+(imag(dd(i-1))*H+real(dd(i-1))*H)*1i, real(dd(i-1))*H-imag(dd(i-1))*(-H)+(imag(dd(i-1))*H+real(dd(i-1))*(-H))*1i, real(dd(i-1))*(-H)-imag(dd(i-1))*H+(imag(dd(i-1))*(-H)+real(dd(i-1))*H)*1i, real(dd(i-1))*(-H)-imag(dd(i-1))*(-H)+(imag(dd(i-1))*(-H)+real(dd(i-1))*(-H))*1i ]; Dist=abs(dd(i)-B); [~,Index]=min(Dist); %找到距离最近的星座图点 x_est(i)=B(Index); %得到发射调制信号的估计 result(sum)=Index; sum=sum+1; bin_res(2*i-1:2*i)=Maptable(Index,:); %得到比特数据的估计 end figure; Len=N; x_t=1:Len; plot(x_t,bin_res(x_t)); title("PI/4QPSK解调后信号"); Ber1=ones(1,20);Ber2=ones(1,20); num=0; for snr=1:20; %信噪比循环,实际使用请改成-10到10 ，SNR_DB改为snr-21 for mt=1:100 %蒙特卡洛循环,因省时实际使用时设置为1000，num也除以1000 SNR_dB=snr; %不同信噪比下生成AWGN信道 SNR=10^(SNR_dB/10); Ps=1; Pn=Ps/SNR; z=sqrt(Pn/2)*(randn(size(dd))+1i*randn(size(dd))); y=dd+z; for i=2:1:length(y) B=[ real(y(i-1))*H-imag(y(i-1))*H+(imag(y(i-1))*H+real(y(i-1))*H)*1i, real(y(i-1))*H-imag(y(i-1))*(-H)+(imag(y(i-1))*H+real(y(i-1))*(-H))*1i, real(y(i-1))*(-H)-imag(y(i-1))*H+(imag(y(i-1))*(-H)+real(y(i-1))*H)*1i, real(y(i-1))*(-H)-imag(y(i-1))*(-H)+(imag(y(i-1))*(-H)+real(y(i-1))*(-H))*1i ]; Dist=abs(y(i)-B); [~,Index]=min(Dist); %找到距离最近的星座图点 result(i)=Index; bin_res(2*i-1:2*i)=Maptable(Index,:); %得到比特数据的估计 end [tureerr1,per2]=symerr(result,s); %per2单次过程中的误符号率 num=num+per2; %求100次循环误码率总和 end Ber1(snr)=num/100; num=0; %100次误码率的平均值 %theoerr1(snr)=2*qfunc(sqrt(2*SNR)*sin(pi/4)); theoerr1(snr)=1/2*erfc(sqrt(SNR/2)); end Ser1=Ber1*log2(2); %%卷积后的Ber/Ser~SNR曲线 num=0;sum=1; code_dd=xkk+1i*ykk; for snr=1:20 %信噪比循环 for mt=1:100 %蒙特卡洛循环，实际使用时设置为1000，num也除以1000 SNR_dB=snr; %不同信噪比下生成AWGN信道 SNR=10^(SNR_dB/10); Ps=1; Pn=Ps/SNR; z=sqrt(Pn/2)*(randn(size(code_dd))+1i*randn(size(code_dd))); y=code_dd+z; for i=2:1:length(y) B=[ real(y(i-1))*H-imag(y(i-1))*H+(imag(y(i-1))*H+real(y(i-1))*H)*1i, real(y(i-1))*H-imag(y(i-1))*(-H)+(imag(y(i-1))*H+real(y(i-1))*(-H))*1i, real(y(i-1))*(-H)-imag(y(i-1))*H+(imag(y(i-1))*(-H)+real(y(i-1))*H)*1i, real(y(i-1))*(-H)-imag(y(i-1))*(-H)+(imag(y(i-1))*(-H)+real(y(i-1))*(-H))*1i ]; Dist=abs(y(i)-B); [~,Index]=min(Dist); code_bin_res(2*i-1:2*i)=Maptable(Index,:); %得到二进制解码信号 end tblen = 18; code_res = vitdec(code_bin_res,Tre,tblen,'trunc','hard'); [tureerr2,per1]=symerr(code_res,bin_s); %per2单次过程中的误符号率 num=num+per1; %求100次循环误码率总和 end Ber2(snr)=num/100; num=0; theoerr1(snr)=2*qfunc(sqrt(2*SNR)*sin(pi/4)); theoerr2(snr)=1/2*erfc(sqrt(SNR/2)); end Ser2=Ber2*log2(2); S=1:20; %信噪比区间 figure;subplot(211); semilogy(S,Ber1,'r',S,theoerr1*2,'b'); xlabel('SNR');ylabel('Ber');title('编码Ber/SNR的关系');legend('实际','理论');grid on; subplot(212); semilogy(S,Ber2,'g',S,theoerr2*2,'b'); xlabel('SNR');ylabel('Ber');title('卷积编码Ber/SNR的关系'); legend('实际','理论');grid on; figure;subplot(211); semilogy(S,Ser1,'r',S,theoerr1,'b'); xlabel('SNR');ylabel('Ser');title('编码Ser/SNR的关系');legend('实际','理论');grid on; subplot(212); semilogy(S,Ser2,'g',S,theoerr2,'b'); xlabel('SNR');ylabel('Ser');title('卷积编码Ser/SNR的关系'); legend('实际','理论');grid on;

采用中频差分解调，下面为基带差分解调代码demo

close all; clc; ps=4*10^6; %码速率 a=0.8; %成形滤波器系数 B=(1+a)*ps; %中频信号处理带宽 Fs=64*10^6; %采样速率 fc=2.4*10^9; %载波频率 N=2000; %仿真数据的长度 t=0:1/Fs:(N*Fs/ps-1)/Fs;%产生长度为N,频率为fs的时间序列 s=(randi(4,1,N)-1); %产生随机四进制数据作为原始数据 A=[1i,2^0.5/2+2^0.5/2*1i,1,2^0.5/2-2^0.5/2*1i,-1i,-2^0.5/2-2^0.5/2*1i,-1,-2^0.5/2+2^0.5/2*1i]; Maptable=[0 0;0 1;1 1;1 0]; % figure ; % scatter(real(A),imag(A),'filled');title('pi/4QPSK星座图'); %将绝对码变换为相对码 xk=ones(1,N); yk=ones(1,N); for i=2:N if s(i)==3 xk(i)=xk(i-1)*cos(pi/4)-yk(i-1)*sin(pi/4); yk(i)=yk(i-1)*cos(pi/4)+xk(i-1)*sin(pi/4); elseif s(i)==2 xk(i)=xk(i-1)*cos(-pi/4)-yk(i-1)*sin(-pi/4); yk(i)=yk(i-1)*cos(-pi/4)+xk(i-1)*sin(-pi/4); elseif s(i)==1 xk(i)=xk(i-1)*cos(3*pi/4)-yk(i-1)*sin(3*pi/4); yk(i)=yk(i-1)*cos(3*pi/4)+xk(i-1)*sin(3*pi/4); elseif s(i)==0 xk(i)=xk(i-1)*cos(-3*pi/4)-yk(i-1)*sin(-3*pi/4); yk(i)=yk(i-1)*cos(-3*pi/4)+xk(i-1)*sin(-3*pi/4); end end %对相对码数据以Fs频率采样 Ads_i=upsample(xk,Fs/ps); Ads_q=upsample(yk,Fs/ps); %加噪声 SNR=50; Ads_i=awgn(Ads_i,SNR); Ads_q=awgn(Ads_q,SNR); %设计平方根升余弦滤波器 n_T=[-2 2]; rate=Fs/ps; T=1; Shape_b = rcosfir(a,n_T,rate,T,'sqrt'); %对采样后的数据进行升余弦滤波; rcos_Ads_i=filter(Shape_b,1,Ads_i); rcos_Ads_q=filter(Shape_b,1,Ads_q); %-------------匹配滤波+降采样------------ rcos_Ads_i=upfirdn(rcos_Ads_i,Shape_b,1,16); rcos_Ads_q=upfirdn(rcos_Ads_q,Shape_b,1,16); [j,k]=size(rcos_Ads_i); %不滤波直接采样 w=zeros(1,k); z=zeros(1,k); i1=zeros(1,k); q1=zeros(1,k); so=zeros(1,k); w=rcos_Ads_i; z=rcos_Ads_q; for i=2:1:k i0(i)=w(i)*w(i-1)+z(i)*z(i-1); q0(i)=z(i)*w(i-1)-w(i)*z(i-1); if i0(i)>0&&q0(i)>0 so(i)=3; elseif i0(i)0 so(i)=1; elseif i0(i)

【本文地址】

公司简介

联系我们