2020年(电子行业企业管理)第一实验室成都信息工程学院电子工程学院

（电子行业企业管理）第一实验室成都信息工程学院
电子工程学院
第一实验室：基础实验篇第Ⅰ部分基本训练题目
第Ⅱ部分简介各题目的原理、程序、效果
第Ⅲ部分基研训练程序软件压缩文件
第Ⅰ部分基本训练题目
1-1-1 序列的图示方法DSP1101
1-1-2 连续信号及采样信号的图示方法DSP 1102
1-1-3单位冲激序列函数impseq
单位冲激序列图示DSP 1103
1-1-4 单位阶跃序列函数stepseq
单位阶跃序列图示DSP1104
1-1-5 矩形序列及图示DSP1105
1-1-6 实指数序列及图示DSP1106
1-1-7 正弦序列及图示DSP1107
1-1-8 复指数序列及图示DSP1108
1-1-9 周期序列及图示DSP 1109
1-1-10 常用5种连续信号及图示DSP1110
1-1-11 离散序列的运算DSP1111
1-1-12 输入序列与系统冲激响应的卷积DSP1112
1-1-13 非零起点时两信号的卷积DSP1113
1-2-1 指数序列的离散时间傅立叶变换DSP 1201
1-2-2 矩形序列的离散时间傅立叶变换DSP1202
1-2-3 离散时间傅立叶变换的性质DSP1203
1-2-4 正弦序列输入，输出为正弦序列，幅度相位因变化DSP1204
1-2-5 模拟信号付氏变换与采样信号的离散时间傅立叶变换
DSP 1205
1-3-1 N点离散傅立叶变换dft(xn,N)
1-3-2 N点离散傅立叶反变换idft(xn,N)
1-3-3 DFT与的Z变换关系DSP1303
1-3-4 DFT与的离散时间傅立叶变换的关系DSP 1304
1-3-5 有限长序列添零填充，得高密度DFT，离散时间付氏频谱不变DSP1305 1-3-6 采样点增多的高分辨率DFT，采样点数少仅添零的高密度DFT DSP1306
1-3-7 DFT的圆周移位函数cirshftt
1-3-8 DFT圆周移位实例DSP1308
1-3-9 圆周卷积DSP1309
1-3-10 复共轭序列的DFT DSP1310
1-3-11 DFT的共轭对称性DSP 1311
1-3-12 补零填充实现线性卷积DSP1312
1-3-13 重迭保留法实现线性卷积DSP1313
1-3-14 重迭保留实现函数ovrlpsav
1-3-15 DET对连续信号作近似谱分析：滤高频，避免混迭频谱；截高时；变有限长序列，避免泄漏频谱DSP1315
1-3-16 采样点为100，进行200点DFT，对进行谱分析DSP 1316 1-3-17 实序列的奇偶分解及DFT的虚实分量DSP1317
1-3-18 实序列的奇偶分解函数DSP1318
1-3-19 用FFT分析信号频率成分DSP1319
1-3-20 用FFT分析语言信号的频谱DSP1320
1-3-21 DCT变换DSP1321
1-3-22 用DCT变换进行语言压缩DSP1322
1-3-23 线性调频Z变换DSP1323
1-3-24 利用CZT计算滤波器100—150HZ频率特性的细节DSP1324
2-1-1 直接型实现系统函数H（Z）的IIR数字滤波器DSP 2101
2-1-2 级联型实现系统H（Z）的IIR数字滤波器DSP2102
2-1-3 级联型实现H（Z）的IIR数字滤波器DSP2103
2-1-4 直接型实现H（Z）的IIR数字滤波器DSP2104
2-1-5 并联型实现H（Z）的IIR数字滤波器DSP 2105
2-1-6 并联型DSP 2106
2-1-7 直接型DSP2107
2-1-8 最终的级联，并联DSP2108
2-1-9 直接型级联型dir2cas(b,a)
2-1-10 级联型直接型cas2par(b0,B,A)
2-1-11 直接型并联型dir2par
2-1-12 并联型直接型par2dire
2-1-13 并联型级联型casfilter
2-1-14 级联型并联型parfilter
2-2-1 FIR直接型滤波器DSP 2201
2-2-2 FIR级联型滤波器DSP2202
2-2-3 FIR的频率取样形式结构DSP 2203
2-2-4 （原例11）由频率样本
求频率采样形式，及单位冲激响应DSP2204
2-2-5 窄带滤波器中的频率采取滤波器是由直接型转换为频率采样型dir2fs(n)
3-1-1 偶对称奇序列的⒈型FIR滤波器的振幅响应hr_type1
3-1-2 偶对称奇序列的及零极点分布DSP3102
3-1-3 偶对称偶序列的Ⅱ型FIR滤波器的振幅响应hr_type2
3-1-4 偶对称偶序列的及零极点分布DSP3104
3-1-5 奇对称奇序列的Ⅲ型FIR滤波器的振幅响应hr_type3
3-1-6 奇对称奇序列的及零极点分布DSP3106
3-1-7奇对称偶序列的Ⅳ型FIR滤波器的振幅响应hr_type4
3-1-8 奇对称偶序列的及零极点分布DSP3108
3-1-9 线性相位FIR滤波器的零点位置有4种可能DSP3109
3-1-10 常用加窗函数DSP3210
3-1-11 对信号用加窗函数的DFT分析频谱DSP3211
3-2-1 计算理想低通滤波器的DSP3201
3-2-2 计算FIR滤波器的绝对和相对的幅度响应DSP3202
3-2-3 提取大于50dB衰减的汉明窗FIR低通滤波器DSP3203
3-2-4 理想高通，偶对称因果序列，N为奇的窗函数，滤波器的单位冲激响应DSP3204
3-2-5 汉宁窗，44dB最小阻带衰减，过度带DSP3205
3-2-6 理想高通，奇对称因果序列，N为偶的窗函数，滤波器的单位冲激响应DSP3206
3-2-7 汉宁窗，44dB最小阻带衰减，过度带DSP3207
3-2-8 理想高通，偶对称因果序列，N为奇的窗函数，滤波器的单位冲激响应DSP3208
3-2-9 设计一个数字FIR带通滤波器DSP3209
3-2-10 理想带通数字滤波器的频率响应DSP3210
3-2-11 设计一个具有相移的数字FIR带通滤波器DSP3211
3-2-12 理想带阻，偶对称因果序列，N为奇的窗函数，滤波器的单位冲激响应ideal-be()
3-2-13 设计一个数字FIR带阻滤波器DSP3213
3-3-1 采样点=0处的频率采样法DSP3301
3-3-2 在过渡带上加两个T1和T2 DSP3302
3-3-3 设计2型FIR低通滤波器DSP3303
3-3-4 设计1型FIR高通滤波器DSP3304
3-3-5 设计4型FIR高通滤波器DSP3305
3-3-6 设计2型FIR带通滤波器DSP3306
3-3-7 设计1型FIR带阻滤波器DSP3307
3-3-8 设计1型FIR低通滤波器DSP3308
3-3-9 设计1型FIR高通滤波器DSP3309
3-3-10 设计4型FIR高通滤波器DSP3310
3-3-11 设计3型FIR带通滤波器DSP3311
3-4-1 用频率响应采样法1设计具有线性相位DSP3401
3-4-2 用窗函数法设计具有线性相位DSP3402
3-4-3 用频率采样法1设计低通滤波器对其进行除噪DSP3403
4-1-1 在MATLAB中用afd_butt(Omegap,Omegar,Ap,Ar)函数来设计巴特沃斯模拟低通滤波器DSP4101
4-1-2 若设计非归一化（≠1）巴特沃斯模拟低通滤波器原型DSP4102
4-1-3 freqs_m(b,a,Omega_max)函数DSP4103
4-1-4 sdir2cas函数DSP4104
4-1-5 设计一个巴特沃斯模拟滤波器DSP4105
4-2-1 用来实现N阶、通带波动为δ的归一化切比学夫1型模拟低通滤波器DSP4201
4-2-2 根据技术指标设计切比学夫1型模拟滤波器DSP4202
4-2-3 设计一个低通切比学夫1型滤波器DSP4203
4-2-4 设计归一化切比学夫2型模拟滤波器DSP4204
4-2-5 根据给定指标设计切比学夫2型模拟滤波器DSP4205
4-2-6 设计一个切比学夫2型低通滤波器DSP4206
4-3-1 用imp_invr函数实现脉冲响应不变法DSP 4301
4-3-2 设计一个巴特沃斯模拟滤波器DSP4302
4-3-3 设计低通数字滤波器DSP4303
4-3-4 设计低通数字滤波器DSP4304
4-4-1 双线性变换法设计低通数字滤波器DSP4401
4-4-2 切比雪夫滤波器原型用双线性变换法设计低通数字滤波器DSP4402 4-5-1
4-5-2
4-5-3
4-5-4 分别设计一个巴特沃斯滤波器和切比雪夫高通滤波器DSP4504
4-5-5 分别设计一个巴特沃斯滤波器和切比雪夫高通滤波器DSP4505
4-5-6 设计一个巴特沃斯带通滤波器DSP4506
4-5-7 设计一个切比雪夫带通滤波器DSP4507
4-5-8 设计一个滤波器DSP4508
4-5-9 设计一个滤波器DSP4509
4-5-10 设计一个滤波器DSP4510
4-6-1 zampping DSP4601
4-6-2用zmapping函数实现例11中的高通滤波器DSP4602
4-6-3切比雪夫1型高通数字滤波器，上述过程由chebhpf函数实现DSP4603 4-6-4用数字频域变换法，设计一个切比雪夫1型高通数字滤波器DSP4604 4-6-5 用双线性变换法设计低通滤波器DSP4605
4-6-6 用脉冲响应不变法设计的低通滤波器对其除噪DSP4606
4-6-7 模拟信号DSP4607
5-1-1下采样DSP 5101
5-1-2 例题DSP5102
5-1-3上采样DSP5103
5-1-4 程序DSP5104
5-1-5 采样率的非整数倍转换DSP5105
5-1-6 程序DSP5106
5-1-7 例题DSP5107
5-1-8 用傅立叶变换对信号进行消噪声处理DSP5108 5-1-9 信号特定频率的提取DSP5109
5-1-10例题DSP5110
5-1-11信号特定频率区间的抑制DSP5111
第Ⅱ部分简介各题目的原理、程序、效果
1-1-1 序列的图示方法DSP1101
原理：数字信号处理中，所有信号都是离散时间信号——序列。

