无论是身处学校还是步入社会,大家都尝试过写作吧,借助写作也可以提高我们的语言组织能力。相信许多人会觉得范文很难写?以下是我为大家搜集的优质范文,仅供参考,一起来看看吧
现代信号处理课后答案篇一
细胞识别
目录
第一部分
1、实验课题名称--------------------3
2、实验目的--------------------------3
第 1 页
3、实验内容概要--------------------3 第二部分
1、建立工程文件--------------------3
2、图像信息获取--------------------4
3、如何建立下拉菜单--------------6
4、标记mark点----------------------6
5、二值化9
6、填洞---9
7、收缩---10
8、获取中心点------------------------11
9、细胞计数---------------------------13
10、all-steps---------------------------13
11、扩展功能-------------------------14 第三部分
12、各步骤结果和错误举例------16 第四部分
13、心得体会--------------------------22
第一部分
1、实验课题:细胞识别
2、实验目的:对血液细胞切片图片进行各种处理,最终得出细胞的数目、面积等信息。
3、实验内容概要:基于vc++6.0软件下的细胞识别,通过细胞的标记、二
第 2 页 值化、填洞、收缩、获取中心点、计数等过程完成实验目的。
第二部分——实验具体步骤
1、建立工程文件
① 新建mfc工程项目:--mfc appwizard、工程名
② 拷贝cdib.h,到工程文件夹,再向工程里添加
③ doc.h添加变量:m_lpdib 和头文件#include”cdib.h”
④ :变量(m_lpdib)的new、delete
第 3 页
⑤ : serialize()
2、图像信息获取
① : ondraw()m_pdib->draw()如果图像不为空的话,那么就执行如下主要代码:
② 点击键,建立类向导,在messages中添加oninitialupdate()函数,添加代码实现对自动打开固定图片。
③ 通过鼠标右击,点击建立类向导,在messages中添加onmousemove()函数,添加代码实现获取所要信息,即实现鼠标在图像任一位置移动时可以直观的读取相对应位置的信息。可以在屏幕上显示鼠标所指点的坐标以及rgb、hsi和灰度值,通过hsi的可以选取合适的阈值来找到细胞以及边界。
第 4 页 ④ 为了rgb图像转化为人眼更容易识别的hsi模型,我们可以通过添加成员函数rgbtohsi来实现这一功能。
hsi模型与rgb模型的转化关系
(添加函数时,可以右击类窗口中的view.h,选中add member function,之后选择函数的返回值类型和函数描述,其它默认不变)
确定后在里面添加实现函数功能的代码。
3、添加下拉菜单
在resourceview那栏的找到菜单按键设置
第 5 页 双击,后在里面添加所需按键
每个按键的id号为注意在填写为idr_加菜单大写。
之后右击按键,建立类向导添加按键所需函数
4、标记mark 分为四步
1.找出mark(red)点和maybemark(blue)点
2.将maybemark(blue)点变成mark(black red)点
3.将mark(black red)点变成edge(yellow(fullred&&fullgreen))点 点滤波
基本思想:mark点指的是我们要寻找的细胞内的点。我们先获取每一个像素点的rgb分量,然后我们将其转化成hsi分量,将h分量进行归一化,因为s的范围是0到1,所以我们要进行尺度的一致,这样才具有可计算性。然后我们通过每个像素点的h分量和s分量的值与细胞内部的h分量和s分量计算欧几里得距离,设定一个mark门限值(我们这里将markdoor设置为0.09,大家可以行设置合适的参数),小于这个门限值我们就当做是细胞的内部,然后对细胞进行标记(red)。还需要设定一个maybe mark门限值(我们这里将maybe markdoor设置为0.15,大家可以行设置合适的参数),我们大于mark门限值小于maybe mark门限值时,我们暂时看成是细胞,我们进行maybe mark的标记(blue)。否则的话,我们需要考虑,一些不是mark和maybe mark点的*lpsrc==0我们区别一下赋值为1,*lpsrc==255
第 6 页 我们区别一下赋值为254,*(lpsrc+1)==255我们区别一下赋值为254.这样的话,我们在后面判断是否为mark点的时候,我们只需要判断*lpsrc是否为0就可以了,判断maybe mark点时只需要判断*(lpsrc)是否为255就可以了。对于边缘的判断只需要判断*(lpsrc+1)是否为255就可以了。
将细胞标记为mark用红色(255,0,0)标记出来,将可能是的细胞标记为maybe mark用蓝色(0,0,255)标记出来。将maybemark to mark的区域用亮红(128,0,0)表示,将不可能是细胞的区域、细胞边界分别用绿色标记出来。操作过程:(1)根据h、s的欧几里得距离sqrt(s2+h2)来大致的确定哪些是细胞(mark)和可能是细胞(maybe mark)的点。
(2)根据maybe mark点周围的情况,如果它的上下左右四个方向有mark点,则将maybe mark点变成mark点。
(3)用sobel算子来做边缘的提取边界(0,255,255)(255,255,0),使用3*3的模板,使用欧几里得距离来判断是否为边缘。
两种sobel算子如下:
第 7 页
主要代码如下:
doubletmp1=pixel[0]+2*pixel[1]+pixel[2]-pixel[6]-2*pixel[7]-pixel[8];doubletmp2=pixel[0]+2*pixel[3]+pixel[6]-pixel[2]-2*pixel[5]-pixel[8];double edge=sqrt(tmp1*tmp1+tmp2*tmp2);
if(edge>edgedoor)*(lpdst+1)=255;//sobel判断该点是否edge//edgedoor=40(4)edge滤波
就是去除全边缘点(四周都是背景或边缘)(强度为5)
5、二值化
基本思想:将原有彩色图像变换为二值图像,其中细胞0x80(128)用gray(灰色)标记出来,边缘0xf0(240)用bright(亮色)标记出来,其他表示为0。主要代码:
第 8 页
6、填洞
将细胞中或者细胞相邻的地方的较小的背景填成细胞的背景,填完的细胞背景的灰度值是129,因为都被访问过了,然后将边缘去掉。
填洞的基本思想:首先将细胞或边缘内的黑点置为vistied = 0x01,以该黑点为中心,在其上下左右侧进行访问是否有未访问的黑点,若有则将上方黑点压栈,且上下左右侧的黑点置已访问。将堆栈顶端的数据弹出,作为新的种子进行扩散,即以该元素为基点,判断其周围是否存在未访问黑点,若有则继续压栈,重复操作。直到找到最后一点,此点四周均不存在未访问黑点,结束访问。若洞像素数小于100大于50,洞内像素数及其初进栈的点(56,(409,222))时,则进行填洞。填洞的过程就是将非mark点转化为mark点。
主要标记访问代码:
填洞函数主要代码分析:
填完洞后,进行下面操作:
如果图像中只有已访问黑点0x01则将其恢复成0;如果图像中只有edge点0xf0则将 edge置为黑点。这样图像中只有黑色的背景以及灰色的细胞mark(0x80)点。
主要代码:
第 9 页
7、收缩
收缩的目的是为了方便计数。通过扫描图像,对图像进行预先的3次腐蚀,判断所生成边界点,然后根据原理判定是否标注该点,存放所标志的中心点,便于统计细胞个数及计算细胞半径。
由mark生成边界,我们有四邻域生成边界和八邻域生成边界。判断该点是否为mark点,如果是mark点的话,我们判断i、j是否是我们选取图片的边界,如果是的话,我们将该点变成边缘点,否则我们判断它的上下左右(周围八个点)是否有非mark点,如有有,则将这边变成边缘点,反之,不变。
8邻域收缩操作代码(4邻域与8邻域思想相同):
第 10 页
8、获取中心点
根据前面所作工作统计获得的中心点个数,去掉一系列不符合要求的点得出最终的细胞个数、细胞的平均半径和平均面积,用对话框输出统计结果。
操作过程如下:
① 首先我们要去除访问标志,是我们先前一次在判断是否需要保存点的时候(markit(int i, int j)),我们将边缘点都标记成访问过了,这时在处理下一次遍历图片发现中心点的时候,我们要进行判断点是否要保存就没有办法做了,所以在没进行一次图片的遍历之前我们都需要去除访问标志。主要代码:*lpsrc&=no_visited;//0xfe// 清除visited标志
最后位 置0操作
② 需要判断是否是边界以外的点,这里我们只处理边界内部的点,对于边界外部的点不加以查找中心点。对于内部的点,我们先要判断是否是孤立的边缘点,即判断该边缘点的上下左右四个点都不是mark点和边缘点我们认为是孤立的边缘点,但是我们在这里也要去除半径不大于2的孤立点,因为我们认为它的半径太小,是噪声。如果是半径大于2的孤立点,我们对他进行标记成中心点,对半径做一点补偿(=k+pre_shrink_count+4,4为补偿)。然后在入队。
主要代码如下:
{
第 11 页 if(k<3)// 如果进行第一次收缩即消失的点则认为该点是噪点,不进行保存直接进行下一次收缩
continue;
// 孤立的点
*lpsrc |=centered;//0x2 对孤立点加上中心点标志0xf2//后面shrink时为0x02
// 保存一下center_point信息(圆心,半径)
pt.x=i;pt.y=j;=k + pre_shrink_count + 4;// +4放大补偿,k为把此圆收缩到一点所经历的收缩次数
_back(pt);continue;}
③ 需要判断是否需要保存该点,我们在判断它的上下左右是否有没有访问过的边缘点,这里我们运用递归函数来找相连通的边缘点,如果是全边缘点的时候我们就需要保存,即将m_bfulledge=1。在保存的函数中我们将该点变成中心点,然后半径补偿然后我们来对该点入队(保存),然后将该点设置成没有访问过的点,因为下面在做该点的上下左右是否为访问过的,访问过才保存,因为这个点已经保存过了,所以要将保存过的点设置成没有访问过,我们必须将这点变成没有访问过。然后判断该点的周围八点是否是访问过的,如果访问过就保存该点,这里也是运用递归函数来实现的。
对中心点处理:
a)独立的中心点直接存储
b)相邻的中心点通过递归求质心作为圆心,最大半径作为新的半径,合并各中心为一点
主要代码如下:
pt.x=tot_x/tot_num;//质心 pt.y=tot_y/tot_num;=max_radius;//取最大半径作为质点中心点半径 *(lpsrc-(pt.y-j)*llinebytes+pt.x-i)|=centered;//质点置为中心点 _back(pt);c)相近但不相邻的点,求质心为圆心,最大半径为半径(直到无相近点)
主要代码如下:
if(abs(x0-x)+abs(y0-y)<10)// 相近
{
/*pt=(j);if((i).radius
pt=(i);*/ (i).x=(x+x0)/2;//取均值,保存最大半径
第 12 页
(i).y=(y+y0)/2;(i).radius=max((i).radius,(j).radius)+4;pt=pt=(i);d)在无相近点的情况下,若半径小于8,则删除。
主要代码如下:
if(bdelete)
{
} e)两圆相交,若其中一圆非相交部分面积小于50%,则删除
主要代码如下:
if(total
9、细胞计数
打开我们处理前的图片,根据前面保存中心点的队列,我们知道中心点的位置和细胞的半径,然后我们重新的导入细胞的图片,在上面画圆,标出细胞。然后我们获取细胞内部的hsi的最大值和最小值,计算出细胞的平均面积和个数。主要代码如下:
(“共有%d个细胞,平均半径%d,平均面积%d : h(%3.1f,%3.1f)s(%3.2f,%3.2f)i(%3.2f,%3.2f)”, (),(int)(totr/()+.5),(int)(tota/()+.5), 360.0*min[0]/255.0,360.0*max[0]/255.0, 1.0*min[1]/255.0,1.0*max[1]/255.0, 1.0*min[2]/255.0,1.0*max[2]/255.0);
10、all-steps
第 13 页 可以一次性实现细胞识别的所有操作步骤
设置控制按键的权限,点击update_command_ui ,键入控制条件
每步操作时给cellprocess设置不同数值,表示那步进行过,只能进行规定的下步操作,从而在运行过程中放置按键误触导致程序崩溃。
11、扩展:
区域选择:
第 14 页
建立类向导:onbuttondown 和onbuttonup 键入代码:
在ondraw中添加下列代码
注意:bool 文件头部,定义在view.h会出现第一次区域选择时出现错误。
添加复位按键:
第 15 页
点击该键后会重新读取图像(和图像自动打开代码一样)
第三部分
12、各步骤结果和错误举例
① 各步骤结果图
(red)& maybe mark(blue)
maybemark to mark(black red)第 16 页
edge information and edge filter
twovalue
fillholes 第 17 页
shrink
findcenter 第 18 页
count
出现的错误举例:
mousemove 程序中出现问题: 1.错误:
没有加#include “mainfrm.h”头文件 2.错误
第 19 页
error c2248: 'm_wndstatusbar' : cannot access protected member declared in class 'cmainframe'
需将protected: // control bar embedded members cstatusbar m_wndstatusbar;ctoolbar
m_wndtoolbar;
protected变为public 供用户操作使用
mousemove函数中((cmainframe*)afxgetmainwnd())->etext(0,str);中内容显示在标准窗口图像的左下方bar
ove的坐标判决放在for循环外,鼠标移动到图像外,程序会崩溃。解决:改变坐标判决代码的位置
通见问题:
在classview中的类视图不见了,文件仍然存在
解决方案:先保存workspace,然后关闭工程,文件,重新打开workspace 原因:classview显示混乱 在类中添加的成员变量和成员函数不能显示出来,即使显示出来了变量或函数,双击后不能跳至正确的位置。
第 20 页
edge information中出现问题
正常
不正常
memcpy(lpnewdibbits,lpsrc,lheight *llinebytes);代码应放在图像处理前,参考图像是初始状态的图像,把第一步的四个小步骤分开写在不同的函数内时,因为每一小步的操作都会改变图像的状态,如果把: memcpy(lpnewdibbits,lpsrc,lheight *llinebytes);写在maybemark_to_mark之后那么参考图像就不是原始图像
发现:做到shrink时,看到收缩后的图像效果很差,和标准收缩图像相差较大,经调试后发现问题(没注意 ppt最后一页有,老师在qq群在中也提到过)。
shrink操作后,关闭图像,出现问题 genedge4()函数中出现问题 for(int j=0;j { lpsrc =(unsigned char *)pdoc->m_lpdib->m_lpimage + llinebytes*(lheight-1-j)-1; for(int i=0;i { lpsrc++;......和 for(int j=0;j for(int i=0;i { 第 21 页 lpsrc =(unsigned char *)pdoc->m_lpdib->m_lpimage + llinebytes*(lheight-1-j)-i;....for循环问题 有差异 后面起始点地址为dot1 前面起始点地址为 dot1+1,不处理四周边界 可以把第二个for循环语句中前面的i=0 改成i=1.观察图像,发现shirnk后未保存填洞点,结果只有 0x00 0x80 0xf0 // *lpsrc &=no_edge_point;//至关重要 少写后,程序一直崩溃。 这段代码在genedge函数中,在shrink步骤中没影响,但在找中心点过程中用到这个函数时,这段代码就十分必要。否则程序会在运行findcenter时直接崩溃。 中心点标志 0x02 不是0xf2 (n).y 等同于 temp[n].y 13、心得体会 通过对本次数字图像处理课程设计的学习,进一步加深了对数图知识的理解,同时也基本掌握了vc++软件的使用方法。从一开始的连图像都无法打开,到最后在老师上课的资料以及同学的帮助和学长的参考程序下终于完成了一幅细胞图像的整个识别过程。自己一向在编程方面有所欠缺,通过咨询老师和同学还有百度,自己也慢慢理解了所写程序代码的含义,间接地提高了自己写代码与识别代码的能力。 通过一个星期多的学习,我对细胞识别的基本思想有了深一步的理解,也让我对c语言相关的知识得到了回顾。此次课程设计给我们提供了一个既能学习又能锻炼的机会,使我们养成了查找资料(主要是在百度上查阅一些代码的含义)的习惯,将理论与实际相结合起来,锻炼了分析问题和实际解决问题的能力。提高了适应能力,为今后的学习和实践打下了基础。 第 22 页 中南大学 课程设计报告 题 目 现代信号处理课程设计 学生姓名 万义武 指导教师 周扬、支国明 学 院 信息科学与工程学院 学 号 0909118219 专业班级 电子信息专业1102班 一、课程设计题目 1、信号发生器 用户根据测试需要,可任选以下两种方式之一生成测试信号:(1)直接输入(或从文件读取)测试序列;(2)输入由多个不同频率正弦信号叠加组合而成的模拟信号公式(如式 1-1 所示)、采样频率(hz)、采样点数,动态生成该信号的采样序列,作为测试信号。12 100sin(2)100sin(2)100sin(2)n f t f t f t (1-1) 2、频谱分析 使用 fft 对产生的测试信号进行频谱分析并展示其幅频特性与相频特性,指定需要滤除 的频带,通过选择滤波器类型(iir / fir),确定对应的滤波器(低通、高通)技术指标。 3、滤波器设计 根据以上技术指标(通带截止频率、通带最大衰减、阻带截止频率、阻带最小衰减),设 计数字滤波器,生成相应的滤波器系数,并画出对应的滤波器幅频特性与相频特性。(1)iir df 设计:可选择滤波器基型(巴特沃斯或切比雪夫型); (2)fir df 设计:使用窗口法(可选择窗口类型,并比较分析基于不同窗口、不同阶数 所设计数字滤波器的特点)。 4、数字滤波 根据设计的滤波器系数,对测试信号进行数字滤波,展示滤波后信号的幅频特性与相频特 性,分析是否满足滤波要求(对同一滤波要求,对比分析各类滤波器的差异)。(1)iir df:要求通过差分方程迭代实现滤波(未知初值置零处理); (2)fir df:要求通过快速卷积实现滤波(对于长序列,可以选择使用重叠相加或重叠 保留法进行卷积运算)。 5、选做内容 将一段语音作为测试信号,通过频谱展示和语音播放,对比分析滤波前后语音信号的变化,进一步加深对数字信号处理的理解。 二、设计过程 《1》、第一、二题:(1).信号发生器。 ①直接输入(或从文件读取)测试序列; ②输入由多个不同频率正弦信号叠加组合而成的模拟信号公式。 ③使用fft对产生的测试信号进行频谱分析并展示其幅频特性与相频特性。(2).源代码 t=(0:0.00001:1);n=0:100;f1=50;y=sin(2*pi*f1*t);f=input('please f=');t=1/f;x=sin(2*pi*f1*n*t);m=fft(x);h=abs(m);figure(1);subplot(321)plot(t,y);subplot(322)stem(n,x,'.');title('xulitu');subplot(323)plot(n,h);title('fupintu');subplot(324)xi=interp1(n,x,t*f1,'linear');plot(t,xi);title('chongjiantu'); (3)结果 (4)分析: 采样原理:对模拟信号进行采样可以看作是一个模拟信号通过一个电子开关s。设电子开关每隔周期t合上一次,每次合上的时间为τ,在电子开关输出端得到其采样信号,一般τ很小, τ越小,采样输出脉冲的幅度越接近输入信号在离散时间点上的瞬时值。 《2》、第三、四题 (1)题目(滤波器设计与数字滤波) 滤波器设计—根据输入的数字滤波器的技术指标,包括通带截止频率,通带最大衰减,阻带截止频率,阻带最小衰减,设计滤波器,生成相应的滤波器系数,并画出对应的滤波器幅频、相频特性。iir df设计:可选择滤波器基型(巴特沃斯或切比雪夫型); (2)源代码 i=input('please input i(choose fuction)=');switch fix(i) case {1}%低通滤波 wp=input('please input wp=');ws=input('please input ws=');ap=input('please input ap=');as=input('please inout as=');fs=1;t=1/fs; wp1=(2/t)*tan(wp/2);ws1=(2/t)*tan(ws/2); [n,wn]=buttord(wp1,ws1,ap,as,'s');[b,a]=butter(n,wn,'s');[bz,az]=bilinear(b,a,fs);w=linspace(0,2*pi,1000);h=freqz(bz,az,w); subplot(311) plot(w(1:500)/pi,abs(h(1:500)));grid; title(['n=',num2str(n)]);text(0.1,0.8,['b=',num2str(bz)]);text(0.1,0.4,['a=',num2str(az)]);xlabel('w/pi');ylabel('幅度(db)'); subplot(312)plot(w/pi,angle(h)); xlabel('w/pi');ylabel('相位');grid; subplot(313)y=real(ifft(h));x=0:999;plot(x,y); title('单位脉冲响应');grid;clear; case{2}%高通滤波 wp=input('please input wp=');ws=input('please input ws=');ap=input('please input ap=');as=input('please inout as=');fs=1;t=1/fs; wp1=(2/t)*tan(wp/2);ws1=(2/t)*tan(ws/2); [n,wn]=buttord(wp1,ws1,ap,as,'s');[b,a]=butter(n,wn,'high','s');[bz,az]=bilinear(b,a,fs);w=linspace(0,2*pi,1000);h=freqz(bz,az,w); subplot(311) plot(w(1:500)/pi,abs(h(1:500)));grid; title(['n=',num2str(n)]);text(0.1,0.9,['b=',num2str(bz)]);text(0.1,0.4,['a=',num2str(az)]);xlabel('w/pi');ylabel('幅度(db)'); subplot(312)plot(w/pi,angle(h)); xlabel('w/pi');ylabel('相位');grid; subplot(313)y=real(ifft(h));x=0:999;plot(x,y); title('单位脉冲响应 ');grid;clear; case{3}%带通滤波 wpl=input('please input wpl=');wph=input('please input wph=');wsl=input('please input wsl=');wsh=input('please input wsh=');ap=input('please input ap=');as=input('please inout as=');wp=[wpl,wph];ws=[wsl,wsh];fs=1;t=1/fs; wp2=(2/t)*tan(wp/2);ws2=(2/t)*tan(ws/2); [n,wn]=buttord(wp2,ws2,ap,as,'s');[b,a]=butter(n,wn,'s');[bz,az]=bilinear(b,a,fs);w=linspace(0,2*pi,1000);h=freqz(bz,az,w); subplot(311) plot(w(1:500)/pi,abs(h(1:500)));grid; title(['n=',num2str(n)]);text(0.1,1.2,['b=',num2str(bz)]);text(0.1,0.4,['a=',num2str(az)]);xlabel('w/pi');ylabel('幅度(db)'); subplot(312)plot(w/pi,angle(h)); xlabel('w/pi');ylabel('相位');grid; subplot(313)y=real(ifft(h));x=0:999;plot(x,y); title('单位脉冲响应 ');grid;clear; end (3)结果 低通滤波 高通滤波 带通滤波(4)分析 用双线性变换法设计无限脉冲响应数字滤波器(iif df)时,先把数字滤波器指标转换成模拟滤波器的指标,然后根据模拟滤波器的指标设计模拟滤波器,再经过线性变换把模拟滤波器转换成数字滤波器。该系统要能够设计巴特沃兹型低通、带通、高通滤波器,并能够输入数字滤波器的性能指标,显示出滤波器的阶数和系数。该系统的关键部分是滤波器的设计部分,按照双线性变换法设计滤波器的步骤进行设计即可。 三、设计总结与心得体会 在课程设计的这段时间,我获益匪浅。不但进一步掌握了数字信号处理的基础知识及matlab的基本操作。虽然在做的过程中遇到了一些问题,但都通过自己的努力解决了它们。这次课程设计对我各方面的综合能力有了很大的提高,对我以后的实践都有很大的帮助。 本次课程设计不但让我又学到了一些知识,而且也提高了我的综合能力。使我在各个方面都得到了锻炼,以后有这样的机会一定会更加的很好利用,它不仅可以提高学习的针对性而且可以很好的锻炼动手能力以及自己的逻辑设计能力和处理问题的能力,希望在以后这方面的能力会很好的加强。 四、课程设计指导书 [1] 《数字信号处理(第二版)》.丁玉美等 西安电子科技大学出版社 [2] 《数字信号处理及其matlab实现》,陈怀琛等译,电子工业出版社; [3] 《matlab及在电子信息课程中的应用》,陈怀琛等,电子工业出版社 五、鸣谢 此次的课程社真心感谢那些为我们提供良好的上机环境已经良好的知道的老师们。同时也感谢中南大学给了我这一次检验自己的动手能力以及发现自己错误的机会! (一).信号分析 1、编制信号生成程序,产生下述各序列,绘出它们的时域波形 1)单位抽样序列 (n) 2)矩形序列 rn(n) 3)三角波序列 n1,0n3x3(n)8n,4n7 0,其它 4)反三角波序列 4n,0n3x4(n)n3,4n7 0,其它 5)gaussian序列 (np) q,0n15x5(n)e 0,其它2 6)正弦序列 x6(n)sin16t 取 fs64hz,n16 7)衰减正弦序列 (t)aesin(2ft)u(t)对连续信号x70进行采样,可得到测试序列 x 7(n)ae antf 0 nt)sin(2u。令(n)a=50,采样周期t=1ms,即fs=1000hz,f0=62.5,a=100。 2.对上述信号完成下列信号分析 1)对三角波序列x3(n)和反三角波序列x4(n),作n=8点的fft,观察比较它们的幅频特 性,说明它们有什么异同?绘出两序列及其它们的幅频特性曲线。at在x3(n)和x4(n)的尾部补零,作n=16点的fft,观察它们的幅频特性发生了什么变化? 分析说明原因。 2)、观察高斯序列x5(n),固定信号x5(n)中的参数p=8,令q分别等于2,4,8,观察它们的时域和幅频特性,了解当q取不同值时,对信号序列的时域幅频特性的影响;固定q=8,令p分别等于8,13,14,观察参数p变化对信号序列的时域及幅频特性的影响,观察p等于多少时,会发生明显的泄漏现象,混叠是否也随之出现?记录实验中观察到的现象,绘出相应的时域序列和幅频特性曲线。 3)对于正弦序列x4(n),取数据长度n分别等于8,16,32,分别作n点fft,观察它们的的时域和幅频特性,说明它们的差别,简要说明原因。 4)、观察衰减正弦序列x7(n)的时域和幅频特性,绘出幅频特性曲线,改变采样频率fs,使 fs=300hz,观察此时的频谱的形状和谱峰出现位置?说明产生现象的原因。 3.设有一连续时间信号s(t),其由20hz、220hz和750hz的正弦信号叠加而成,分析确定采样频率及数据分析长度,计算并绘出信号的频谱,指出各个频率份量。 你们先自己看一下matlab的书,对照书上的例题仿真一下,多练习。 先给出信号分析部分的题目给你们,你们可以先做做,最好使用gui,将所有的部分集成在一起。滤波器部分的题目开学后再给你们,如果matlab熟练了,那部分做起来很快的。 如果题目中的公式看不到的话,可能是公式编辑器的版本问题,我采用的是公式编辑器5.2 追求完美。他还告诫在场的师生:“每个清华人都负有责任,建设这个国家。为学,要扎扎实实,不可沽名钓誉。做事,要公正廉洁,不要落身后骂名。” 信息科学与工程学院 数字信号处理课程设计实验报告 课题名称: 简单信号滤波演示系统 学生姓名: 学 号: 专业班级: 指导老师: 实验时间: 2014.10.8 目 录 第一章 概述.................................3 第二章 第三章 第四章第五章第六章 1.1 fir、iir概述.................................3 1.2题目要求......................................3 设计分析.............................5 2.1算法分析......................................5 2.2 在matlab中实现的分析........................6 程序实现.......................8 3.1 程序主体介绍..................................8 3.2 子程序........................................9 3.3 程序调试及运行结果............................9 3.4 结果分析及问题分析...........................16 心得体会............................17 参考文献............................18 源代码..............................19 matlab 第一章 概述 1.1 fir、iir概述 数字滤波器是指输入输出均为数字信号,通过一定的运算关系改变输入信号所含频率成分的相对比例或者滤除某些频率成分的器件。数字滤波器与模拟滤波器相比数字滤波器具有精度高、稳定、体积小、重量小、灵活等特点。主要分为两种:有限脉冲响应fir和无限脉冲响应iir。设计滤波器的主要要求有两种,一是幅频特性,一是相频特性。一般的滤波器主要是对幅频特性作出要求,如果对输出相频特性也有要求,就需要用到线性相位滤波器。iir滤波器的设计主要有两类,一是借助于模拟滤波器设计进行,二是直接在频域或时域中进行设计。fir滤波器的设计不能借助于模拟滤波器,也有两类设计方法,一是窗函数法,二是频率采样法。还有一种比较有效的方法是切比雪夫等波纹逼近法,需通过计算机辅助进行。 1.2题目要求 设计一个工作流程如图所示的信号滤波演示系统: 图2 滤波演示图 ⑴ 信号发生器—根据信号选择分为两大类: ① 静态型:直接输入(或从文件读取)测试信号序列; ② 动态型:输入由多个不同频率正弦信号叠加组合而成的模拟信号公式 100sin(2f1t)100sin(2f2t)...100sin(2fnt)采样频率(hz)以及采样点数,动态生成该信号的采样序列,作为测试信号。 ⑵ 频谱分析—使用fft 对产生的测试信号进行频谱分析并展示其幅频、相频特 性,指定需要滤除或保留的频带,通过选择滤波器类型(iir/fir),确定对 应的滤波器(低通、高通、带通、带阻)技术指标。 ⑶ 滤波器设计—根据iir/fir 数字滤波器技术指标设计滤波器,生成相应的滤 波器系数,并展示对应的滤波器幅频(衰减)、相频特性。 ① iir df 设计:使用双线性变换法,可选择滤波器基型(巴特沃斯或切比雪夫型); ② fir df 设计:使用窗口法,可选择窗口类型。 ⑷ 数字滤波—根据设计的滤波器系数,对测试信号进行滤波。 ① iir df:要求通过差分方程迭代实现滤波,未知初值置零处理; ② fir df:要求通过快速卷积实现滤波。可以选择使用重叠相加或重叠保留法进行卷积运算,并动态展示卷积运算的详细过程。 ⑸ 输出信号分析—展示滤波后信号的幅频与相频特性,分析是否满足滤波要求。 对同一滤波要求,根据输出信号频谱,对比分析各类滤波器的差异。选做部分 将一段语音作为测试信号,通过频谱展示和语音播放,对比分析滤波前后语音信号的变化,进一步加深对数字信号处理的理解。 第二章 设计分析 2.1算法分析 此题目的实现可分为三个某块的设计实现:输入信号模块,设计滤波器模块,滤波模块。 首先明确各模块间的数据依赖关系:在输入信号模块得到信号后,对信号 进行频域分析,从而确定滤波器的相关技术指标,根据滤波器的技术指标与类型,在滤波器设计模块完成滤波器的设计,然后将滤波器的设计结果传递给滤波模块,滤波模块同时接收输入信号,从而通过运算,实现信号的滤波处理。从数据传递关系上分析,滤波模块的实现依赖于其他模块的数据输出,因此放在最后设计。先设计输入模块。 因为此设计相对复杂,分模块设计,通过参数传递和接口实现分模块设计即检验,提高程序的稳定性与健壮性。 输入的实现可以有两种方式:静态输入和动态输入。静态输入选择从文本输入数据,将信号取样值以矩阵的形式存放在文本中。采用文件的读取就可以实现,比较容易。动态输入,即输入由一系列频率的正弦信号相加的组成的信号,需要经过采样的,注意在设置采样频率时一定要符合奈奎斯特准则,提高采样点数,增加频谱分辨率。最后输出一采样信号向量,传递给其余两模块。 滤波器的设计,通过输入信号的频谱分析,设置滤波器的参数,然后才可以设计滤波器。第一步需要总结设计滤波器需要哪些参数,通过学习可以总结,所有滤波器的参数有四个:通带截止频率、阻带截止频率、通带最大衰减、阻带最小衰减。 对滤波器的设计分两类:fir和iir,二则所需的函数及设计方法不同。iir采用借助于模拟滤波器的方式,包括巴特沃斯滤波器和切比雪夫滤波器两种类型。fir采用窗函数方式,有矩形窗、三角窗、汉明窗、汉宁窗、布莱克曼和凯塞窗。通过调用不同的函数来实现滤波器的设计。特别在实现窗函数滤波器时,各个函数的主要区别是不同的频率采样,可以通过选择结构实现,简化程序。通过滤波器的设计最后可以得到滤波器的系统函数h(z)的系数。分析滤波器的幅频特性和相频特性,如果不符合要求重新设定滤波器参数或者换成其他滤波器类型。如果性能符合要求,则将系数传递给滤波模块。滤波模块不调用滤波函数,实现滤波功能根据滤波器类型的不同,有两种方式可以选择,一种是通过差分方程运算,一种是通过线性卷积运算。前者适合对iir滤波器进行滤波,后者适合对fir滤波器进行滤波。且线性卷积为实现对长序列的卷积运算,采用重叠相加法,采用动态变化展示运算过程。 2.2 在matlab中实现的分析 输入模块通过读取文件和直接输入数据运算可以很容易实现。在输入信号确定后,要对其进行频谱分析,从而确定滤波器的参数和类型(低通、高通、带通、带阻),此模块作用也就完成,将数据分别用全局变量传递给下面的模块。 设置模拟信号: 我采用的测试信号是两个正弦信号叠加而成的信号,信号为 y=sin(2*pi*f1*n*ts)+sin(2*pi*f2*n*ts)其中频率f1=30;频率f2=50;频率f3=200;采样频率fs=3000;采样间隔 ts=1/fs;采样点数 n=1024;n=1:n 采集模拟信号的程序代码: f1=30;% 频率1 f2=50;% 频率2 f3=200;% 频率3 fs=3000;% 采样频率 ts=1/fs;% 采样间隔 n=1024;% 采样点数 n=1:n; y=sin(2*pi*f1*n*ts)+sin(2*pi*f2*n*ts)+sin(2*pi*f1*n*ts);% 正弦波混合 频谱分析: 使用fft对产生的测试信号进行频谱分析并展示其幅频特性与相频特性,指定需要滤除的频带,通过选择滤波器类型iir,确定对应的滤波器(低通、高通)技术指标fp、fc、rp、rs。 滤波器的设计: 根据以上技术指标(通带截止频率、通带最大衰减、阻带截止频率、阻带最 小衰减),设计数字滤波器,生成相应的滤波器系数,并画出对应的滤波器 幅频特性与相频特性。分别设计巴特沃斯滤波器、切比雪夫i型滤波器、切 比雪夫ii型滤波器、椭圆滤波器、yulewalk滤波器。 巴特沃斯和切比雪夫的滤波函数调用为: [n,wc]=buttord(wp,ws,rp,rs);[n,wc]=cheb1ord(wp,ws,rp,rs);[b,a]=butter(n,wc,’property’);[b,a]=cheby1(n,rp,wc,’property’); property对于低通和高通为’’,带通’high’,带阻’stop’;窗函数滤波器设计的调用函数: 求窗函数的阶数: n=ceil((h*pi)/wdel);%wdel为窗函数的过渡带宽,h对应不同窗函数的值 wn=boxcar(n+1);wn=triang(n+1);wn=hanning(n+1);wn=hamming(n+1);wn=blackman(n+1);wn=kaiser(n+1,beta);%bata为kaiser的a参数 b=fir1(n,ws,'property',wn);property对于低通和高通为’’,带通’high’,带阻’stop’; 数字滤波: 根据设计的滤波器系数,对测试信号通过差分方程迭代实现滤波数字滤波,展示滤波后信号的幅频特性与相频特性,分析是否满足滤波要求(对同一滤波要求,对比分析各类滤波器的差异)。 需要注意的是,窗函数对滤波参数的使用,需要通过运算得到一窗函数阶数,对通带最大衰减无使用。需要注意fir和iir对参数的不同利用。 滤波模块:运用差分方程运算的滤波函数,可以用循环调用简单实现。动态展示线性卷积的重叠相加法可以用流程图来说明: 第三章 matlab程序实现 3.1 程序主体介绍 此程序界面也是采用gui,不过是采用图形用户界面开发环境。和第一个对比,布局更加容易,很容易调整界面。但是在一些函数操作行会受一些限制,比如函数参数的传递,目前来说还是采用全局变量才能完成。 pushbutton1_callback(hobject, eventdata, handles)%获得输入信号 uipanel11_selectionchangefcn(hobject, eventdata, handles)%获取滤波器的设计参数 uipanel15_selectionchangefcn(hobject, eventdata, handles)%设计滤波器 pushbutton4_callback(hobject, eventdata, handles)%进行滤波 以上是一些功能控件,需要注意的是他们的view ballcak属性,对于groupbutton一组按钮,需要选用selectionchangefcn,其余选择callback属性即可。另外,handles函数作用为传递的句柄,是一结构体,调用时注意此处。其余和直接用函数创建gui界面无太大区别。 3.2 子程序 3.2.1 iir滤波器 fir滤波器 巴特沃斯滤波器 3.2.2 矩形窗 三角窗 汉宁窗 哈明窗 3.3 程序调试及运行结果 原始信号幅频特性及相频特性 巴特沃斯的幅频特性及相频特性 信号滤波后的幅频特性及相频特性 矩形窗的幅频特性及相频特性 信号滤波后的幅频特性及相频特性 三角窗的幅频特性及相频特性 信号滤波后的幅频特性及相频特性 汉宁窗的幅频特性及相频特性 信号滤波后的幅频特性及相频特性 哈明窗的幅频特性及相频特性 信号滤波后的幅频特性及相频特性 3.4 结果分析及问题分析 用双线性变换法设计无限脉冲响应数字滤波器(iif df)时,先把数字滤波器指标转换成模拟滤波器的指标,然后根据模拟滤波器的指标设计模拟滤波器,再经过线性变换把模拟滤波器转换成数字滤波器。该系统要能够设计巴特沃兹型低通、带通、高通滤波器,并能够输入数字滤波器的性能指标,显示出滤波器的阶数和系数。该系统的关键部分是滤波器的设计部分,按照双线性变换法设计滤波器的步骤进行设计即可。 第四章 心得体会 本次课程设计可以称为matlab和数字信号处理的综合设计,因为这次有一半的时间在研究matlab,有一半的时间在思考数字信号问题的解决。几天下来收获是很多的,对数字信号处理有了一次全面的回顾,而且也看到了matlab的神秘面容。只要是自己真正的努力去做了,就绝对不会后悔在课程设计上花费的时间。现在发现自己确实会的太少,而这次又学了太多,对自己以后的学习还是有鞭策意义的。 这次课程设计也使我认识到了matlab的强大,或者进一步说工具软件有着你意想不到的功效,能节省你的时间,同时又能够让你从实践上更深层次的认识到所学知识的奇妙,而且可以同时明白了学习与实践的相辅相成。对课程设计也是一直保持很高兴趣的,虽然有时为它焦头烂额,但是也会因此有柳暗花明的喜悦,任何事情都折射着一个普通的道理——付出才有回报。浅显的道理,却是需要我们用毅力来坚持见证的。 具体来说,通过本次课程设计,我掌握了matlab的基本运用,尤其是其中gui可视化图形用户界面的设计,包括函数调用与在图形用户界面开发环境下的调用。函数的调用与参数的传递是两个简单却很容易出错的地方,自由细心才可以保证程序的健壮性。 对数字信号本身内容的理解来说,全书的内容确实是可以融合在这两个课程设计题目中的,一个是dft运用,另一个是滤波器的设计和利用。对全书的内容可以做出最好的概括。其中滤波器的设计中,调用了许多滤波函数,这是本次实验有点欠缺的地方,但是仍然从整体上把握了整个滤波器的设计过程。 对于课程设计中出现的问题,解决的方式有两种,一是自己解决。可以上网,查阅图书和matlab的本身的帮助文件,不断尝试,自己修改错误,可以更好的反省。二是与同学相商。一加一不是等于二那么简单的,相互交流才是进步的最好方式。其中解决问题最重要的应该是坚持不懈,在遇到问题时,不放弃才有可能称为最后的胜利者。 每次课程设计都收获颇多,而且每次都对自己的学习态度做一次调整,这个的作用确实是很大的,希望以后更加珍惜的课程设计的时间。 第五章 参考文献 [1] 丁玉美等.数字信号处理 [m].西安:西安电子科技大学出版社,2002 [2] 程佩青.数字信号处理教程,第二版[m].北京:清华大学出版社,2001 [3] 赵树杰等.数字信号处理[m].西安:西安电子科技大学出版社,1997 [4] 及在电子信息课程中的应用[m],北京:电子工业出版社出 版,2002 [5] 余成波.数字信号处理及matlab实现,第一版[m].北京:清华大学出版社,2008 [6] l signal processing: a computer based approach, 3rdedition[m],new york, usa:mcgraw-hill,2000 [7] tanding digital signal processing,2nd edition[m].new jersey, usa:prentice hall,2005 第六章 源代码 function varargout = df(varargin) gui_singleton = 1; gui_state = struct('gui_name', mfilename,...'gui_singleton', gui_singleton,...'gui_openingfcn', @df_openingfcn,...'gui_outputfcn', @df_outputfcn,...'gui_layoutfcn', [] ,...'gui_callback', []);if nargin && ischar(varargin{1}) _callback = str2func(varargin{1});end if nargout [varargout{1:nargout}] = gui_mainfcn(gui_state, varargin{:});else gui_mainfcn(gui_state, varargin{:});end function df_openingfcn(hobject, eventdata, handles, varargin) = hobject; guidata(hobject, handles); function varargout = df_outputfcn(hobject, eventdata, handles) varargout{1} = ; function edit1_callback(hobject, eventdata, handles) function edit1_createfcn(hobject, eventdata, handles) if ispc set(hobject,'backgroundcolor','white');else set(hobject,'backgroundcolor',get(0,'defaultuicontrolbackgroundcolor'));end function edit2_callback(hobject, eventdata, handles) function edit2_createfcn(hobject, eventdata, handles) if ispc set(hobject,'backgroundcolor','white');else set(hobject,'backgroundcolor',get(0,'defaultuicontrolbackgroundcolor'));end function edit3_callback(hobject, eventdata, handles) function edit3_createfcn(hobject, eventdata, handles) if ispc set(hobject,'backgroundcolor','white');else set(hobject,'backgroundcolor',get(0,'defaultuicontrolbackgroundcolor'));end function edit4_callback(hobject, eventdata, handles) function edit4_createfcn(hobject, eventdata, handles)if ispc set(hobject,'backgroundcolor','white');else set(hobject,'backgroundcolor',get(0,'defaultuicontrolbackgroundcolor'));end function edit5_callback(hobject, eventdata, handles) function edit5_createfcn(hobject, eventdata, handles)if ispc set(hobject,'backgroundcolor','white');else set(hobject,'backgroundcolor',get(0,'defaultuicontrolbackgroundcolor'));end function pushbutton1_callback(hobject, eventdata, handles)f=zeros(1,4); f(1)=str2num(get(2,'string'));f(2)=str2num(get(1,'string'));f(3)=str2num(get(3,'string'));f(4)=str2num(get(4,'string')); t=str2num(get(12,'string'));%采样周期% fre=str2num(get(5,'string'));%采样频率% t=0:1/fre:t; xn=10*sin(2*pi*t*f(1))+10*sin(2*pi*t*f(2))+10*sin(2*pi*t*f(3))+10*sin(2*pi*t*f(4)); set(15,'userdata',xn);%将xn放在用户数据userdata% yn3=abs(fft(xn));%快速傅立叶变换(符频特性)% n1=[0:length(yn3)-1]/length(yn3)*fre;%横坐标% axes(12);%坐标系编号% stem(n1,yn3,'.'); axis([0,fre/2,0,max(yn3)]);%坐标轴单位控制% title('信号的幅频特性');xlabel('频率(hz)');ylabel('|x(ejw)|'); yn4=angle(fft(xn));%相频特性% n4=[0:length(yn4)-1]/length(yn4)*fre;axes(16);stem(n4,yn4,'.'); axis([0,fre/2,min(yn4),max(yn4)]);title('信号的相频特性');xlabel('频率(hz)');ylabel('相位'); function edit6_callback(hobject, eventdata, handles) function edit6_createfcn(hobject, eventdata, handles) if ispc set(hobject,'backgroundcolor','white');else set(hobject,'backgroundcolor',get(0,'defaultuicontrolbackgroundcolor'));end function edit7_callback(hobject, eventdata, handles) function edit7_createfcn(hobject, eventdata, handles) if ispc set(hobject,'backgroundcolor','white');else set(hobject,'backgroundcolor',get(0,'defaultuicontrolbackgroundcolor'));end function edit8_callback(hobject, eventdata, handles) function edit8_createfcn(hobject, eventdata, handles) if ispc set(hobject,'backgroundcolor','white');else set(hobject,'backgroundcolor',get(0,'defaultuicontrolbackgroundcolor'));end function edit9_callback(hobject, eventdata, handles) function edit9_createfcn(hobject, eventdata, handles) if ispc set(hobject,'backgroundcolor','white');else set(hobject,'backgroundcolor',get(0,'defaultuicontrolbackgroundcolor'));end function pushbutton2_callback(hobject, eventdata, handles)xn=get(15,'userdata');fre=str2num(get(5,'string')); fp=str2num(get(6,'string'));%通带最大频率 % fs=str2num(get(7,'string'));%阻带最小频率% rp=str2num(get(8,'string'));rs=str2num(get(9,'string'));wp=2*fp/fre;ws=2*fs/fre; [n,wc]=buttord(wp,ws,rp,rs,'s');%求阶数,截止频率% [b,a]=butter(n,wc,'s');%巴特沃兹模拟低通滤波器系数% [bz,az]=bilinear(b,a,1); [h,w]=freqz(bz,az);%分析数字滤波器% axes(14);plot(w/pi,abs(h)); axis([0,1,0,max(abs(h))]);title('巴特沃兹的幅频特性');xlabel('w/π'); ylabel('|x(ejw)|');axes(17);plot(w/pi,angle(h));title('巴特沃兹的相频特性');xlabel('频率(hz)');ylabel('相位'); yn2=filter(bz,az,xn);%迭代法求解滤波信号% yn=abs(fft(yn2)); n1=[0:length(yn)-1]/length(yn)*fre;axes(15);stem(n1,yn,'.'); title('信号滤波后的幅频特性');xlabel('频率(hz)');ylabel('|x(ejw)|');axis([0,fre/2,0,max(yn)]);yn=angle(fft(yn2)); n1=[0:length(yn)-1]/length(yn)*fre;axes(18);stem(n1,yn,'.'); title('信号滤波后的相频特性');xlabel('频率(hz)'); ylabel('相位'); axis([0,fre/2,min(yn),max(yn)]); pushbutton1_callback(hobject, eventdata, handles); function edit10_callback(hobject, eventdata, handles)function edit10_createfcn(hobject, eventdata, handles) if ispc set(hobject,'backgroundcolor','white');else set(hobject,'backgroundcolor',get(0,'defaultuicontrolbackgroundcolor'));end function edit11_callback(hobject, eventdata, handles) function edit11_createfcn(hobject, eventdata, handles)if ispc set(hobject,'backgroundcolor','white');else set(hobject,'backgroundcolor',get(0,'defaultuicontrolbackgroundcolor'));end function pushbutton5_callback(hobject, eventdata, handles) xn=get(15,'userdata');fp=str2num(get(10,'string'));fs=str2num(get(11,'string'));fre=str2num(get(5,'string'));wp=2*fp/fre*pi;ws=2*fs/fre*pi; bt=ws-wp; n0=ceil(1.8*pi/bt);n=n0+mod(n0+1,2);wc=(wp+ws)/2/pi; hn=fir1(n-1,wc,boxcar(n));%单位序列相应% yn=abs(fft(hn));yn=20*log10(yn); n1=[0:length(yn)-1]/length(yn)*2;axes(14);plot(n1,yn); title('矩形窗的损耗函数');xlabel('w/π'); ylabel('20log|hg(w)|');axis([0,1,min(yn),max(yn)]);axes(17);yn=angle(hn); n1=[0:length(yn)-1]/length(yn)*2;plot(n1,yn); title('矩形窗的相频特性');xlabel('w/π'); ylabel('相位'); axis([0,1,min(yn),max(yn)]); yn2=fftfilt(hn,xn,512);%重叠相加法求滤波后序列% yn=abs(fft(yn2));n1=[0:length(yn)-1]/length(yn)*fre;axes(15);stem(n1,yn,'.'); title('信号滤波后的幅频特性');xlabel('频率(hz)');ylabel('|x(ejw)|');axis([0,fre/2,0,max(yn)]);yn=angle(fft(yn2)); n1=[0:length(yn)-1]/length(yn)*fre;axes(18);stem(n1,yn,'.'); title('信号滤波后的相频特性');xlabel('频率(hz)'); ylabel('相位'); axis([0,fre/2,min(yn),max(yn)]); pushbutton1_callback(hobject, eventdata, handles); function pushbutton6_callback(hobject, eventdata, handles) xn=get(15,'userdata');fp=str2num(get(10,'string'));fs=str2num(get(11,'string'));fre=str2num(get(5,'string'));wp=2*fp/fre*pi;ws=2*fs/fre*pi; bt=ws-wp; n0=ceil(6.1*pi/bt);n=n0+mod(n0+1,2);wc=(wp+ws)/2/pi; hn=fir1(n-1,wc,bartlett(n)); yn=abs(fft(hn));yn=20*log10(yn); n1=[0:length(yn)-1]/length(yn)*2;axes(14);plot(n1,yn); title('三角窗的损耗函数');xlabel('w/π'); ylabel('20log|hg(w)|');axis([0,1,min(yn),max(yn)]);axes(17);yn=angle(hn);n1=[0:length(yn)-1]/length(yn)*2;plot(n1,yn); title('三角窗的相频特性');xlabel('w/π'); ylabel('相位'); axis([0,1,min(yn),max(yn)]); yn2=fftfilt(hn,xn,512);yn=abs(fft(yn2)); n1=[0:length(yn)-1]/length(yn)*fre;axes(15);stem(n1,yn,'.'); title('信号滤波后的幅频特性');xlabel('频率(hz)');ylabel('|x(ejw)|');axis([0,fre/2,0,max(yn)]);yn=angle(fft(yn2)); n1=[0:length(yn)-1]/length(yn)*fre;axes(18);stem(n1,yn,'.'); title('信号滤波后的相频特性');xlabel('频率(hz)'); ylabel('相位'); axis([0,fre/2,min(yn),max(yn)]); pushbutton1_callback(hobject, eventdata, handles); function pushbutton7_callback(hobject, eventdata, handles) xn=get(15,'userdata');fp=str2num(get(10,'string'));fs=str2num(get(11,'string'));fre=str2num(get(5,'string'));wp=2*fp/fre*pi;ws=2*fs/fre*pi;if(fp n0=ceil(6.2*pi/bt);n=n0+mod(n0+1,2);wc=(wp+ws)/2/pi; hn=fir1(n-1,wc,hanning(n));else bt=wp-ws; n0=ceil(6.2*pi/bt);n=n0+mod(n0+1,2);wc=(wp+ws)/2/pi; hn=fir1(n-1,wc,'high',hanning(n));end yn=abs(fft(hn));yn=20*log10(yn); n1=[0:length(yn)-1]/length(yn)*2;axes(14);plot(n1,yn); title('汉宁窗的损耗函数');xlabel('w/π'); ylabel('20log|hg(w)|');axis([0,1,min(yn),max(yn)]);axes(17);yn=angle(hn); n1=[0:length(yn)-1]/length(yn)*2;plot(n1,yn); title('汉宁窗的相频特性');xlabel('w/π'); ylabel('相位'); axis([0,1,min(yn),max(yn)]); yn2=fftfilt(hn,xn,512);yn=abs(fft(yn2)); n1=[0:length(yn)-1]/length(yn)*fre;axes(15);stem(n1,yn,'.'); title('信号滤波后的幅频特性');xlabel('频率(hz)');ylabel('|x(ejw)|');axis([0,fre/2,0,max(yn)]);yn=angle(fft(yn2)); n1=[0:length(yn)-1]/length(yn)*fre;axes(18);stem(n1,yn,'.'); title('信号滤波后的相频特性');xlabel('频率(hz)'); ylabel('相位'); axis([0,fre/2,min(yn),max(yn)]); pushbutton1_callback(hobject, eventdata, handles); function pushbutton8_callback(hobject, eventdata, handles) xn=get(15,'userdata');fp=str2num(get(10,'string'));fs=str2num(get(11,'string'));fre=str2num(get(5,'string'));wp=2*fp/fre*pi;ws=2*fs/fre*pi;if(fp n0=ceil(6.2*pi/bt);n=n0+mod(n0+1,2);wc=(wp+ws)/2/pi; hn=fir1(n-1,wc,hamming(n));else bt=wp-ws; n0=ceil(6.2*pi/bt);n=n0+mod(n0+1,2);wc=(wp+ws)/2/pi; hn=fir1(n-1,wc,'high',hamming(n));end yn=abs(fft(hn));yn=20*log10(yn); n1=[0:length(yn)-1]/length(yn)*2;axes(14);plot(n1,yn); title('哈明窗的损耗函数');xlabel('w/π'); ylabel('20log|hg(w)|');axis([0,1,min(yn),max(yn)]);axes(17);yn=angle(hn); n1=[0:length(yn)-1]/length(yn)*2;plot(n1,yn); title('哈明窗的相频特性');xlabel('w/π'); ylabel('相位'); axis([0,1,min(yn),max(yn)]); yn2=fftfilt(hn,xn,512);yn=abs(fft(yn2)); n1=[0:length(yn)-1]/length(yn)*fre;axes(15);stem(n1,yn,'.'); title('信号滤波后的幅频特性');xlabel('频率(hz)');ylabel('|x(ejw)|');axis([0,fre/2,0,max(yn)]);yn=angle(fft(yn2)); n1=[0:length(yn)-1]/length(yn)*fre;axes(18);stem(n1,yn,'.'); title('信号滤波后的相频特性');xlabel('频率(hz)'); ylabel('相位'); axis([0,fre/2,min(yn),max(yn)]); pushbutton1_callback(hobject, eventdata, handles); %---executes during object creation, after setting all on axes12_createfcn(hobject, eventdata, handles)% hobject handle to axes12(see gcbo) % eventdata reservedhandles not created until after all createfcns called % hint: place code in openingfcn to populate axes12 %---executes on mouse press over axes on axes12_buttondownfcn(hobject, eventdata, handles)% hobject handle to axes12(see gcbo) % eventdata reservedto be defined in a future version of matlab % handles structure with handles and user data(see guidata) % hints: get(hobject,'string')returns contents of edit12 as text % str2double(get(hobject,'string'))returns contents of edit12 as a double %---executes during object creation, after setting all on edit12_createfcn(hobject, eventdata, handles)% hobject handle to edit12(see gcbo) % eventdata reservedhandles not created until after all createfcns called % hint: edit controls usually have a white background on windows.% see ispc and ispc && isequal(get(hobject,'backgroundcolor'), get(0,'defaultuicontrolbackgroundcolor'))set(hobject,'backgroundcolor','white');end %---if enable == 'on', executes on mouse press in 5 pixel border.%---otherwise, executes on mouse press in 5 pixel border or over on edit5_buttondownfcn(hobject, eventdata, handles)% hobject handle to edit5(see gcbo) % eventdata reservedto be defined in a future version of matlab % handles structure with handles and user data(see guidata) %---executes during object deletion, before destroying on edit12_deletefcn(hobject, eventdata, handles)% hobject handle to edit12(see gcbo) % eventdata reserved-to be defined in a future version of matlab % handles structure with handles and user data(see guidata)现代信号处理课后答案篇二
现代信号处理课后答案篇三
现代信号处理课后答案篇四
目录
一、课程设计目的要求 ………………………………………… 1
二、课程设计选题 ……………………………………………… 1
三、程序分析及运行果 …………………………………………
1、第一题 ……………………………………………………
2、第三题 ……………………………………………………
3、第五题 ……………………………………………………
4、第七题 ……………………………………………………
四、课程设计总结 ………………………………………………
参考文献 …………………………………………………………
附.程序源代码 ……………………………………………………
2 4 6 9 12 13 14
一、课程设计目的要求
1.全面复习课程所学理论知识,巩固所学知识重点和难点,将理论与实践很好地结合起来。
2.提高综合运用所学知识独立分析和解决问题的能力; 3.熟练使用一种高级语言进行编程实现。
二、课程设计选题
本次课程设计,按任务书要求选择了1、3、5、7题,第8题为选作,因为能力有限,所以并未选择,现将题目及要求摘录如下:
1.给定模拟信号:xa(t)aeatsin(0t)u(t),式中a444.128,502,0502 rad/s。对xa(t)进行采样,可得采样序列
x(n)xa(nt)aeantsin(0nt)u(nt)1)选择采样频率fs=1 khz,观测时间tp50ms,观测所得序列x(n)及其幅频特性|x(ejw)|
2)改变采样频率fs=300hz,观测此时|x(ejw)|的变化 3)令采样频率fs=200hz,观测此时|x(ejw)|的变化
要求分析说明原理,绘出相应的序列及其它们对应的幅频特性曲线,指出|x(ejw)|的变化,说明为什么?
3.一个连续信号含两个频率分量,经采样得
x(n)=sin2π*0.125n+cos2π*(0.125+δf)n n=0,1……,n-1 已知n=16,δf分别为1/16和1/64,观察其幅频特性;当n=128时,δf不变,其结果有何不同,为什么?分析说明原因,并打印出相应的幅频特性曲线
119n)cos(n),使用fft分析其频谱: 5.一个序列为x(n)0.5cos(20201)使用不同宽度的矩形窗截短该序列为m点长度,取m分别为: a)m=20 b)m=40 c)m=160 ;观察不同长度的窗对谱分析结果的影响; 2)使用汉宁窗、哈明窗重做1)
3)对三种窗的结果进行理论分析及比较。并绘出相应的幅频特性曲线
df的设计
分别利用矩形窗、汉宁窗、哈明窗设计一个n=11的线性相位fir 低通和高通数字滤波器,截止频率c3rad,要求:求出各滤波器的单位脉冲响应h(n);绘出各滤波器的幅频及相频响应曲线;观察各滤波器的通带波纹和阻带波纹;比较不同窗函数对滤波特性的影响。
三、程序分析及运行结果
程序综述及gui界面:程序全部利用matlab语言写成,通过与gui模块的按钮控件pushbutton链接,使之能够通过按动按钮跳出显示窗口来实现程序的运行以及最后结果的呈现,简洁明了。程序最终的gui界面如下图所示:
1、第一题:
设计思想:此程序按第一题题目要求,预先给定了模拟信号at(t)aesin(0t)u(t),同时给出了相关参数,a444.128,502,xa0502按照题目要求直接在程序中定义相关变量并赋值。题目还rad/s。要求采取三种不同的采样频率分别为1khz、300hz和200hz,所以t1=0.001,t2=1/300,t3=0.005。fft函数对序列进行傅立叶变换,得到幅频特性曲线。程序通过subplot、stem、plot等函数将结果显示在弹出的窗口中。这样第一题的要求就基本实现了。程序代码:
a=444.128;%设定参数 w=50*sqrt(2)*pi;a=50*sqrt(2)*pi;t1=0.001;n1=50;n=0:49;x1=a*exp(-a*n*t1).*sin(w*n*t1);%所给信号 figure;subplot(3,2,1);stem(x1,'.');%坐标函数及现实函数 xlabel('n');ylabel('x(n)');title('x(n)序列1');y1=fft(x1,n1);%傅立叶变换 subplot(3,2,2);plot(abs(y1));xlabel('(n=50 wk=2pik/n)k');ylabel('|x1(jw)|');title('幅频特性|x1(jw)|');
t2=1/300;n2=15;n=0:14 x2=a*exp(-a*n*t2).*sin(w*n*t2);subplot(3,2,3);stem(x2,'.');xlabel('n');ylabel('x(n)');title('x(n)序列2');y2=fft(x2,n2);subplot(3,2,4);plot(abs(y2));xlabel('(n=15 wk=2pik/n)k');ylabel('|x2(jw)|');title('幅频特性|x2(jw)|');t3=0.005;n3=10;n=0:9;x3=a*exp(-a*n*t3).*sin(w*n*t3);subplot(3,2,5);stem(x3,'.');xlabel('n');ylabel('x(n)');title('x(n)序列3');y3=fft(x3,n3);subplot(3,2,6);plot(abs(y3));xlabel('(n=10 wk=2pik/n)k');ylabel('|x3(jw)|');title('幅频特性|x3(jw)|');
运行结果:
结果分析:从运行结果看,随着采样频率的减少,幅频特性的失真越明显。所以采样频率越低,越容易失真。
2、第三题:
设计思想:用fft进行谱分析时,首先必须对信号进行采样,使之变成离散信号,然后就可按照前面的方法用fft来对连续信号进行谱分析。按采样定理,采样频率fs应大于2倍信号的最高频率,为了满足采样定理,一般在采样之前要设置一个抗混叠低通滤波器。
按题目要求,讲程序分为n=16和n=128两部分。首先定义出n、f1、f2等参数,接着列出题目所要求的信号,然后利用傅立叶变换函数,最后通过坐标以及输出函数得出相应幅频特性曲线。程序代码:
%第三题(n=16)n=16;n=0:15;f1=1/16;f2=1/64;
x1=sin(2*pi*0.125*n)+cos(2*pi*(0.125+f1)*n);y1=fft(x1,16);%傅立叶变换函数 figure;
subplot(2,2,1);stem(n,abs(y1),'.');
axis([0 15 0 8]);%标注坐标长度
xlabel('(n=16 wk=2pik/n)k');ylabel('|x1(k)|');
title('幅频特性|x1(k)|');
x2=sin(2*pi*0.125*n)+cos(2*pi*(0.125+f2)*n);y2=fft(x2,16);
subplot(2,2,2);stem(n,abs(y2),'.');axis([0 15 0 8]);
xlabel('(n=16 wk=2pik/n)k');ylabel('|x2(k)|');title('幅频特性|x2(k)|');%第三题(n=128)n=128;n=0:127;f1=1/16;f2=1/64;
x1=sin(2*pi*0.125*n)+cos(2*pi*(0.125+f1)*n);y3=fft(x1,128);figure;
subplot(2,2,1);stem(n,abs(y3),'.');axis([0 127 0 65]);
xlabel('(n=128 wk=2pik/n)k');ylabel('|x3(k)|');title('幅频特性|x3(k)|');
x2=sin(2*pi*0.125*n)+cos(2*pi*(0.125+f2)*n);y4=fft(x2,128);
subplot(2,2,2);stem(n,abs(y4),'.');axis([0 127 0 65]);
xlabel('(n=128 wk=2pik/n)k');ylabel('|x4(k)|');title('幅频特性|x4(k)|');
运行结果: n=16:
n=128:
结果分析:由实验运行结果可得出当采样点数n不同时,信号的幅频特性是不同的,这是因为x(k)是x(n)的傅里叶变换x(ejw)在频率区间[0,2π]上的n点等间隔采样,这就是dft的物理意义,由此可见,dft的变换长度n不同,表示对x(ejw)在[0,2π]区间上的采样间隔和采样点数不同,所以dft的变换结果不同。
3、第五题:
设计思想:本题给定了一个序列,要求用窗函数截短。通过调整窗口长度可以有效的控制过渡段的宽度。其中题目中要求的汉宁窗是升余弦窗,能够使能量更集中在主瓣中;哈明窗是改进的升余弦窗,能使能量更加集中在主瓣中。本程序定义了窗长度,给出相应的序列之后,用窗函数去截短相应的长度,最后通过窗口输出截短后的序列。程序代码:
%第五题(矩形窗)
m1=20;m2=40;m3=160;%定义窗长度 n1=0:19;n2=0:39;n3=0:159;
x1=0.5*cos(11*pi*n1/20)+cos(9*pi*n1/20);%相应序列 x2=0.5*cos(11*pi*n2/20)+cos(9*pi*n2/20);x3=0.5*cos(11*pi*n3/20)+cos(9*pi*n3/20);b1=boxcar(m1)'.*x1;%用窗函数截短 b2=boxcar(m2)'.*x2;b3=boxcar(m3)'.*x3;figure;
subplot(3,1,1);stem(abs(fft(b1)),'.');subplot(3,1,2);stem(abs(fft(b2)),'.');subplot(3,1,3);stem(abs(fft(b3)),'.');%第五题(汉宁窗)
m1=20;m2=40;m3=160;%定义窗长度 n1=0:19;n2=0:39;n3=0:159;
x1=0.5*cos(11*pi*n1/20)+cos(9*pi*n1/20);
x2=0.5*cos(11*pi*n2/20)+cos(9*pi*n2/20);x3=0.5*cos(11*pi*n3/20)+cos(9*pi*n3/20);y1=hanning(m1)'.*x1;%窗函数截短 y2=hanning(m2)'.*x2;y3=hanning(m3)'.*x3;figure;
subplot(3,1,1);stem(abs(fft(y1)),'.');subplot(3,1,2);stem(abs(fft(y2)),'.');subplot(3,1,3);stem(abs(fft(y3)),'.');%第五题(哈明窗)
m1=20;m2=40;m3=160;%定义窗长度 n1=0:19;n2=0:39;n3=0:159;
x1=0.5*cos(11*pi*n1/20)+cos(9*pi*n1/20);x2=0.5*cos(11*pi*n2/20)+cos(9*pi*n2/20);x3=0.5*cos(11*pi*n3/20)+cos(9*pi*n3/20);h1=hamming(m1)'.*x1;%窗函数截短 h2=hamming(m2)'.*x2;h3=hamming(m3)'.*x3;figure;
subplot(3,1,1);stem(abs(fft(h1)),'.');subplot(3,1,2);stem(abs(fft(h2)),'.');subplot(3,1,3);stem(abs(fft(h3)),'.');
运行结果: 矩形窗:
汉宁窗:
哈明窗:
结果分析:信号的截短产生了能量泄漏,不同的窗函数,产生泄漏的大小不一样,频率分辨能力也不一样,用fft算法计算频谱又产生了栅栏效应,所以不同的窗函数对信号频谱的影响是不一样的。从运行结果中可以看出矩形窗主瓣比较集中,旁瓣较高,并有负旁瓣;汉宁窗和哈明窗结果从图上粗看基本一样,但哈明窗加权的系数能使旁瓣达到更小。
4、第七题: 设计思想:
题目要求分别利用矩形窗、汉宁窗、哈明窗设计一个n=11的线性相位fir 低通和高通数字滤波器。fir滤波器在保证幅度特性满足技术要求的同时,很容易做到有严格的线性相位特性。fir滤波器设计的任务是选择有限长度的h(n),使传输函数h(ejw)满足技术要求。本题中要求设计的线性相位fir滤波器的,可利y用fir1函数设计出滤波器,再利用freqz函数得出其频率响应。程序代码:
%第七题(矩形窗)n=10;wn=1/3;b1lp=fir1(n,wn,boxcar(n+1));%用fir1函数设计低通滤波器
b1hp=fir1(n,wn,'high',boxcar(n+1));%用fir1函数设计高通滤波器 m=1024;w=0:(m-1);t=0:10;
h1lp=freqz(b1lp,1,w/m*2*pi);%得出频率响应 h1hp=freqz(b1hp,1,w/m*2*pi);figure;
subplot(3,2,1);stem(t,b1lp,'.');subplot(3,2,2);stem(t,b1hp,'.');subplot(3,2,3);
plot(w(1:(m/2))/m*2,abs(h1lp(1:(m/2))));subplot(3,2,4);
plot(w(1:(m/2))/m*2,abs(h1hp(1:(m/2))));subplot(3,2,5);
plot(w/m*2,angle(h1lp));subplot(3,2,6);
plot(w/m*2,angle(h1hp));%第七题(汉宁窗)n=10;wn=1/3;
b2lp=fir1(n,wn,hanning(n+1));
b2hp=fir1(n,wn,'high',hanning(n+1));m=1024;
w=0:(m-1);t=0:10;
h2lp=freqz(b2lp,1,w/m*2*pi);h2hp=freqz(b2hp,1,w/m*2*pi);figure;
subplot(3,2,1);stem(t,b2lp,'.');subplot(3,2,2);stem(t,b2hp,'.');subplot(3,2,3);
plot(w(1:(m/2))/m*2,abs(h2lp(1:(m/2))));subplot(3,2,4);
plot(w(1:(m/2))/m*2,abs(h2hp(1:(m/2))));subplot(3,2,5);
plot(w/m*2,angle(h2lp));subplot(3,2,6);
plot(w/m*2,angle(h2hp));
%第七题(哈明窗)n=10;wn=1/3;
b3lp=fir1(n,wn,hamming(n+1));
b3hp=fir1(n,wn,'high',hamming(n+1));m=1024;w=0:(m-1);t=0:10;
h3lp=freqz(b3lp,1,w/m*2*pi);h3hp=freqz(b3hp,1,w/m*2*pi);figure;
subplot(3,2,1);stem(t,b3lp,'.');subplot(3,2,2);stem(t,b3hp,'.');subplot(3,2,3);
plot(w(1:(m/2))/m*2,abs(h3lp(1:(m/2))));subplot(3,2,4);
plot(w(1:(m/2))/m*2,abs(h3hp(1:(m/2))));subplot(3,2,5);
plot(w/m*2,angle(h3lp));subplot(3,2,6);
plot(w/m*2,angle(h3hp));
运行结果:
矩形窗:
汉宁窗:
哈明窗:
结果分析:如运行结果所示,通过三个窗口设计出的滤波器,得到的频率响应特性是不同的。此题再次验证了第五题的结论:矩形窗主瓣比较集中,旁瓣较高,并有负旁瓣;汉宁窗和哈明窗都有使能量集中在主瓣中的作用,但哈明窗加权的系数能使旁瓣达到更小。
四、课程设计总结:
随着课程设计报告的基本完成,本次课程设计终于接近了尾声。本次课程设计要求我们利用上学期所学的现代信号处理课程的知识结合matlab编程工具,完成几道信号处理以及滤波器设计的题目,通过实际操作,回顾所学内容,夯实基础,强化理论知识,并体验理论与实际相结合的过程。
设计过程中遇到的第一个问题便是对于matlab语言的不熟悉,其实现在想想这个本不应该成为问题。虽然本专业并没有开设matlab程序设计这门课,但是在上学期上现代信号处理课时老师已经向我们强调过matlab语言的重要性,并且也在课堂上通过ppt的形式向我们展示了其部分功能。同时,老师还在上学期结课时特别提醒过课程设计会用到matlab语言,要求我们在暑假里去熟悉它。可惜的是,老师的的话并没有引起我的足够重视,才导致等到要真正设计的时候一头雾水,无从下手。通过这件事,我明白了对于一件事情,想要做得很好,提前做功课和准备是必不可少的,机会只垂青有准备的人。
当然,通过本次课程设计,我还是基本熟悉了一些matlab模块以及与本课程有关的一些函数的用法。但是由于上学期先代信号处理课程学得不是很好,所以也就是在懂的同学的指导下按部就班的写了一些代码,也没有什么创新,但还是很好的回顾的所学的频谱分析以、傅立叶变换以及滤波器设计的内容。也实际设计了一个滤波器,虽然借助与matlab函数,但也算一个不错的结果。老师在任务书中提到的目的要求,我自认为多多少少完成到达了一些。但不可否认的是还有很多没有达到。总之,很多东西还是需要自己在不断的摸索中找到答案。
另外,通过这次课程设计让我认识了一门新的语言———matlab,它在以后的学习中还会用到,相信在以后学习的过程的中,我会比较好的掌握并运用它。
参考文献
1、《数字信号处理》第二版,丁玉美 高西全等,西安电子科技大学出版社,2001.1
2、《matlab及在电子信息课程中的应用》,陈怀琛等,电子工业出版社出版,2002.4
3、《matlab数字信号处理与应用》,张德丰等,清华大学出版社,2010
4、《matlab数字信号处理》,王彬等,机械工业出版社,2010
5、《信号与信息处理基础》,彭军 李宏等,中国铁道出版社,2009.2
附.程序源代码
function varargout = untitled(varargin)%untitled m-file for
% untitled, by itself, creates a new untitled or raises the existing % singleton*.%
% h = untitled returns the handle to a new untitled or the handle to % the existing singleton*.%
% untitled('property','value',...)creates a new untitled using the % given property value gnized properties are passed via % varargin to calling syntax produces a % warning when there is an existing singleton*.%
% untitled('callback')and untitled('callback',hobject,...)call the % local function named callback in untitled.m with the given input % arguments.%
% *see gui options on guide's tools “gui allows only one % instance to run(singleton)”.%
% see also: guide, guidata, guihandles
% edit the above text to modify the response to help untitled
% last modified by guide v2.5 21-oct-2011 14:29:04
% begin initialization codedo not edit
%---executes just before untitled is made on untitled_openingfcn(hobject, eventdata, handles, varargin)% this function has no output args, see outputfcn.% hobject handle to figure
% eventdata reservedto be defined in a future version of matlab % handles structure with handles and user data(see guidata)
% get default command line output from handles structure varargout{1} = ;
%---executes on button press in on pushbutton1_callback(hobject, eventdata, handles)% hobject handle to pushbutton1(see gcbo)
% eventdata reservedto be defined in a future version of matlab % handles structure with handles and user data(see guidata)
%第三题(n=16)n=16;n=0:15;f1=1/16;f2=1/64;
x1=sin(2*pi*0.125*n)+cos(2*pi*(0.125+f1)*n);y1=fft(x1,16);%傅立叶变换函数 figure;
subplot(2,2,1);stem(n,abs(y1),'.');axis([0 15 0 8]);%标注坐标长度
xlabel('(n=16 wk=2pik/n)k');ylabel('|x1(k)|');title('幅频特性|x1(k)|');
x2=sin(2*pi*0.125*n)+cos(2*pi*(0.125+f2)*n);y2=fft(x2,16);
subplot(2,2,2);stem(n,abs(y2),'.');axis([0 15 0 8]);
xlabel('(n=16 wk=2pik/n)k');ylabel('|x2(k)|');title('幅频特性|x2(k)|');
%---executes on button press in on pushbutton11_callback(hobject, eventdata, handles)% hobject handle to pushbutton11(see gcbo)
% eventdata reservedto be defined in a future version of matlab
% handles structure with handles and user data(see guidata)
%第五题(矩形窗)
m1=20;m2=40;m3=160;%定义窗长度 n1=0:19;n2=0:39;n3=0:159;
x1=0.5*cos(11*pi*n1/20)+cos(9*pi*n1/20);%相应序列 x2=0.5*cos(11*pi*n2/20)+cos(9*pi*n2/20);x3=0.5*cos(11*pi*n3/20)+cos(9*pi*n3/20);b1=boxcar(m1)'.*x1;%用窗函数截短 b2=boxcar(m2)'.*x2;b3=boxcar(m3)'.*x3;figure;
subplot(3,1,1);stem(abs(fft(b1)),'.');subplot(3,1,2);stem(abs(fft(b2)),'.');subplot(3,1,3);stem(abs(fft(b3)),'.');
%---executes on button press in on pushbutton5_callback(hobject, eventdata, handles)% hobject handle to pushbutton5(see gcbo)
% eventdata reservedto be defined in a future version of matlab% handles structure with handles and user data(see guidata)
%第五题(哈明窗)
m1=20;m2=40;m3=160;%定义窗长度 n1=0:19;n2=0:39;n3=0:159;
x1=0.5*cos(11*pi*n1/20)+cos(9*pi*n1/20);x2=0.5*cos(11*pi*n2/20)+cos(9*pi*n2/20);x3=0.5*cos(11*pi*n3/20)+cos(9*pi*n3/20);h1=hamming(m1)'.*x1;%窗函数截短 h2=hamming(m2)'.*x2;h3=hamming(m3)'.*x3;figure;
subplot(3,1,1);stem(abs(fft(h1)),'.');subplot(3,1,2);stem(abs(fft(h2)),'.');subplot(3,1,3);stem(abs(fft(h3)),'.');
%---executes on button press in on pushbutton7_callback(hobject, eventdata, handles)
% hobject handle to pushbutton7(see gcbo)
% eventdata reservedto be defined in a future version of matlab % handles structure with handles and user data(see guidata)
%第七题(矩形窗)n=10;wn=1/3;
b1lp=fir1(n,wn,boxcar(n+1));%用fir1函数设计低通滤波器
b1hp=fir1(n,wn,'high',boxcar(n+1));%用fir1函数设计高通滤波器 m=1024;w=0:(m-1);t=0:10;
h1lp=freqz(b1lp,1,w/m*2*pi);%得出频率响应 h1hp=freqz(b1hp,1,w/m*2*pi);figure;
subplot(3,2,1);stem(t,b1lp,'.');subplot(3,2,2);stem(t,b1hp,'.');subplot(3,2,3);
plot(w(1:(m/2))/m*2,abs(h1lp(1:(m/2))));subplot(3,2,4);
plot(w(1:(m/2))/m*2,abs(h1hp(1:(m/2))));subplot(3,2,5);
plot(w/m*2,angle(h1lp));subplot(3,2,6);
plot(w/m*2,angle(h1hp));
%---executes on button press in on pushbutton9_callback(hobject, eventdata, handles)% hobject handle to pushbutton9(see gcbo)
% eventdata reservedto be defined in a future version of matlab % handles structure with handles and user data(see guidata)
%第七题(哈明窗)n=10;wn=1/3;
b3lp=fir1(n,wn,hamming(n+1));
b3hp=fir1(n,wn,'high',hamming(n+1));m=1024;w=0:(m-1);t=0:10;
h3lp=freqz(b3lp,1,w/m*2*pi);h3hp=freqz(b3hp,1,w/m*2*pi);figure;
subplot(3,2,1);stem(t,b3lp,'.');subplot(3,2,2);stem(t,b3hp,'.');subplot(3,2,3);
plot(w(1:(m/2))/m*2,abs(h3lp(1:(m/2))));subplot(3,2,4);
plot(w(1:(m/2))/m*2,abs(h3hp(1:(m/2))));subplot(3,2,5);
plot(w/m*2,angle(h3lp));subplot(3,2,6);
plot(w/m*2,angle(h3hp));
10>3>
一键复制