数字信号处理（二）

%解差分方程  y[n]-0.95y[n-1]+0.9025y[n-2]=1/3[x[n]+x[n-1]+x[n-2]], n>=0
% x[n]=cos(pin/3), y[-1]=-2, y[-2]=-3, x[-1]=1, x[-2]=1
b=[1,1,1]/3; a=[1,-0.95, 0.9025];
Y=[-2,-3]; X=[1,1];
xic=filtic(b,a,Y,X)
bxplus=[1,-0.5]; axplus=[1,-1,1];
ayplus=conv(a,axplus)
byplus=conv(b,bxplus)+conv(xic,axplus)
[R,p,C]=residuez(byplus,ayplus)
Mp=abs(p),Ap=angle(p)/pi
n=[0:50];
x=cos(pi*n/3);
y=filter(b,a,x,xic);
plot(n,y),title(‘系统响应曲线‘);

%解差分方程  y[n]-0.95y[n-1]+0.9025y[n-2]=1/3[x[n]+x[n-1]+x[n-2]], n>=0
% x[n]为采集的信号, y[-1]=-2, y[-2]=-3, x[-1]=1, x[-2]=1
b=[1,1,1]/3; a=[1,-0.95, 0.9025];
Y=[-2,-3]; X=[1,1];
xic=filtic(b,a,Y,X)
bxplus=[1,-0.5]; axplus=[1,-1,1];
ayplus=conv(a,axplus)
byplus=conv(b,bxplus)+conv(xic,axplus)
[R,p,C]=residuez(byplus,ayplus)
Mp=abs(p),Ap=angle(p)/pi
n=[0:50];
      [x,Fs]=audioread(‘E:test1.wma‘);  %自己录入的wma文件
X=x(:,1);  %单声道？
y=filter(b,a,X,xic);

figure
subplot(211),plot(X), title(‘输入信号曲线‘);
subplot(212),plot(y),title(‘差分方程运算处理后信号曲线‘)

sound(X,Fs);
sound(y,Fs)

%解差分方程  y[n]-0.95y[n-1]+0.9025y[n-2]=1/3[x[n]+x[n-1]+x[n-2]], n>=0
% x[n]为采集的信号, y[-1]=-2, y[-2]=-3, x[-1]=1, x[-2]=1
b=[1,1,1]/3; a=[1,-0.95, 0.9025];
Y=[-2,-3]; X=[1,1];
xic=filtic(b,a,Y,X)
bxplus=[1,-0.5]; axplus=[1,-1,1];
ayplus=conv(a,axplus)
byplus=conv(b,bxplus)+conv(xic,axplus)
[R,p,C]=residuez(byplus,ayplus)
Mp=abs(p),Ap=angle(p)/pi
n=[0:50];
      [x,fs,bits]=wavread(‘D:liulili.wav‘);  %自己录入的wav文件
X=x(:,1);  %单声道？
y=filter(b,a,X,xic);

figure
subplot(211),plot(X), title(‘输入信号曲线‘);
subplot(212),plot(y),title(‘差分方程运算处理后信号曲线‘)

sound(X,fs);
sound(y,fs)