利用MATLAB进行信号分析
? 数字滤波器的频率响应:
0.2?0.3z?1?z?2例:对一数字滤波器H(z)=,求其幅频特性和相频特性: ?1?21?0.4z?z a=[1 0.4 1]; b=[0.2 0.3 1]; freqz(b,a,128)
? 离散系统零极点图: b =[0.2 0.1 0.3 0.1 0.2]; a=[1.0 -1.1 1.5 -0.7 0.3]; zplane(b,a)
? 数字滤波器的冲激响应: b=[0.2 0.1 0.3 0.1 0.2];
a=[1.0 -1.1 1.5 -0.7 0.3]; impz(b,a,50)
? 计算离散卷积:
x=[1 1 1 1 0 0]; y=[2 2 3 4]; z=conv(x,y);
stem(0:length(x)+length(y)-2,z)
? 系统函数转换:
H(z)=
(z?0.1)(z?0.5)(z?1)
(z?0.3)(z?2)(z?5)y=[0.1 0.5 1]; x=[0.3 2 5]; k=1;
b=k*poly(y); a=poly(x);
--------------------------------------------------- b b =
1.0000 -1.6000 0.6500 -0.0500
a a =
1.0000 -7.3000 12.1000 -3.0000
1?1.6z?1?0.65z?2?0.05z?3则H(z)= ?1?2?31?7.3z?12.1z?3.0z ------------------------------------------------------
a=[1 –7.3 12.1 –3]; b=[ –1 6 0.65 –0.05]; y=roots(b); x=roots(a); d=b(1)/a(1);
--------------------------------------------- y y =
1.0000 0.5000 0.1000
x x =
5.0000 2.0000
0. 3000
-----------------------------------------------
b=[8 -4 11 -2];
a=[1 -1.25 0.75 -0.125]; y=roots(b); x=roots(a)
---------------------------------------------- y =
0.1550 + 1.1367i 0.1550 - 1.1367i 0.1900
x =
0.5000 + 0.5000i 0.5000 - 0.5000i
0.2500
? 计算离散信号的能量和功率:
序列的能量:E=
?x(n)n?0N?1n?0?2=
?x(n)x*(n)n?0?
1序列的功率:P=
N?x(n)21=N?x(n)x*(n)n?0?
例:x(n)=0.8nu(n),计算前10点的能量占总能量的百分比。 前十点能量:
N=10;
n=0:N-1; x=(0.8).^n;
e=sum(abs(x).^2)
E=1/(1-0.8.^2) e/E
? 用FFT计算周期序列的频谱:X(ej)=
ω
1X(k) , N为序列x(n)的周期。 Nx(n)=cos(
??n+),用FFT计算其频谱。 83N=16;n=0:N-1;
x=cos(pi/8*n+pi/3); X=1/N*fft(x,N);
omega=2*pi/N*(n-N/2);
subplot(2,1,1);stem(omega,abs(fftshift(X)));axis([-pi,pi,0,1]); ylabel('Magnitude');xlabel('Frequency(rad)');
subplot(2,1,2);stem(omega,angle(fftshift(X)));axis([-pi,pi,-4,4]); ylabel('Phase');xlabel('Frequency(rad)'); ? 利用DFT实现线性卷积:
已知x(n)=[1 2 3 4],y(n)=[5 6 7 8],求z(n)=x(n)*y(n)。
x=[1 2 3 4];y=[5 6 7 8]; z=conv(x,y)
x=[1 2 3 4];y=[5 6 7 8]; X=fft(x);Y=fft(y); Z=X.*Y;z=ifft(Z)
x=[1 2 3 4];y=[5 6 7 8];L=7; X=fft(x,L);Y=fft(y,L); Z=X.*Y;z=ifft(Z,L)
x=[1 2 3 4];y=[5 6 7 8];L=8; X=fft(x,L);Y=fft(y,L); Z=X.*Y;z=ifft(Z,L)
? 求逆Z变换:
求X(z)=
1 z﹥0.9的逆Z变换。
(1?0.9z?1)2(1?0.9z) b=1;
a=poly([0.9 0.9 -0.9]); [R P C]=residuez(b,a) 执行结果为: R =
0.2500 0.5000 0.2500 P =
0.9000 0.9000 -0.9000 C =
[] 则:X(z)=
0.250.5+
1?0.9z?1(1?0.9z?1)20.25
1?0.9z?1相应的逆Z变换为:x(n)=0.25(0.9)nu(n)+
5(n+1)(0.9)n+1u(n+1)+ 0.25(-0.9)nu(n) 9 =0.75(0.9)nu(n)+ 0.5n(0.9)n+1u(n)+ 0.25(-0.9)nu(n)
FFT的应用
? 试用FFT计算其频谱:
x(n)=cos(
7?9?n)+cos(n) 1716% pp.m
M=80;n=0:M-1;
x=(cos(7*pi/17*n)+cos(9*pi/16*n)); N=256; X=fft(x,N);
omega=2*pi/N*[(0:N-1)-N/2]; figure(1);
subplot(2,1,1);stem(n,x);xlabel('Time index n');ylabel('x(n)'); subplot(2,1,2);plot(omega,abs(X));axis([-pi,pi,0,50]); ylabel('X');xlabel('omega (rad)');
? 用FFT分析信号频率成分:
一被噪声污染的信号,很难看出它所包含的频率分量,如一个由50Hz和120Hz正弦信号构成的信号,受到均值随机噪声的干扰,数据采样为1000Hz。通过FFT来分析其信号频率成分。 ?t1.m
t=0:0.001:0.6;
x=sin(2*pi*50*t)+sin(2*pi*120*t); y=x+randn(1,length(t)); Y=fft(y,512);
subplot(2,1,1);plot(y);title('受噪声污染的信号'); k=0:511;f=1000*k/512;% w=2*pi*f/fs= 2*pi*k/N subplot(2,1,2);plot(f,abs(Y));title('FFT');
? 用FFT分析语音信号的频谱:
?t2.m load mtlb;
subplot(2,2,1);plot(mtlb);title('原始语音信号'); y=fft(mtlb);
subplot(2,2,2);plot(abs(y));title('FFT变换'); y(abs(y)<1)=0;x=ifft(y);
subplot(2,2,3);plot(abs(y));title('去掉幅值小于1的变换值'); subplot(2,2,4);plot(real(x));title('重构语音信号');