matlab:qam调制和解调, rcosdesign和upfirdn的理解以及运用

在QAM信号生成中，为了导入VPI软件或者示波器，需要约束输出的波形采样点数，才能正确的显示波形。以下是一些代码和理解。

比如我们需要生成一个QAM的波形，要求如下：

采样频率(带宽):fs = 16GHz

波特率（信号所占带宽）:baudrate = 1GHz

QAM的符号数: n_symbols = 256
每个符号采样次数:n_Sample=fs/baud = 16 
采样间隔:t = 1/fs
时间窗长度:T=(n_qam*n_Sample)*t=256 *16*(1/16e9)=256ns

代码中出现的图片如下：

发射的QAM载波信号频谱图如下

时域图在VPI里显示，这里主要强调时间窗长度和采样点数

rcosdesign的成型滤波器

传输和接收的qam信号星座图如下

可以设置为4QAM,16QAM,64QAM等改一下代码就可以

16QAM

64QAM

matlab代码 ：

pfft.m

function[fscale,fre] = pfft(fs,wave); 
N = length(wave);
s = wave;
Px = (abs(fft(s)).^2/(N^2))./50;
ef = 10*log10(1000*Px);
fre = fftshift(ef);
fscale = fs*(1:N)/N-fs/2;
end

main.m

clear all
close all
clc
%% ==================参数输入=====================
n_symbol = 256;% QAM符号数
n_QAM = 64; % QAM调制进制数
fs = 16e9; % 采样率
fc = 2e9; % 载波频率
baud = 1e9; % 波特率--传输的带宽
desireWavePoint = 8192 % 示波器采样点数
angle = 0; 
%% ==================产生随机序列=====================
bit_per_symbol = log2(n_QAM); % 每个符号比特数
bit_length = n_symbol*bit_per_symbol;
bit_sequence = round(rand(1,bit_length))'; % 列向量
%% ===================bit to symbol 1===================
dataInMatrix = reshape(bit_sequence,length(bit_sequence)/bit_per_symbol,bit_per_symbol); % 重组，使其每行只有4列
Symbols = bi2de(dataInMatrix);%将每行转化为十进制数，[0,15]
symb = qammod(Symbols,n_QAM,'UnitAveragePower',true); % mean(abs(real(symb))+abs(imag(symb)))
X_blocks = reshape(symb,1,length(symb)/1);
x_ifft = ifft(X_blocks); %将QAM的复序列(频域)转换为时域，为生成基带波形做准备
x_s   = x_ifft(:);%串并转换
%% =================升余弦滤波器设置===================
rolloff = 0.5; % 滚降因子
span = 10 ; % 截断长度
sps = fs/baud;  % Samples per symbol
%%  =================基带&载波信号生成===================
% ifft后的序列时域表示为矩形电平脉冲，导致频谱很宽，需要脉冲成型控制传输带宽baud=1GHz
% rcosdesign 参数：span 代表滤波器符号数，sps 代表每个符号采样数
rrcFilter=rcosdesign(rolloff,span,sps,'sqrt'); 
figure(1);stem(rrcFilter);title('系统函数');

% upfirdn参数:对信号X_s进行sps个上采样后与rrcFilter进行卷积，低通滤波
% 带宽为baud。注意：sps=fs/baud暗示了滤波器带宽 baud = 1GHz
waveform_B_0 = upfirdn(x_s, rrcFilter, sps); %升余弦滤波，基带信号

ts=1/fs;%时间间隔
L=length(waveform_B_0);
T_total = L/fs;
n = (0:ts:(L-1)*ts)';
% 载波信号
waveform_C = real(waveform_B_0).*cos(2.*pi.*fc.*(n)) + imag(waveform_B_0).*sin(2.*pi.*fc.*(n));%OFDM射频信号

% 为了满足示波器采样的输入点数 8192，需要重采样
% 注意!! tx_sam 就是传入示波器的信号，他是实信号，注意看tx_sam的类型
tx_sam = resample(waveform_C,desireWavePoint,length(waveform_C))

% 载波信号频谱
[fscale1_temp,fre_temp] = pfft(fs,waveform_C );
figure(4);
plot (fscale1_temp,fre_temp);

%%  =================信号解调===================
% 降采样，恢复waveform_C的采样数
rx_sam = resample(tx_sam,length(waveform_C),length(tx_sam))
receive_w = rx_sam;
ts=1/fs;%时间间隔
% 下变频，生成基带信号
L=length(receive_w);
n = (0:ts:(L-1)*ts)';
I=2*receive_w.*cos(2.*pi.*fc.*(n)+angle);
Q=2*receive_w.*sin(2.*pi.*fc.*(n)+angle);
IQ=I+1j*Q;

% 滤波器降采样
RxSig_0=upfirdn(IQ, rrcFilter, 1, sps);
tx_length = length(x_s)/2
RxSig = RxSig_0(span+1:end-span) %去首尾，保证QAM符号数对应
X_re_blocks = fft(RxSig(:)); % 时域变频域

% 扩展因子，保证我们解调时幅度的正确
factor = mean(abs(real(symb))+abs(imag(symb)))./mean(abs(real(X_re_blocks))+abs(imag(X_re_blocks)));
final_rx = factor*X_re_blocks;

%%  =================星座图解调===================
symb_demod = qamdemod(final_rx,n_QAM,'UnitAveragePower',true);
scatterplot(symb) % 传输星座图
title('Tx')

scatterplot(final_rx) % 接收的星座图
title('Rx')


% ======================EVM===========================
evm = comm.EVM;
evm.ReferenceSignalSource = 'Estimated from reference constellation';
evm.ReferenceConstellation = qammod(0:n_QAM-1,n_QAM,'UnitAveragePower',true);
% test1 = qammod(0:n_QAM-1,n_QAM);
%归一化yMod
% A=abs(mean(syncSignal));
rmsEVM = evm(final_rx);

% ======================误码率1=========================
symb_demod = qamdemod(final_rx,n_QAM,'UnitAveragePower',true);
error_symb = sum(symb_demod~=Symbols)/length(symb_demod);
% =====================误比特率1=========================
bit_demod = de2bi(symb_demod); %将每行转化为十进制数，[0,15]
bit_demod =bit_demod(:);
error_bit = sum(bit_demod~=bit_sequence)/length(bit_demod);
% =====================误比特率2=========================
% bit_demod =  qamdemod(syncSignal,n_QAM,'OutputType','bit','UnitAveragePower',true);

这里需要注意，upfirdn()上采样后，采样点数是比较难设置的。对生成波形waveform_c进行重采样，是为了适应Matlab和VPI联合仿真的输入，VPI要求输入的采样点数为2的N次幂。

在这个过程中可能会有人对rcosdesign() upfirdn()不理解，以下是一个例子

clc;
clear all;
close all;

x=ones(1,20);        
figure(1);
subplot(4,1,1)
stem(x);
title('sig_in');

span = 4
sps = 4
h=rcosdesign(0.3,span,sps,'sqrt');
subplot(4,1,2)
stem(h);
title('h(x)');


y=upfirdn(x,h,4);   
% fvtool(y,'Analysis','impulse');
subplot(4,1,3)
stem(y);title('upsample');


z=upfirdn(y,h,1,4); 
z1 = z(span+1:end-span)
subplot(4,1,4)
stem(z1);title('sig_out');