x(n)={…，x(-1),x(0),x(1),…}
如果x(n)={0,5,7,9,6,3,2,1}，-1<=n<=6。

程序：
n=-1:6; %序列的序号依兴次取1至6的各整数
x=[0,5,7,9,6,3,2,1]; %对应序号的序列各值
stem(n,x); %调绘离散序列函数
图形如下：
1-1-2 连续信号及采样信号的图示方法DSP1102
原理：例,当f1=50HZ, f2 =120HZ,
fs =1000HZ时的信号为。

程序：
f1=50;f2=120;fs=1000; %fs表采样频率
t=0:1/fs:1;n=t*fs; %时刻t从0至1，步长为1/fs
y=sin(2*pi*f1*t)+2*sin(1*pi*f2*t);
subplot(211);plot(t(1:50),y(1:50)); title('y(t)');
subplot(212);stem(n(1:50),y(1:50)); title('y(n)');
图形如下：
1-1-3 单位冲激序列函数impseq
单位冲激序列图示DSP 1103
原理：产生,可以用函数impseq
function[x,n]=impseq(n0,n1,n2)
n=[n1:n2]; %n取从n1至n2的各整数
x=[(n-n0)==0]; %n仅当n=n0时x值为1，其它x值为0 单位冲激序列图示DSP 1103
n1=-4; %指定参数
n2=6;
n0=2;
impseq(n0,n1,n2); %调用函数
stem(n,x)
1-1-4单位阶跃序列stepseq
单位阶跃序列图示DSP 1104
function[x,n]=stepseq(n0,n1,n2)
n=[n1:n2]; %n取从n1至n2的各整数
x=[(n-n0)>=0]; %n当n>=n0时x值为1，其它x值为0
单位阶跃序列图示DSP 114
n1=-4;
n2=10;
n0=2;
stepseq (n0,n1,n2);
stem(n,x)
1-1-5矩形序列R N(n)及图示DSP1105
n=[n1:n2]; %n取从n1至n2的各整数
x=[n>=0&n<=N-1]; %当n在0至N-1时x值为1，其它x值为0 n1=-5;
n2=20;
N=8;
stem(n,x)
1-1-6单边实指数序列及图示DSP 1106
N=30;
n=0:N-1;
x=a.^n;
stem(n,x);
1-1-7 单边正弦序列及图示DSP 1107
N=50;
n=0:N-1;
x=sin(*n);;
stem(n,x);
1-1-8 单边复指数序列及图示DSP 1108
N=50;
;
n=0:N-1;
x=exp(*n);;
stem(n,x);
1-1-9 周期序列及图示DSP 1109
n=-1:6;
x1=[0,5,7,9,6,3,2,1];
x=[x1 x1 x1]; %x1是x中的一个周期，要产生3个周期的x序列
stem(n,x);
1-1-10 常用5种连续信号的图示DSP 1110
t=0:0.0001:0.2;
x=sawtooth(2*pi*50*t,1); %调用锯齿波函数
subplot(3,2,1),plot(t,x); %调绘连续曲线函数
x=sawtooth(2*pi*50*t,0.5); %调用三角波函数，与锯齿波差异仅参数由1变0.5
subplot(3,2,2),plot(t,x);
x=square(2*pi*50*t); %调用方波函数
subplot(3,2,3),plot(t,x); axis([0,0.2,-1.5, 1.5]); %后者指定x、y轴取值范围x=tripuls(t,0.1); %调用非周期三角波函数
subplot(3,2,4),plot(t,x); axis([0,0.2,-0.1, 1.1]); %后者指定宽度0.1, 为正轴值二倍
x=rectpuls(t,0.1); %调用非周期方波函数
subplot(3,2,5),plot(t,x); axis([0,0.2,-0.1, 1.1]);
t=-5:0.1:5;x=sinc(t);
subplot(3,2,6),plot(t,x); axis([-5,5-0.4,1.1]);
结果图待补
1-1-11 离散序列的运算DSP1111
n=n1:n2;
X a(n)=X1(n)+X2(n); %信号x1(n)与x2相加
; %信号x1(n)与x2相乘
X c(n)=a*X1(n);
X d(n)=fliplr(X1(n));
X e(n)=sum(x(n1:n2));
X f(n)=proa(X(n1:n2));
X g(n)=sum(abs(X)^2);
X h(n)=X((n1-m):(n2-m));
subplot(2,5,1),stem(n, X1);
subplot(2,5,2),stem(n, X2);
subplot(2,5,3),stem(n, X a);
subplot(2,5,4),stem(n, X b);
subplot(2,5,5),stem(n, X c);
subplot(2,5,6),stem(n, X d);
subplot(2,5,7),stem(n, X e);
subplot(2,5,8),stem(n, X f);
subplot(2,5,9),stem(n, X g);
subplot(2,5,10),stem(n, X h);
1-1-12 输入序列与冲激响应卷积DSP1112
x=[ones(1,10)]; %输入为矩形脉冲序列，为一行10列向量N1=length(x);n1=0:N1-1; %序列长度为N1
N2=20;n2=0:N2-1;
h=0.8.^n2; %冲激响应
N=N1+N2-1;n=0:N-1;
y=conv(x,h); %调用卷积函数，x、h是参数
subplot(311);stem(n1,x);
subplot(312);stem(n2,h);
subplot(313);stem(n,y);
结果待补
1-1-13 非零起点时两信号的卷积
非零起点时两信号的卷积图示DSP 1113
若x,h的起点不为0，则用conv-m计算卷积。

