医学信号处理实验报告.docx

资源描述

医学信号处理实验报告.docx

《医学信号处理实验报告.docx》由会员分享，可在线阅读，更多相关《医学信号处理实验报告.docx（5页珍藏版）》请在冰点文库上搜索。

医学信号处理实验报告.docx

医学信号处理实验报告

　　　　　　医学信号处理实验报告　　　　　　　　　　班级：

生医1201　　姓名：

葛然　　　　学号：

12282002　　　　实验一：

　　基本函数图形　　　　　　1　　　　　　　　2　　　　　　　　　　第一题　　　　3　　　　　　　　第二题第一小题　　　　源代码：

　　n=-10:

10;　　x1=n.*n.*（[n>=-5]-[n>=6]）;subplot（2,2,1）;stem（n,x1）;　　title（‘n*n*[u（n+5￡?

-u（n-6）]’）;xlabel（‘n’）;ylabel（‘x（n）’）;　　x2=10.*[n==0];subplot（2,2,2）;stem（n,x2）;　　title（‘10*[n==0]’）;xlabel（‘n’）;ylabel（‘x（n）’）;　　x3=20.*

.*（[n>=4]-[n>=10]）;subplot（2,2,3）;　　4　　　　

　　　　　　stem（n,x3）;　　title（‘20.*

.*（[n>=4]-[n>=10]）’）;xlabel（‘n’）;ylabel（‘x（n）’）;　　　　x4=n.*n.*（[n>=-5]-[n>=6]）+10.*[n==0]+20.*

.*（[n>=4]-[n>=10]）;subplot（2,2,4）;stem（n,x4）;title（‘x1+x2+x3’）;xlabel（‘n’）;ylabel（‘x（n）’）;　　运行结果：

　　　　　　5　　　　第二题第二小题　　　　源代码：

　　n=-20:

20;x1=cos（*pi.*n）;x2=cos（*pi.*n）;subplot（3,1,1）;stem（n,x1）;　　title（‘cos（*pi.*n）’）;xlabel（‘n’）;ylabel（‘x1（n）’）;gridon;subplot（3,1,2）;stem（n,x2）;　　title（‘cos（*pi.*n）’）;xlabel（‘n’）;ylabel（‘x2（n）’）;gridon;　　n=-40:

40;x1=cos（*pi.*n）;subplot（3,1,3）;stem（n,x1）;　　title（‘cos（*pi.*n）’）;xlabel（‘n’）;ylabel（‘x1（n）’）;gridon;　　　　　　　　6　　运行结果：

　　　　　　A小问：

cos（）（-20　　　　答：

上图可以看出cos（*pi*n）是周期序列，基本的周期是20.序列cos（*pi*n）是周期为5的周期序列。

对比两者可得N实际上是基波的周期，基波的频率为2*pi/N而K是一个任意的正整数与N互为素数。

　　　　B小问：

cos（）（-40　　答：

此序列是非周期的，对比a组中的结论我们可得到一下结论：

1）当2*pi/w0是整数或者有理数时，K取适当的整数值可以得到最小周期如前两个图所示；2）当2*pi/w0是一个无理数时，任何的K值都不能使得N为正整数，所以余弦序列是非周期的，如第三个图所示。

　　　　7　　　　实验二：

　　　　第一题　　运行结果：

　　　　第二题：

　　源代码：

　　p=[,,,];d=[1];n=-5:

20;x=[n==0];x=double（x）;h=filter（p,d,x）;subplot（2,1,1）;stem（n,h）;　　title（‘impulseresponse’）;xlabel（‘time’）;ylabel（‘amplitude’）;x=[n>=0];　　8　　　　x=double（x）;s=conv（x,h）;subplot（2,1,2）;n=-10:

40;stem（n,s）;　　title（‘stepresponse’）;xlabel（‘time’）;ylabel（‘amplitude’）;　　　　　　运行结果：

　　　　　　理论计算结果及对比分析：

　　　　答：

通过理论计算得到n=0,1,2,3,4五个点的结果分别为：

　　n=0:

h（n）=0　　n=1,2,3,4时h（n）=　　利用得到的单位冲击响应函数求得卷积的结果分别为0，，0.5，，1　　对比图像，理论结果与实验结果一致　　9　　

　　　　　　实验三：

　　　　第一题　　运行结果：

　　　　　　第二题：

　　源代码：

　　num=[-1];den=[-1];[z,p,k]=tf2zp（num,den）;n=0:

20;　　x=2.*

.*[n>=0];y=filter（num,den,x）;　　subplot（2,2,1）;zplane（z,p）;w=[0:

512]*pi/512;H=freqz（num,den,w）;magH=abs（H）;phaH=angle（H）;　　10　　　　　　subplot（2,2,2）;plot（w/pi,magH）;grid;xlabel（‘?

μ?

êμ￥?

pi’）;ylabel（‘·ù?

è’）;title（‘?

μ?

ê?

ìó|’）;　　subplot（2,2,3）;plot（w/pi,phaH）;grid;　　xlabel（‘?

μ?

êμ￥?

pi’）;ylabel（‘?

à?

’）;title（‘?

à?

