《数字调制解调技术的MATLAB与FPGA实现》-matlab ASK编码debug记录

本文是针对"第5.2 ASK调制信号的matlab仿真"的源码的debug记录，因为原书中有一些代码在2016后的matlab版本下有bug，所以这里自己用matlab2020版本跑完后发一下debug记录

1、对rectpulse的理解：

原文链接：https://blog.csdn.net/Explore_OuO/article/details/108947411

y = rectpulse(x, nsamp);

对x应用矩形脉冲整形，即x中的每个符号重复nsamp次，形成输出y。

函数是对列向量进行处理的，把列向量的每个元素在列上重复nsamp次。如果输入有多列，那么就对每列分别处理。
nsamp = 4; % number of samples per symbol
nsymb = 3; % number of symbol
s = RandSteam('mt19937ar', 'Seed', 0); % 只是给随机数一个新的随机方式
ch1 = randi(s, [0 1], nsymb, 1);
ch2 = [1:nsymb]';
x = [ch1 ch2]; % 把两个列向量ch1 ch2拼起来
y = rectpulse(x, nsamp);

原文为CSDN博主「西土城山羊卷」的原创文章，遵循CC 4.0 BY-SA版权协议

2、对rcosflt函数的修改

这里报错是历史版本用的升余弦函数rcosflt已不再使用，改用rcosdesign，rcosflt 使用升余弦滤波器对输入信号进行滤波，rcosdesign函数Matlab推荐用来做升余弦成型滤波，但是两者的函数定义不一致，打"help rcosflt"之后，

Y = rcosflt(X, Fd, Fs, TYPE_FLAG, R, DELAY) filters the input signal X
using a raised cosine FIR filter. The sample frequency for the input, X,
is Fd (Hz). The sample frequency for the output, Y, is Fs (Hz). Fs must be
an integer multiple of Fd. （这句话看出Fd为Fs的倍数，为采样率）

The TYPE_FLAG gives specific filter design
or filtering options. The roll-off factor, R, determines the width of the
transition band of the filter. DELAY is the time delay from the beginning
of the filter to the peak of the impulse response.

R, the roll-off factor specifies the excess bandwidth of the filter. R must
be in the range [0, 1]. For example, R = .5 means that the bandwidth of the
filter is 1.5 times the input sampling frequency, Fd. This also means that
the transition band of the filter extends from .5 * Fd and 1.5 * Fd. Since
R is normalized to the input sampling frequency, Fd, it has no units.
Typical values for R are between 0.2 to 0.5.

如果rcosflt只有三个参数的话，传输参数是X\Fd\Fs：

Y = rcosflt(X, Fd, Fs) filters the input signal X using a raised cosine
filter and the following default parameters
DELAY = 3
R = .5
TYPE_FLAG = ‘fir/normal’

这里我认为初始变量里面
%设置函数的默认参数值
if nargin < 1
Len=1000; %数据长度为1000
IsPlot=0; %不绘图
IsOutput=0; %不将数据写入文本文件中
end

code_rate=1*10^6; %码元速率 Fs=8*code_rate; %采样频率
LenData=Len*Fs/code_rate; %数据长度=8*len
Fc=70*10^6; %载波频率
Qn=8; %量化位数
a=0.8; %成形滤波器滚降因子
码元速率应该是输入采样率Fd，8倍的频率是Fs，在rcosflt的后面两个参数Fd和Fs选取了归一化到Fd的倍数，所以分别为1和8

所以这行代码为Y = rcosflt(X, Fd, Fs)

code_2ask_filter=rcosdl(code_2ask,1,Fs/code_rate);

为了方便理解rcosflt和rcosdesign的应用，我找到了一篇讲得比较好的文章

---------------------以下片段引用自脉冲成型滤波器仿真 - fangying - 博客园 (cnblogs.com)]

例如：设计一个16阶升余弦滤波器，载波频率Fc = 1KHz，滚降系数0.25，采样率为8KHz。

N = 16; 
Fc = 1000; 
R = 0.25; 
Fs = 8000; 
h = firrcos(N, Fc, R, Fs, 'rolloff', 'normal');
figure();plot(h)

如果用rcosdesign函数来设计这个滤波器，那么要用下面的调用方法来实现。

h1 = rcosdesign(R, N/(Fs/Fc/2), Fs/Fc/2,‘normal’);