function[y,ny]=convm(x,nx,h,nh)
nyb=nx(1)+nh(1); %两信号起始序号相加，作为输出的起始序号
bye=nx(length(x))+nh(length(h)); %两信号终止序号相加，作输出终止序号ny=[nyb,bye]; %卷积的序号各值
y=conv(x,h);
程序DSP 1113
例
运行结果为
1-2-1指数序列的离散时间傅立叶变换DSP 1201
研究序列的离散傅立叶变换。

解：x(n)是绝对可和的，因此它的DTFT存在。

X(e jw)= =e jw/(e jw-0.8)
流程图如图示
程序实现如下：
n=0:50;x=(0.8).^n;
subplot(221);stem(n,x);title('输入序列');
w=[0:1:500]*2*pi/500;
X=exp(j*w)./(exp(j*w)-0.8*ones(1,501));
magx=abs(X);angx=angle(X);
subplot(223);plot(w/pi,magx);
xlabel('以pi为单位的频率');title('离散时间傅立叶变换幅度'); subplot(224);plot(w/pi,angx);
xlabel('以pi为单位的频率');title('离散时间傅立叶变换相位');
1-2-2矩形序列的离散时间傅立叶变换DSP 1202
原理：，设N-=7
程序：
N=7;n=0;x=[ones(1,N);
k=0:199;w=(pi/100)*k; %将0至轴分为200点
X=x*(exp(-j*pi/100).^(n’*k); %用矩阵向量乘法求DTFT
MagX=abs(X);angX=angle(X);
subplot(3,1,1),stem(n,x);
subplot(3,1,2),plot(w/pi,magX);
subplot(3,1,3), plot(w/pi,angX/pi);
结果待补
1-2-3 离散时间傅立叶变换的性质DSP1203
暂缺
1-2-4 正弦序列输入，输出为正弦序列，幅度相位因变化DSP1204
线性时不变系统，当输入为正弦序列时，则输出也为同频正弦序列，其幅度和相位受H(e jw)影响。

流程图
图1.
程序实现如下：
b=[1,0.5];a=[1.-0.5];
d=impseq(0,0,30);
n=0:30;%在0<=n<=30之间，h(n)截取有限长度
x=cos(0.2*pi*m+pi/4);
h=filter(b,a,d);
y=filter(b,a,x);w=[0:500]*2*pi/500;
w=[0:500]*2*pi/500;
H=freqz(b,a,w);
M=abs(H);A=angle(H);
subplot(231);stem(n,d);title('单位脉冲响应');
subplot(234);stem(n,h);title('单位脉冲响应');
subplot(233);stem(n,x);title('输入信号');
subplot(236);stem(n,y);title('输出信号');
subplot(232);plot(w/pi,M);title('幅度响应');
subplot(235);plot(w/pi,A/pi);title('相位响应');
1-2-5 模拟信号付氏变换与采样信号的离散时间傅立叶变换DSP1205
令x a(t)=，求出并绘制其傅立叶变换Xa(jΩ)。

用f s=5kHz进行采样，求出并画出离散时间傅立叶变换X(e jw)。

程序实现如下：
Dt=0.00005;t=-0.005:Dt:0.005;xa=exp(-1000*abs(t));%模拟信号
Wmax=2*pi*2000;K=500;k=0:1:K;W=k*Wmax/K;
Xa=xa*exp(-j*t'*W)*Dt;Xa=real(Xa);%连续时间傅立叶变换
W=[-fliplr(W),W(2:501)];%频率从-Wmax to Wmax
Xa=[fliplr(Xa),Xa(2:501)];%Xa介于-Wmax和Wmax间
subplot(221);plot(t*1000,xa);xlabel('时间(毫秒)');
ylabel('xa(t)');title('模拟信号')
subplot(222);plot(W/(2*pi*1000),Xa*1000);xlabel('频率(kHz)');
ylabel('Xa(jw)');title('连续时间傅立叶变换')
Ts=0.0002;n=-25:1:25;x=exp(-1000*abs(n*Ts));%离散信号
K=500;k=0:1:K;w=pi*k/K;
X=x*exp(-j*n'*w);X=real(X);%离散时间傅立叶变换
w=[-fliplr(w),w(2:K+1)];
X=[fliplr(X),X(2:K+1)];
subplot(223);stem(n*Ts*1000,x);
xlabel('时间(毫秒)');gtext('Ts=0.2毫秒');
ylabel('x1(n)');title('离散信号')
subplot(224);plot(w/pi,X);xlabel('频率(弧度)');
ylabel('X1(w)');title('离散傅立叶变换')
1-3-1 N点离散傅立叶变换dft(xn,N)
设x(n)是一个长度为N的有限长序列，定义x(n)的N点离散傅立叶变换为X(k)= k=0,1,…,N-1
function[Xk]=dft(xn,N)
n=[0:1:N-1];%n的行向量
k=[0:1:N-1];%k的行向量
WN=exp(-j*2*pi/N);%旋转因子
nk=n'*k;%产生一个含nk值的N乘N维矩阵
WNnk=WN.^nk;%DFT矩阵
Xk=xn*WNnk;%DFT系数的行向量
1-3-2 N点离散傅立叶反变换idft(Xk,N)
function[xn]=idft(Xk,N)
n=[0:1:N-1];%n的行向量
k=[0:1:N-1];%k的行向量
WN=exp(-j*2*pi/N);%旋转因子
nk=n'*k;
WNnk=WN.^(-nk);%DFT矩阵
xn=(Xk*WNnk)/N;%DFT系数的行向量
1-3-3DFT与的Z变换关系DSP1303
X(k)=X（Z）
序列DFT的物理意义：序列x(n)的N点DFT是x(n)的z变换在单位圆上的N点等间隔采样；X(k)为x(n)的离散傅立叶变换X(e jw)在区间[0:2π]上的N点等间隔采样。

暂缺待补充
1-3-4 DFT与的离散时间傅立叶变换的关系DSP1304
后者在区间[0，2pi]上的N个等间隔采样
暂缺待补充
1-3-5 有限长序列添零填充，得高密度DFT，离散时间付氏频谱不变DSP1305 例x(n)=R5(n),求X(e jw)及N分别取10，20的X(k)。

解：设N=10，则
k=0,1,…,9
设N＝20，则
k=0,1, (19)
流程图：
图1-3-5
程序实现如下：
n=0:4;x=[ones(1,5)];
N1=10;n1=0:1:N1-1;
N2=20;n2=0:1:N2-1;
x1=[ones(1,5),zeros(1,N1-5)];
X1=fft(x1,N1);%N=10点离散傅立叶变换
magX1=abs(X1);
k1=(0:length(magX1)'-1)*N1/length(magX1);
x2=[ones(1,5),zeros(1,N2-5)];
X2=fft(x2,N2); magX2=abs(X2);
k2=(0:length(magX2)'-1)*N2/length(magX2);
subplot(321);stem(n,x);ylabel('x(n)');
subplot(323);stem(n1,x1);ylabel('x(n)');
subplot(324);stem(k1,magX1);ylabel('/X(k)/');
subplot(325);stem(n2,x2);ylabel('x(n)');
subplot(326);stem(k2,magX2);ylabel('/X(k)/');
结论：
①填零是给原序列填零的运算，会给原始序列的离散时间傅立叶变换提供间隔较密的样本。

②为画出X(e jw)，只需要5点的X(k)用内插公式即可得到X(e jw)。

但实际上是用10或20点的X(k)来填充X(e jw)的值。

③填零运算提供了较密的频谱，而没有增加任何新的信息，因此它不能提供高分辨率的频谱。

④为得到高分辨率的频谱，需从实验或观察中取得更多的数据。

1-3-6 采样点增多的高分辨率DFT，采样点数少仅添零的高密度DFT DSP1306 为了说明高密度和高分辨率之间的区别，考察序列
x(n)=2cos(0.35πn)+cos(0.5πn)
①当0≤n＜10时，确定并画出x(n)的；离散傅立叶变换。

②当x(n)=｛x(n) 0≤n＜10
时，确定并画出x(n)的离散傅立叶变换。

0 10≤n＜40
③当0≤n＜40时，确定并画出x(n)的离散傅立叶变换。

流程图：
图1-3-6
程序实现如下：
N1=10;N2=40;n1=0:N1-1;n2=0:N2-1;
x=2*cos(0.35*pi*n)+cos(0.5*pi*n);
x1=x(1:N1);
Y1=dft(x1,N1);magY1=abs(Y1);
k1=0:N1-1;w1=2*pi/N1*k1;
x2=[x1 zeros(1,N2-N1)];
Y2=dft(x2,N2);magY2=abs(Y2);
k2=0:N2-1;w2=2*pi*k2/N2;
Y3=dft(x,N2);magY3=abs(Y3);
k3=0:N2-1;w3=2*pi/N2*k;
subplot(231);stem(n1,x1);title('没有足够点的采样信号');
subplot(234);stem(w1/pi,magY1);title('信号的频谱');
subplot(232);stem(n2,x2);title('添零信号');
subplot(235);stem(w2/pi,magY2);title('高密度频谱');
subplot(233);stem(n2,x);title('有足够采样点的信号');
subplot(236);stem(w3/pi,magY3);title('高分辨率频谱');
图1-3-6高密度与高分辨率频谱
结论：
①当0≤n＜10时的序列x(n)与X(k)，从X(k)几乎无法看出有关信号的频谱的信息。

②将x(n)补30个零时的x(n)y与X(k),这时的频谱相当密，但从中很难看出信号的频谱成分，故成为高密度频谱。

③将x(n)的长度加长到40时的x(n)与X(k),这是可以清晰的看出信号的频谱成分（w1=0.35π,w2=0.5π）,故成为高分辩率频谱。

1-3-7圆周移位性质cirshftt(x,m,N)
设x(n)是一个长度为N的有限长序列，圆周移位定义为
y(n)=x((n+m))N R N(n)
（1－7）
将x(n)以N为周期进行周期延拓得到～
x(n)=x((n))N,再将～
x(n)
左移m位得到
～
x(n+m)，最后取～
x(n+m)
的主值序列，则得到有限长序列x(n)的周期移位序列y(n)，y(n)
仍为长度为N的有限长序列
((n))N表示n对N求余，即如果n=MN+n1, 0≤n1＜N-1,M为整数，则((n))N = n1 。

圆周移位用cirshftt实现如下：
function y=cirshftt(x,m,N)
if length(x)>N
error('N 必须>=x的长度')
end
x=[x zeros(1,N-length(x))];
n=0:N-1;
n=mod(n-m,N);
y=x(n+1);
1-3-8 DFT园周移位实例DSP1308
例序列x(n)={9,8,7,6,5,4,3,2,1},求分别移位1，3，5，7，9位的圆周移位。

程序实现如下：
n=0:8;x=[9,8,7,6,5,4,3,2,1];
y1=cirshftt(x,1,9);
y2=cirshftt(x,3,9);
y3=cirshftt(x,5,9);
y4=cirshftt(x,7,9);
y5=cirshftt(x,9,9);
subplot(611);stem(n,x);ylabel('x(n)');
subplot(612);stem(n,y1);ylabel('y1(n)');
subplot(613);stem(n,y2);ylabel('y2(n)');
subplot(614);stem(n,y3);ylabel('y3(n)');
subplot(615);stem(n,y4);ylabel('y4(n)');
subplot(616);stem(n,y5);ylabel('y5(n)');
图1-3-8 序列圆周移位
1-3-9 圆周卷积DSP1309
例计算两序列x1(n)={1,2,2,3}；x2(n)={1,2,3,4,3,2}的圆周卷积。

流程图：
图1-3-9
程序实现如下：
x1=[1,2,2,3];x2=[1,2,3,4,3,2];
N=length(x1)+length(x2)-1;n=0:N-1;n1=0:N-2;n2=0:N+1;
y1=circonvt(x1,x2,N-1);
y2=circonvt(x1,x2,N);
y3=circonvt(x1,x2,N+2);
y4=conv(x1,x2);
M=N+2;m=0:M-1;
x1=[x1 zeros(1,M-length(x1))];
x2=[x2 zeros(1,M-length(x2))];
X1=dft(x1,M);
X2=dft(x2,M);
X=X1.*X2;
x=idft(X,M);x=real(x);
subplot(241);stem(m,x1);title('x1(n)');
subplot(242);stem(m,x2);title('x2(n)');
subplot(243);stem(n1,y1);title('N-1点圆周卷积');
subplot(244);stem(n,y2);title('N点圆周卷积');
subplot(245);stem(n2,y3);title('N+2点圆周卷积');
subplot(246);stem(n,y4);title('一般卷积运算');
subplot(247);stem(m,x);title('x(n)=IDFT[X(k)]');
图1-1-9 圆周卷积
结论：两序列，若x1的长度为N，x2的长度为M。

L≥N＋M－1时，循环卷积等于线性卷积。

L＝N＋M－1时，不管时循环卷积也好，还是线性卷积也好，可以用一般卷积公式进行计算，因为三者的结果时一样的。

1-3-10 计算共轭序列x(n)={1-j,2+2j,3-3j,-4+4j,5-5j}的DFT和x*(n)的
DFT DSP1310
xn=[1-j,2+2j,3-3j,-4+4j,5-5j];
Xk=dft(xn,5);x1=(xn').';X=dft(x1,5);n=0:4;
subplot(221);stem(n,abs(xn));title('/x(n)/');
subplot(222);stem(n,abs(Xk));title('/x(k)/');
subplot(223);stem(n,abs(x1));title('/x*(n)/');
subplot(224);stem(n,abs(X));title('/X*(N-k)/');
图1-3-10 复共轭序列的DFT
1-3-11DFT共轭对称性DSP1311
暂缺待补
1-3-12补零填充实现线性卷积DSP1312
暂缺待补
1-3-13设x(n)={10,9,8,7,6,5,4,3,2,1},h(n)={1,1,-1}按N=6用重叠保留方法计算y(n)=x(n)*h(n) DSP1313
x=[10,9,8,7,6,5,4,3,2,1];h=[1,1,-1];
x1=[0,0,10,9,8,7];x2=[8,7,6,5,4,3];x3=[4,3,2,1,0,0];
y1=circonvt(x1,h,6);y2=circonvt(x2,h,6);
y3=circonvt(x3,h,6);y=ovrlpsav(x,h,6);;
n=0:5;N=length(x)+length(h)-1;n1=0:N-1;n2=0:9;
subplot(241);stem(n,x1);title('x1');axis([0,6,0,10]);
subplot(245);stem(n,y1);title('y1');axis([0,6,-10,20]);
subplot(242);stem(n,x2);title('x');axis([0,6,0,10]); subplot(246);stem(n,y2);title('y2');axis([0,6,-10,20]); subplot(243);stem(n,x3);title('x3');axis([0,6,0,10]); subplot(247);stem(n,y3);title('y3');axis([0,6,-10,20]); subplot(244);stem(n2,x);title('x');axis([0,11,0,10]); subplot(248);stem(n1,y);title('y');axis([0,11,-10,20]);
1-3-14重叠保留法实现函数ovrlpsav
function[y]=ovrlpsav(x,h,N)
Lenx=length(x);M=length(h);%X输入序列，h脉冲响应M1=M-1;L=N-M1;%N段长
h=[h zeros(1,N-M)];
x=[zeros(1,M1),x,zeros(1,N-1)];%予置M-1个零
k=floor((Lenx+M1-1)/(L));% 段数
Y=zeros(k+1,N);
for k=0:k %各段园卷积
xk=x(k*L+1:k*L+n);
Y(k+1,:)=circonvt(xk,h,N);
end
Y=Y(:,M:N)'; %去掉前M-1个值
y=(Y(:))': %装成输出
1-3-15DFT对连续信号作近似谱分析DSP1315
例
x a(t)幅度的估计对模拟信号x a(t)＝2sin（4πt）＋5cos（8πt）以时间x a(t)间隔T对其采样，得到N点序列x（n），用N点DFT得到对x a(t)幅度的估计。

（1）T＝0.01s。

N＝40或N＝50，一个能提精确x a(t) 的幅度谱，画出DFT 的幅度普。

（2）T＝0.005s，N＝40或N＝50，画出DFT的幅度谱。

流程图
图1.8
程序实现如下：
T=0.01;N=40;n=0:N-1;t=n*T;
xn=2*sin(4*pi*t)+5*cos(8*pi*t);
Xk=dft(xn,N);
magXk=abs(Xk);
k=(0:length(magXk)'-1)*N/length(magXk);
subplot(241);plot(t,xn);axis([0,0.4,-7.5,7]);
title('T=0.01s,t=0.4s');ylabel('x(t)');
subplot(245);stem(k,magXk);title('T=0.01s,N=40');ylabel('X(k)');
图1－12 用DFT进行频谱分析
从图上可以看出，采样间隔T＝0.01s，采样点数N＝50是的幅度频谱是最精确。

T＝0.005两种情况都存在频谱泄漏。

1-3-16 采样间隔、采样点数变化时频谱样值比较DSP1316
例已知一模拟信号x a(t)＝e-t u(t),现以采样率fs＝20Hz进行采样。

用DFT计算当序列长度①L＝100，②L＝20时，N＝200点地幅度频谱样值并通过作图与理论上准确地频谱样值进行比较。

解：原信号的傅立叶变换
其幅度为| Xa(jΩ)|=1/(1+Ω2)
流程图
图1.9
程序实现如下：
fs=20;
L1=100;N=200;n1=0:L1-1;t1=n1/fs;L2=20,n2=0:L2-1;t2=n2/fs;
xn1=exp(-t1);xn=[xn1,zeros(1,N-L1)];
Xk1=dft(xn,N);magXk1=abs(Xk1);
k1=(0:length(magXk1)'-1)*N/length(magXk1);
xn2=exp(-t2);xn=[xn2,zeros(1,N-L2)];
Xk2=dft(xn,N);magXk2=abs(Xk2);
k2=(0:length(magXk2)'-1)*N/length(magXk2);
Omeger=0:0.1:20*pi;
Xa=1./(1+Omeger.^2);
subplot(231);plot(t1,xn1);title('xa(t) t=5s');
subplot(232);plot(t2,xn2);title('xa(t) t=1s');
subplot(233);plot(k1,magXk1);title('X(k) L1=100 N=200');
subplot(234);plot(k2,magXk2);title('X(k) L2=20 N=200');
subplot(235);plot(Omeger/pi,Xa);title('/Xa(j\Omega)/');
图1－13 用DFT计算的频谱
结论：
①当序列长度为100，进行200点DFT计算的结果混叠与泄漏的影响比较小，基本上接近原信号的频谱。

因为按给定的fs＝20Hz，相当于取信号的最好频率fh＝10Hz，故在[0,fh]
| Xa(jΩ)|2dΩ=0.495
频率范围内的信号能量为E h＝1/2π∫20π
-20π
| Xa(jΩ)|2dΩ=0.5
信号的总能量为E x＝1/2π∫∞
-∞
E h/ E x = 99%，基本上满足频谱不混叠的要求。

②当序列长度为20，进行200点DFT计算，由于截取x(n)长度太短x(t)|t=LT=e-LT=1/e=0.3079>>0
所以频谱泄漏出现较大的波动，以致与原信号频谱有较大差别。

㈡用DFT对离散信号进行频谱分析序列x(n)在单位圆上的z变换就是傅立叶变换X（e jw），即X（e jw）＝X(z)|z= e jw
对序列x(n)进行N点DFT得到X（k），X（k）是X（e jw）在区间[0,2π]上的N点等间隔采样，因此序列的傅立叶变换可利用DFT来计算。

1-3-17实序列奇偶分解及DFT的虚实分解DSP1317
例：设x（n）＝0.5（0.8）n 0≤n≤20
（1）分解x（n）成x ec（n）和x oc（n）；（奇偶部分）
（2）检验序列的性质。

DFT[x ec（n）]=Re[X(k)]
DFT[x oc（n）]=Im[X(k)]3
程序实现如下：
N=20;n=0:N-1;x=5*(0.8).^n;
[xec,xoc]=circevod(x);
X=dft(x,N);Xec=dft(xec,N);Xoc=dft(xoc,N);
subplot(241);stem(n,x);title('x(n)');
subplot(242);stem(n,abs(X));title('abs[X(k)]');
subplot(243);stem(n,real(X));title('Re[X(k)]');
subplot(244);stem(n,imag(X));title('Im[X(k)]');
subplot(245);stem(n,xec);title('xec(n)');
subplot(246);stem(n,xoc);title('xoc(n)');
subplot(247);stem(n,real(Xec));title('DFT[xec(n)]');
subplot(248);stem(n,imag(Xoc));title('DFT[xoc(n)]');
图1－15 DFT的实部和虚部
结论：
实序列的偶分量关于N/2点对称，奇分量关于N/2点反对称，偶分量的
DFT等于实序列的DFT的实部，奇分量的DFT等于实序列的DFT的虚部。

1-3-18 实序列奇偶分解函数circevod(x)
function [xec,xoc]=circevod(x)
if any(imag(x)~=0)
error('x非实数序列')
end
N=length(x);n=0:N-1;
xec=0.5*(x+x(mod(-n,N)+1));
xoc=0.5*(x-x(mod(-n,N)+1));
1-3-19用FFT分析信号频率成分DSP1319
一被噪声污染的信号，很难看出它所包含的频率分量，如一个由50Hz和120Hz正弦信号构成的信号，受到均值随机噪声的干扰，数据采样率为1000Hz。

通过FFT来分析其信号频率成分，用MATLAB实现如下：
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(311);plot(x);title('原噪声信号');
subplot(312);plot(y);title('受噪声污染的信号');
n=0:511;f=1000*n/512;
subplot(313);plot(f,abs(Y));title('FFT');
图1－16
1-3-20 用FFT分析信号频率成分DSP1320
（1）用FFT分析语音信号的频谱
load mtlb;
subplot(221);plot(mtlb);title('原语音信号');
y=fft(mtlb);
subplot(222);plot(abs(y));title('FFT变换');
y(abs(y)<1)=0;x=ifft(y);
subplot(223);plot(abs(y));title('去掉幅值小于1的FFT变换值'); subplot(224);plot(real(x));title('重构语音信号');
图1－17 语音信号的频谱
1-3-21 计箅信号的DCT变换DSP1321
N=100;n=1:N;
x=n+50*cos(2*pi*n/40);
y=dct(x);
subplot(121);stem(x);title('ԭʼÐźÅ')
subplot(122);stem(y);title('DCT±ä»¯')
1-3-22 用DCT变换进行语音压缩DSP1322
load mtlb
subplot(221);plot(mtlb);title('ԭʼÓïÒôÐźÅ');
y=dct(mtlb);
subplot(222);plot(y);title('DCT񄯯');
y(abs(y)<1)=0;x=idct(y);
subplot(223);plot(y);title('È¥µô·ùֵСÓÚ1µÄDCT±ä»»Öµ'); subplot(224);plot(y);title('Öع¹ÓïÒôÐźÅ');
1-3-23 线形调频Z变换DSP1323
已知信号, 用Z变换求其前30点的复频请谱，给定A0=0.8, 程序如下
subplot(1,3,1),zplane(zk.’);
;
subplot(1,3,2),stem(k,abs(Z));title(‘CZT’);
y=fft(X,M);
subplo t(1,3,3),stem(k,abs(y));title(‘FFT’);
1-3-24利用czt计算滤波器x在100-150hz频率特性的细节DSP1324
X=firl(30,125/500,boxcar(31);
Fs=1000;
f1=100;f2=150;
n=1024;
w=wxp(-j*2*pi*(f2-1)/(m*Fs);
a=exp(j*2*pi*f1/Fs)
y=fft(X,1000);
z=czt(X,m,w,a);
fy=0:(length(y)-1)’*1000/:leng th(y);
fz=((0:length(z)-1)’*(f2-f1)/length(z))+f1;
subplot(1,2,1),plot(fy(1:500),abs(y(1:500)));title(‘FFT’);
subplot(1,2,2),plot(fz,abs(z));title(‘CZT’);
2-1-1 直接型系统实现函数filter
直接型系统的IIR数字滤波器DSP2101
例用直接型实现系统函数H（z）＝的IIR数字滤波器，求单位脉冲响应和单位阶跃信号的输出。

图2.1
程序实现如下：
b=[1,-3,11,27,18];a=[16,12,2,-4,-1];
N=25;delta=impseq(0,0,N);
x=[ones(1,5),zeros(1,N-5)];%单位阶跃信号
h=filter(b,a,delta);%直接型单位脉冲响应
y=filter(b,a,x);
subplot(121);stem(h);title('直接型h(n)');
subplot(122);stem(y);title('直接型y(n)');
效果图：
图2-1-1 直接型单位脉冲响应和输出信号
2-1-2 级联型实现函数casfiltr
级联型系统函数可表示为
H（z）＝b0
在MATLAB中给定级联型系统函数，由扩展函数casfiltr实现IIR的级联形式。

function y=casfiltr(b0,B,A,x);
[K,L]=size(B);
N=length(x);
w=zeros(K+1,N);
for i=1:K
w(i+1,:)=filter(B(i,:),A(i,:),w(i,:));
end
y=b0*w(K+1,:);
2-1-3 级联型系统的IIR数字滤波器DSP2103
例系统函数H（z）＝，用级联型结构实现。

流程图：
图2.2
级联型实现系统DSP2102 程序实现如下:
b0=4;B=[1,1,0;1,-1.4142136,1];A=[1,-0.5,0;1,0.9,0.81]; N=60;delta=impseq(0,0,N);
h=casfiltr(b0,B,A,delta);
x=[ones(1,5),zeros(1,N-5)];
y=casfiltr(b0,B,A,x);
subplot(121);stem(h);title('级联型h(n)');
subplot(122);stem(y);title('级联型y(n)');
图2－2级联型单位脉冲响应和输出信号
2-1-4 并联型实现函数parfilter
并联型系统函数可表示为
H（z）＝
在MATLAB用给定并联系统函数，由扩展函数parfiltr实现IIR的并联形式。

function y=parfiltr(C,B,A,x);
[K,L]=size(B);
N=length(x);
w=zeros(K+1,N);
w(1,:)=filter(C,1,x);
for i=1:K
w(i+1,:)=filter(B(i,:),A(i,:),x);
end
y=sum(w);
2-1-5 并联型系统的IIR数字滤波器DSP2105
例用并联型实现系统函数为（z）＝的滤波器。

流程图：
图2.5
程序实现如下：
C=0;B=[-14.75,-12.90;24.50,26.82];A=[1,-7/8,3/32;1,-1,0.5];
N=60;delta=impseq(0,0,N);
h=parfiltr(C,B,A,delta);
x=[ones(1,5),zeros(1,N-5)];
y=parfiltr(C,B,A,x);
subplot(121);stem(h);title('并联型h(n)');
subplot(122);stem(y);title('并联型y(n)');
图2－5并联型单位脉冲想响应和输出信号
2-1-6 直接型级联型函数dir2cas
编写函数dir2cas来实现IIR直接型转换为级联型。

其在matalb中实现程序为：
function[b,a]=dircas(b,a); %直接型转换成级连型);
b0=b(1);b=b/b0; %计算系统增益系数a0=a(1);a=a/a0;
b0=b0/a0;
M=length(b);N=length(a); %分子、分母长度补齐if N>M
b=[b zeros(1,N-M)];
elseif M>N
a=[a zeros(1,M-N)];N=M;
else
NM=0;
end
k=floor(N/2);B=zeros(k,3);A=zeros(k,3); %初始化系数矩阵
if k*2==N;
b=[b 0];
a=[a 0];
end
broots=cplxpair(roots(b)); %实数从小到大成对排列
aroots=cplxpair(roots(a));
for i=1:2:2*k %变换多项式系数
Brow=broots(i:1:i+1,:);
Brow=real(poly(Brow));
B(fix(i+1)/2,:)=Brow;
Arow=aroots(i:1:i+1,:);
Arow=real(poly(Arow));
A(fix(i+1)/2,:)=Arow;
end
程序框图如下：
2-1-7 直接型级联型的IIR滤波器实现DSP2107
用级联型实现系统函数为的IIR数字滤波器，求出单位脉冲响应和单位阶跃信号的输出。

程序原代码如下：
N=25;delta=impseq(0,0,N);
b=[1 -3 11 27 18];a=[16 12 2 -4 -1];
[b0,B,A]=dir2cas(b,a);
h=casfiltr(b0,B,A,delta);
x=[ones(1,5),zeros(1,N-5)];
y=casfiltr(b0,B,A,x);
subplot(121);stem(h);title('级连型h(n)');
subplot(122);stem(y);title('级连型y(n)')
波形展示：
2-1-8 级联型直接型函数cas2dir
编写函数cas2dir来实现IIR级联型转换为直接型。

其在matalb中实现程序为：
function[b,a]=cas2dir(b0,B,A); %级连型转换成直接型
[K,L]=size(B);
b=[1];
a=[1];
for i=1:K %将多项式相乘
b=conv(b,B(i,:));
a=conv(a,A(i,:));
end
b=b*b0;
程序框图如下：
2-1-9 级联型直接型的IIR滤波器实现DSP2109
用直接型实现系统函数为的IIR数字滤波器，求出单位脉冲响应和单位阶跃信号的输出。

程序原代码如下：
b0=4;B=[1,1,0;1,-1.4142136,1];
A=[1,-0.5,0;1,0.9,0.81];
N=60;
delta=impseq(0,0,N);
x=[ones(1,5),zeros(1,N-5)];
[b,a]=cas2dir(b0,B,A);
h=filter(b,a,delta);
y=filter(b,a,x);
subplot(121);stem(h);title('直接型h(n)');
subplot(122);stem(y);title('直接型y(n)');
波形展示：
2-1-10 直接型并联型函数dir2par
编写函数dir2par来实现IIR直接型转换为并联型。

其在matalb中实现程序为：
function[C,B,A]=dir2par(b,a); %直接型转换成并联型
M=length(b);N=length(a);
[r1,p1,C]=residue(b,a); %多项式分解
p=cplxpair(p1,10000000*eps);
I=cplxcomp(p1,p); %复共轭比较
r=r1(I);
K=floor(N/2);B=zeros(K,2);A=zeros(K,3); %初始化矩阵
if K*2==N
for i=1:2:N-2 %合并因子为二阶多项式分式Brow=r(i:1:i+1,:);
Arow=p(i:1:i+1,:);
[Brow, Arow]=residuez(Brow,Arow,[]);
B(fix((i+1)/2),:)=real(Brow);
A(fix((i+1)/2),:)=real(Arow);
end
[Brow,Arow]=residuez(r(N-1),p(N-1),[]);
B(K,:)=[real(Brow') 0];
A(K,:)=[real(Arow') 0];
else
for i=1:2:N-1
Brow=r(i:1:i+1,:);
Arow=p(i:1:i+1,:);
[Brow,Arow]=residuez(Brow,Arow,[]);
B(fix((i+1)/2),:)=real(Brow);
A(fix((i+1)/2),:)=real(Arow);
end
end
程序框图如下：
2-1-11 直接型系统用并联型实现DSP2111
用并联型实现系统函数为的IIR数字滤波器，求出单位脉冲响应和单位阶跃信号的输出。

程序原代码如下：
b=[1 -3 11 27 18];a=[16 12 2 -4 -1];
N=25;delta=impseq(0,0,N);
[C,B,A]=dir2par(b,a);
h=parfiltr(C,B,A,delta);
x=[ones(1,5),zeros(1,N-5)];
y=parfiltr(C,B,A,x);
subplot(121);stem(h);title('并联型h(n)');
subplot(122);stem(y);title('并联型y(n)');
波形展示：
2-1-12 并联型直接型函数par2dir
若给定系统函数并联形式，也可由扩展函数par2dir将其转换成直接形式，再由filter 直接形式实现。

function[b,a]=par2dir(C,B,A);
[K,L]=size(A);R=[];P=[];
for i=1:K
[r,p,k]=residuez(B(i,:),A(i,:));
R=[R;r];P=[P;p];
end
[b,a]=residuez(R,P,C);
b=b(:)';a=a(:)';
2-1-13 并联型系统用直接型实现DSP2113
用直接型实现系统函数为的IIR数字滤波器，求出单位脉冲响应和单位阶跃信号的输出。

程序原代码如下：
C=0;B=[-14.75,-12.90;24.50,26.82];
A=[1,-7/8,3/32;1,-1,0.5];
N=60;
delta=impseq(0,0,N);
[b,a]=par2dir(C,B,A);
h=filter(b,a,delta);
x=[ones(1,5),zeros(1,N-5)];
y=filter(b,a,x);
subplot(121);stem(h);title('直接型h(n)');
subplot(122);stem(y);title('直接型y(n)');
波形展示：
2-1-14 并联型级联型DSP2114
IIR数字滤波器级联型与并联型网络结构的相互转换可通过直接型来完成
并联型级联型，可先将并联型转换为直接型，再将其转换为级联型。

例如：
用级联型实现系统函数为的IIR数字滤波器，求出单位脉冲响应和单位阶跃信号的输出。

程序原代码如下：
C=0;B=[-14.75,-12.90;24.50,26.82];
A=[1,-7/8,3/32;1,-1,0.5];
N=60;
delta=impseq(0,0,N);
[b,a]=par2dir(C,B,A);
[b0,B,A]=dir2cas(b,a);
h=casfiltr(b0,B,A,delta);
x=[ones(1,5),zeros(1,N-5)];
y=casfiltr(b0,B,A,x);
subplot(121);stem(h);title('级连型h(n)');
subplot(122);stem(y);title('级连型y(n)');
2-1-15 级联型并联型DSP2115
IIR数字滤波器级联型与并联型网络结构的相互转换可通过直接型来完成
级联型并联型，可先将级联型转换为直接型，再将其转换为并联型。

例如：用并联型实现系统函数为的IIR数字滤波器，求出单位脉冲响应和单位阶跃信号的输出。

程序原代码如下：
b0=4;B=[1,1,0;1,-1.4142136,1];
A=[1,-0.5,0;1,0.9,0.81];
N=60;
delta=impseq(0,0,N);
x=[ones(1,5),zeros(1,N-5)];
[b,a]=cas2dir(b0,B,A);
[C,B,A]=dir2par(b,a);
h=parfiltr(C,B,A,delta);
x=[ones(1,5),zeros(1,N-5)];
y=parfiltr(C,B,A,x);
subplot(121);stem(h);title('并联型h(n)');
subplot(122);stem(y);title('并联型y(n)');
波形展示：。