第一篇:电子科技大学 数字信号处理 第二次编程作业
第二次编程作业
加载信号(EEGdata.txt);分析信号的幅度谱,确定低通FIR数字滤波器的指标;分别利用各种窗函数(Rectangular, Hamming, Kaiser)设计此FIR滤波器;对信号加随机噪声,并用设计的FIR滤波器对含噪声信号进行滤波。要求:
1)给出程序,画出滤波器幅度谱及其增益图(dB);分析设计的滤波器是否达到指标。
2)利用设计的三种滤波器对加载的信号分别进行滤波,对比分析滤波结果。程序:
data = load('EEGdata.txt');n=1:length(data);
xn=0.1*rand(length(data),1)-0.05;x=xn+data;
[h0,w0]=freqz(data);[h1,w1]=freqz(x);figure(1)
plot(w0/pi,abs(h0),w1/pi,abs(h1),'r');legend('Ô-dataÐźÅƵÆ×ͼ','¼ÓÔëÐźÅƵÆ×ͼ');
wp = 0.15*pi;ws = 0.2*pi;DB = ws-wp;
N = ceil(2*3.11*pi/DB+1);wc =(wp+ws)/2;wc = wc/pi;
win = rectwin(N+1);
b = fir1(N, wc, 'low', win);figure(2)freqz(b,1);
wp = 0.15*pi;ws = 0.2*pi;DB = ws-wp;
N = ceil(2*3.11*pi/DB+1);wc =(wp+ws)/2;wc = wc/pi;
win = hamming(N+1);
b1 = fir1(N, wc, 'low', win);figure(3)freqz(b1,1);
fpts = [0.15, 0.2];mag = [1, 0];dev = [0.01, 0.01];
[N, Wn, beta, ftype] = kaiserord(fpts, mag, dev)Kwin = kaiser(N+1, beta);b2 = fir1(N, Wn, Kwin);[h,omega] = freqz(b2,1,512);figure(4)freqz(b2,1);
y=filter(b, 1, x);figure(5)subplot(2,1,1);
plot(n,data,n,x,'g',n,y,'r');title('¾ØÐδ°Â˲¨·ù¶ÈÆ×Ч¹û¶Ô±Èͼ');legend('Ô-ÐźÅdata','¼ÓÔëÐźÅx','Â˲¨ºóÐźÅy');[H,W]=freqz(y);subplot(2,1,2);
plot(w0/pi,abs(h0),W/pi,abs(H),'r');title('¾ØÐδ°Â˲¨ÆµÆ×Ч¹û¶Ô±Èͼ');legend('Ô-ÐźÅdata','Â˲¨ºóÐźÅy');
y1=filter(b1, 1, x);figure(6)subplot(2,1,1);
plot(n,data,n,x,'g',n,y1,'r');title('Hamming´°Â˲¨·ù¶ÈÆ×Ч¹û¶Ô±Èͼ');legend('Ô-ÐźÅdata','¼ÓÔëÐźÅx','Â˲¨ºóÐźÅy1');[H,W]=freqz(y1);subplot(2,1,2);
plot(w0/pi,abs(h0),W/pi,abs(H),'r');title('Hamming´°Â˲¨ÆµÆ×Ч¹û¶Ô±Èͼ');legend('Ô-ÐźÅdata','Â˲¨ºóÐźÅy1');
y2=filtfilt(b2, 1, x);figure(7)subplot(2,1,1);
plot(n,data,n,x,'g',n,y2,'r');title('Kaiser´°Â˲¨·ù¶ÈÆ×Ч¹û¶Ô±Èͼ');legend('Ô-ÐźÅdata','¼ÓÔëÐźÅx','Â˲¨ºóÐźÅy2');[H,W]=freqz(y2);subplot(2,1,2);
plot(w0/pi,abs(h0),W/pi,abs(H),'r');title('Kaiser´°Â˲¨ÆµÆ×Ч¹û¶Ô±Èͼ');legend('Ô-ÐźÅdata','Â˲¨ºóÐźÅy2');
figure(8)
plot(n,data,n,y,'c',n,y1,'r',n,y2,'g');title('ÈýÖÖ´°º¯ÊýÂ˲¨Ð§¹û¶Ô±Èͼ');legend('Ô-dataÐźÅ','¾ØÐδ°Â˲¨','Hamming´°Â˲¨','Kaiser´°Â˲¨');
***6420 0原data信号频谱图加噪信号频谱图 0.10.20.30.40.50.60.70.80.91 50Magnitude(dB)0-50-10000.10.20.30.40.50.60.70.8Normalized Frequency( rad/sample)0.910Phase(degrees)-1000-202_-300000.10.20.30.40.50.60.70.8Normalized Frequency( rad/sample)0.91
50Magnitude(dB)0-50-100-15000.10.20.30.40.50.60.70.8Normalized Frequency( rad/sample)0.910Phase(degrees)-1000-202_-300000.10.20.30.40.50.60.70.8Normalized Frequency( rad/sample)0.91 50Magnitude(dB)0-50-100-15000.10.20.30.40.50.60.70.8Normalized Frequency( rad/sample)0.910Phase(degrees)-500-1000-1500-200000.10.20.30.40.50.60.70.8Normalized Frequency( rad/sample)矩形窗滤波幅度谱效果对比图0.20.10-0.1-0.2 ***原信号data加噪信号x滤波后信号y300350 0.91矩形窗滤波频谱效果对比图1510原信号data滤波后信号y 50 00.10.20.30.40.50.60.70.80.91Hamming窗滤波幅度谱效果对比图0.20.10-0.1-0.2 050100150200原信号data加噪信号x滤波后信号y1250300350 Hamming窗滤波频谱效果对比图1510原信号data滤波后信号y1 50 00.10.20.30.40.50.60.70.80.91Kaiser窗滤波幅度谱效果对比图0.20.10-0.1-0.2 050100150200原信号data加噪信号x滤波后信号y2250300350 Kaiser窗滤波频谱效果对比图20151050 00.10.20.30.40.50.60.70.80.91原信号data滤波后信号y2 三种窗函数滤波效果对比图0.2原data信号矩形窗滤波Hamming窗滤波Kaiser窗滤波 0.150.10.050-0.05-0.1-0.***0200250300350
结果分析:滤波效果皆不尽如人意,原因是加载的噪声信号是随机信号,各频率皆有,滤波器通带部分的噪声无法滤掉。尤其是矩形窗函数和Hamming窗函数构造的滤波器长度大于data的三分之一,无法使用filtfilt函数滤波,只能用filter函数,致使滤波结果有相位平移。
第二篇:数字信号处理实验教学-电子教案
数字信号处理实验 电子教案
讲义1 Matlab简介及其安装使用说明..................................................2 讲义2 Matlab基本语句..........................................................................8 讲义3 Matlab基本数值运算................................................................13 讲义4 Matlab函数、及其调用方法....................................................16 实验1 常见离散信号产生和实现.........................................................20 实验2 离散系统的时域分析.................................................................22 实验3 FFT算法的应用..........................................................................24 实验4 离散系统的变换域分析.............................................................27 实验5 有限冲激响应数字滤波器设计.................................................32 实验6 无限冲激响应数字滤波器设计.................................................36 实验7 设计性和研究性实验.................................................................41讲义1 Matlab简介及其安装使用说明
一.MATLAB程序设计语言简介
MATLAB,Matrix Laboratory的缩写,是由Mathworks公司开发的一套用于科学工程计算的可视化高性能语言,具有强大的矩阵运算能力。与大家常用的Fortran和C等高级语言相比,MATLAB的语法规则更简单,更贴近人的思维方式,被称之为“草稿纸式的语言”。截至目前,MATLAB已经发展到7.x版,适用于所有32位的Windows操作系统,按NTFS(NT文件系统)格式下完全安装约需 850 MB。MATLAB软件主要由主包、仿真系统和工具箱三大部分组成。
二.MATLAB应用入门
1.MATLAB的安装与卸载
MATLAB软件在用户接口设计上具有较强的亲和力,其安装过程比较典型,直接运行光盘中的安装向导支撑程序SETUP.exe,按其提示一步步选择即可。MATLAB自身带有卸载程序,在其安装目录下有uninstall子目录,运行该目录下的uninstall.exe即可; 也可以通过Windows系统的安装卸载程序进行卸载。2.MATLAB的启动与退出
MATLAB安装完成后,会自动在Windows桌面上生成一个快捷方式,它是指向安装目录下binwin32matlab.exe的链接,双击它即可来到MATLAB集成环境的基本窗口,通常称之为命令窗口。MATLAB的退出与普通WIN32的程序一样,值得一提的是它有一个自身专有的快捷键Ctrl+Q。3.MATLAB界面简介
图 1 MATLAB基本界面——命令窗口
图2
图3
图4
图5
图6 指令历史
图7 1)菜单栏
菜单栏中包括File、Edit、View、Web、Window和Help六个菜单项。这里着重介绍File、help项。
File项:数据输入/输出的接口,包括10个子项,这里重点介绍其中的5个子项: New:新建文件项。有四个选择:M File(*.M,文本格式的MATLAB程序文件,可以直接通过文件名的方式在MATLAB环境下解释运行;Figure(图形);Model(仿真模型文件)和GUI(可视化界面文件)。
Open:打开所有MATLAB支持的文件格式,系统将自动识别并采用相应的程序对文件进行处理。例如, 打开一个.m文件,系统将自动打开M文件编辑器对它进行编辑。
Set Path...:设置工作路径。可以打开路径设置(Set Path)对话框(图2),将用户自己建立的目录加入MATLAB的目录系统中,以便所编制的文件能够在MATLAB环境中直接调用。
图8 路径设置对话框
单击Add Folder...按钮可以将你的一个文件夹加入到系统路径中;Add with Subfolders...允许把一个文件夹包括其所有的子文件夹加入到系统路径中。这两种操作均可以直观地在右侧的路径栏内看到结果。
选中一个加入的文件夹,你可以利用Move to Top(移至所有路径的最前面),Move Up(上移一个),Move Down(下移一个),Move to Bottom(移至所有路径的最后面)等四个按钮将改变文件在系统路径中的排列位置以利于对文件的搜索使用,也可以利用Remove按钮将其删除。
对路径操作完毕后,按Save按钮予以保存;Help 项:
Matlab Help:打开以html超文本形式存储的的帮助文件主页; Demos:打开matlab演示窗主页;
About Matlab:Matlab注册图标、版本、制造商和用户信息;
图9 Help选项
图10 Help窗口
2)命令行区
进行通用操作,数值计算,编程和数据类型,输入输出,绘图,用户界面等命令,例如,命令:help函数名(*.m文件); edit编辑函数、文件
对输入命令的解释MATLAB按以下顺序进行:
① 检查它是否是工作空间中的变量,是则显示变量内容。② 检查它是否是嵌入函数,是则运行之。③ 检查它是否是子函数。④ 检查它是否是私有函数。
⑤ 检查它是否是位于MATLAB搜索路径范围内的函数文件或脚本文件。
请注意,如果有两个以上的方案与输入的命令相匹配,MATLAB将只执行第一个匹配。
讲义2 Matlab基本语句
一.程序控制语句
1.循环语句
MATLAB的循环语句包括for循环和while循环两种类型。1)for循环 语法格式:
for 循环变量 = 起始值:步长:终止值 循环体 end 起始值和终止值为一整形数,步长可以为整数或小数,省略步长时,默认步长为1。执行for循环时,判定循环变量的值是否大于(步长为负时则判定是否小于)终止值,不大于(步长为负时则小于)则执行循环体,执行完毕后加上步长,大于(步长为负时则小于)终止值后退出循环。例1 给矩阵A、B赋值。MATLAB 语句及运行结果如下: k=5;a=zeros(k, k)%矩阵赋零初值 for m=1 : k for n=1: k a(m,n)=1/(m+n-1);end end for i=m :-1 : 1 b(i)=i;end 运行结果: a= 1.0000 0.5000 0.3333 0.2500?0.202_ 0.5000 0.3333 0.2500 0.202_ 0.1667 0.3333 0.2500 0.202_ 0.1667 0.1429 0.2500 0.202_ 0.1667 0.1429 0.1250 0.202_ 0.1667 0.1429 0.1250 0.1111 b= 1 2?3 4 5 2)while循环 语法格式:
while 表达式 循环体 end 其执行方式为:若表达式为真(运算值非0),则执行循环体;结果为0),则退出循环体,执行end后的语句。
例2 a=3;while a a=a-1 end 输出: a=2 a=1 a=0 2.条件转移语句
条件转移语句有if和switch两种。
若表达式为假(运算 1)if语句
MATLAB中if语句的用法与其他高级语言相类似,其基本语法格式有以下几种: 格式一: if 逻辑表达式 执行语句 end 格式二: if 逻辑表达式 执行语句1 else 执行语句2 end 格式三:
if 逻辑表达式1 执行语句1
else? if 逻辑表达式2 执行语句2 end 2)switch语句
switch语句的用法与其他高级语言相类似,其基本语法格式为:
switch表达式(标量或字符串)case 值1 语句1 case 值2 语句2 … otherwise 语句n end
二.绘图语句 常用的MATLAB绘图语句有figure、plot、subplot、stem等,图形修饰语句有title、axis、text等。
1.figure figure有两种用法,只用一句figure命令,会创建一个新的图形窗口,并返回一个整数型的窗口编号。figure(n)表示将第n号图形窗口作为当前的图形窗口,并将其显示在所有窗口的最前面;如果该图形窗口不存在,则新建一个窗口,并赋以编号n。2.plot 线型绘图函数。用法为plot(x,y,'s')。参数x为横轴变量,y为纵轴变量,s用以控制图形的基本特征如颜色、粗细等,通常可以省略,常用方法如表1所示。
表1
3.Stem 绘制离散序列图,常用格式stem(y)和stem(x,y)分别和相应的plot函数的绘图规则相同,只是用stem命令绘制的是离散序列图。4.subplot subplot(m,n,i)图形显示时分割窗口命令,把一个图形窗口分为m行,n列,m×n个小窗口,并指定第i个小窗口为当前窗口。5.绘图修饰命令
在绘制图形时,我们通常需要为图形添加各种注记以增加可读性。在plot语句后使用title('标题')可以在图形上方添加标题,使用xlabel('标记')或ylabel('标记')为X轴或Y轴添加说明,使用text(X值、Y值、'想加的标示')可以在图形中任意位置添加标示。例3 画图基本语句如图1所示。
MATLAB 语句及运行结果如下: x=0:0.1*pi:2*pi;%定义x向量
figure(1);%创建一个新的图形窗口,编号为1 subplot(2,2,1);%将窗口划分为2行,2列,在第1个窗口中作图 plot(x,sin(x));%画图
title('正弦线');%给图形加标题 subplot(2,2,2);%在第2个窗口中作图 plot(x,sin(x),'r');%画一正弦波,红色 xlabel('X');%给x轴加说明 ylabel('SIN(X)');%给y轴加说明 subplot(2,2,3);%在第3个窗口中作图 plot(x,sin(x),'--');%画一正弦波,破折线 subplot(2,2,4);%在第4个窗口中作图 plot(x,sin(x),'r+');%画一正弦波,红色+线 text(4,0,'注记');
图1 讲义3 Matlab基本数值运算
一.MATLAB内部特殊变量和常数
MATLAB内部有很多变量和常数,用以表达特殊含义。常用的有:
1.变量ans: 指示当前未定义变量名的答案;
2.常数eps:表示浮点相对精度,其值是从1.0到下一个最大浮点数之间的差值。该变量值作为一些MATLAB函数计算的相对浮点精度,按IEEE标准,如:,近似为2.2204e-016;
3.常数Inf:表示无穷大。当输入或计算中有除以0时产生Inf;
4.虚数单位i,j:表示复数虚部单位,相当于5.NaN:表示不定型值,是由0/0运算产生的。
6.常数pi:表示圆周率π,其值为3.14***…;
;
MATLAB中可表示的数字的近似范围从,1.有效数字表示的典型例子如下:1234.56789,123456.789E-2,1.23456789e3(format指令可以控制显示格式)
2.复数形式:3.5+4*j,-2.1-7.4*j(i 也可以)
取绝对值:abs()语法格式:abs(x)。当x为实数时计算x的绝对值;x为复数时得到的是复数的模值;x为字符串时得到各字符的ASCII码。取相角:angle()语法格式:angle(z)。求复矢量或复矩阵的相角,结果为一个以弧度为单位介于-π和+π之间的值。
二.变量
1.变量命名规则
MATLAB中对变量的命名应遵循以下规则:
(1)变量名可以由字母、数字和下划线混合组成,但必须以字母开头;(2)字符长度不能大于31;(3)变量命名区分大小写。2.局部变量和全局变量 局部变量是指那些每个函数体内自己定义的,不能从其他函数和MATLAB工作空间访问的变量。全局变量是指用关键字“global”声明的变量。全局变量名应尽量大写,并能反映它本身的含义。如果需要在工作空间和几个函数中都能访问一个全局变量,必须在工作空间和这几个函数中都声明该变量是全局的。
三.矩阵及其运算
MATLAB具有强大的矩阵运算和数据处理功能,对矩阵的处理必须遵从代数规则。1.矩阵生成1)一般矩阵的生成
对于一般的矩阵MATLAB的生成方法有多种。最简单的方法是从键盘直接输入矩阵元素。直接输入矩阵元素时应注意:各元素之间用空格或逗号隔开,用分号或回车结束矩阵行,用中括号把矩阵所有元素括起来。
例1 在工作空间产生一个3×3矩阵A可用MATLAB语言描述如下: A=[1 2 3;4 5 6;7 8 9] 或 A=[1 2 3 4 5 6 7 8 9] 运行结果: A= 1 2 3 4 5 6 7 8 9 Size(A)得到矩阵的大小,ans = 3 3 2)特殊矩阵的生成
对于特殊的矩阵可直接调用MATLAB的函数生成。
用函数zeros生成全0矩阵:格式 B=zeros(m,n)生成m×n的全0阵。用函数ones生成全1矩阵:格式 B=ones(m,n)生成m×n的全1阵。
用函数eye生成单位阵:格式 B=eye(m,n)生成m×n矩阵,其中对角线元素全为1,其他元素为0。2.矩阵的运算
矩阵的运算有基本运算和函数运算两种类型。基本运算包括矩阵的加、减、乘、除、乘方、求转置、求逆等,其主要特点是通过MATLAB提供的基本运算符+、-、*、/()、^等即可完成。函数运算主要是通过调用MATLAB系统内置的运算函数来求取矩阵的行列式(det(A)),求秩(rank(A)),求逆(inv(A)),求特征值和特征向量([V,D]=eig(A)), 求Jordan标准形(jordan(A))和矩阵分解等。需要用时可以参阅联机帮助和相关参考书。(.*,.,表示逐个元素的乘积和相除;矩阵X/Y 相当于X*inv(Y), XY相当于inv(Y)*X)
例2 矩阵的基本运算。A=[1, 2, 3;4, 5, 6];B =[6, 5, 4;3, 2, 1];C =A+B %计算两个矩阵的和 D =B' %计算矩阵B的转置
E=A*D %做矩阵乘法,必须要满足矩阵乘法的基本要求E应该是2阶方阵 F=det(E)%求E的行列式值 G=E^(-1)%求E的逆
输出结果: C= 7 7 7 7 7 7 D= 6 3 5 2 4 1 E= 28 10 73 28 F=54 G= 0.5185-0.1852-1.3519 0.5185 讲义4 Matlab函数、及其调用方法
在MATLAB语言中,M文件有两种形式:脚本和函数。
脚本没有输入/输出参数,只是一些函数和命令的组合。它可以在MATLAB环境下直接执行,也可以访问存在于整个工作空间内的数据。由脚本建立的变量在脚本执行完后仍将保留在工作空间中可以继续对其进行操作,直到使用clear命令对其清除为止。
函数是MATLAB语言的重要组成部分。MATLAB提供的各种工具箱中的M文件几乎都是以函数的形式给出的。函数接收输入参数,返回输出参数,且只能访问该函数本身工作空间中的变量,从命令窗或其他函数中不能对其工作空间的变量进行访问。
1.函数结构
MATLAB语言中提供的函数通常由以下五个部分组成:
(1)函数定义行: 以function开头,函数名(必须与文件名相同)及函数输入输出参数在此定义;(2)H1行:第一注释行,供lookfor和help在线帮助使用;
(3)函数帮助文件;通常包括函数输入输出参数的含义,调用格式说明;
(4)函数体:它包括进行运算和赋值的所有MATLAB程序代码。函数体中可以包括流程控制、输入/输出、计算、赋值、注释以及函数调用和脚本文件调用等。在函数体中完成对输出参数的计算;
(5)注释。
这五个部分中最重要的是函数定义行和函数体。
函数定义行是一个MATLAB函数所必需的,其他各部分的内容可以没有,这种函数称为空函数。
例1 function sa=circle(r,s)% Circle plot a circle of radii r in the line specified by s % r raddi % s line color % sa area of the circle % % circle(r)use blue line to draw a circle of radii r % circle(r,s)use 's' to draw circle % sa=circle(r)compute the area of the circle and draw it in blue % sa=circle(r,s)compute the area of the circle and draw it in color 's' if nargin>2 error('to many input ');end clf;t=0:pi/100:2*pi;x=r*exp(i*t);if nargout==0 plot(x,s);else sa=pi*r*r;fill(real(x),imag(x),s)end axis('square');%makes the current axis box square in size.例2 function y = imp_fun(n,n0)% IMP_FUN Unit impulse function.% % IMP_FUN(N,N0), where N is a vector of sequential integers, % returns a vector the same length as N with zeros everywhere except % N = N0.[m1,n1] = size(n);if ~(m1 == 1 | n1 == 1)error('The sample vector must be one-dimensional.');end y = zeros(m1,n1);i = find(n >= n0);if isempty(i)return end y(i(1))= 1;% where n = n0, set output to 1 2.函数调用 函数调用的过程实际上就是参数传递的过程。例如,在一个脚本文件里调用函数“max”可采用如下方式: n=1:20;
a=sin(2*pi*n/20); [Y,I]=max(a);该调用过程把变量“a”传给了函数中的输入参数“x”,然后把函数运算的返回值传给输出参数“Y”和“I”。其中,Y是a序列的最大值,I是最大值Y对应的坐标值。
例3 构造:y[n] = δ[n-3]: 调用函数: n = 0:6;y = imp_fun(n,3);stem(n,y)
图1 例4 构造:y[n] = 5δ[n]imp_fun(n,2);stem(n,y);
图2 实验1 常见离散信号产生和实现
一、实验目的
1、加深对常用离散信号的理解;
2、熟悉使用MATLAB在时域中产生一些基本的离散时间信号。
二、实验原理
1、单位抽样序列
在MATLAB中可以利用
函数实现。
2、单位阶越序列
在MATLAB中可以利用
函数实现:
3、正弦序列
在MATLAB中实现过程如下:
4、复指数序列
在MATLAB中实现过程如下:
5、指数序列
在MATLAB中实现过程如下:
三、预习要求
1、预先阅读实验讲义(MATLAB基础介绍);
2、讨论正弦序列、复指数序列的性质。
A.绘出信号,当、时、、时的信号实部和虚部图;当B.绘出信号
时呢?此时信号周期为多少? 的频率是多少?周期是多少?产生一个数字频率为0.9的正弦序列,并显示该信号,说明其周期。
3、使用帮助功能学习square(方波),sawtooth(锯齿波)和sinc函数,并绘图。
四、实验内容
编制程序产生上述5种信号,长度可输入确定,函数需要的参数可输入确定,并绘出其图形。
实验2 离散系统的时域分析
一、实验目的
1、熟悉并掌握离散系统的差分方程表示法;
2、加深对冲激响应和卷积分析方法的理解。
二、实验原理
在时域中,离散时间系统对输入信号或者延迟信号进行运算处理,生成具有所需特性的输出信号,具体框图如下:
其输入、输出关系可用以下差分方程描述:
输入信号分解为冲激信号,记系统单位冲激响应,则系统响应为如下的卷积计算式:
当称系统为IIR系统。
时,h[n]是有限长度的(),称系统为FIR系统;反之,三、预习要求
1、在MATLAB中,熟悉利用函数
2、在MATLAB中,熟悉用函数 响应的过程。
实现差分方程的仿真;
计算卷积,用求系统冲激
四、实验内容
1、以下程序中分别使用conv和filter函数计算h和x的卷积y和y1,运行程序,并分析y和y1是否有差别,为什么要使用x[n]补零后的x1来产生y1;具体分析当h[n]有i个值,x[n]有j个值,使用filter完成卷积功能,需要如何补零? % Program P2_7 clf;h = [3 2 1-2 1 0-4 0 3];%impulse response x = [1-2 3-4 3 2 1];%input sequence y = conv(h,x);n = 0:14;subplot(2,1,1);stem(n,y);xlabel('Time index n');ylabel('Amplitude');title('Output Obtained by Convolution');grid;x1 = [x zeros(1,8)];y1 = filter(h,1,x1);subplot(2,1,2);stem(n,y1);xlabel('Time index n');ylabel('Amplitude');title('Output Generated by Filtering');grid;
2、编制程序求解下列两个系统的单位冲激响应和阶跃响应,并绘出其图形。要求分别用 filter、conv、impz三种函数完成。,,给出理论计算结果和程序计算结果并讨论。实验3 FFT算法的应用
一、实验目的
1、加深对离散信号的DFT的理解;
2、在MATLAB中实现FFT算法。
二、实验原理
N点序列的DFT和IDFT变换定义式如下: , , 利用旋转因子具有周期性,可以得到快速算法(FFT)。
在MATLAB中,可以用函数反变换。
和计算N点序列的DFT正、三、预习要求
1、在MATLAB中,熟悉函数fft、ifft的使用;
2、阅读扩展练习中的实例,学习在MATLAB中的实现FFT算法的实现;
3、利用MATLAB编程完成计算,绘出相应图形。并与理论计算相比较,说明实验结果的原因。
四、实验内容 1、2N点实数序列
N=64。用一个64点的复数FFT程序,一次算出,形。
并绘出 的图
2、已知某序列在单位圆上的N=64等分样点的Z变换为:。
用N点IFFT程序计算出和。
五、扩展练习
例1:对连续的单一频率周期信号 按采样频率和N =16,观察其DFT结果的幅度谱。
采样,截取长度N分别选N =20解:此时离散序列算并作图,函数fft用于计算离散傅里叶变换DFT,程序如下: k=8;n1=[0:1:19];xa1=sin(2*pi*n1/k);subplot(2,2,1)plot(n1,xa1)xlabel('t/T');ylabel('x(n)');xk1=fft(xa1);xk1=abs(xk1);subplot(2,2,2)stem(n1,xk1)xlabel('k');ylabel('X(k)');n2=[0:1:15];,即k=8。用MATLAB计xa2=sin(2*pi*n2/k);subplot(2,2,3)plot(n2,xa2)xlabel('t/T');ylabel('x(n)');xk2=fft(xa2);xk2=abs(xk2);subplot(2,2,4)stem(n2,xk2)xlabel('k');ylabel('X(k)');
图 不同的截取长度的正弦信号及其DFT结果
计算结果示于图,(a)和(b)分别是N=20时的截取信号和DFT结果,由于截取了两个半周期,频谱出现泄漏;(c)和(d)分别是N=16时的截取信号和DFT结果,由于截取了两个整周期,得到单一谱线的频谱。上述频谱的误差主要是由于时域中对信号的非整周期截断产生的频谱泄漏。
实验4 离散系统的变换域分析
一、实验目的
1、熟悉对离散系统的频率响应分析方法;
2、加深对零、极点分布的概念理解。
二、实验原理
离散系统的时域方程为
其变换域分析方法如下: 频域:
系统的频率响应为:
Z域:
系统的转移函数为:
分解因式:,其中和称为零、极点。
三、预习要求
1.在MATLAB中,熟悉函数tf2zp、zplane、freqz、residuez、zp2sos的使用,其中:[z,p,K]=tf2zp(num,den)求得有理分式形式的系统转移函数的零、极点;zplane(z,p)绘制零、极点分布图;h=freqz(num,den,w)求系统的单位频率响应;[r,p,k]=residuez(num,den)完成部分分式展开计算;sos=zp2sos(z,p,K)完成将高阶系统分解为2阶系统的串联。
2.阅读扩展练习中的实例,学习频率分析法在MATLAB中的实现;
3.编程实现系统参数输入,绘出幅度频率响应和相位响应曲线和零、极点分布图。
四、实验内容
求系统 的零、极点和幅度频率响应和相位响应。
五、扩展练习
例1: 求下列直接型系统函数的零、极点,并将它转换成二阶节形式
解:用MATLAB计算程序如下:
num=[1-0.1-0.3-0.3-0.2];den=[1 0.1 0.2 0.2 0.5];[z,p,k]=tf2zp(num,den);m=abs(p);disp('零点');disp(z);disp('极点');disp(p);disp('增益系数');disp(k);sos=zp2sos(z,p,k);disp('二阶节');disp(real(sos));zplane(num,den)输入到“num”和“den”的分别为分子和分母多项式的系数。计算求得零、极点增益系数和二阶节的系数: 零点: 0.9615-0.5730-0.1443 + 0.5850i-0.1443-0.5850i 极点: 0.5276+0.6997i 0.5276-0.6997i-0.5776+0.5635i-0.5776-0.5635i 增益系数: 1 二阶节: 1.0000-0.3885-0.5509 1.0000 1.15520 0.6511 1.0000 0.28850 0.36300 1.0000-1.0552 0.7679 系统函数的二阶节形式为:
极点图见图:
图 系统函数的零、极点图
例2: 差分方程
所对应的系统的频率响应。
解:差分方程所对应的系统函数为:
用MATLAB计算的程序如下: k=256;num=[0.8-0.44 0.36 0.02];den=[1 0.7-0.45-0.6];w=0:pi/k:pi;h=freqz(num,den,w);subplot(2,2,1);plot(w/pi,real(h));grid title('实部')xlabel('omega/pi');ylabel('幅度')subplot(2,2,2);plot(w/pi,imag(h));grid title('虚部')xlabel('omega/pi');ylabel('Amplitude')subplot(2,2,3);plot(w/pi,abs(h));grid title('幅度谱')xlabel('omega/pi');ylabel('幅值')subplot(2,2,4);plot(w/pi,angle(h));grid title('相位谱')xlabel('omega/pi');ylabel('弧度')
图
实验5 有限冲激响应数字滤波器设计
一、实验目的:
1、加深对数字滤波器的常用指标理解。
2、学习数字滤波器的设计方法。
二、实验原理:
低通滤波器的常用指标:
(1)通带边缘频率;
(2)阻带边缘频率;
(3)通带起伏;
(4)通带峰值起伏,(5)阻带起伏,最小阻带衰减。
三、预习要求
1.在MATLAB中,熟悉函数fir1、kaiserord、remezord、remez的使用;
B = fir1(n,Wn,'high','noscale')设计滤波器;
[n,Wn,beta,ftype] = kaiserord(f,a,dev)估计滤波器阶数;
[n,fo,ao,w] = remezord(f,a,dev,fs)计算等波纹滤波器阶数n和加权函数w(ω);
B=remez(n,f,a)进行等波纹滤波器的设计。
2.阅读扩展练习中的实例,学习FIR滤波器的设计方法及其在MATLAB中的实现;
3.给出FIR数字滤波器的冲激响应,绘出它们的幅度和相位频响曲线,讨论它们各自的实现形式和特点。
数字滤波器有IIR和FIR两种类型,它们的特点和设计方法不同。
四、实验内容:
利用MATLAB编程,分别用窗函数法和等波纹滤波器法设计两种FIR数字滤波器,指标要求如下:
通带边缘频率:,通带峰值起伏:。
阻带边缘频率:,最小阻带衰减:。
五、扩展练习
例1: 用凯塞窗设计一FIR低通滤波器,通带边界频率,阻带衰减
不小于50dB。,阻带边界频率解: 首先由过渡带宽和阻带衰减
来决定凯塞窗的N和
上图给出了以上设计的频率特性,(a)为N=30直接截取的频率特性(b)为凯塞窗设计的频率特性。凯塞窗设计对应的MATLAB程序为: wn=kaiser(30,4.55);nn=[0:1:29];alfa=(30-1)/2;hd=sin(0.4*pi*(nn-alfa))./(pi*(nn-alfa));h=hd.*wn';[h1,w1]=freqz(h,1);或者:b = fir1(29,0.4,kaiser(30,4.55));[h1,w1]=freqz(b,1);plot(w1/pi,20*log10(abs(h1)));axis([0,1,-80,10]);grid;xlabel('归一化频率/p');ylabel('幅度/dB');还可以使用[n,Wn,beta,ftype]=kaiserord(f,a,dev)函数来估计滤波器阶数等,得到凯塞窗滤波器:
fcuts = [0.3 0.5];%归一化频率omega/pi mags = [1 0];devs = [0.05 10^(-2.5)];[n,Wn,beta,ftype] = kaiserord(fcuts,mags,devs);%计算出凯塞窗N,beta的值 hh = fir1(n,Wn,ftype,kaiser(n+1,beta),'noscale');freqz(hh);实际中,一般调用MATLAB信号处理工具箱函数remezord来计算等波纹滤波器阶数N和加权函数W(ω),调用函数remez可进行等波纹滤波器的设计,直接求出滤波器系数。函数remezord中的数组fedge为通带和阻带边界 例2:利用雷米兹交替算法设计等波纹滤波器,设计一个线性相位低通FIR数字滤波器,其指标为:通带边界频率fc=800Hz,阻带边界fr=1000Hz,通带波动At=40dB,采样频率fs=4000Hz。
阻带最小衰减解:
在MATLAB中可以用remezord 和remez两个函数设计,其结果如图2,MATLAB程序如下: fedge=[800 1000];mval=[1 0];dev=[0.0559 0.01];fs=4000;[N,fpts,mag,wt]=remezord(fedge,mval,dev,fs);b=remez(N,fpts,mag,wt);[h,w]=freqz(b,1,256);plot(w*202_/pi,20*log10(abs(h)));grid;xlabel('频率/Hz');ylabel('幅度/dB');所得图像如下所示:
实验6 无限冲激响应数字滤波器设计
一、实验目的
1、掌握双线性变换法及脉冲相应不变法设计IIR数字滤波器的具体设计方法;
2、熟悉用双线性变换法及脉冲响应不变法设计低通、高通和带通IIR数字滤波器的计算机编程。
二、实验原理
在MATLAB中,可以用下列函数辅助设计IIR数字滤波器:
1)利用buttord和cheb1ord可以确定低通原型巴特沃斯和切比雪夫滤波器的阶数和截止频率; 2)[num,den]=butter(N,Wn)(巴特沃斯)和[num,den]=cheby1(N,Wn),[num,den]=cheby2(N,Wn)(切比雪夫1型和2型)可以进行滤波器的设计;
3)lp2hp,lp2bp,lp2bs可以完成低通滤波器到高通、带通、带阻滤波器的转换; 4)使用bilinear可以对模拟滤波器进行双线性变换,求得数字滤波器的传输函数系数; 5)利用impinvar可以完成脉冲响应不变法的模拟滤波器到数字滤波器的转换。
三、预习要求
1.在MATLAB中,熟悉函数butter、cheby1、cheby2的使用,其中:
[num,den]=butter(N,Wn)巴特沃斯滤波器设计; [num,den]=cheby1(N,Wn)切比雪夫1型滤波器设计; [num,den]=cheby2(N,Wn)切比雪夫2型滤波器设计。
2.阅读扩展练习中的实例,学习在MATLAB中进行数字滤波器的设计;
3.给出IIR数字滤波器参数和滤波器的冲激响应,绘出它们的幅度和相位频响曲线,讨论它们各自的实现形式和特点。
四、实验内容
利用MATLAB编程,用脉冲响应不变法和双线性变换法设计一个数字带通滤波器,指标要求如下: 通带边缘频率:阻带边缘频率:,,通带峰值起伏:,最小阻带衰减:
。
五、扩展练习例1:设采样周期T=250μs(采样频率fs =4kHz),用脉冲响应不变法和双线性变换法设计一个三阶巴特沃兹滤波器,其3dB边界频率为fc =1kHz。[B,A]=butter(3,2*pi*1000,'s');[num1,den1]=impinvar(B,A,4000);[h1,w]=freqz(num1,den1);[B,A]=butter(3,2/0.00025,'s');[num2,den2]=bilinear(B,A,4000);[h2,w]=freqz(num2,den2);f=w/pi*202_;plot(f,abs(h1),'-.',f,abs(h2),'-');grid;xlabel('频率/Hz ')ylabel('幅值/dB')程序中第一个butter的边界频率2π×1000,为脉冲响应不变法原型低通滤波器的边界频率;第二个butter的边界频率2/T=2/0.00025,为双线性变换法原型低通滤波器的边界频率.图1给出了这两种设计方法所得到的频响,虚线为脉冲响应不变法的结果;实线为双线性变换法的结果。脉冲响应不变法由于混叠效应,使得过渡带和阻带的衰减特性变差,并且不存在传输零点。同时,也看到双线性变换法,在z=-1即Ω=π或f=2000Hz处有一个三阶传输零点,这个三阶零点正是模拟滤波器在ω=∞处的三阶传输零点通过映射形成的。
下图给出了MATLAB计算的结果。
例2: 设计一数字高通滤波器,它的通带为400~500Hz,通带内容许有0.5dB的波动,阻带内衰减在小于317Hz的频带内至少为19dB,采样频率为1,000Hz。wc=2*1000*tan(2*pi*400/(2*1000));wt=2*1000*tan(2*pi*317/(2*1000));[N,wn]=cheb1ord(wc,wt,0.5,19,'s');[B,A]=cheby1(N,0.5,wn,'high','s');[num,den]=bilinear(B,A,1000);[h,w]=freqz(num,den);f=w/pi*500;plot(f,20*log10(abs(h)));axis([0,500,-80,10]);grid;xlabel('')ylabel('幅度/dB')下图给出了MATLAB计算的结果。
例3: 设计一巴特沃兹带通滤波器,其3dB边界频率分别为f2=110kHz和f1=90kHz,在阻带f3 = 120kHz处的最小衰减大于10dB,采样频率fs=400kHz。w1=2*400*tan(2*pi*90/(2*400));w2=2*400*tan(2*pi*110/(2*400));wr=2*400*tan(2*pi*120/(2*400));[N,wn]=buttord([w1 w2],[0 wr],3,10,'s');[B,A]=butter(N,wn,'s');[num,den]=bilinear(B,A,400);[h,w]=freqz(num,den);f=w/pi*200;plot(f,20*log10(abs(h)));axis([40,160,-30,10]);grid;xlabel('频率/kHz')ylabel('幅度/dB')下图给出了MATLAB计算的结果。
例4: 一数字滤波器采样频率fs=1kHz,要求滤除100Hz的干扰,其3dB的边界频率为95Hz和105Hz,原型归一化低通滤波器为: w1=95/500;w2=105/500;[B,A]=butter(1,[w1, w2],'stop');[h,w]=freqz(B,A);f=w/pi*500;plot(f,20*log10(abs(h)));axis([50,150,-30,10]);grid;xlabel('频率/Hz')ylabel('幅度/dB')下图为MATLAB的计算结果。
实验7 设计性和研究性实验
设计性实验1 图像信号的抽取与插值
实验目的
1、熟悉图像处理常用函数和方法;
2、培养通过查阅文献解决问题的能力。实验要求
给出一个二维灰度图像,3、编程实现对该图像的任意比例的放大及缩小;
4、编程实现对该图像的任意角度旋转;
5、解决缩放及旋转时产生的锯齿等不图像不平滑问题。实验提示
6、利用上采样、下采样等方法对信号进行缩放变换;
7、观察对图像进行缩放或旋转时,图像是否会出现锯齿等不平滑现象?
8、分析产生锯齿现象的原因;
9、查阅文献了解解决锯齿现象的方法。(例如平滑滤波、双线性插值、双立方插值等处理)
设计性实验2 语音及音乐信号的采样、滤波
实验目的
1、理解采样率和量化级数对语音信号的影响;
2、设计滤波器解决实际问题。实验要求
利用电脑的声卡录一段语音信号及音乐信号,(1)观察使用不同采样率及量化级数所得到的信号的听觉效果,从而确定对不同信号的最佳的采样率;
(2)分析音乐信号的采样率为什么要比语音的采样率高才能得到较好的听觉效果;(3)注意观察信号中的噪声(特别是50hz交流电信号对录音的干扰,设计一个滤波器去除该噪声。实验提示
(1)推荐录音及播放软件:CoolEdit;
(2)分析语音及音乐信号的频谱,根据信号的频率特性理解采样定律对信号数字化的工程指导意义;
(3)可用带阻滤波器对50Hz交流电噪声进行去噪处理;
(4)也可研究设计自适应滤波器对50Hz噪声及其它随机环境噪声进行滤波处理。设计性实验3 双音多频(DTMF)信号的合成和识别
二、实验目的
1、了解电话按键音形成的原理,理解DTMF音频产生软件和DTMF解码算法;
2、利用FFT算法识别按键音;
三、实验要求
(1)设计音频产生函数,音频信号见下图,每个数据信号持续半秒;(2)实现解码函数:接受(1)产生的DTMF信号,识别信号的频率,并生成包含拨号数字的序列;
四、实验提示
(1)DTFT音频可以用两个正弦波按比例叠加产生;
(2)检测算法可以用FFT算法的DFT,或是用一组滤波器实现;
(3)Goertzel算法可以实现调谐滤波器;
设计性实验4 音乐信号处理
五、实验目的
1、了解回声的产生和梳妆滤波器;
2、混音效果的原理和均衡器的设计;
六、实验要求
(1)设计函数实现一段语音或音乐的回声产生;
(2)设计均衡器,使得得不同频率的混合音频信号,通过一个均衡器后,增强或削减某些频率区域,以便修正低频和高频信号之间的关系;
七、实验提示
(1)回声产生可以使用梳妆滤波器,y(n)=x(n)+ax(n-R), a<1(回声
zRH(z),1R1z衰减系数);或者传输函数为的全通滤波器实现;比较这两种实现方式的区别,分析为什么会有这样的区别;
(2)可以用许多一阶和二阶参数可调的滤波器级联来实现均衡器的功能,滤波器的结构选择结构要求是调整方便,最好调一个参数只影响一个应用指标,且可调参数少;
第三篇:电子科技大学成都学院数字信号处理复习重点
数字信号处理复习重点
第1章
(1)求序列的周期、Z变换和逆Z变换,会求简单序列的DTFT。掌握DTFT与Z变换的关系。
(2)判断系统是否为线性时不变系统,判断系统的稳定性与因果性。
第2章
(1)计算简单序列的DFT。掌握DFT的一些重要性质及应用(圆周共轭对称性,时域、频域循环移位性质,循环卷积)。
(2)DFT对连续信号进行频谱分析时可能出现的几个问题。
(3)掌握基-2 DFT—FFT的基本思想及特点(算法思想,运算量,运算流图,结构规则等),会计算FFT算法的运算量(加法和乘法次数),会画运算流图。
(4)计算两个有限长序列的线性卷积和循环卷积,掌握线性卷积与循环卷积的关系,线性卷积等于循环卷积的条件。
第3章
(1)映射变换的两个基本原则。
(2)脉冲响应不变法和双线性变换法的原理与特点。(优点与缺点,适合设计哪些滤波器)
(3)掌握由模拟低通原型到数字各型滤波器的设计步骤。(重点是低通、高通)
第4章
(1)FIR数字滤波器的线性相位特性。(线性相位满足的条件)
(2)窗口设计法设计FIR滤波器时怎样减少振荡,怎样使过渡带变窄。熟悉四种窗函数的特点,掌握窗长对频谱的影响。
(3)理解频率采样设计法的概念及理论依据,掌握设计步骤及要点。
(4)IIR与FIR数字滤器的比较。(掌握两类滤波器的特点)
第5章
(1)IIR滤波器的各种结构以及各种结构的特点(直接型、级联、并联)
(2)已知系统的差分方程或系统函数,用各种结构来实现滤波器。(掌握IIR直接型、级联型和并联型)
第四篇:数字信号处理课后习题Matlab作业
数字信号处理MATLAB
第1页
习题数字信号处理MATLAB习题
M1-1 已知g1(t)cos(6t),g2(t)cos(14t),g3(t)cos(26t),以抽样频率fsam10Hz对上述三个信号进行抽样。在同一张图上画出g1(t),g2(t)和g3(t)及抽样点,对所得结果进行讨论。
解:
第2页
从以上两幅图中均可看出,三个余弦函数的周期虽然不同,但它们抽样后相应抽样点所对应的值都相同。那么这样还原回原先的函数就变成相同的,实际上是不一样的。这是抽样频率太小的原因,我们应该增大抽样频率才能真实还原。如下图:f=50Hz
第3页
程序代码
f=10;
t=-0.2:0.001:0.2;g1=cos(6.*pi.*t);g2=cos(14.*pi.*t);g3=cos(26.*pi.*t);k=-0.2:1/f:0.2;h1=cos(6.*pi.*k);h2=cos(14.*pi.*k);h3=cos(26.*pi.*k);% subplot(3,1,1);
% plot(k,h1,'r.',t,g1,'r');% xlabel('t');% ylabel('g1(t)');% subplot(3,1,2);
% plot(k,h2,'g.',t,g2,'g');% xlabel('t');% ylabel('g2(t)');% subplot(3,1,3);
% plot(k,h3,'b.',t,g3,'b');% xlabel('t');% ylabel('g3(t)');
plot(t,g1,'r',t,g2,'g',t,g3,'b',k,h1,'r.',k,h2,'g.',k,h3,'b.')
第4页
xlabel('t');ylabel('g(t)');
legend('g1(t)','g2(t)','g3(t)');
M2-1 利用DFT的性质,编写一MATLAB程序,计算下列序列的循环卷积。
(1)g[k]={1,-3,4,2,0,-2,},h[k]={3,0,1,-1,2,1};(2)x[k]=cos(k/2),y[k]=3k,k=0,1,2,3,4,5。解:(1)循环卷积结果
6.0000-3.0000 17.0000-2.0000 7.0000-13.0000
程序代码
第5页
g=[1-3 4 2 0-2];h=[3 0 1-1 2 1];l=length(g);L=2*l-1;GE=fft(g,L);HE=fft(h,L);y1=ifft(GE.*HE);for n=1:l
if n+l<=L
y2(n)=y1(n)+y1(n+l);else
y2(n)=y1(n);
end end y2
stem(0:l-1,y2)xlabel('k')ylabel('y(k)')title('循环卷积')
(2)循环卷积结果
-71.0000-213.0000 89.0000 267.0000 73.0000 219.0000
第6页
程序代码
k=0:5;
x=cos(pi.*k./2);y=3.^k;l=length(x);L=2*l-1;GE=fft(x,L);HE=fft(y,L);y1=ifft(GE.*HE);for n=1:l
if n+l<=L
y2(n)=y1(n)+y1(n+l);
else
y2(n)=y1(n);
end end y2
stem(0:l-1,y2)xlabel('k')ylabel('y’(k)')title('循环卷积')
第7页
M2-2 已知序列x[k]cos(k/2N),|k|N
0,其他(1)计算序列DTFT的表达式X(ej),并画出N=10时,X(ej)的曲线。
(2)编写一MATLAB程序,利用fft函数,计算N=10时,序列x[k]的DTFT在m2m/N的抽样值。利用hold函数,将抽样点画在X(ej)的曲线上。
解:
(1)X(e)DTFT{x[k]}jkx[k]ejkkNcos(k/2N)eNjk
程序代码
N=10;k=-N:N;
x=cos(k.*pi./(2*N));W=linspace(-pi,pi,512);
第8页
X=zeros(1,length(W));for k=-N:N
X1=x(k+N+1).*exp(-j.*W.*k);X=X+X1;end
plot(W,abs(X))xlabel('W');ylabel('abs(X)');
(2)
程序代码
N=10;k=-N:N;
x=cos(k.*pi./(2*N));X_21=fft(x,21);L=-10:10;
W=linspace(-pi,pi,1024);X=zeros(1,length(W));for k=-N:N
X1=x(k+N+1).*exp(-j.*W.*k);X=X+X1;end
第9页
plot(W,abs(X));hold on;
plot(2*pi*L/21,fftshift(abs(X_21)),'o');xlabel('W');ylabel('abs(X)');
M2-3 已知一离散序列为x[k]Acos0kBcos[(0)k]。用长度N=64的Hamming窗对信号截短后近似计算其频谱。试用不同的A和B的取值,确定用Hamming窗能分辨的最小的谱峰间隔wc的值。
解:f1=100Hz f2=120Hz时
2中cN
f2=140Hz时
第10页
f2=160Hz时
第11页
由以上三幅图可见
f2=140Hz时,各谱峰可分辨。则f又
wc2N
40Hz
且
wT2fT2401 800所以c=3.2(近似值)
程序代码
N=64;L=1024;
f1=100;f2=160;;fs=800;
A=1;B1=1;B2=0.5;B3=0.25;B4=0.05;T=1/fs;ws=2*pi*fs;k=0:N-1;
x1=A*cos(2*pi*f1*T*k)+B1*cos(2*pi*f2*T*k);x2=A*cos(2*pi*f1*T*k)+B2*cos(2*pi*f2*T*k);x3=A*cos(2*pi*f1*T*k)+B3*cos(2*pi*f2*T*k);x4=A*cos(2*pi*f1*T*k)+B4*cos(2*pi*f2*T*k);hf=(hamming(N))';x1=x1.*hf;x2=x2.*hf;x3=x3.*hf;x4=x4.*hf;
X1=fftshift(fft(x1,L));X2=fftshift(fft(x2,L));X3=fftshift(fft(x3,L));X4=fftshift(fft(x4,L));
W=T*(-ws/2+(0:L-1)*ws/L)/(2*pi);subplot(2,2,1);plot(W,abs(X1));title('A=1,B=1');xlabel('W');ylabel('X1');subplot(2,2,2);
第12页
plot(W,abs(X2));title('A=1,B=0.5');xlabel('W');ylabel('X2');subplot(2,2,3);plot(W,abs(X3));title('A=1,B=0.25');xlabel('W');ylabel('X3');subplot(2,2,4);plot(W,abs(X4));title('A=1,B=0.05');xlabel('W');ylabel('X4');
M2-4 已知一离散序列为x[k]cos0k0.75cos1k,0k63。其中, 02/15,12.3/15。
(1)对x[k]做64点FFT, 画出此时信号的谱。
(2)如果(1)中显示的谱不能分辨两个谱峰,是否可对(1)中的64点信号补0而分辨出两个谱峰。通过编程进行证实,并解释其原因。
解:(1)
第13页
程序代码
W0=2*pi/15;W1=2.3*pi/15;N=64;k=0:N-1;
x=cos(W0*k)+0.75*cos(W1*k);X=fft(x);
plot(k/N,abs(X));grid on;
title('64点FFT');
(2)
第14页
第15页
由以上三幅图看出:不能对(1)中的64点信号补零而分辨出两个谱峰,这样的方法只能改变屏幕分辨率,但可以通过加hamming窗来实现对谱峰的分辨。程序代码
W0=2*pi/15;W1=2.3*pi/15;N=64;L=1024;k=0:N-1;
x=cos(W0*k)+0.75*cos(W1*k);X=fft(x,L);
plot((0:L-1)/N,abs(X));grid on;
title('1024点FFT');
M2-5 已知一连续信号为x(t)=exp(-3t)u(t),试利用DFT近似分析
第16页
其频谱。若要求频率分辨率为1Hz,试确定抽样频率fsam、抽样点数N以及持续时间Tp。
解:
本题使用矩形窗,则Nfsamfsam1fsam,Tp1 f1f
第17页
由以上三幅图可以看出当fsam越来越大时,近似值越来越接近
第18页
于实际值。即fsam越大拟合效果越好,造成的混叠也是在可以允许的范围内。程序代码
fs=100;ws=2*pi*fs;Ts=1/fs;N=fs;
x=exp(-3*Ts*(0:N-1));y=fft(x,N);l=length(y);
k=linspace(-ws/2,ws/2,l);
plot(k,Ts*fftshift(abs(y)),'b:');hold on;
w=linspace(-ws/2,ws/2,1024);y1=sqrt(1./(9+w.^2));plot(w,y1,'r')
title('fs=100Hz时的频谱')legend('近似值','实际值);
M2-6 试用DFT近似计算高斯信号g(t)exp(dt2)的频谱抽样值。
π2通过和频谱的理论值G(j)exp()比较,讨论如何根据时域的信
d4d号来恰当地选取截短长度和抽样频率使计算误差能满足精度要求。
解:
第19页
第20页
由以上三幅图可以看出:
当时域截取长度相同时,抽样间隔越小时误差越小,当抽样间隔一定时,时域截取长度越长,误差越小。当取抽样间隔为1S,时域截取长度为2S时,误差较大,绝对误差在0.5左右;当抽样间隔为0,5S,时域截取长度为2S时,误差比间隔为1S时小,绝对误差不大于0.2;当抽样间隔为0.5S时域截取长度为4S时,误差更小,绝对误差不大于0.04。因为时域截取长度越长,保留下来的原信号中的信息越多,抽样间隔越小,频谱越不容易发生混叠,所以所得频谱与理论值相比,误差更小。
程序代码
Ts=0.5;N=4;N0=64;
k=(-N/2:(N/2))*Ts;
第21页
x=exp(-pi*(k).^2);X=Ts*fftshift(fft(x,N0));
w=-pi/Ts:2*pi/N0/Ts:(pi-2*pi/N0)/Ts;XT=(pi/pi)^0.5*exp(-w.^2/4/pi);subplot(2,1,1)
plot(w/pi,abs(X),'-o',w/pi,XT);xlabel('omega/pi');ylabel('X(jomega)');
legend('试验值','理论值');
title(['Ts=',num2str(Ts)subplot(2,1,2)plot(w/pi,abs(X)-XT)ylabel('实验误差')
xlabel('omega/pi');
'N=',num2str(N)]);第22页
' '
第五篇:数字信号处理课程设计
目 录
摘要...........................................................................................................................................1 1 绪论..............................................................................................................................................2
1.1 DSP系统特点和设计基本原则......................................................................................2 1.2 国内外研究动态.............................................................................................................2 2系统设计........................................................................................................................................3 3硬件设计........................................................................................................................................5
3.1 硬件结构...........................................................................................................................5 3.2 硬件电路设计...................................................................................................................7
3.2.1 总输入电路...........................................................................................................7 3.2.2 总输出电路...........................................................................................................7 3.2.3 语音输入电路.......................................................................................................9 3.2.4 语音输出电路.......................................................................................................9 实验结果及分析.........................................................................................................................10 4.1 实验结果.........................................................................................................................10 4.2 实验分析.........................................................................................................................12 5 总结与心得体会.........................................................................................................................13 参考文献.........................................................................................................................................14 致谢................................................................................................................................................15
摘要
基于DSP的语音信号处理系统,该系统采用TMS320VC5509作为主处理器,TLV320AIC23B作为音频芯片,在此基础上完成系统硬件平台的搭建和软件设计,从而实现对语音信号的采集、滤波和回放功能,它可作为语音信号处理的通用平台。
语音是人类相互之间进行交流时使用最多、最自然、最基本也是最重要的信息载体。在高度信息化的今天,语音信号处理是信息高速公路、多媒体技术、办公自动化、现代通信及智能系统等新兴领域应用的核心技术之一。通常这些信号处理的过程要满足实时且快速高效的要求,随着DSP技术的发展,以DSP为内核的设备越来越多,为语音信号的处理提供了良好的平台。本文设计了一个基于TMS320VC5509定点的语音信号处理系统,实现对语音信号的采集、处理与回放等功能,为今后复杂的语音信号处理算法的研究和实时实现提供一个通用平台。
关键词:语音处理;DSP;TMS320VC5509;TLV320AIC23B
1 绪论
语音是人类相互间所进行的通信的最自然和最简洁方便的形式,语音通信是一种理想的人机通信方式。语音通信的研究涉及到人工智能、数字信号处理、微型计算机技术、语言声学、语言学等许多领域,所以说语音的通信是一个多学科的综合研究领域,其研究成果具有重要的学术价值。另外通过语音来传递信息是人类最重要的、最有效、最常用的交换信息的形式。语言是人类特有的功能,声音是人类常用的工具,是相互传递信息的主要手段。同时也是众构成思想交流和感情沟通的最主要的途径。
1.1 DSP系统特点和设计基本原则
DSP(digital signal processor)是一种独特的微处理器,是以数字信号来处理大量信息的器件。其工作原理是接收模拟信号,转换为0或1的数字信号。再对数字信号进行修改、删除、强化,并在其他系统芯片中把数字数据解译回模拟数据或实际环境格式。它不仅具有可编程性,而且其实时运行速度可达每秒数以千万条复杂指令程序,远远超过通用微处理器,是数字化电子世界中日益重要的电脑芯片。它的强大数据处理能力和高运行速度,是最值得称道的两大特色。
1.2 国内外研究动态
语音信号处理作为一个重要的研究领域,已经有很长的研究历史。但是它的快速发展可以说是从1940年前后Dudley的声码器和Potter等人的可见语音开始的;20世纪60年代中期形成的一系列数字信号处理的理念和技术基础;到了80年代,由于矢量量化、隐马尔可夫模型和人工神经网络等相继被应用于语音信号处理,并经过不断改进与完善,使得语音信号处理技术产生了突破性的进展。一方面,对声学语音学统计模型的研究逐渐深入,鲁棒的语音识别、基于语音段的建模方法及隐马尔可夫模型与人工神经网络的结合成为研究的热点。另一方面,为了语音识别实用化的需要,讲者自适应、听觉模型、快速搜索识别算法以及进一步的语言模型的研究等课题倍受关注。
在通信越来越发达的当今世界,尤其最近几十年,语音压缩编码技术在移动 通信、IP电话通信、保密通信、卫星通信以及语音存储等很多方面得到了广泛的应用。因此,语音编码一直是通信和信号处理的研究热点,并其取得了惊人的进展,目前在PC机上的语音编码已经趋于成熟,而如何在嵌入式系统中实时实现语音压缩编码则是近些年来语音信号处理领域的研究热点之一。
2系统设计
在实际生活中,当声源遇到物体时会发生反射,反射的声波和声源声波一起传输,听者会发现反射声波部分比声源声波慢一些,类似人们面对山体高声呼喊后可以在过一会儿听到回声的现象。声音遇到较远物体产生的反射会比遇到较近的反射波晚些到达声源位置,所以回声和原声的延迟随反射物体的距离大小改变。同时,反射声音的物体对声波的反射能力,决定了听到的回声的强弱和质量。另外,生活中的回声的成分比较复杂,有反射、漫反射、折射,还有回声的多次反射、折射效果。
当已知一个数字音源后,可以利用计算机的处理能力,用数字的方式通过计算模拟回声效应。简单的讲,可以在原声音流中叠加延迟一段时间后的声流,实现回声效果。当然通过复杂运算,可以计算各种效应的混响效果。如此产生的回声,我们称之为数字回声。
本次实验的程序流程图如下:
图2.1 程序流程图
本次实验的系统框图如下:
图2.2 系统框图
3硬件设计
3.1 硬件结构
图3.1是系统的硬件结构框图, 系统主要包括VC5509和A IC23 两个模块。
图3.1系统硬件结构框图
利用VC5509 的片上外设I2C(Inter-Integrated Circuit, 内部集成电路)模块配置AIC23 的内部寄存器;通过VC5509 的McBSP(Multi channel Buffered Serial Ports, 多通道缓存串口)接收和发送采样的音频数据。控制通道只在配置AIC23 的内部寄存器时工作, 而当传输音频数据时则处于闲置状态。
AIC23通过麦克风输入或者立体声音频输入采集模拟信号, 并把模拟信号转化为数字信号, 存储到DSP的内部RAM中,以便DSP处理。
当DSP完成对音频数据的处理以后, AIC23再把数字信号转化为模拟信号, 这样就能够在立体声输出端或者耳机输出端听到声音。
AIC23能够实现与VC5509 DSP的McBSP端口的无缝连接, 使系统设计更加简单。接口的原理框图, 如下图所示。
图3.2 AIC23与VC5509接口原理图
系统中A IC23的主时钟12 MHz直接由外部的晶振提供。MODE接数字地, 表示利用I2 C控制接口对AIC23传输控制数据。CS接数字地, 定义了I2 C总线上AIC23的外设地址, 通过将CS接到高电平或低电平, 可以选择A IC23作为从设备在I2 C总线上的地址。SCLK和SDIN是AIC23控制端口的移位时钟和数据输入端,分别与VC5509的I2C模块端口SCL和SDA相连。
收发时钟信号CLKX1和CLKR1由A IC23的串行数据输入时钟BCLK提供, 并由A IC23的帧同步信号LRCIN、LRCOUT启动串口数据传输。DX1和DR1分别与A IC23 的D IN 和DOUT 相连, 从而完成VC5509与AIC23间的数字信号通信。
3.2 硬件电路设计
3.2.1 总输入电路
图3.3 总输入电路
从左到右各部分电路为:
话筒,开关,语音输入电路,UA741高增益放大电路,有源二阶带 通滤波器。
3.2.2 总输出电路
图3.4 总输出电路
从左到右各部分电路为:
LM386高频功率放大器及其外围器件连接电路,语音输出电路,开关,扬声器。
3.2.3 语音输入电路
图3.5语音输入电路
3.2.4 语音输出电路
图3.6 语音输出电路
语音信号通道包括模拟输入和模拟输出两个部分。模拟信号的输入输出电路如图所示。上图中MICBIAS 为提供的麦克风偏压,通常是3/4 AVDD,MICIN为麦克风输入,可以根据需要调整输入增益。下图中LLINEOUT 为左声道输出,RLINEOUT为右声道输出。用户可以根据电阻阻值调节增益的大小,使语音输入输出达到最佳效果。从而实现良好的模拟语音信号输入与模拟信号的输出。4 实验结果及分析
4.1 实验结果
按“F5”键运行,注意观察窗口中的bEcho=0,表示数字回声功能没有激活。这时从耳机中能听到麦克风中的输入语音放送。将观察窗口中bEcho的取值改成非0值。这时可从耳机中听到带数字回声道语音放送。
分别调整uDelay和uEffect的取值,使他们保持在0-1023范围内,同时听听耳机中的输出有何变化。
当uDelay和uEffect的数值增大时,数字回声的效果就会越加的明显。
图4.1 修改前程序图
图4.2 修改前程序图
图4.3 频谱分析
图4.4 左声道及右声道波形 4.2 实验分析
所以,从本实验可知当已知一个数字音源后,可以利用计算机的处理能力,用数字的方式通过计算模拟回声效应。简单的讲,可以在原声音流中叠加延迟一段时间后的声流,实现回声效果。当然通过复杂运算,可以计算各种效应的混响效果。
声音放送可以加入数字回声,数字回声的强弱和与原声的延迟均可在程序中设定和调整。5 总结与心得体会
通过本次课程设计,我明白了细节决定成败这句话的道理,在实验中,有很多注意的地方,都被忽视了,导致再花费更多的时间去修改,这严重影响了试验的进度。同时,在本次实验中我了解了ICETEK – VC5509 – A板上语音codec芯片TLV320AIC23的设计和程序控制原理,并进一步掌握了数字回声产生原理、编程及其参数选择、控制,以及了解了VC5509DSP扩展存储器的编程使用方法。
这一学期的理论知识学习加上这次课程设计,使我对DSP有了更加深刻的了解,对数字信号的处理功能,软硬件相结合,语音信号的采集与放送等等方面都有了很深的了解,相信本次课程设计,无论是对我以后的学习,还是工作等方面都有一个很大的帮助。因此,本次课程设计让我受益匪浅。
参考文献
[1]李利.DSP原理及应用[M].北京:中国水利水电出版社,202_.[2]王安民,陈明欣,朱明.TMS320C54xxDSP实用技术[M].北京:清华大学出版社,202_ [3]彭启琮,李玉柏.DSP技术[M].成都:电子科技大学出版社,1997 [4]李宏伟,等.基于帧间重叠谱减法的语音增强方法[J].解放军理工大学学报,202_(1):41~44 [5]TexasInstrumentsIncorporated.TMS320C54x系列DSP的CPU与外设[M].梁晓雯,裴小平,李玉虎,译.北京:清华大学出版社,202_ [6]赵力.语音信号处理[M].北京:机械工业出版社,202_比较图4和图5,可以看到1200Hz以上的频谱明显得到了抑制。
[7]江涛,朱光喜.基于TMS320VC5402的音频信号采集与系统处理[J].电子技术用,202_,28(7):70~72[8]TexasInstrumentsIncorporated:TMS320VC5402Datasheet,202_
致谢
在本次课程设计的即将完成之际,笔者的心情无法平静,本文的完成既是笔者孜孜不倦努力的结果,更是指导老师樊洪斌老师亲切关怀和悉心指导的结果。在整个课程设计的选题、研究和撰写过程中,老师都给了我精心的指导、热忱的鼓励和支持,他的精心点拨为我开拓了研究视野,修正了写作思路,对课程设计的完善和质量的提高起到了关键性的作用。另外,导师严谨求实的治学态度、一丝不苟的工作作风和高尚的人格魅力,都给了学生很大感触,使学生终生受益。在此,学生谨向老师致以最真挚的感激和最崇高的敬佩之情。
另外,还要感谢这段时间来陪我一起努力同学,感谢我们这个小团队,感谢每一个在学习和生活中所有给予我关心、支持和帮助的老师和同学们,几年来我们一起学习、一起玩耍,共同度过了太多的美好时光。我们始终是一个团结、友爱、积极向上的集体。