数字信号处理与MATLAB 实现1． n1=[ns:nf]；x1=[zeros(1,n0-ns),1,zeros (1,nf-n0)]； %单位抽样序列的产生2. subplot(2,2,4) 画2行2列的第4个图3. stem(n,x) %输出离散序列，(plot 连续)4. 编写子程序可调用4.1 单位抽样序列)(0n n -δ生成函数impseq.m[x,m]=impseq(n0,ns,nf); %序列的起点为ns ，终点为nf ，在n=n0点处生成一个单位脉冲n=[-5:5];x1=3*impseq(2,-5,5)-impseq(-4,-5,5)x1 =0 -1 0 0 0 0 0 3 0 0 0n=[-5:5];x1=3*impseq(2,-4,5)-impseq(-4,-5,4) %起点到终点长度要一致x1 =0 -1 0 0 0 0 3 0 0 04.2 单位阶跃序列)(0n n u -生成函数stepseq.m[x,n]=stepseq(no,ns,nf) %序列的起点为ns ，终点为nf ，在n=n0点处生成一个单位阶跃4.3 两个信号相加的生成函数sigadd.m[y,n]=sigadd(x1,n1,x2,n2)4.4 两个信号相乘的生成函数sigmult.m[y,n]=sigmult(x1,n1,x2,n2)4.5 序列移位y(n)=x(n-n0)的生成函数sigshift.m[y,n]=sigshift(x,m,n0)4.6 序列翻褶y(n)=x(-n)的生成函数sigfold.m[y,n]=sigfold(x,n)4.7 evenodd.m 函数可以将任一给定的序列x(n)分解为xe(n)和xo(n)两部分[xe,xo,m]=evenodd(x,n)4.8 序列从负值开始的卷积conv_m, conv 默认从0开始function [y,ny]=conv_m(x,nx,h,nh)有{x(n):nx1≤n ≤nx2},{h(n):nh1≤n ≤nh2}, 卷积结果序列为{y(n):nx1+nh1≤n ≤nx2+nh2}例. 设1132)(-++=z z z X ,1225342)(-+++=z z z z X ,求)()()(21z X z X z Y += 程序：x1=[1,2,3];n1=-1:1;x2=[2,4,3,5];n2=-2:1;[y,n]=conv_m(x1,n1,x2,n2)结果：y =2 8 17 23 19 15n =-3 -2 -1 0 1 2因此21231519231782)(--+++++=z zz z z z Y 计算1121))9.01()9.01(()(---+-=z z z X ，9.0>z 得Z 反变换b=1;a=poly([0.9 0.9 -0.9]);[r,p,k]=residuez(b,a)结果：r =0.25000.50000.2500p =0.90000.9000-0.9000k =[] 因此得到,9.0125.0)9.01(5.09.0125.0)(1211---++-+-=z z z z X 9.0>z 相应的)()9.0(25.0)1()9.0)(1(95)()9.0(25.0)(1n u n u n n u n X n n n -++++=+ 5. [H,w]=freqz[B,A,M]计算出M 个频率点上的频率响应，存放于H 向量中，M 个频率存放在w 向量中，freqz 函数自动将这M 个频率点均匀设置在频率围[0，π]之间。
若缺省w 和M 时，函数自动选取512个频率点计算。
不带输出向量的freqz 函数将自动绘制幅频和相频曲线。
也可[H,w]=freqz(B,A);plot(w/pi,abs(H)); 绘出幅频特性6. zplane(z,p)绘制出列向量z 中的零点（以符号o 表示）和列向量p 中的极点（以符号x 表示）以及参考单位圆。
7. 傅里叶变换信号时域和频域变换连续对应着非周期，离散对应着周期。
一个域的离散必然导致另一个域的周期延拓8. fft 和ifft ：一维快速傅里叶变换和逆傅里叶变换X=fft(x,N)采用FFT 算法计算序列向量X 的N 点DFT ，缺省N 时，fft 函数自动按X 的长度计算DFT 。
当N 为2的整数次幂时，fft 按基2算法计算，否则用混合基算法。
ifft 的调用格式类似。
9. fft2和ifft2: 二维快速正傅里叶变换和逆傅里叶变换（1）Y=fft2(x)数据二维傅里叶变换参数X 是向量fft2(x)相当于fft(fft(x)’)’,即先对X 的列做一维傅里叶变换，然后再对变换结果的行作一维傅里叶变换；若X 是向量，则此傅里叶变换即变成一维傅里叶变换fft; 若X 是矩阵，则是计算该矩阵的二位傅里叶变换。
（2）Y=fft2(X, M, N)通过对X 进行补零或截断，使得X 成为（M*N ）的矩阵。
函数iff2的参数应用与函数fft2完全相同。
10. czt: 线形调频Z 变换Y=czt(x,m,w,a)此函数计算由z=a*w.^(-(0:m-1))定义的z 平面螺旋线上各点的z 变换，a 规定了起点，w 规定了相邻点的比例，m 规定了变换的长度，后三个变量默认值为a=1,w=exp(j*2*pi/m)及m=length(x)因此y=czt(x)就等于y=fft(x).11. dct 和idct ：离散余弦正变换和离散余弦逆变换Y=dct(x,N)完成如下变换，N 的默认值为length(x).∑=+=N n n k n n x K Y 1))12(2cos()(2)(π k=0, 1, …，N-112. fftshiftY=fftshift(x)用来重新排列X=fft(x)的输出，当X 为向量时，把X 的左右两半进行交换，从而将零频分量移至频谱的中心；如果X 为二维傅里叶变换的结果，它同时将X 的左右和上下部分进行交换。
13. fftfiltY=fftfilt(b,x)采用重叠相加法FFT 对信号向量x 快速滤波，得到输出序列向量y ，向量b 为FIR 滤波的单位脉冲相应，h(n)=b(n+1),n=0,1,…,length(b)-1.Y=fftfilt(b,x,N)自动选取FFT 长度NF=2^nextpow2(N), 输入数据x 分段长度M=NF-length(b)+1, 其中nextpow2(N)函数求的一个整数，满足2^(nextpow2(N)-1)<=N<=2^nextpow2(N)N 缺省时，fftfilt 自动选择合适的FFT 长度NF 和对x 的分段长度M 。
p = nextpow2(A)， p 满足 2^p >= abs(A)L=pow2(nextpow2(M+N-1)) M, N 分别为序列延拓周期14. 用DFT 进行谱分析时，必须将序列阶段为长度为N 的有限长序列。
造成频谱泄露和谱间干扰。
泄露使频谱变得模糊，分辨率降低；旁瓣引起不同分量间的干扰。
截断效应无法完全消除，可以加宽窗和缓慢截断。
矩形窗比海明窗的频率分辨率高（泄漏小），但谱间干扰大，因此海明窗是以牺牲分辨率来换取谱间干扰的降低。
栅栏效应是只能在离散点的地方看到真实的像，其余频谱被遮挡。
为减少栅栏效应，可以在时域数据末端增加一些零点，使周期点数增加，但不改变原有数据，即增加频域抽样点数N ，频域抽样为k Nπ2，这样必然使谱线更密，这样原来看不到的谱分量就可能看到了。
15. x=[one(1,5),zero(1,N-5)]; %单位阶跃信号y = filter(b,a,x) %直接型输出信号16. N=25;h=impz(b,a,N) %N 次采样，b,a 为零极点多项式系数，a 不算1如：543214321421216118271131)(-----------++++++-=zz z z z z z z z z H %直接型 b=[1，-3，11，27，18]；a=[16，12，2，-4，-1])81.09.01)(5.01()4142136.11)(1(4)(211211------++-+-+=z z z z z z z H %级联型 b0=4;B=[1,1,0;1,-1.4142136,1];A=[1.-0.5.0;1.0.9.0.81]211211126243211214)(------+--++---=z z z z z z z H %并联型 C=0;B=[-14,-12;24,26];A=[1,-2,3;1,-1,1]delta=impseq(0,0,N);x=[one(1,5),zero(1,N-5)];h=filter(b,a,delta); %直接型单位脉冲响应y = filter(b,a,x); %直接型输出信号h=casfilter(b0,B,A,delta); %级联型单位脉冲响应 h=casfilter(b0,B,A,x); %级联型输出响应h=parfiltr(C,B,A,delta); %并联型单位脉冲响应h=parfiltr(C,B,A,x); %并联型输出响应[b0,B,A]=dir2cas(b,a); %直接型转换成级联型[b,a]=cas2dir(b0,B,A); %级联型转换成直接型[C,B,A]=dir2par(b,a); %直接型转换成并联型[b,a]=par2dir(C,B,A); %并联型转换成直接型17. buttord.m用来确定数字低通或模拟低通滤波器的阶次，其调用格式分别是(1)[N,Wn]=buttord(Wp,Ws,Rp,Rs)(2) [N,Wn]=buttord(Wp,Ws,Rp,Rs,’s ’)格式(1)对应数字滤波器, 式中Wp,Ws 分别是通带和阻带的截止频率，实际上它们是归一化频率，其值在0~1之间，1对应抽样频率的一半。
对低通和高通滤波器，Wp,Ws 都是标量，对带通和带阻滤波器，Wp,Ws 都是1⨯2的向量。
Rp,Rs 分别是通带和阻带的衰减，单位为dB 。
N 是求出的相应低通滤波器的阶次，Wn 是求出的3dB 频率，它和Wp 稍有不同。
格式（2）对应模拟滤波器，式中各个变量的含义和格式（1）相同，但Wp,Ws 及Wn 的单位为rad/s, 因此，他们实际上是频率。
18. buttap.m用来设计模拟低通原型滤波器G(p), 其调用格式是[z,p,k]=buttap(N)N 是与设计的低通原型滤波器的阶次，z,p 和k 分别是设计出的G(p)的极点、零点及增益。
19. lp2lp.m, lp2hp.m, lp2bp.m, lp2bs.m以上4个文件的功能是将模拟低通原型滤波器G(p)分别转换为实际的低通、高通、带通及带阻滤波器。