信号课设报告解析.docx
《信号课设报告解析.docx》由会员分享，可在线阅读，更多相关《信号课设报告解析.docx（29页珍藏版）》请在冰豆网上搜索。
信号课设报告解析
课程设计报告
（2013--2014年度第1学期）
名称：
信号分析与处理课程设计
题目：
使用Matlab工具进行数字信号处理
院系：
控制与计算机工程学院
班级：
测控1104
学号：
1111160410
学生姓名：
李伯尧
指导教师：
高峰
设计周数：
一周
成绩：
日期：
2014年1月2日
《信号分析与处理（自）》课程设计
任务书
一、目的与要求
是使学生通过上机使用Matlab工具进行数字信号处理技术的仿真练习，加深对《信号分析与处理（自）》课程所学基本理论和概念的理解，培养学生应用Matlab等工具进行数字信号处理的基本技能和实践能力，为工程应用打下良好基础。
二、主要内容
1．了解Matlab基本使用方法，掌握Matlab数字信号处理的基本编程技术。
掌握数字信号的基本概念。
1）使用Matlab（生成几种典型数字信号（正弦信号、周期信号、高斯随机信号等），2）编程计算离散信号的特征值（均值、方差等）。
3）进行信号加减运算。
2．Matlab编程实现典型离散信号（正弦信号、周期信号）的离散傅立叶变换，显示时域信号和频谱图形（幅值谱和相位谱）；以正弦周期信号为例，观察讨论基本概念（混叠、泄漏、整周期截取、频率分辨率等）。
3．设计任意数字滤波器（分别采用IIR和FIR方法），并对某类型信号进行滤波，并对结果进行显示和分析。
三、进度计划
序号
设计（实验）内容
完成时间
备注
1
布置课程设计任务、示范讲解Matlab程序使用方法
0.5天
2
上机进行课程设计
3天
3
整理课程设计报告
1天
4
演示课程设计内容并答辩
0.5天
四、设计成果要求
1．提交完成设计内容的程序
2．提交设计报告
五、考核方式
课程设计报告、设计内容演示和答辩相结合。
考核内容：
考勤、纪律、课程设计报告、实际编程能力和基本概念掌握程度等。
学生姓名：
李伯尧
指导教师：
高峰
2014年1月2日
一、课程设计的目的与要求
通过上机使用Matlab工具进行数字信号处理技术的仿真练习，加深对《数字信号处理（自）》课程所学基本理论和概念的理解，培养学生应用Matlab等工具进行数字信号处理的基本技能和实践能力，为工程应用打下良好基础。
二、设计正文
1.了解Matlab基本使用方法，掌握Matlab数字信号处理的基本编程技术。
掌握数字信号的基本概念。
1.1使用Matlab（生成几种典型数字信号（正弦信号、周期信号、高斯随机信号等））。
%生成正弦信号
omega=pi/8;%设置数字角频率
ns1=0;
nf1=32;
n1=[ns1:
nf1];
xn1=sin（omega*n1）;
subplot（3,1,1）;
stem（n1,xn1,'.'）;
axis（[0,35,-1.2,1.2]）;
xlabel（'n'）;ylabel（'xn'）;title（'（a）正弦信号'）;grid
%生成周期信号
x=[1234567];%周期为7
xn2=x'*ones（1,3）;%产生3个周期长度的序列
xn2=xn2（:
）;
xn2=xn2';
n2=0:
length（xn2）-1;
subplot（3,1,2）;
stem（n2,xn2,'.'）;
axis（[0,20,-8,8]）;
xlabel（'n'）;ylabel（'xn'）;title（'（b）周期信号'）;grid
%生成高斯随机信号
n3=40;%设置序列长度为40
xn3=randn（1,n3）;%产生均值为0方差为1的高斯随机序列
subplot（3,1,3）;
stem（xn3）;
xlabel（'n'）;ylabel（'xn'）;title（'（c）高斯随机信号'）;grid
1.2编程计算离散信号的特征值（均值、方差等）。
%生成正弦信号
omega=pi/8;%设置数字角频率
ns1=0;
nf1=32;
n1=[ns1:
nf1];
xn1=sin（omega*n1）;
%生成周期信号
x=[1234567];%周期为7
xn2=x'*ones（1,3）;%产生3个周期长度序列
xn2=xn2（:
）;
xn2=xn2';
n2=0:
length（xn2）-1;
%生成高斯随机信号
n3=40;%设置序列长度为40
xn3=randn（1,n3）;%产生均值为0方差为1的高斯随机序列
fprintf（'正弦信号的均值为%.4f方差分别为%.4f\n',mean（xn1）,var（xn1,1））;
fprintf（'周期信号的均值为%.4f方差分别为%.4f\n',mean（xn2）,var（xn2,1））;
fprintf（'高斯随机信号均值为%.4f方差分别为%.4f\n',mean（xn3）,var（xn3,1））;
1.3信号的加减运算
%生成正弦信号
omega=pi/8;
ns1=0;
nf1=32;
n1=[ns1:
nf1];
xn1=sin（omega*n1）;
%生成周期信号
x=[1234567];
xn2=x'*ones（1,3）;
xn2=xn2（:
）;
xn2=xn2';
ns2=0;
nf2=length（xn2）-1;
n2=0:
nf2;
ny=0:
max（nf1,nf2）;%y的向量位置
y1=zeros（1,length（ny））;y2=y1;%y1,y2序列初始化
y1（find（ny<=nf1））=xn1;
y2（find（ny<=nf2））=xn2;
ya=y1+y2;%两序列相加
ys=y1-y2;%两序列相减
%画图
subplot（4,2,1）;stem（n1,xn1,'.'）;
xlabel（'n'）;ylabel（'xn'）;title（'正弦信号'）;grid
subplot（4,2,2）;stem（ny,y1,'.'）;
xlabel（'n'）;ylabel（'xn'）;title（'修正后的正弦信号'）;grid
subplot（4,2,3）;stem（n2,xn2,'.'）;
xlabel（'n'）;ylabel（'xn'）;title（'周期信号'）;grid
subplot（4,2,4）;stem（ny,y2,'.'）;
xlabel（'n'）;ylabel（'xn'）;title（'修正后的周期信号'）;grid
subplot（4,2,6）;stem（ny,ya,'.'）;
xlabel（'n'）;ylabel（'xn'）;title（'序列相加'）;grid
subplot（4,2,8）;stem（ny,ys,'.'）;
xlabel（'n'）;ylabel（'xn'）;title（'序列相减'）;grid
2．Matlab编程实现典型离散信号（正弦信号、周期信号）的离散傅立叶变换，显示时域信号和频谱图形（幅值谱和相位谱）；以正弦周期信号为例，观察讨论基本概念（混叠、泄漏、整周期截取、频率分辨率等）。
2.1Matlab编程实现典型离散信号的离散傅立叶变换。
%生成正弦信号
omega=pi/8;f=omega/（2*pi）;
ns1=0;
nf1=32;
n1=[ns1:
nf1];
xn1=sin（omega*n1）;
xk164=fft（xn1,256）;%计算xn1的256点dft
%以下为正弦信号的绘图
k=0:
255;wk=2*k/256;fk=wk/（2*pi）;%产生256点dft对应采样频率（关于pi的归一化值）
subplot（3,2,1）;stem（n1,xn1,'.'）;
xlabel（'n'）;ylabel（'xn'）;title（'正弦信号'）;grid
subplot（3,2,2）;stem（wk,abs（xk164）,'.'）;%绘制256点dft的幅频特性图
xlabel（'\omega/\pi'）;ylabel（'幅度'）;title（'256点dft的幅频特性图'）;grid
subplot（3,2,4）;stem（fk,abs（xk164）,'.'）;
xlabel（'f'）;ylabel（'幅度'）;title（'256点dft的幅频特性图'）;grid
subplot（3,2,6）;stem（k,abs（xk164）,'.'）;
xlabel（'k'）;ylabel（'幅度'）;title（'256点dft的幅频特性图'）;grid
subplot（3,2,5）;stem（wk,angle（xk164）,'.'）;%绘制256点dft的相频特性图
xlabel（'\omega/\pi'）;ylabel（'相位'）;title（'256点dft的相频特性图'）;grid
%周期信号
x=[1234567];
xn=x'*ones（1,10）;
xn=xn（:
）;
xn=xn';
n=0:
length（xn）-1;
xk=fft（xn,256）;%计算xn的256点dft
%以下为周期信号的绘图
k=0:
255;wk=2*k/256;fk=wk/（2*pi）;%产生256点dft对应采样频率（关于pi的归一化值）
subplot（3,2,1）;stem（n,xn,'.'）;
xlabel（'n'）;ylabel（'xn'）;title（'周期信号'）;grid
subplot（3,2,2）;stem（wk,abs（xk）,'.'）;%绘制256点dft的幅频特性图
xlabel（'\omega/\pi'）;ylabel（'幅度'）;title（'256点dft的幅频特性图'）;grid
subplot（3,2,4）;stem（fk,abs（xk）,'.'）;
xlabel（'f'）;ylabel（'幅度'）;title（'256点dft的幅频特性图'）;grid
subplot（3,2,6）;stem（k,abs（xk）,'.'）;
xlabel（'k'）;ylabel（'幅度'）;title（'256点dft的幅频特性图'）;grid
subplot（3,2,5）;stem（wk,angle（xk）,'.'）;%绘制256点dft的相频特性图
xlabel（'\omega/\pi'）;ylabel（'相位'）;title（'256点dft的相频特性图'）;grid
%生成高斯随机信号
n3=40;
xn3=randn（1,n3）;
xk364=fft（xn3,128）;%计算xn3的128点dft
%以下为高斯随机信号的绘图
k=0:
127;wk=2*k/128;fk=wk/（2*pi）;%产生128点dft对应采样频率（关于pi的归一化值）
subplot（3,2,1）;stem（xn3,'.'）;
xlabel（'n'）;ylabel（'xn'）;title（'高斯随机信号'）;grid
subplot（3,2,2）;stem（wk,abs（xk364）,'.'）;%绘制128点dft的幅频特性图
xlabel（'\omega/\pi'）;ylabel（'幅度'）;title（'128的dft的幅频特性图'）;grid
subplot（3,2,4）;stem（fk,abs（xk364）,'.'）;
xlabel（'f'）;ylabel（'幅度'）;title（'128点dft的幅频特性图'）;grid
subplot（3,2,6）;stem（k,abs（xk364）,'.'）;
xlabel（'k'）;ylabel（'幅度'）;title（'128点dft的幅频特性图'）;grid
subplot（3,2,5）;stem（wk,angle（xk364）,'.'）;%绘制128点dft的相频特性图
xlabel（'\omega/\pi'）;ylabel（'相位'）;title（'128点dft的相频特性图'）;grid
2.2时域混叠的分析
以正弦周期信号为例，一个周期信号长度为7，截取10个周期。
对其进行128点DFT，如图所示（a）、（b），然后抽取其中的64点进行IDFT还原，如图（c）。
因为采样点数N小于序列长度M，故无法由频域采样恢复原信号，在两周期处发生时域混叠。
对128点DFT进行IDFT还原，可以无失真的恢复原信号，如图（d）所示。
上述混叠实质上是频域采样定理所描述的情形。
%生成正弦周期信号
omega=pi/8;
ns1=0;
nf1=32;
n1=[ns1:
nf1];
xn=sin（omega*n1）;
xk=fft（xn,128）;%计算xn1的128点dft
xn1=ifft（xk,128）;
xk1=xk（1:
2:
128）;%抽取64点dft信号
xn2=ifft（xk1,128）;
%以下为正弦周期信号的绘图
k=0:
127;wk=2*k/128;%产生128点dft对应的采样频率（关于pi的归一化值）
subplot（4,1,1）;stem（xn,'.'）;%绘制原始信号
xlabel（'n'）;ylabel（'xn'）;title（'（a）正弦周期信号'）;grid
subplot（4,1,2）;stem（wk,abs（xk）,'.'）;%绘制128点dft的幅频特性图
xlabel（'\omega/\pi'）;ylabel（'幅度'）;title（'（b）128点dft的幅频特性图'）;grid
subplot（4,1,3）;stem（xn1,'.'）;
xlabel（'n'）;ylabel（'xn'）;title（'（c）抽取128点dft还原的信号'）;grid
subplot（4,1,4）;stem（xn2,'.'）;
xlabel（'n'）;ylabel（'xn'）;title（'（d）抽取64点dft还原的信号'）;grid
2.3频谱泄露的分析
正弦周期信号是无限长的，要用DFT分析其频谱，必须对周期信号进行截断。
用矩形窗截断后信号的频谱向两边展宽，使频谱变模糊，谱分辨率降低。
这种现象称为频谱泄露。
一般矩形窗的长度越长，展宽就越窄。
%生成正弦周期信号
omega=pi/8;
ns1=0;
nf1=16;
nf2=64;
fori=1:
2;
switchi;
case1
n1=[ns1:
nf1];
case2
n1=[ns1:
nf2];
end
figure（i）;
xn=sin（omega*n1）;
xk=fft（xn,128）;
%以下为正弦周期信号的绘图
k=0:
127;wk=2*k/128;
subplot（2,1,1）;stem（xn,'.'）;
xlabel（'n'）;ylabel（'xn'）;title（'正弦周期信号'）;grid
subplot（2,1,2）;stem（wk,abs（xk）,'.'）;
xlabel（'\omega/\pi'）;ylabel（'幅度'）;title（'128点幅频特性图'）;grid
end
前两图是矩形窗为0:
16个周期长度截取的信号及其频谱，后两图是矩形窗为0:
64个周期长度截取的信号及其频谱。
观察两者频谱，很容易理解频谱泄露的含义。
2.4整周期截取的分析
%生成正弦周期信号
omega=pi/8;
ns1=0;
nf1=64;
n1=[ns1:
nf1];
xn=sin（omega*n1）;
xk=fft（xn,128）;
xn1=xn（1:
40）;
xk1=fft（xn1,128）;
%以下为正弦周期信号的绘图
k=0:
127;wk=2*k/128;
subplot（4,1,1）;stem（xn,'.'）;%绘制原始信号
xlabel（'n'）;ylabel（'xn'）;title（'（a）整周期截取的正弦周期信号'）;
subplot（4,1,2）;stem（wk,abs（xk）,'.'）;
xlabel（'\omega/\pi'）;ylabel（'幅度'）;title（'（b）整周期截取128点dft的幅频特性图'）;grid
subplot（4,1,3）;stem（xn1,'.'）;%绘制非整周期截取的正弦周期信号
xlabel（'n'）;ylabel（'xn'）;title（'（c）非整周期截取的正弦周期信号'）;grid
subplot（4,1,4）;stem（wk,abs（xk1）,'.'）;
xlabel（'\omega/\pi'）;ylabel（'幅度'）;title（'（d）非整周期截取128点dft的幅频特性图'）;grid
很明显，非正整周期截取时它的最大数字频率不在0~0.2
之间，这与理论值相差是比较大的。
在截取长度不是很长时必须要整周期截取，这样才能保证截断误差不是很大。
2.5频率分辨率的分析
率分辨率，顾名思义，就是将信号中两个靠的很近的频谱分开的能力。
信号x（t）长度为Ts，通过傅氏变换后得到X，其频率分辨率为Δf=1/T（Hz），若经过采样后，假设采样频率为fs=1/Ts，而进行频谱分析时要将这个无穷长的序列使用窗函数截断处理，假设使用矩形窗，我们知道，矩形窗的频谱为sinc函数，主瓣宽度可以定义为2*pi/M，M为窗宽，那么，时域相乘相当于频域卷积，频域内，这一窗函数能够分辨出的最近频率肯定不可能小于2*pi/M了，也就是如果数据长度不能满足2*pi/M<|w2-w1|（w2,w1为两个靠的很近的频率），那么在频谱分析时，频谱上将不能分辨出这两个谱，由于w2-w1=2*pi（f2-f1）/fs=2*pi*Δf/fs也就是2*pi/M<2*piΔf/fs，得到Δf的限制为fs/M，这就是窗函数宽度的最小选择，就是说，根据Shannon采样定理确定了采样频率后，要根据靠的最近的谱峰来确定最小的采样长度，这样，所作出来的频谱才能分辨出那两个谱峰，也就是拥有了相应的频率分辨率。
几个例子：
考虑双正弦信号：
x=sin（2*pi*10*n）+sin（2*pi*9.8*n）;根据Shannon采样定理，采样频率要大于截止频率的两倍，这里选采样频率为80，那么，我们可以看到，Δf为0.2Hz，那么，最小的数据长度为0.2/80=400，但是对正弦信号的频谱分析经验告诉我们，在截断时截断时的数据要包含整周期，并且后面不宜补零以避免频谱泄露（这一点见胡广书《数字信号处理导论》，清华大学出版社），那么，我们要选择至少980个点，才能保含到一个整周期，另外，FFT的经验告诉我们作分析时最好选择2的整数次幂，我们选择靠的最近的1024点。
分析结束。
Fs=80
n=0:
1/Fs:
1023*1/Fs;
x=sin（2*pi*10*n）+sin（2*pi*9.8*n）;
N=length（n）;figure
（1）;
X=fftshift（fft（x））;
plot（（-N/2:
N/2-1）*Fs/N,abs（X）*2/N）;
gridon;
axis（[01501]）;
可以看出这两个谱峰很好的被分辨开来，9.8Hz不在谱线上，所以幅值不为1，以下是一些对比：
例1
Fs=80;
n=0:
1/Fs:
1023*1/Fs;
x=sin（2*pi*10*n）+sin（2*pi*9.8*n）;
N=length（n）;
X=fftshift（fft（x））;
figure
（1）;
subplot（211）
plot（（-N/2:
N/2-1）*Fs/N,abs（X）*2/N）;
gridon;
axis（[01501]）;
title（'窗长1024，包含整周期'）;
例2
n=0:
1/Fs:
979*1/Fs;
x=sin（2*pi*10*n）+sin（2*pi*9.8*n）;
N=length（n）;
X=fftshift（fft（x））;
subplot（212）
plot（（-N/2:
N/2-1）*Fs/N,abs（X）*2/N）;
gridon;
axis（[01501]）;
title（'窗长980，包含整周期'）;
例3
n=0:
1/Fs:
399*1/Fs;
x=sin（2*pi*10*n）+sin（2*pi*9.8*n）;
N=length（n）;
X=fftshift（fft（x））;
figure
（2）;
subplot（211）
plot（（-N/2:
N/2-1）*Fs/N,abs（X）*2/N）;
gridon;
axis（[01501]）;
title（'窗长400，不包含整周期'）;
例4
Fs=20;
n=0:
1/Fs:
1024*1/Fs;
x=sin（2*pi*10*n）+sin（2*pi*9.8*n）;
N=length（n）;
X=fftshift（fft（x））;
subplot（212）
plot（（-N/2:
N/2-1）*Fs/N,abs（X）*2/N）;
gridon;
axis（[01501]）;
title（'窗长1024，采样率过小'）;
3．设计数字滤波器，并对某类型信号进行滤波，并对结果进行显示和分析。
3.1设计IIR数字滤波器
%巴特沃兹滤波器
clc
wp=input（'请输入通带边界频率wp='）;
ws=input（'请输入阻带边界频率ws='）;
rp=input（'请输入通带最大衰减rp='）;
rs=input（'请输入阻带最小衰减rs='）;
[N,wc]=buttord（wp,ws,rp,rs,'s'）%计算相应的模拟低通滤波器阶数N和3dB截止频率
[B,A]=butter（N,wc,'s'）;%计算相应的模拟滤波器的系统函数
tf（B,A）
[Bz,Az]=impinvar（B,A,8000）;%用脉冲响应不变法将模拟滤波器转换成数字滤波器
[H,w]=freqz（Bz,Az）;%计算数字滤波器的频率响应
db=20*log10（abs（H））;
figure
（1）;holdon
plot（w/pi,db）;
grid
xlabel（'\omega/\pi'）;ylabel（'幅度（db）'）;
title（sprintf（'%d阶数字低通滤波器的损耗函数',N））;
%以下滤波过程
x=0:
0.1:
40;
y1=10*sin（pi*x/8）;
figure
（2）;
subplot（2,2,1）;
plot（x,y1）;grid
title（'正弦信号'）
y2=randn（1,length（x））;
subplot（2,2,2）;
plot（x,y2）;grid
title（'高斯随机信号'）
y3=y1+y2;
subplot（2,2,3）;
plot（x,y3）;grid
title（'滤波前信号'）;
y=filter（Bz,Az,y3）;
subplot（2,2,4）;
plot（x,y）;grid
title（'滤波后信号'）;
3.2FIR数字滤波器
t=0:
0.02/32:
0.06;
f1=50;
u1=cos（2*pi*f1*t）+0.8
升级会员 