h1 = h1 / max(h1) / (Fs/Fc/2);　　%重新量化一下系数

figure;plot(h1)

下面重点介绍一下 rcosdesign函数，这个函数数Matlab推荐的用来做升余弦成型滤波的函数。

[语法]

b = rcosdesign(beta,span,sps)

b = rcosdesign(beta,span,sps,shape)

[描述]

b = rcosdesign(beta,span,sps) 函数返回一个滚降系数为beta的均方根升余弦函数。函数被截断为span个符号并且每个符号周期有sps个采样点。滤波器的阶数为 span*sps并且必须为偶数。滤波器的能量为1。。

b = rcosdesign(beta,span,sps,shape)中，最后一个参数shape，当shape设置为’sqrt’时返回均方根升余弦滤波器系数，当shape被设置为’normal’时，返回一个升余弦滤波器的系数。

[例如]

设计一个滚降系数为0.25，符号截断数为6，每符号采样点为4的滤波器。

h = rcosdesign(0.25, 6, 4);

mx = max(abs(h - rcosdesign(0.25, 6, 4)));

fvtool(h,‘Analysis’,‘impluse’);

笔记：引用的文章里面Fc当作Fd来用，但是这里注意到一点就是rcosdesign跟rcosflt在和输入码的交互上不相同

---

下面的程序片段做了一个成型滤波的过程。

设计一个48阶的升余弦FIR滤波器，滚降因子为0.5，发送速率为1被速率，过采样率为8倍过采样。

%% filter: order = 48, rolloff factor: alpha = 0.5, sps = 8
　　% b = firros(n,Fc,df,Fs) Fc:cutoff frequency,df: transmition bandwidth, Fs:oversampling frequency
　　h = firrcos(48, 0.5, 0.5, 8); %Fc = 1/2, h = rcosdesign(0.5, 6, 8);也可以
　　figure(1);
　　plot(h);
　　grid on;
　　xlabel(‘Time’);
　　ylabel(‘Amplitude’);
　　title(’ raised cosine rolloff filter’);

tx_bits1x = randint(100,1)*2 -1; %%100个0-1的数据乘2-1，一列
　　tx_bits8x = zeros(1,800); %%一行
　　tx_bits8x(1:8:end) = tx_bits1x; %%%从1到end，以8为步进

tx_shaped = filter(h,1,tx_bits8x);%%y = filter(b,a,x) 使用由分子和分母系数 b 和 a 定义的有理传递函数 对输入数据 x 进行滤波。如果 a(1) 不等于 1，则 filter 按 a(1) 对滤波器系数进行归一化。因此，a(1) 必须是非零值。

tx_sampled = tx_shaped(1:8:end);

figure(2);
　　stem(tx_bits1x(1:40));
　　title(‘orginal bitstream’);

figure(3);
　　plot(tx_shaped(1:100));
　　title(‘output waveform’);
　　grid on;

figure(4);
　　stem(tx_sampled(1:40));
　　grid on;
　　title(‘sampled output’);

----------------------引用结束线----------------------------------

代码分析：

randint（n,m)产生的是一个n*m维的矩阵，矩阵的元素或者是0或者是1，是随机的。如果想产生一个范围的数，可以设置一个区间，如randint(2,3,[1 6])，就是产生一个2x3随机矩阵，这个矩阵的元素是区间[1 6]的随机数。

所以将代码改成randi，randi(10,3,4) 返回一个由介于 1 和 10 之间的伪随机整数组成的 3×4 数组，或者是X = randi([imin,imax],___) 使用以上任何语法返回一个数组，其中包含从区间 [imin,imax] 的均匀离散分布中得到的整数。

3、升余弦函数的时域和频域响应理解

---------------------以下片段引用自https://www.zhihu.com/question/37222536/answer/360630367，作者：云飞

画出不同滚降系数α升余弦滤波器的时域脉冲响应和频域幅频响应就好理解了。

时域上，容易看出α越大其衰减速度就越快，α=0也就是此时对应于理想矩形窗的sinc脉冲其衰减最慢，α=1时衰减最快。为了减弱定时同步偏差的影响，滚降系数α越大越好。但通过增加滚降系数α加快衰减速度是以牺牲其频谱占用为代价的。

