数字信号处理课程设计报告.docx
《数字信号处理课程设计报告.docx》由会员分享，可在线阅读，更多相关《数字信号处理课程设计报告.docx（17页珍藏版）》请在冰豆网上搜索。
数字信号处理课程设计报告
山东工商学院
课程设计报告
课程名称：
数字信号处理A
班级：
电子131
姓名:
XXXX
学号：
XXXXXX
指导教师：
XXXXX
时间：
2016年1月1日
一、课程设计题目
题目1：
（1）、已知Xa（t）=e^-1000|t|,求其傅立叶变换Xa（jΩ），画出模拟信号及其连续时间傅里叶变换的曲线图。
（2）、以Xa（t）为例，说明采样频率对频率响应的影响，分别采用fs=1000Hz和fs=5000Hz，绘出X（e^jw）曲线。
（1）
代码:
closeall
clear;clc;
W=10;f=1000;n=-10:
W-1;t=n/f;
X=exp（-1000*abs（t））;
subplot（1,2,1）;plot（t,X）;%画模拟信号曲线
xlabel（'t/s'）;ylabel（'xa（n）'）;
title（'模拟信号'）;%标题模拟信号
tf=10;N=100;dt=10/N;t=（1:
N）*dt;
wf=25;Nf=50;
w1=linspace（0,wf,Nf）;%0-25之间分成50点
dw=wf/（Nf-1）;W1=-50:
50;
Xat=exp（-1000*abs（t））;%表达式
F1=Xat*exp（-1i*t'*w1）*dt;%傅立叶变换
w=[-fliplr（w1）,w1（2:
Nf）];%负频率的频谱
Y1=（exp
（2）-1）./（exp
（2）-exp（1-1i*W1）-exp（1+1i*W1）+1）;
F=[fliplr（F1）,F1（2:
Nf）];t=[-fliplr（t）,t];
subplot（1,2,2）;
plot（w,F,'linewidth',1）;%画傅立叶变换曲线
xlabel（'w/pi'）;ylabel（'Xa（jΩ）'）;
title（'傅里叶变换'）;%标题傅立叶变换
结果：
分析：
模拟信号在[-0.01，0.01]区间为连续信号,其傅立叶变换曲线在[-10，10]内为连续曲线。
（2）
代码：
closeall
clear
clc
Dt=0.00005;%步长为0.00005s
t=-0.005:
Dt:
0.005;
xa=exp（-1000*abs（t））;%取时间从-0.005s到0.005s这段模拟信号
Ts1=0.001;Ts2=0.0002;%周期
n=-25:
1:
25;
x1=exp（-1000*abs（n*Ts1））;
x2=exp（-1000*abs（n*Ts2））;
K=100;k=0:
1:
K;w=pi*k/K;%求模拟角频率
X1=x1*exp（-j*n'*w）;%求其傅立叶变换
X2=x2*exp（-j*n'*w）;%求其傅立叶变换
X11=real（X1）;X12=real（X2）;
w=[-fliplr（w）,w（2:
101）];%将角频率范围扩展为从-到+
X11=[fliplr（X11）,X11（2:
101）];
X12=[fliplr（X12）,X12（2:
101）];
subplot（2,1,1）;plot（w/pi,X11）;%画出fs=1000Hz的频率响应
xlabel（'w/pi'）;ylabel（'X1（jw）'）;
title（'fs=1000Hz的DTFT'）;%标题fs=1000Hz的DTFT
subplot（2,1,2）;plot（w/pi,X12）;%画出fs=5000Hz的频率响应
xlabel（'w/pi'）;ylabel（'X2（jw）'）;
title（'fs=5000Hz的DTFT'）;%标题fs=5000Hz的DTFT
结果：
分析：
当采样频率越大的时候，采样信号频谱越陡峭，而其失真情况也越来越小。
题目2：
已知时域信号x（n）=cos（0.48πn）+cos（0.52πn），求下面5种情况的X（ejω）和X（k）。
（1）取x（n）的前10点数据，求N=10点的X（ejω）和X（k），并作图。
（2）将
（1）中的x（n）补零至100点，求N=100点的X（ejω）和X（k），并作图。
（3）取x（n）的前100点数据，求N=100点的X（ejω）和X（k），并作图。
（4）取x（n）的前128点数据，求N=128点的X（ejω）和X（k），并作图。
（5）取x（n）的前50点数据，求N=50点的X（ejω）和X（k），并作图。
讨论以上5种情况的区别。
（1）
代码：
closeall
clear
clc
n=（0:
1:
9）;
y=cos（0.48*pi*n）+cos（0.52*pi*n）;
n1=（0:
1:
9）;x=y（1:
1:
10）;
subplot（3,1,1）;
stem（n1,x）;%画出x（n）曲线
title（'x（n）（0<=n<=9）'）;%标题0<=n<=9）
xlabel（'n'）;ylabel（'x（n）'）;
axis（[0,10,-2.5,2.5]）;%axis（[xminxmaxyminymax]）
w=linspace（0,2*pi,length（x））;%0-2*pi区域分为10点
xw=x*exp（-j*[1:
length（x）]'*w）;
magx=abs（xw）;%对xw取绝对值
subplot（3,1,2）;
plot（w,magx）;%画出x（jw）曲线
title（'DTFT'）;%标题DTFT
xlabel（'w'）;ylabel（'x（jw）'）;
axis（[0,2*pi,0,10]）;%axis（[xminxmaxyminymax]）
subplot（3,1,3）;
x1=fft（x）;%对x进行傅立叶变换
magx1=abs（x1）;%对x1取绝对值
stem（n1,abs（magx1））;%画出x（k）曲线
title（'DFT'）;%标题DFT
xlabel（'k'）;ylabel（'x（k）'）;
axis（[0,10,0,10]）;%axis（[xminxmaxyminymax]）
结果：
分析：
由图可见，由于截断函数的频谱混叠作用，X（K）不能正确分辨w1=0.48*pi,w2=0.52*pi这两个频率分量。
（2）
代码：
closeall
clear
clc
n=（0:
1:
9）;
y=cos（0.48*pi*n）+cos（0.52*pi*n）;
n1=（0:
1:
99）;
x=[y（1:
1:
10）,zeros（1,90）];%第10位到100位的数据都为0
subplot（3,1,1）;
stem（n1,x）;%画出x（n）曲线
title（'x（n）（0<=n<=9+90zeros）'）;%标题0<=n<=9+90zeros）
xlabel（'n'）;ylabel（'x（n）'）;
axis（[0,100,-2.5,2.5]）;%axis（[xminxmaxyminymax]）
w=linspace（0,2*pi,length（x））;%0-2*pi区域分为100点
xw=x*exp（-j*[1:
length（x）]'*w）;
magx=abs（xw）;%对xw取绝对值
subplot（3,1,2）;
plot（w,magx）;%画出x（jw）曲线
title（'DTFT'）;%标题DTFT
xlabel（'w'）;ylabel（'x（jw）'）;
axis（[0,2*pi,0,10]）;%axis（[xminxmaxyminymax]）
subplot（3,1,3）;
x1=fft（x）;%对x进行傅立叶变换
magx1=abs（x1）;%对x1取绝对值
stem（n1,abs（magx1））;%画出x（k）曲线
title（'DFT'）;%标题DFT
xlabel（'k'）;ylabel（'x（k）'）;
axis（[0,100,0,10]）;%axis（[xminxmaxyminymax]）
结果：
分析：
由图可见，虽然x（n）补零至100点，X（K）的密度，截断函数的频谱混叠作用没有改变，这时的物理分辨率使X（K）仍不能正确分辨w1=0.48*pi,w2=0.52*pi这两个频率分量。
（3）
代码：
closeall
clear
clc
n=（0:
1:
99）;
y=cos（0.48*pi*n）+cos（0.52*pi*n）;
n1=（0:
1:
99）;x=y（1:
1:
100）;
subplot（3,1,1）;
stem（n1,x）;%画出x（n）曲线
title（'x（n）（0<=n<=99）'）;%标题0<=n<=99）
xlabel（'n'）;ylabel（'x（n）'）;
axis（[0,100,-2.5,2.5]）;%axis（[xminxmaxyminymax]）
w=linspace（0,2*pi,length（x））;%0-2*pi区域分为100点
xw=x*exp（-j*[1:
length（x）]'*w）;
magx=abs（xw）;%对xw取绝对值
subplot（3,1,2）;
plot（w,magx）;%画出x（jw）曲线
title（'DTFT'）;%标题DTFT
xlabel（'w'）;ylabel（'x（jw）'）;
axis（[0,2*pi,0,54]）;%axis（[xminxmaxyminymax]）
subplot（3,1,3）;
x1=fft（x）;%对x进行傅立叶变换
magx1=abs（x1）;%对x1取绝对值
stem（n1,abs（magx1））;%画出x（k）曲线
title（'DFT'）;%标题DFT
xlabel（'k'）;ylabel（'x（k）'）;
axis（[0,100,0,54]）;%axis（[xminxmaxyminymax]）
结果：
分析：
由图可见，截断函数的加宽且为周期序列的整数倍，改变了频谱混叠作用，提高了“物理”分辨率使X（K）能正确分辨w1=0.48*pi,w2=0.52*pi这两个频率分量。
（4）
代码：
closeall
clear
clc
n=（0:
1:
127）;
y=cos（0.48*pi*n）+cos（0.52*pi*n）;
n1=（0:
1:
127）;x=y（1:
1:
128）;
subplot（3,1,1）;
stem（n1,x）;%画出x（n）曲线
title（'x（n）（0<=n<=127）'）;%标题0<=n<=127）
xlabel（'n'）;ylabel（'x（n）'）;
axis（[0,128,-2.5,2.5]）;%axis（[xminxmaxyminymax]）
w=linspace（0,2*pi,length（x））;%0-2*pi区域分为128点
xw=x*exp（-j*[1:
length（x）]'*w）;
magx=abs（xw）;%对xw取绝对值
subplot（3,1,2）;
plot（w,magx）;%画出x（jw）曲线
title（'DTFT'）;%标题DTFT
xlabel（'w'）;ylabel（'x（jw）'）;
axis（[0,2*pi,0,65]）;%axis（[xminxmaxyminymax]）
subplot（3,1,3）;
x1=fft（x）;%对x进行傅立叶变换
magx1=abs（x1）;%对x1取绝对值
stem（n1,abs（magx1））;%画出x（k）曲线
title（'DFT'）;%标题DFT
xlabel（'k'）;ylabel（'x（k）'）;
axis（[0,128,0,65]）;%axis（[xminxmaxyminymax]）
结果：
分析：
由图可见，截断函数虽进一步加宽，但不是周期序列的整数倍，所以尽管X（K）能正确分辨w1=0.48*pi,w2=0.52*pi这两个频率分量，但还是呈现除了频谱泄露。
（5）
代码：
closeall
clearclc
n=（0:
1:
49）;
y=cos（0.48*pi*n）+cos（0.52*pi*n）;
n1=（0:
1:
49）;x=y（1:
1:
50）;
subplot（3,1,1）;stem（n1,x）;%画出x（n）曲线
title（'x（n）（0<=n<=49）'）;%标题0<=n<=49）
xlabel（'n'）;ylabel（'x（n）'）;
axis（[0,50,-2.5,2.5]）;%axis（[xminxmaxyminymax]）
w=linspace（0,2*pi,length（x））;%0-2*pi区域分为50点
xw=x*exp（-j*[1:
length（x）]'*w）;
magx=abs（xw）;%对xw取绝对值
subplot（3,1,2）;plot（w,magx）;%画出x（jw）曲线
title（'DTFT'）;%标题DTFT
xlabel（'w'）;ylabel（'x（jw）'）;
axis（[0,2*pi,0,30]）;%axis（[xminxmaxyminymax]）
subplot（3,1,3）;
x1=fft（x）;%对x进行傅立叶变换
magx1=abs（x1）;%对x1取绝对值
stem（n1,abs（magx1））;%画出x（k）曲线
title（'DFT'）;%标题DFT
xlabel（'k'）;ylabel（'x（k）'）;
axis（[0,50,0,30]）;%axis（[xminxmaxyminymax]）
结果：
分析：
由图可见，截断函数的宽度正好为序列的周期，即这时的“物理”分辨率使X（K）正好能正确分辨w1=0.48*pi,w2=0.52*pi这两个频率分量。
题目3：
利用MATLAB编程设计一个数字带通滤波器，指标要求如下：
通带边缘频率：
，
，通带峰值起伏：
阻带边缘频率：
，
最小阻带衰减：
分别用IIR和FIR两种数字滤波器类型进行设计。
实验要求：
给出IIR数字滤波器参数和FIR数字滤波器的冲激响应，绘出它们的幅度和相位频响曲线，讨论它们各自的实现形式和特点。
（1）
代码：
closeall
clear
clc
ws=[0.3,0.75];
wp=[0.45,0.65];
Rp=1;Rs=40;
[N,wc]=buttord（wp,ws,Rp,Rs）;%求阶数N及频率参数wc
[B,A]=butter（N,wc）;%设计巴特沃斯滤带通波器
[H,w]=freqz（B,A）;%求频响
subplot（2,1,1）;
plot（w/pi,abs（H））;%画幅频曲线
title（'IIR幅度频响曲线'）;%标题IIR幅度频响曲线
xlabel（'w/pi'）;
ylabel（'幅度'）;
subplot（2,1,2）;
plot（w/pi,angle（H））;%画相频曲线
title（'IIR相位频响曲线'）;%标题IIR相位频响曲线
xlabel（'w/pi'）;
ylabel（'相位'）;
IIR数字滤波器参数：
N=7
wc=0.42840.6693
B=0.00030-0.001900.00570-0.0095
00.00950-0.005700.00190
-0.0003
A=1.00001.74514.92826.11959.81349.224510.43237.51546.40913.45952.26010.84700.41670.08560.0299
IIR幅度和相位频响曲线：
（2）
FIR设计：
由于
，查表可选hamming窗，其阻带最小衰减-53dB满足要求。
代码：
closeall
clear
clc
ws1=0.3*pi;ws2=0.75*pi;
wp1=0.45*pi;wp2=0.65*pi;
wb=（wp1-ws1+ws2-wp2）/2;%计算过渡带宽
wc=[wp1/pi-（wp1-ws1）/（2*pi）;%设置理想带通截止频率
wp2/pi+（ws2-wp2）/（2*pi）];
N=ceil（3.3*2*pi/wb）;%计算窗口长度
b=fir1（N,wc）;%设计滤波器系数
n=0:
N;
hn=b（n+1）;%求冲激响应
[H,w]=freqz（b,1）;%求频率响应
subplot（2,1,1）;
plot（w/pi,abs（H））;%画幅频曲线
title（'FIR幅度频响曲线'）;%标题FIR幅度频响曲线
xlabel（'w/pi'）;ylabel（'幅度'）;
subplot（2,1,2）;
plot（w/pi,angle（H））;%画相频曲线
title（'FIR相位频响曲线'）;%标题FIR相位频响曲线
xlabel（'w/pi'）;ylabel（'相位'）;
冲激响应hn=b（n+1）：
0.00110.00050.00010.0006-0.0031-0.00080.0051-0.0002
-0.0002-0.0016-0.01090.00870.0149-0.0098-0.0018-0.0119
-0.01450.04630.0118-0.04450.0003-0.03400.02500.1591
-0.1126-0.23970.21350.2135-0.2397-0.11260.15910.0250
-0.03400.0003-0.04450.01180.0463-0.0145-0.0119-0.0018
-0.00980.01490.0087-0.0109-0.0016-0.0002-0.00020.0051
-0.0008-0.00310.00060.00010.00050.0011
FIR幅度和相位频响曲线：
分析IIR与FIR的实现形式与特点：
1、IIR滤波器阶数比FIR少。
IIR滤波器存在着输出对输入的反馈，因此可以用比FIR滤波器少的阶数来满足技术指标。
2、FIR滤波器可得到严格的线性相位，而IIR滤波器则做不到这一点。
IIR滤波器的选频特性越好，则相位的非线性就越严重，在需要严格线性相位的情况下应该选择FIR滤波器。
3、IIR滤波器必须采用递归结构实现，FIR滤波器主要采用非递归结构。
IIR滤波器必须采用递归结构实现，只有当所有极点都在单位圆内时滤波器才是稳定的。
但实际中由于存在有限字长效应，滤波器有可能变得不稳定。
而FIR滤波器主要采用非递归结构，因而从理论上以及从实际的有限精度的运算中，都是稳定的。
4、IIR滤波器有现成的设计公式、数据和表格，FIR滤波器没有现成的设计公式。
IIR滤波器可利用模拟滤波器现成的设计公式、数据和表格，因而计算工作量较小，对计算工具要求不高。
FIR滤波器没有现成的设计公式，窗函数法只给出窗函数的计算公式，但计算通带和阻带衰减仍无显式表达式。
5、IIR滤波器主要是设计规格化的、频率特性为分段常数的标准低通、高通、带通、带阻和全通滤波器，而FIR滤波器可设计出理想正交变换器、理想微分器、线性调频器等各种网络，适应性较广。
6、FIR滤波器可以采用快速傅立叶变换（FFT）来实现，在相同阶数下，运算速度可以快得多
三、自我评价和总结
在这将近三天的课设里，我学习了一些与MATLAB相关的基础知识，增强的自己的知识，这段日子，时间过的真的好快，解决一个问题要花费很长时间，这是因为自己的能力缘故，但是每当解决一个问题之后，很开心，真的很开心，收获得很大，刚开始的时候，没有头绪，不知道从何做起，后来自己翻课本，在网上查阅相关的资料，再加上老师的讲解，慢慢的有了头绪，然后自己还挺高兴的，就这样到了最后一天下午，老师看到截图后说周期有问题，我们的缺失了一个周期，必须要调试，于是我们走上了漫漫的调试之旅，很长时间后，问题都还没有解决，然后老师就开始慢慢的指导我们了，真的很佩服老师，感觉老师真的很厉害，不到5分钟就调试结束了，期间还给我们讲解了过程，删了一些代码，又加了一些代码，运行后，非常成功，然后我们自己对老师的这种方法理解了大部分，剩下的一些题，跟老师讲的这个题是同类型的，只是参数不同，然后我们又开始调试了，这一次比较快，一个多小时，初步完成了，然后又接着写代码注释，截图，分析，在宿舍写了一晚上，最后终于基本完成了本次课程设计，我相信这次经历对以后的学习与工作都有很大的帮助，另外，衷心的感谢老师的指导，如果没有老师的指导，课设不会这么快的完成，老师，您辛苦了。
升级会员 