ìó|’）;　　subplot（2,2,4）;plot（n,y）;grid;　　xlabel（‘n’）;ylabel（‘y’）;title（‘y（n）’）;　　　　　　运行结果：

　　　　11　　　　　　实验四：

tude’）;　　subplot（4,1,3）;stem（w/pi,angle（y））;title（‘phase’）;　　　　13　　　　运行结果：

　　　　　　源代码：

　　N=10;n=0:

9;　　x=n.*（1./n）.*[n　　y=fft（x,N）;x1=ifft（y,N）;　　subplot（4,1,1）;stem（n,x）;　　title（‘originalsignal’）;　　subplot（4,1,2）;stem（n,abs（y））;title（‘magnitude’）;　　subplot（4,1,3）;　　14　　

　　　　　　stem（n,angle（y））;title（‘phase’）;　　subplot（4,1,4）;stem（n,x1）;　　title（‘resumedsignal’）;　　　　运行结果：

　　　　　　15　　　　实验五：

　　第一题运行结果：

　　　　　　16　　　　第二题　　源代码：

　　N=64;k=0:

N-1;w=2*pi*k/N;z=exp（j*w）;Xk=（z）./（）;xn=real（ifft（Xk,N））;xtilde=xn’*ones（1,2）;xtilde=（xtilde（:

））’;subplot（2,1,1）;stem（0:

127,xtilde）;title（‘64μ?

μ?

IFFT·ù?

è?

×’）;xlabel（‘n’）;　　ylabel（‘xtilde（n）’）;N=16;k=0:

N-1;w=2*pi*k/N;z=exp（j*w）;Xk=（z）./（）;xn=real（ifft（Xk,N））;xtilde=xn’*ones（1,8）;xtilde=（xtilde（:

））’;subplot（2,1,2）;stem（0:

127,xtilde）;title（‘16μ?

μ?

IFFT·ù?

è?

×’）;xlabel（‘n’）;　　ylabel（‘xtilde（n）’）;　　　　　　　　17　　运行结果：

　　　　　　　　18　　　　实验七：

　　第一题：

源代码：

n=50;　　Wn=;　　b1=fir1（n,Wn,boxcar（n+1））;[H,w]=freqz（b1,1,512,2）;subplot（4,1,1）;plot（w,abs（H））;grid;　　title（‘矩形窗低通FIR数字滤波器’）;b2=fir1（n,Wn,hanning（n+1））;[H,w]=freqz（b2,1,512,2）;subplot（4,1,2）;　　plot（w,abs（H））;grid;　　title（‘汉宁窗低通FIR数字滤波器’）;b3=fir1（n,Wn,hamming（n+1））;[H,w]=freqz（b3,1,512,2）;subplot（4,1,3）;　　plot（w,abs（H））;grid;　　title（‘海明窗低通FIR数字滤波器’）;b4=fir1（n,Wn,blackman（n+1））;[H,w]=freqz（b4,1,512,2）;subplot（4,1,4）;　　plot（w,abs（H））;grid;　　title（‘希拉克曼窗低通FIR数字滤波器’）;　　运行结果：

　　矩形窗低通FIR数字滤波器希拉克曼窗低通FIR数字滤波器海明窗低通FIR数字滤波器汉宁窗低通FIR数字滤波器　　19　　

　　　　　　第二题：

　　　　n=50;Wn=,　　b1=fir1（n,Wn,’high’,boxcar（n+1））;[H,w]=freqz（b1,1,512,2）;subplot（4,1,1）;　　plot（w,abs（H））;grid;　　title（‘矩形窗高通FIR数字滤波器’）;b2=fir1（n,Wn,’high’,hanning（n+1））;[H,w]=freqz（b2,1,512,2）;subplot（4,1,2）;　　plot（w,abs（H））;grid;　　title（‘汉宁窗高通FIR数字滤波器’）;b3=fir1（n,Wn,’high’,hamming（n+1））;[H,w]=freqz（b3,1,512,2）;　　subplot（4,1,3）;　　plot（w,abs（H））;grid;　　title（‘海明窗高通FIR数字滤波器’）;b4=fir1（n,Wn,’high’,blackman（n+1））;[H,w]=freqz（b4,1,512,2）;　　subplot（4,1,4）;　　plot（w,abs（H））;grid;　　title（‘希拉克曼窗高通FIR数字滤波器’）;　　运行结果：

　　矩形窗高通FIR数字滤波器希拉克曼窗高通FIR数字滤波器海明窗高通FIR数字滤波器汉宁窗高通FIR数字滤波器　　　　　　第三题：

　　　　n=50;　　Wn=[,];　　b1=fir1（n,Wn,’bandpass’,boxcar（n+1））;[H,w]=freqz（b1,1,512,2）;subplot（4,1,1）;　　plot（w,abs（H））;grid;　　title（‘矩形窗带通FIR数字滤波器’）;b2=fir1（n,Wn,’bandpass’,hanning（n+1））;[H,w]=freqz（b2,1,512,2）;subplot（4,1,2）;　　plot（w,abs（H））;grid;　　title（‘汉宁窗带通FIR数字滤波器’）;b3=fir1（n,Wn,’bandpass’,hamming（n+1））;[H,w]=freqz（b3,1,512,2）;　　subplot（4,1,3）;　　plot（w,abs（H））;grid;　　title（‘海明窗带通FIR数字滤波器’）;b4=fir1（n,Wn,’bandpass’,blackman（n+1））;[H,w]=freqz（b4,1,512,2）;　　subplot（4,1,4）;　　plot（w,abs（H））;grid;　　title（‘希拉克曼窗带通FIR数字滤波器’）;　　运行结果：

　　矩形窗带通FIR数字滤波器希拉克曼窗带通FIR数字滤波器海明窗带通FIR数字滤波器汉宁窗带通FIR数字滤波器　　　　　　第四题：

　　　　n=50;　　Wn=[,];　　b1=fir1（n,Wn,’stop’,boxcar（n+1））;[H,w]=freqz（b1,1,512,2）;subplot（4,1,1）;　　plot（w,abs（H））;grid;　　title（‘矩形窗带阻FIR数字滤波器’）;b2=fir1（n,Wn,’stop’,hanning（n+1））;[H,w]=freqz（b2,1,512,2）;subplot（4,1,2）;　　plot（w,abs（H））;grid;　　title（‘汉宁窗带阻FIR数字滤波器’）;b3=fir1（n,Wn,’stop’,hamming（n+1））;[H,w]=freqz（b3,1,512,2）;　　subplot（4,1,3）;　　plot（w,abs（H））;grid;　　title（‘海明窗带阻FIR数字滤波器’）;b4=fir1（n,Wn,’stop’,blackman（n+1））;[H,w]=freqz（b4,1,512,2）;　　subplot（4,1,4）;　　plot（w,abs（H））;grid;　　title（‘希拉克曼窗带阻FIR数字滤波器’）;　　运行结果：

　　矩形窗带阻FIR数字滤波器希拉克曼窗带阻FIR数字滤波器海明窗带阻FIR数字滤波器汉宁窗带阻FIR数字滤波器　　22　　　　答：

从低通、高通、带通、带阻的设计所得到的波形可以看出矩形窗的过渡带最窄，但是它的余振最大，布拉克曼窗的过渡带最长，几乎没有余振。

　　于实际的DF在理想特性不连续的边缘处，H（w）产生一个过渡带。

过渡带的宽度取决于窗口频谱的主瓣宽度.而余振取决于窗口频谱的旁瓣，旁瓣愈多，余振也越多。

依据能量守恒，主瓣越宽必然导致旁瓣越小，主瓣越窄必然导致旁瓣越大。

所以导致过渡带越宽余振越小，过渡带越窄余振越大。

　　　　　　　　23　　　　实验八：

　　实验　　源代码：

　　wp1=*pi;wp2=*pi;ws1=*pi;ws2=*pi;rp=1;rs=40;wp=[wp1,wp2];ws=[ws1,ws2];　　[n,wn]=buttord（wp/pi,ws/pi,rp,rs）;[b,a]=butter（n,wn,’bandpass’）;freqz（b,a）;　　得到图形为：

　　200Magnitude（dB）Frequency（?

rad/sample）（degrees）-Frequency（?

rad/sample）　　源代码2：

　　i=0:

;x1=sin（600*pi*i）;x2=sin（200*pi*i）;　　x3=sin（900*pi*i）;　　x=x1+x2+x3;　　y=filter（b,a,x）;　　　　24　　

　　　　　　figure　　subplot（5,1,1）;plot（i,x1）;subplot（5,1,2）;plot（i,x2）;subplot（5,1,3）;plot（i,x3）;subplot（5,1,4）;plot（i,x）;subplot（5,1,5）;plot（i,y）;　　得到的图形为：

　　　　得到的参数如下>>bb=　　Columns1through11　　0-　　0　　0-　　0　　0-Columns12through15　　0　　0->>aa=　　Columns1through11　　　　　　　　Columns12through15　　　　25　　　　>>nn=　　7>>wnwn=　　　　所以计算得到的传递函数为：

　　h（n）?

10?

12　　　　?

14]/[1?

10?

11?

12?

13?

15]　　其单位响应及零级点分布：

　　num=[0-00-00-00-]　　den=[　　　　][z,p,k]=tf2zp（num,den）;subplot（2,1,1）;　　zplane（z,p）;grid;w=[0:

512]*pi/512;H=freqz（num,den,w）;　　n=0:

50;　00-20幅度-频率　　答：

FIR系统只有零点，因此FIR系统不像IIR系统那样易于取得比较好的通带与阻带衰减特性，要取得好的衰减特性，滤波器的阶次要求更高，既N更大。

本实验中设计的等波纹最优FIR带通滤波器要符合相应的滤波特性，最少需n=31阶。

　　FIR数字滤波器可以很容易实现严格的线性相位特性。

实际应用中，像语音处理、图像处理以及数据传输要求线性相位，任意幅度，FIR滤波器满足要求。

同时FIR数字滤波器的单位抽样响应是有限长的，因而滤波器一定是稳定的，只须经过一定的延时，任何非因果的有限长序列都变成因果的有限序列。

展开阅读全文