如下图，不同滚降系数的幅频响应，可以看出α=0时频谱占用最小，这时对应的就是奈奎斯特第一定律下的奈奎斯特带宽，随着α增大其频谱占用增加。

---------------------以下片段引用自脉冲成型滤波器仿真 - fangying - 博客园 (cnblogs.com)]

脉冲型滤波器用成型脉冲即数字1用矩形脉冲表示用升余弦脉冲或高斯脉冲表示主要用于基带数据处理。

在数字通信系统中，基带信号进入调制器前，波形是矩形脉冲，突变的上升沿和下降沿包含高频成分较丰富，信号的频谱一般比较宽。从本质上说，脉冲成形就是一种滤波。数字通信系统的信号都必须在一定的频带内，但是基带脉冲信号的频谱是一个Sa函数，在频带上是无限宽的，单个符号的脉冲将会延伸到相邻符号码元内产生码间串扰，这样就会干扰到其他信号，这是不允许的。为了消除干扰，信号在发射之前要进行脉冲成形滤波，把信号的频率约束在带内。因此在信道带宽有限的条件下，要降低误码率，提升信道频带利用率，需要在信号传输前，对其进行脉冲成形处理，改善其频谱特，产生适合信道传输的波形。符号/秒代表单位波特（Baud），波特率是符号信息的比特率。一般的脉冲成型是要过采样的，不然没有意义，因为成型滤波会扩展带宽，过采样是为了减少频谱混叠。

常用的脉冲成型滤波器有RC成型（升余弦）、Gaussian成型等。

----------------------引用结束线----------------------------------

修正原代码中的初始序列矩阵，从2列恢复到1列

%产生随机分布的二进制数据
code_2ask=randi([0,1],Len,1);

修正后的源码为：
clc
clear
Len=10; %数据长度为10
IsPlot=0; %不绘图
IsOutput=0; %不将数据写入文本文件中

code_rate=1*10^6; %码元速率 Fs=8*code_rate; %采样频率
LenData=Len*Fs/code_rate; %数据长度 8*len
Fc=50*10^6; %载波频率 Qn=8; %量化位数 R=0.8; %成形滤波器滚降因子 N=16; %产生载波信号 t=0:1/Fs:Len/code_rate; carrier=cos(2*pi*Fc*t);
carrier=carrier(1:LenData);%取len*8时间长度的载波 figure set (gcf,'Position',[500,100,700,600], 'color','w') subplot(3,3,1) %产生随机分布的二进制数据 code_2ask=randi([0,1],Len,1); stem(code_2ask) title('原始随机码'); %对基带数据以Fs速率采样 code_2ask_upsamp=rectpulse(code_2ask,Fs/code_rate); subplot(3,3,2) stem(code_2ask_upsamp) title('采样后的随机码'); %设计成形滤波器 code_2ask_filter=rcosdesign(R,N/(Fs/code_rate/2),Fs/code_rate/2,'normal'); code_2ask_filter=code_2ask_filter/max(code_2ask_filter)/(Fs/code_rate/2); subplot(3,3,3) plot(code_2ask_filter) title('滤波器自身'); code_2ask_filter = filter(code_2ask_filter,1,code_2ask_upsamp); subplot(3,3,4) plot(code_2ask_filter) title('滤波后'); subplot(3,3,5) plot(carrier) title('载波'); %产生未进行成形滤波的2ASK信号 ASK2=carrier.*code_2ask_upsamp(1:LenData)';
subplot(3,3,6)
plot(ASK2)
title('未进行成形滤波的2ASK信号');
%产生成形滤波后的2ASK信号
ASK2_filter=carrier.*code_2ask_filter(1:LenData)';
subplot(3,3,7)
plot(ASK2_filter)
title('成形滤波后的2ASK信号');
结果产出

原书上

所以，对2ASK的仿真达到书中要求，debug done~~

4、4ASK调制改动

4ASK的改动点在初始序列里面，
原文代码：code_4ask=randint(1,Len,4);
改成：code_4ask=randi([0,3],Len,1);
产出图：

原文中图：

所以，对4ASK的仿真达到书中要求，debug done~~

5、总结

总结：本例程，集中精神一个一个函数去查询攻克，这样从原理、时域、频域三方面都理解清楚了ASK的编码方式、成形滤波器的作用。