matlab数字信号处理

1.解差分方程

差分方程如下：（k0、k1、k2为常数）

u(k)=u(k-1)+k0*e(k)+k1*e(k-1)+k2*e(k-2)

matlab代码如下：

clc
clear all
N = 6
k0=2.01
k1=-3
k2=1
x = [k0 ,k1, k2 , 0,0,0]
y = [1 , -1, 0, 0, 0,0]
e = [1 zeros(1,N-1)]
k=0:1:N-1
u = filter(x,y,e)
stem(k,u,'k','filled')
xlabel('n')
ylabel('U')
grid on

2.低通滤波器的逼近

clc;clear all;
syms w n
wc=0.5*pi
X=exp(-1i*w*n)
u(n)=int(X,w,-wc,wc)/(2*pi)
lim=limit(u,n,0)
n=-10:10
zero=find(n==0)
n(zero)=[]
figure(1)
stem(n,u(n),'k','filled')
hold on
stem(0,lim,'k','filled')
xlabel('n')
ylabel('hd')
grid on

syms  w H
M=50
H(w)=0
N=-M:M
for n=N
H(w)=H(w)+sinc(n/2)*exp(-1i*w*n)/2
end
x=-10:pi/100:10
figure(2)
plot(x,abs(H(x)),'k')
xlabel('w')
ylabel('|HM(exp(jw))|')
grid on

当然也可以用symsum函数：

clc ;clear all;
syms  w  n
M=50
f=sinc(n/2)*exp(-1i*w*n)/2
H=symsum(f,n,-M,-1)+1+symsum(f,n,1,M)
fplot(abs(H))

结果同上

3.截断函数的影响

Dg(exp(jw))=sin(wN/2)/sin(w/2)
原理为：

matlab代码如下：

clear all ;clc;Omega_step=0.001;Omega1=0.4*pi;Omega2=0.6*pi;N=32;
OmegaX=[-pi:Omega_step:pi];X=zeros(1,length(OmegaX));
X(fix((abs(OmegaX(1))-Omega1)/Omega_step):fix((abs(OmegaX(1))+Omega1)/Omega_step))=31;

 subplot(2,2,1) ;plot(OmegaX,X,'k');ax=gca;ax.FontSize=24;xlim([-pi pi]);ylim([-1,32]);xlabel('Omega');ylabel('X');
grid on
OmegaDg=[-pi:Omega_step:pi];      Dg=sin(OmegaDg/2*N)./sin(OmegaDg/2);
 subplot(2,2,2);plot(OmegaDg,Dg,'k');ax=gca;ax.FontSize=24;xlim([-pi pi]);xlabel('Omega');ylabel('Dg');
grid on
XDg=conv(X/Omega_step,Dg)*Omega_step/2/pi;  OmegaXDg=OmegaX(1)+OmegaDg(1):Omega_step:OmegaX(end)+OmegaDg(end);
 subplot(2,2,[3 4]);plot(OmegaXDg,XDg,'k');ax=gca;ax.FontSize=24;xlim([-pi pi]);xlabel('Omega');ylabel('X*Dg');
grid on
结果如下：