communicationsystemof16qam.m

一个16QAM数字调制后的通信系统,包括信号发生器,调制器,信道传输,接收器,均衡器等部分

clear all;
clc;
close all;
%参数初始化
N=10000;        % 二进制码流比特数
T=0.001;        % 16QAM符号长度(s)
Fd=1/T;         % 符号率
Fs=2000;        % 过采样率(Hz)
r=Fs*T;         % 过采样因子，r=Fs*T
alpha=0.25;     % 滚降因子
Fc=500;         % 信道截止频率(Hz)，200 – 500
SNR=15;         % 信道的信噪比

%16QAM星座图
qam16=[-3 -1 1 3];
A=repmat(qam16,4,1);
B=flipud(A');
konst_qam=A+j*B; 
N=10000;
qam16=konst_qam(randint(N,1,16)+1); 
figure(1);
plot(qam16,'*'); 
title('16QAM 星座图');
axis([-4,4,-4,4]); 

%平方根升余弦滤波器
len=length(qam16);                                                            
qam=zeros(size(1:r*len));
qam(1:r:r*len)=qam16;
dt=1/Fs;
t=-5*T:dt:5*T-1/Fs;
t=t+0.00000001;
Tx=rcosine(Fd,Fs,'sqrt',alpha);
[trash,pos] = max(Tx); 
figure(2);
hold on;
plot(Tx);
stem(Tx);
title('平方根升余弦滤波器');
xlabel('Time [s]');
ylabel('Amplitude');
hold off;
grid;

%发送信号的眼图和发送信号星座图
TxOut=filter(Tx,1,qam);
figure(3);
hold on;
for i=pos:r*2:length(qam16)-r*2;
    plot(real(real(TxOut(i:i+r*2))));
end;
title('发送信号的眼图');
figure(4);
plot(TxOut,'*');
title('发送信号星座图');
axis([-4,4,-4,4]);

%原始QAM信号频谱和Tx滤波器频谱的对比
f=-Fs/2:Fs/2;    
fre_qam=freqz(qam,1,f,Fs);   
fre_Tx=freqz(Tx,1,f,Fs);   
fre_TxOut=freqz(TxOut,1,f,Fs);   
figure(5);
subplot(2,1,1);
plot(f,20*log10(real(fre_qam)),f,20*log10(real(fre_Tx)),'g');
title('BLUE：原始QAM信号频谱 GREEN：Tx滤波器频谱');
grid on;           
subplot(2,1,2);
plot(f,20*log10(real(fre_TxOut)));
title('通过Tx滤波器信号的频谱');
grid on;

%3阶巴特沃斯滤波器
[ButterB,ButterA]=butter(3,2*Fc/Fs);                               
figure(6);
freqz(ButterB,ButterA,512); 
title('3阶巴特沃斯滤波器');

%信号经过信道后的频率响应
ButterOut=filter(ButterB,ButterA,TxOut);
AddNoise=awgn(ButterOut,SNR); 
fre_AddNoise=freqz(AddNoise,1,f,Fs);
figure(7);
plot(f,20*log10(real(fre_AddNoise)));
title('信号经过信道后的频率响应');
grid on;

%通过RX滤波器的接收信号的眼图
Rx=Tx;
RxOut=filter(Rx,1,AddNoise);
figure(8);
hold on;
for i=100:r*2:9000-r*2;
    plot(real(RxOut(i:i+r*2)));
end;
title('通过RX滤波器的接收信号的眼图');

%经过匹配滤波和采样的星座图
Impulse=[1 zeros(1,100)];
TxImpOut=filter(Tx,1,Impulse);
ChannelImpOut=filter(ButterB,ButterA,TxImpOut);
RxImpOut=filter(Rx,1,ChannelImpOut);
[Trash,Pos]=max(abs(RxImpOut)); 
RxSampOut=RxOut(Pos:Fs*T:end);

figure(9);
plot(real(RxSampOut),imag(RxSampOut),'*');
title('经过匹配滤波和采样的星座图');

%通过均衡器后的星座图
RxSampOut_x=real(RxSampOut);
qam16_x=real(qam16);
beta=0.01; 
c=zeros(size(1:31)).'; 
for i=1:9000;
   rk=flipud(RxSampOut_x(i:i+30).');                                        
   Ek(i)=qam16_x(i+15)-c.'*rk; 
   c=c+beta*Ek(i)*rk; 
end;

EqOut=filter(c,1,RxSampOut);

figure(10);
plot(real(EqOut),imag(EqOut),'*');
title('通过均衡器后的星座图');

[trash,pos_c] = max(c);
SampOut = EqOut(pos_c:end); 

for i=1:length(SampOut),;
    dist=abs(konst_qam(:)-SampOut(i));
    [y,ind_2(i)]=min(dist); 
end

%结果输出
qam_det=konst_qam(ind_2).';                                                
temp=qam16(1:1:length(qam_det));
Err_num=sum(qam_det~=temp) %误码数
Err_ratio=Err_num/N %误码率