数字信号处理实验参考程序
数字信号处理实验参考程序
此程序仅供参考非标准答案。
实验一Q1-4clear, close all;N = 40;num = [0.9 -0.45 0.35 0.002];den = [1 0.71 -0.46 -0.62];y = impz(num,den,N);% Plot the impulse responsestem(y);xlabel('Time index n'); ylabel('Amplitude');title('Impulse Response'); grid;Q1-5clear,close all,n=0:39;x=delta(n);num = [0.9 -0.45 0.35 0.002];den = [1 0.71 -0.46 -0.62];y=filter(num,den,x);stem(n,y);xlabel('Time index n');ylabel('Amplitute');title('Impulse Response');grid;Q1-6clear, close all;N = 40;num = [0.9 -0.45 0.35 0.002];den = [1 0.7 -0.46 -0.62];y = stepz(num,den,N);% Plot the impulse responsestem(y);xlabel('Time index n'); ylabel('Amplitude');title('Impulse Response'); grid;Q1-8% Program_P1_10clear, close all;h = [3 2 1 -2 1 0 -4 0 3 3 2 1 -2 0 -4]; % impulse response x = [1 -2 3 -4 3 2 1 1 -2 3]; % input sequencey = conv(h,x); % convolutionn = 0:23;subplot(2,1,1);stem(n,y, '. ');title('Output Obtained by Convolution'); grid;x1 = [x zeros(1,14)];y1 = filter(h,1,x1);subplot(2,1,2);stem(n,y1, '. ');xlabel('Time index n'); ylabel('Amplitude');title('Output Generated by Filtering'); grid;% Program_P1_4% Stability test based on the sum of the absolute values of the impulse response samples clear,close all;num = [1 -4 3]; den = [1 -1.7 1]; % Define a systemN = 200;h = impz(num,den,N+1);parsum = 0;for k = 1:N+1;parsum = parsum + abs(h(k));if abs(h(k)) < 10^(-6), break, endend% Plot the impulse responsen = 0:N;stem(n,h)xlabel('Time index n'); ylabel('Amplitude');% Print the value of abs(h(k))disp('V alue =');disp(abs(h(k)));实验二% Program P2_3% Evaluation of the DTFTclear, close all;% Compute the frequency samples of the DTFTw = -4*pi:8*pi/511:4*pi;num = [0.7 -0.5 0.3 1]; den = [1 0.3 -0.5 0.7];h = freqz(num, den, w);% Plot the DTFTsubplot(2,2,1)plot(w/pi,real(h)); grid; title('Real part of H(e^{j\omega})')xlabel('\omega /\pi'); ylabel('Amplitude');subplot(2,2,2)plot(w/pi,imag(h)); grid;title('Imaginary part of H(e^{j\omega})')xlabel('\omega /\pi'); ylabel('Amplitude');subplot(2,2,3)plot(w/pi,abs(h)); grid;title('Magnitude Spectrum |H(e^{j\omega})|')subplot(2,2,4)plot(w/pi,angle(h)); grid;title('Phase Spectrum arg[H(e^{j\omega})]')xlabel('\omega /\pi'); ylabel('Phase in radians');% Program P2_4% Evaluation of the DTFTclear, close all;% Compute the frequency samples of the DTFTw = -4*pi:8*pi/511:4*pi;num = [1 3 5 7 9 11 13 15 17]; den = [1 0 0 0 0 0 0 0 0];h = freqz(num, den, w);% Plot the DTFTsubplot(2,2,1)plot(w/pi,real(h)); grid; title('Real part of H(e^{j\omega})')xlabel('\omega /\pi'); ylabel('Amplitude');subplot(2,2,2)plot(w/pi,imag(h)); grid;title('Imaginary part of H(e^{j\omega})')xlabel('\omega /\pi'); ylabel('Amplitude');subplot(2,2,3)plot(w/pi,abs(h)); grid;title('Magnitude Spectrum |H(e^{j\omega})|')xlabel('\omega /\pi'); ylabel('Amplitude');subplot(2,2,4)plot(w/pi,angle(h)); grid;title('Phase Spectrum arg[H(e^{j\omega})]')xlabel('\omega /\pi'); ylabel('Phase in radians');2.5% Program2.5clear, close all;w = -4*pi:0.02:4*pi;n = 0:8;x = (n+1);X = freqz(x,1,w);subplot(221)stem(n,x,'.'), title('Sequence x[n]'), xlabel('Time index n')subplot(222)plot(w/pi,abs(X)), title('M agnitude of DTFT of x[n]'), axis([-4,4,0,max(abs(X))]) subplot(224)plot(w/pi,angle(X)), title('Phase of DTFT of x[n]'), xlabel('Frequency*pi in ratians') axis([-4,4,min(angle(X)),max(angle(X))])% Program_2_12clear, close all;L = input('L=');k = 0:L-1;x = input('x=');X = fft(x,L);y=ifft(X,L);subplot(221)stem(x,'.'), title('Sequence x[n]'), xlabel('Time index n'), axis([0,15,min(x),max(x)]) subplot(222)stem(k,abs(X),'.'), title('M agnitude of DFT of x[n]'), axis([0,15,0,max(abs(X))]) subplot(224)stem(k,angle(X),'.'), title('Phase of DFT of x[n]'), xlabel('Frequency index k')axis([0,15,min(angle(X)),max(angle(X))])subplot(223)stem(x,'.'), title('Sequence y[n]'), xlabel('Time index n'), axis([0,15,min(x),max(x)])% Program_2_15clear,close all;g=input('g=');h=input('h=');N=length(g);x=g+j*h;X=fft(x,N);for k=0:N-1G(k+1)=(X(k+1)+conj(X(mod(-k,N)+1)))/2,H(k+1)=(X(k+1)-conj(X(mod(-k,N)+1)))./2jendsubplot(221)stem(g,'.'), title('Sequence g[n]'), xlabel('Time index n'), axis([0,N,0,max(g)]) subplot(222)stem(h,'.'), title('Sequence h[n]'), xlabel('Time index n'), axis([0,N,0,max(h)]) subplot(223)stem(abs(G),'.'), title('M agnitude of DFT of g[n]'), axis([0,N,0,max(abs(G))]) subplot(224)stem(abs(H),'.'), title('M agnitude of DFT of h[n]'), axis([0,N,0,max(abs(H))])% Program_2_17clear,close all;x=0:9;M=input('M=');y=circshift(x,M);subplot(221)stem(x,'.'), title('Sequence x[n]'), xlabel('Time index n'), axis([0,9,0,max(x)]) subplot(222)stem(y,'.'), title('Sequence y[n]'), xlabel('Time index n'), axis([0,9,0,max(y)])% Program_2_18clear,close all;g=[1 2 0 1];h=[2 2 1 1];y=circonv(g,h);subplot(221)stem(g,'.'), title('Sequence g[n]'), xlabel('Time index n'), axis([0,4,0,3])subplot(222)stem(h,'.'), title('Sequence h[n]'), xlabel('Time index n'), axis([0,4,0,3])subplot(223)stem(y,'.'), title('Sequence y[n]'), xlabel('Time index n'), axis([0,4,0,max(y)])% Program_2_19clear, close all;g=input('g=');h=input('h=');M=length(h);N=length(g);ge=[g zeros(1,M-1)];he=[h zeros(1,N-1)];Ge=fft(ge,M+N-1);He=fft(he,M+N-1);y=ifft(Ge.*He,M+N-1);subplot(221)stem(g,'.'), title('Sequence g[n]'), xlabel('Time index n'), axis([0,M+N-1,0,max(g)]) subplot(222)stem(h,'.'), title('Sequence h[n]'), xlabel('Time index n'), axis([0,M+N-1,0,max(h)]) subplot(223)stem(y,'.'), title('Sequence y[n]'), xlabel('Time index n'), axis([0,M+N-1,0,max(y)])% Program_2_20clear,close all,N=input('Type in the DFT length N=');n=0:N-1;x=input('Type in the sequence x(n)=');k=0:N-1;X=zeros(1,N);for m=0:N-1;for j=0:N-1;X(m+1)=X(m+1)+x(j+1)*exp(-i*2*pi*m*j/N);endendsubplot(211)stem(n,x,'.')axis([0,N,1.2*min(x),1.2*max(x)])title('The original sequence x[n]')xlabel('Index n')subplot(212)stem(k,abs(X),'r.')axis([0,N,-0.2*max(abs(X)),1.2*max(abs(X))]) title('The DFT of the seqence x[n]')xlabel('Index k')% Program_2_21clear,close all;N=input('Type in the DFT length N=');M=log2(N);n=0:N-1;x=input('Type in the sequence x(n)=');x1=zeros(1,N);x1(1)=x(1);x1(N)=x(N);J=0;for I=1:N-2;k=1;c=0;while k~=M+1&c~=1;if J>=N/(2.^k);J=J-N/(2.^k);c=0;else J=J+N/(2.^k);c=1;endk=k+1;endx1(I+1)=x(J+1);endm = 0:N/2-1;W = exp(-i*2*pi*m/N);for r=1:M;B=2^(r-1);for J=0:B-1;p=2^(M-r)*J;for k=J+1:2^r:N-1;q=x1(k);x1(k)=q+x1(k+B)*W(p+1);x1(k+B)=q-x1(k+B)*W(p+1);endendendk=0:N-1;subplot(211)stem(k,x,'.'),title('The sequence x[n]')xlabel('Index k'),axis([0,max(k),1.1*min(x),1.1*max(x)])subplot(212)stem(k,abs(x1),'r.'),title('The DFT of the sequence x[n]')xlabel('Index k'),axis([0,max(k),-0.2*max(abs(x1)),1.2*max(abs(x1))])实验四%Q4_1clear,close all;n=0:40;wp = input('wp=');ws = input('ws=');Rp = input('Rp=');Rs = input('Rs=');[N,wn] = buttord(wp,ws,Rp,Rs,'s');[b,a]= butter(N,wn,'s');[H,w] = freqs(b,a);subplot(221)plot(w/pi,abs(H));grid on;title('Magnitude response')subplot(222)plot(w/pi,angle(H));grid on;title('Phase response')subplot(223)plot(w/(2*pi),20*log(abs(H)));axis([0,6000,-150,10])grid ontitle('The gain of the butterworth filter')xlabel('Frequency index w/2*pi');subplot(224)h = impz(b,a,12);n = 0:11;stem(n,h/max(abs(h)),'.');grid ontitle('The impulse response--h=h/max(abs(h))')%Q4_3clear,close all;wp = input('wp=');ws = input('ws=');Wp = tan(wp/2);Ws = tan(ws/2);Rp = input('Rp=');Rs = input('Rs=');[N,wn] = buttord(Wp,Ws,Rp,Rs,'s');[b,a]= butter(N,wn,'s');[b,a] = bilinear(b,a,0.5);[H,w] = freqz(b,a);subplot(221)plot(w/pi,abs(H));grid on;title('Magnitude response')subplot(222)plot(w/pi,angle(H));grid on;title('Phase response')subplot(223)plot(w/(pi),20*log(abs(H)));axis([0,1,-750,5])grid ontitle('The gain of the IIR filter')xlabel('Frequency index w/pi')subplot(224)h = impz(b,a,50);n = 0:49;stem(n,h,'.');grid ontitle('The impulse response')%Q4_5clear,close all;x = [-4, -2, 0, -4, -6, -4, -2, -4, -6, -6, -4, -4, -6, -6, -2, 6, 12, 8, 0, -16, -38, -60, -84, -90, -66, -32, -4, -2, -4, 8, 12, 12, 10, 6, 6, 6, 4, 0, 0, 0, 0, 0, -2, -4, 0, 0, 0, -2, -2, 0, 0, -2, -2, -2, -2, 0];dw = 0.01;w = -pi:dw:pi;L = length(x);n = 0:L-1;subplot(222)stem(n,x,'.')title('The original signal')gridwp = 0.2*pi;ws = 0.3*pi;Wp = tan(wp/2);Ws = tan(ws/2);Rp = 1;Rs = 30;[N,wn] = buttord(Wp,Ws,Rp,Rs,'s');[b,a]= butter(N,wn,'s');[b,a] = bilinear(b,a,0.5);H = freqz(b,a,w);subplot(221)plot(w/(pi),20*log(abs(H)));axis([0,0.5,-150,0])grid ontitle('The gain of the IIR filter')xlabel('Frequency index w/pi')subplot(223)h = impz(b,a,L);stem(n,h,'.')title('The impulse response')gridxlabel('Time index n')X = freqz(x,1,w);Y = X.*H;y = dw*Y*exp(j*w'*n)/(2*pi);subplot(224)stem(n,y,'.')title('The output signal of the filter')gridxlabel('Time index n')% Q4_6clear,close all;wp=input('Type in the passband edge frequency (in Radian) wp='); ws=input('T ype in the passband edge frequency (in Radian) ws='); deltaw = ws-wp;wc=(wp+ws)/2;L=ceil(8*pi/(ws-wp));if L/2==ceil(L/2);L=L+1;elseL=L;enda=(L-1)/2;n=0:L-1;w=0.5*(1-cos(2*pi*n/(L-1))).*(u(n)-u(n-L));h=(wc/pi)*sin(wc*(n-a+eps))./(wc*(n-a+eps)).*w;N=2000;H=fft(h,N);H=H/max(abs(H));W=0:2*pi/N:2*pi-2*pi/N;subplot(222)plot(W/pi,abs(H)),title('The magnitude response of the FIR filter') axis([0,1,0,max(abs(H))]),grid on;subplot(223)plot(W/pi,angle(H))axis([0,1,-pi,pi]),grid on;title('The phase response of the FIR filter'),xlabel('Frequency (*pi)')subplot(224)plot(W/pi,20*log10(abs(H)));axis([0,1,-100,0]),grid on;title('The frequency response in dB'),xlabel('Frequency (*pi)')subplot(221)stem(n,h,'.'),axis([0,L-1,min(h),max(h)]),grid on;title('The impulse response of the FIR filter')xlabel('Index n')% Q4_8clear,close all;x = [-4, -2, 0, -4, -6, -4, -2, -4, -6, -6, -4, -4, -6, -6, -2, 6, 12, 8, 0, -16, -38, -60, -84, -90, -66, -32, -4, -2, -4, 8, 12, 12, 10, 6, 6, 6, 4, 0, 0, 0, 0, 0, -2, -4, 0, 0, 0, -2, -2, 0, 0, -2, -2, -2, -2, 0];wp=input('Type in the passband edge frequency (in Radian) wp=');ws=input('T ype in the passband edge frequency (in Radian) ws=');deltaw = ws-wp;wc=(wp+ws)/2;L=ceil(8*pi/(ws-wp));if L/2==ceil(L/2);L=L+1;elseL=L;enda=(L-1)/2;n=0:L-1;w=0.5*(1-cos(2*pi*n/(L-1))).*(u(n)-u(n-L));h=(wc/pi)*sin(wc*(n-a+eps))./(wc*(n-a+eps)).*w;y = conv(x,h);subplot(211)n = 0:length(x)-1;stem(n,x,'.')axis([0,100,-100,50])gridtitle('The original signal')subplot(212)n = 0:length(y)-1;stem(n,y,'.')gridaxis([0,100,-100,50])title('The output signal of the filter')xlabel('Time index n')。
数字信号处理实验报告
实验一 信号、系统及系统响应一、实验目的1、熟悉理想采样的性质，了解信号采样前后的频谱变化，加深对时域采样定理的理解。
2、熟悉离散信号和系统的时域特性。
3、熟悉线性卷积的计算编程方法：利用卷积的方法，观察、分析系统响应的时域特性。
4、掌握序列傅里叶变换的计算机实现方法，利用序列的傅里叶变换对离散信号、系统及其系统响应进行频域分析。
二、 实验原理1.理想采样序列：对信号x a (t)=A e −αt sin(Ω0t )u(t)进行理想采样，可以得到一个理想的采样信号序列x a (t)=A e −αt sin(Ω0nT ),0≤n ≤50，其中A 为幅度因子，α是衰减因子，Ω0是频率，T 是采样周期。
2.对一个连续时间信号x a (t)进行理想采样可以表示为该信号与一个周期冲激脉冲的乘积，即x ̂a (t)= x a (t)M(t),其中x ̂a (t)是连续信号x a (t)的理想采样；M(t)是周期冲激M(t)=∑δ+∞−∞(t-nT)=1T ∑e jm Ωs t +∞−∞,其中T 为采样周期，Ωs =2π/T 是采样角频率。
信号理想采样的傅里叶变换为X ̂a (j Ω)=1T ∑X a +∞−∞[j(Ω−k Ωs )],由此式可知：信号理想采样后的频谱是原信号频谱的周期延拓，其延拓周期为Ωs =2π/T 。
根据时域采样定理，如果原信号是带限信号，且采样频率高于原信号最高频率分量的2倍，则采样以后不会发生频率混叠现象。
三、简明步骤产生理想采样信号序列x a (n),使A=444.128，α=50√2π，Ω0=50√2π。
（1） 首先选用采样频率为1000HZ ，T=1/1000，观察所得理想采样信号的幅频特性，在折叠频率以内和给定的理想幅频特性无明显差异，并做记录；（2） 改变采样频率为300HZ ，T=1/300，观察所得到的频谱特性曲线的变化，并做记录；（3） 进一步减小采样频率为200HZ ，T=1/200，观察频谱混淆现象是否明显存在，说明原因，并记录这时候的幅频特性曲线。
数字信号处理实验指导书(M)
数字信号处理实验电子信息科学与技术实验室2007年7月目录实验一离散时间信号的时域表示 (3)实验二离散信号的卷积和 (6)实验三离散傅立叶变换及其特性验证 (8)实验四信号处理中FFT的应用 (11)实验五离散系统的Z域分析 (15)实验六无限冲激响应(IIR)数字滤波器的三种结构 (19)实验七冲激响应不变法IIR数字滤波器设计 (23)实验八双线性变换法IIR数字滤波器设计 (26)实验一 离散时间信号的时域表示一、实验目的1、熟悉Matlab 命令，掌握离散时间信号－序列的时域表示方法。
2、掌握用Matlab 描绘二维图像的方法。
3、掌握用Matlab 对序列进行基本的运算和时域变换的方法。
二、实验原理与计算方法（一）序列的表示方法 序列的表示方法有列举法、解析法和图形法，相应的用Matlab 也可以有这样几种表示方法，分别介绍如下：1、列举法 在Matlab 中，用一个列向量来表示一个有限长序列，由于一个列向量并不包含位置信息，因此需要用表示位置的n 和表示量值的x 两个向量来表示任意一个序列，如：例1.1：>>n=[-3,-2,-1,0,1,2,3,4]; >>x=[2,1,-1,0,1,4,3,7];如果不对向量的位置进行定义，则Matlab 默认该序列的起始位置为n=0。
由于内存有限，Matlab 不能表示一个无限序列。
2、解析法对于有解析表达式的确定信号，首先定义序列的范围即n 的值，然后直接写出该序列的表达式，如：例1.2：实现实指数序列nn x )9.0()(=，100≤≤n 的Matlab 程序为：>>n=[0:10]; >>x=(0.9).^n;例 1.3：实现正余弦序列)5.0sin(2)31.0cos(3)(n n n x πππ++=，155≤≤n 的Matlab 程序为：>>n=[5:15];>>x=3*cos(0.1*pi*n+pi/3)+2*sin(0.5*pi*n); 3、图形法在Matlab 中用图形法表示一个序列，是在前两种表示方法的基础上将序列的各个量值描绘出来，即首先对序列进行定义，然后用相应的画图语句画图，如：例1.4：绘制在例1.1中用列举法表示的序列的图形，则在向量定义之后加如下相应的绘图语句：>>stem(n,x);此时得到的图形的横坐标范围由向量n 的值决定，为-3到4，纵坐标的范围由向量x 的值决定，为-1到7。
数字信号处理实验全部程序MATLAB
实验一熟悉MATLAB环境一、实验目的(1)熟悉MATLAB的主要操作命令。
(2)学会简单的矩阵输入和数据读写。
(3)掌握简单的绘图命令。
(4)用MATLAB编程并学会创建函数。
(5)观察离散系统的频率响应。
二、实验内容认真阅读本章附录，在MATLAB环境下重新做一遍附录中的例子，体会各条命令的含义。
在熟悉了MATLAB基本命令的基础上，完成以下实验。
上机实验内容：(1)数组的加、减、乘、除和乘方运算。
输入A=[1 2 3 4]，B=[3 4 5 6]，求C=A+B，D=A-B，E=A.*B，F=A./B，G=A.^B并用stem语句画出A、B、C、D、E、F、G。
实验程序：A=[1 2 3 4];B=[3 4 5 6];n=1:4;C=A+B;D=A-B;E=A.*B;F=A./B;G=A.^B;subplot(4,2,1);stem(n,A,'fill');xlabel ('时间序列n');ylabel('A');subplot(4,2,2);stem(n,B,'fill');xlabel ('时间序列n ');ylabel('B');subplot(4,2,3);stem(n,C,'fill');xlabel ('时间序列n ');ylabel('A+B');subplot(4,2,4);stem(n,D,'fill');xlabel ('时间序列n ');ylabel('A-B');subplot(4,2,5);stem(n,E,'fill');xlabel ('时间序列n ');ylabel('A.*B');subplot(4,2,6);stem(n,F,'fill');xlabel ('时间序列n ');ylabel('A./B');subplot(4,2,7);stem(n,G,'fill');xlabel ('时间序列n ');ylabel('A.^B');运行结果：(2)用MATLAB实现以下序列。
数字信号处理实验报告_五个实验
实验一 信号、系统及系统响应一、 实验目的1、熟悉连续信号经理想采样前后的频谱变化关系，加深对时域采样定理的理解；2、熟悉时域离散系统的时域特性；3、利用卷积方法观察分析系统的时域特性；4、掌握序列傅立叶变换的计算机实现方法，利用序列的傅立叶变换对连续信号、离散信号及系统响应进行频域分析。
二、 实验原理及方法采样是连续信号数字处理的第一个关键环节。
对采样过程的研究不仅可以了解采样前后信号时域和频域特性发生变化以及信号信息不丢失的条件，而且可以加深对傅立叶变换、Z 变换和序列傅立叶变换之间关系式的理解。
对一个连续信号)(t x a 进行理想采样的过程可用下式表示：)()()(^t p t t xx aa=其中)(^t x a 为)(t x a 的理想采样，p(t)为周期脉冲，即∑∞-∞=-=m nT t t p )()(δ)(^t x a的傅立叶变换为)]([1)(^s m a m j X T j a XΩ-Ω=Ω∑∞-∞=上式表明^)(Ωj Xa为)(Ωj Xa的周期延拓。
其延拓周期为采样角频率（T /2π=Ω）。
只有满足采样定理时，才不会发生频率混叠失真。
在实验时可以用序列的傅立叶变换来计算^)(Ωj X a 。
公式如下：Tw jw ae X j X Ω==Ω|)()(^离散信号和系统在时域均可用序列来表示。
为了在实验中观察分析各种序列的频域特性，通常对)(jw e X 在[0，2π]上进行M 点采样来观察分析。
对长度为N 的有限长序列x(n)，有：n jw N n jw k ke m x eX--=∑=)()(1其中，k Mk πω2=，k=0,1,……M-1 时域离散线性非移变系统的输入/输出关系为 ∑∞-∞=-==m m n h m x n h n x n y )()()(*)()(上述卷积运算也可在频域实现)()()(ωωωj j j e H e X eY =三、 实验程序s=yesinput(Please Select The Step Of Experiment:\n 一.(1时域采样序列分析 s=str2num(s); close all;Xb=impseq(0,0,1); Ha=stepseq(1,1,10);Hb=impseq(0,0,3)+2.5*impseq(1,0,3)+2.2*impseq(2,0,3)+impseq(3,0,3); i=0;while(s);%时域采样序列分析 if(s==1) l=1; k=0;while(1)if(k==0)A=yesinput('please input the Amplitude:\n',...444.128,[100,1000]); a=yesinput('please input the Attenuation Coefficient:\n',...222.144,[100,600]); w=yesinput('please input the Angle Frequence(rad/s):\n',...222.144,[100,600]); end k=k+1;fs=yesinput('please input the sample frequence:\n',...1000,[100,1200]); Xa=FF(A,a,w,fs); i=i+1;string+['fs=',num2str(fs)]; figure(i)DFT(Xa,50,string); 1=yesinput 1=str2num(1); end%系统和响应分析else if(s==2)kk=str2num(kk);while(kk)if(kk==1)m=conv(Xb,Hb);N=5;i=i+1;figure(i)string=('hb(n)');Hs=DFT(Hb,4,string);i=i+1;figure(i)string('xb(n)');DFT(Xb,2,string);string=('y(n)=xb(n)*hb(n)');else if (kk==2)m=conv(Ha,Ha);N=19;string=('y(n)=ha(n)*(ha(n)');else if (kk==3)Xc=stepseq(1,1,5);m=conv(Xc,Ha);N=14;string=('y(n)=xc(n)*ha(n)');endendendi=i+1;figure(i)DFT(m,N,string);kk=yesinputkk=str2num(kk);end卷积定理的验证else if(s==3)A=1;a=0.5;w=2,0734;fs=1;Xal=FF(A,a,w,fs);i=i+1;figure(i)string=('The xal(n)(A=1,a=0.4,T=1)'); [Xa,w]DFT(Xal,50,string);i=i+1;figure(i)string =('hb(n)');Hs=DFT(Hb,4,string);Ys=Xs.*Hs;y=conv(Xal,Hb);N=53;i=i+1;figure(i)string=('y(n)=xa(n)*hb(n)');[yy,w]=DFT(y,N,string);i=i+1;figure(i)subplot(2,2,1)plot(w/pi,abs(yy));axis([-2 2 0 2]);xlabel('w/pi');ylabel('|Ys(jw)|');title(FT[x(n)*h(n)]');subplot(2,2,3)plot(w/pi,abs(Ys));axis([-2 2 0 2]);xlabel('w/pi');ylabel('|Ys(jw)|');title('FT[xs(n)].FT[h(n)]');endendend子函数:离散傅立叶变换及X(n),FT[x(n)]的绘图函数function[c,l]=DFT(x,N,str)n=0:N-1;k=-200:200;w=(pi/100)*k;l=w;c=x*Xc=stepseq(1,1,5);子函数：产生信号function c=FF(A,a,w,fs)n=o:50-1;c=A*exp((-a)*n/fs).*sin(w*n/fs).*stepseq(0,0,49); 子函数：产生脉冲信号function [x,n]=impseq(n0,n1,n2)n=[n1:n2];x=[(n-n0)==0];子函数：产生矩形框信号function [x,n]=stepseq(n0,n1,n2) n=[n1:n2];x=[(n-n0>=0)];四、 实验内容及步骤1、认真复习采样理论，离散信号与系统，线性卷积，序列的傅立叶变换及性质等有关内容，阅读本实验原理与方法。
数字信号处理上机实验及参考程序
数字信号处理上机实验及参考程序数字信号处理实验实验⼀离散时间信号与系统及MA TLAB实现1.单位冲激信号：n = -5:5;x = (n==0);subplot(122);stem(n, x);2.单位阶跃信号：x=zeros(1,11);n0=0;n1=-5;n2=5;n = n1:n2;x(:,n+6) = ((n-n0)>=0);stem(n,x);3.正弦序列：n = 0:1/3200:1/100;x=3*sin(200*pi*n+1.2);stem(n,x);4.指数序列n = 0:1/2:10;x1= 3*(0.7.^n);x2=3*exp((0.7+j*314)*n);subplot(221);stem(n,x1);subplot(222);stem(n,x2);5.信号延迟n=0:20;Y1=sin(100*n);Y2=sin(100*(n-3));subplot(221);stem(n,Y1);subplot(222);stem(n,Y2);6.信号相加X1=[2 0.5 0.9 1 0 0 0 0];X2=[0 0.1 0.2 0.3 0.4 0.5 0.6 0.7];X=X1+X2;stem(X);7.信号翻转X1=[2 0.5 0.9 1];n=1:4;X2=X1(5-n);subplot(221);stem(n,X1);subplot(222);stem(n,X2);8.⽤MATLAB计算序列{-2 0 1 –1 3}和序列{1 2 0 -1}的离散卷积。
a=[-2 0 1 -1 3]; b=[1 2 0 -1];c=conv(a,b);M=length(c)-1;n=0:1:M;stem(n,c);xlabel('n');ylabel('幅度');9．⽤MA TLAB计算差分⽅程当输⼊序列为时的输出结果。
N=41;a=[0.8 -0.44 0.36 0.22];b=[1 0.7 -0.45 -0.6];x=[1 zeros(1,N-1)];k=0:1:N-1;y=filter(a,b,x);stem(k,y)xlabel('n');ylabel('幅度')10.冲激响应impzN=64;a=[0.8 -0.44 0.36 0.22];b=[1 0.7 -0.45 -0.6];x=[1 zeros(1,N-1)];k=0:1:N-1;y=impz(a,b,N);stem(k,y)xlabel('n');ylabel('幅度')11.传递函数频率响应a=[0.8 -0.44 0.36 0.22];%分⼦的系数数组b=[1 0.7 -0.45 -0.6];%分母的系数数组n=(0:500)*pi/500%在pi范围内取501个采样点[h,f]=freqz(a,b,n);%求系统的频率响应subplot(2,1,1),plot(n/pi,abs(h));grid%作系统的幅度频响图axis([0,1,1.1*min(abs(h)),1.1*max(abs(h))]);ylabel('幅度');subplot(2,1,2),plot(n/pi,angle(h));grid %作系统的相位频响图axis([0,1,1.1*min(angle(h)),1.1*max(angle(h))]);ylabel('相位');xlabel('以pi为单位的频率');12.系统零极点图a=[0.8 -0.44 0.36 0.22];b=[1 0.7 -0.45 -0.6];h=zplane(a,b);实验⼆离散信号变换1.解⽅程y(n)-2y(n-1)+3y(n-2)=4x(n)-5x(n-1)+6x(n-2)-7x(n-3)a=[4,-5,6,-7];b=[1,-2,3];n=[0:7]; x=ones(length(n));Y=[-1,1];X=[1,-1];xic=filtic(b,a,Y,X);y1=filter(b,a,x,xic)stem(n,y1);xlabel('n');ylabel('y(n)');2．对连续的单⼀频率周期信号按采样频率采样，截取长度N分别选N =20和N =16，观察其DFT结果的幅度谱。
数字信号处理实验指导书
《数字信号处理》实验指导书信息与机电工程学院实验中心2017-11-20实验一 常见离散信号的MATLAB 产生和图形显示一、实验目的：加深对常用离散信号的理解； 二、实验原理：1、基础知识：R1.1 单位样本序列10[]0n n n δ=⎧=⎨≠⎩如果()n δ在时间轴上延迟了k 个单位，得到()n k δ-，即：1[]0n k n k n kδ=⎧-=⎨≠⎩R1.2 单位阶跃序列10[]0n u n n ≥⎧=⎨<⎩ R1.3 指数序列[]n x n A α=，其中()00j e σωα+=，j A A e φ=，则前式化为()000000[]cos()sin()n j n n n x n A eA e n j A e n σωφσσωφωφ++==+++R1.4 正弦序列0[]cos()x n A n ωφ=+，其中A ，0ω，φ是实数，分别称为正弦序列的振幅、角频率和初始相位。
00/2f ωπ=称为频率。
2、用到的MATLAB 命令 运算符和特殊符号 ： . + -* / .^ ; %基本矩阵和矩阵控制 i ones pirand randnzeros基本函数 cos sin exp imag real二维图形 axis gird legendplotstem title xlabel ylabelstairs 通用图形函数 clf subplot三、实验内容及要求：编制程序产生信号，并绘出其图形。
例1.1单位样本和单位阶跃序列% 程序 P1.1% 一个单位样本序列的产生clf;% 产生一个从－10到20的向量n = -10:20;% 产生单位样本序列u = [zeros(1,10) 1 zeros(1,20)];% 绘制单位样本序列stem(n,u);xlabel('时间序号 n');ylabel('振幅');title('单位样本序列');axis([-10 20 0 1.2]);习题：Q1.1 运行程序P1.1，以产生单位样本序列u[n]并记录它。
数字信号处理实验一
数字信号处理实验一实验目的：掌握利用Matlab产生各种离散时间信号，实现信号的相加、相乘及卷积运算实验函数：参考课本77-19页，注意式(2.11.1)的表达与各matlab子函数间的关系。
1、stem(x,y) % 绘制以x为横轴，y为纵轴的离散序列图形2、[h ,t] = impz(b, a) % 求解数字系统的冲激响应h，取样点数为缺省值[h, t] = impz(b, a, n) % 求解数字系统的冲激响应h，取样点数为nimpz(b, a) % 在当前窗口用stem(t, h)函数出图3、[h ,t] = dstep(b, a) % 求解数字系统的阶跃响应h，取样点数为缺省值[h, t] = dstep (b, a, n) % 求解数字系统的阶跃响应h，取样点数为ndstep (b, a) % 在当前窗口用stairs(t, h)函数出图4、y = filter(b,a,x) % 在已知系统差分方程或转移函数的情况下求系统输出实验原理：一、常用的时域离散信号及其程序1、产生单位抽样函数δ(n)n1 = -5;n2 = 5;n0 = 0;n = n1:n2;x = [n==n0]; % x在n=n0时为1，其余为0stem(n,x,'filled'); %filled:序列圆心处用实心圆表示axis([n1,n2,0,1.1*max(x)])title('单位抽样序列')xlabel('time(n)')ylabel('Amplitude:x(n)')2、产生单位阶跃序列u(n)n1 = -2;n2 = 8;n0 = 0;n = n1:n2;x = [n>=n0]; % x在n>=n0时为1，其余为0stem(n,x,'filled');axis([n1,n2,0,1.1*max(x)])title('单位阶跃序列')xlabel('time(n)')ylabel('Amplitude:x(n)')3、复指数序列复指数序列的表示式为()(),00,0j n e n x n n σω+⎧≥⎪=⎨<⎪⎩，当0ω=时，()x n 为实指数序列；当0σ=时，()x n 为虚指数序列，即()()cos sin j n e n j n ωωω=+，即其实部为余弦序列，虚部为正弦序列。
数字信号处理实验报告MATLAB
数字信号处理实验报告姓名：班级：09电信一班学号：2)]得下图二，图二图一3．将如下文件另存为：sigadd.m文件function [y,n] = sigadd(x1,n1,x2,n2)% 实现y(n) = x1(n)+x2(n)% -----------------------------% [y,n] = sigadd(x1,n1,x2,n2)% y = 在包含n1 和n2 的n点上求序列和,% x1 = 在n1上的第一序列% x2 = 在n2上的第二序列(n2可与n1不等)n = min(min(n1),min(n2)):max(max(n1),max(n2)); % y(n)的长度y1 = zeros(1,length(n)); y2 = y1; % 初始化y1(find((n>=min(n1))&(n<=max(n1))==1))=x1; % 具有y的长度的x1y2(find((n>=min(n2))&(n<=max(n2))==1))=x2; % 具有y的长度的x2y = y1+y2;在命令窗口输入：x1=[1,0.5,0.3,0.4];n1=-1:2;x2=[0.2,0.3,0.4,0.5,0.8,1];n2=-2:3; [y,n] = sigadd(x1,n1,x2,n2)得：y =n=-1:10;x=sin(0.4*pi*n);y=fliplr(x);n1=-fliplr(n);subplot(2,1,1),stem(n,x) subplot(2,1,2),stem(n1,y在命令窗口键入:n=-1:10; x=sin(0.4*pi*n);n (samples)实验结果：1.（1）在命令窗口输入：tic; [am,pha]=dft1(x)N=length(x);w=exp(-j*2*pi/N);for k=1:Nsum=0;for n=1:Nsum=sum+x(n)*w^((k-1)*(n-1));endam(k)=abs(sum);pha(k)=angle(sum);end;toc得到如下结果：am =Columns 1 through 11120.0000 41.0066 20.9050 14.3996 11.3137 9.6215 8.6591 8.1567 8.0000 8.1567 8.6591Columns 12 through 169.6215 11.3137 14.3996 20.9050 41.0066pha =Columns 1 through 110 1.7671 1.9635 2.1598 2.3562 2.5525 2.7489 2.9452 3.1416 -2.9452 -2.7489Columns 12 through 16-2.5525 -2.3562 -2.1598 -1.9635 -1.7671Elapsed time is 0.047000 seconds.（2）在命令窗口输入：tic;[am,pha]=dft2(x)N=length(x);n=[0:N-1];k=[0:N-1];w=exp(-j*2*pi/N);nk=n’*k;wnk=w.^(nk); Xk=x*wnk; am= abs(Xk); pha=angle(Xk); toc得到下图：figure(1)00.10.20.30.40.50.60.70.80.91signal x(n), 0 <= n <= 99（2）在命令窗口键入：n3=[0:1:99];y3=[x(1:1:10) zeros(1,90)]; ％添90个零。
《数字信号处理》实验及代码
实验一 用Matlab 进行数字信号处理一、实验目的1．掌握 MATLAB 基本的操作；2．学习典型的离散信号的Matlab 实现方法。
3．学习离散时间序列的基本运算：相加、相乘、移位等； 二、实验内容1. 练习把y=sin(x), z=cos(x),u=2sin(x),v=sin(x)/cos(x)在[0，2pi]区间内的4个子图分别用不同的颜色、点型和线型绘制在同一个窗口中，并加上纵坐标、标题、图例和网格线。
2. 利用MATLAB 编程产生和绘制下列有限长序列： (1)单位脉冲序列()n δ (2)单位阶跃序列()u n (3)矩形序列8()R n(4)正弦序列()sin()53x n n ππ=+(5)实指数序列0.9()n u n3．上机调试并打印或记录实验结果。
4．完成实验报告。
三、实验结果1. 实验程序如下：x=0:pi/10:2*pi;y=sin(x); plot(x,y,'r*:');grid on %绘制网格线 hold on z=cos(x); plot(x,z,'y+-'); hold on u=2*sin(x); plot(x,u,'bx-'); hold on v=sin(x)/cos(x); plot(x,v,'ko-');hold onxlabel('x 轴'); %x轴注释ylabel('y 轴'); %y轴注释legend({'y=sin(x)','z=cos(x)', 'u=2sin(x) ', 'v=sin(x)/cos(x) '}); %图形注解2. 程序如下：n=-20:20;n0=0;n1=8;w0=pi/5; w1=pi/3;x1=[(n-n0)==0];x2=[(n-n0)>=0];x3=[(n-n0)>=0& (n-n1)<=0];x4=sin(w0*n+w1);x5=0.9.^n.*x2;subplot(511);stem(n,x1);axis([ -20 20 0 2]);ylabel('\sigma(n)'); subplot(512);stem(n,x2);axis([ -20 20 0 2]);ylabel('u(n)');subplot(513);stem(n,x3);axis([ -20 20 0 2]);ylabel('R8(n)'); subplot(514);stem(n,x4);axis([ -20 20 -2 2]);ylabel(' sin(w0n+w1) '); subplot(515);stem(n,x5);axis([ -20 20 0 2]);ylabel('0.9nu(n)'); xlabel('n');实验二离散信号与系统一、实验目的1．掌握卷积定理、熟悉离散信号和系统的时域特性；2．学习Matlab进行卷积计算方法；3．学习Matlab求解差分方程二、实验内容1．离散信号的基本运算：对序列x(n)={2,3,4,1,2,5} ,n=0,1,2,3,4,5,的移位、乘法、加法、翻转及尺度变换。
数字信号处理实验程序清单
实验一: 系统响应及系统稳定性%实验1：系统响应及系统稳定性close all;clear all%======内容1：调用filter解差分方程，由系统对u(n)的响应判断稳定性==A=[1,-0.9];B=[0.05,0.05]; %系统差分方程系数向量B和Ax1n=[1 1 1 1 1 1 1 1 zeros(1,50)]; %产生信号x1(n)=R8(n)x2n=ones(1,128); %产生信号x2(n)=u(n)hn=impz(B,A,58); %求系统单位脉冲响应h(n)subplot(2,2,1);y='h(n)';tstem(hn,y); %调用函数tstem绘图title('(a) 系统单位脉冲响应h(n)');box ony1n=filter(B,A,x1n); %求系统对x1(n)的响应y1(n)subplot(2,2,2);y='y1(n)';tstem(y1n,y);title('(b) 系统对R8(n)的响应y1(n)');box ony2n=filter(B,A,x2n); %求系统对x2(n)的响应y2(n)subplot(2,2,4);y='y2(n)';tstem(y2n,y);title('(c) 系统对u(n)的响应y2(n)');box on%===内容2：调用conv函数计算卷积========x1n=[1 1 1 1 1 1 1 1 ]; %产生信号x1(n)=R8(n)h1n=[ones(1,10) zeros(1,10)];h2n=[1 2.5 2.5 1 zeros(1,10)];y21n=conv(h1n,x1n);y22n=conv(h2n,x1n);figure(2)subplot(2,2,1);y='h1(n)';tstem(h1n,y); %调用函数tstem绘图title('(d) 系统单位脉冲响应h1(n)');box onsubplot(2,2,2);y='y21(n)';tstem(y21n,y);title('(e) h1(n)与R8(n)的卷积y21(n)');box onsubplot(2,2,3);y='h2(n)';tstem(h2n,y); %调用函数tstem绘图title('(f) 系统单位脉冲响应h2(n)');box onsubplot(2,2,4);y='y22(n)';tstem(y22n,y);title('(g) h2(n)与R8(n)的卷积y22(n)');box on%=========内容3：谐振器分析================un=ones(1,256); %产生信号u(n)n=0:255;xsin=sin(0.014*n)+sin(0.4*n); %产生正弦信号A=[1,-1.8237,0.9801];B=[1/100.49,0,-1/100.49]; %系统差分方程系数向量B和Ay31n=filter(B,A,un); %谐振器对u(n)的响应y31(n)y32n=filter(B,A,xsin); %谐振器对u(n)的响应y31(n)figure(3)subplot(2,1,1);y='y31(n)';tstem(y31n,y);title('(h) 谐振器对u(n)的响应y31(n)');box onsubplot(2,1,2);y='y32(n)';tstem(y32n,y);title('(i) 谐振器对正弦信号的响应y32(n)');box ontstem 程序清单function tstem(xn,yn)%时域序列绘图函数% xn:信号数据序列，yn:绘图信号的纵坐标名称（字符串）n=0:length(xn)-1;stem(n,xn,'.');box onxlabel('n');ylabel(yn);axis([0,n(end),min(xn),1.2*max(xn)])实验二时域采样与频域采样1 时域采样理论的验证程序清单% 时域采样理论验证程序exp2a.mTp=64/1000; %观察时间Tp=64微秒%产生M长采样序列x(n)% Fs=1000;T=1/Fs;Fs=1000;T=1/Fs;M=Tp*Fs;n=0:M-1;A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5;xnt=A*exp(-alph*n*T).*sin(omega*n*T);Xk=T*fft(xnt,M); %M点FFT[xnt)]yn='xa(nT)';subplot(3,2,1);tstem(xnt,yn); %调用自编绘图函数tstem绘制序列图box on;title('(a) Fs=1000Hz');k=0:M-1;fk=k/Tp;subplot(3,2,2);plot(fk,abs(Xk));title('(a) T*FT[xa(nT)],Fs=1000Hz');xlabel('f(Hz)');ylabel('幅度');axis([0,Fs,0,1.2*max(abs(Xk))])%========================================% Fs=300Hz和Fs=200Hz的程序与上面Fs=1000Hz完全相同。
数字信号处理实验文档
实验一：离散信号的MATLAB实现一、实验目的：1、掌握离散时间信号的一般表示方法。
2、熟悉连续信号经理想采样后的频谱变化关系，加深对时域采样定理的理解。
3、掌握离散信号序列的操作。
二、实验内容：M1-1 已知g1(t)=cos(6*pi*t), g1(t)=co 14*pi*t), g1(t)=cos(26*pi*t),以抽样频率fsam=10Hz对上述三个信号进行抽样。
在同一张图上画出g1(t),g2(t)和g3(t)及其抽样点，对所得结果进行讨论。
解：代码如下：100:100)*1/100;g1t=cos(6*pi*t);g2t=cos(14*pi*t);g3t=cos(26*pi*t);subplot(3,1,1);plot(g1t);subplot(3,1,2);plot(g2t);subplot(3,1,3);plot(g3t);绘出的图形如图1_1所示:图1_1采样频率为fsam=10Hz，采样时间为0.1s，而f1=3Hz，f2=7Hz，f3=13Hz，使得三个信号的采样图形相似，这样不能很好还原原来的信号图像。
所以对信号的采样频率应足够大，应满足fsam>=2fm.M1-2利用MATLAB的filter函数，求出下列系统的单位脉冲响应，并判断系统是否稳定。
讨论题所获得的结果。
代码1：k=1:300;x=zeros(1,300);x(1)=1;b1=[1];a1=[1,-1.845,0.850586];y1=filter(b1,a1,x);subplot(1,2,1);plot(k,y1);xlabel('k');ylabel('幅度y1');b2=[1];a2=[1,-1.85,0.85];y2=filter(b2,a2,x);subplot(1,2,2);plot(k,y2);xlabel('k');ylabel('幅度y2');图1_2_1代码2：x=zeros(1,500);x(1)=1;b1=[1];a1=[1,-1.845,0.850586];y1=filter(b1,a1,x);plot(k,y1);b2=[1];a2=[1,-1.85,0.85];y2=filter(b2,a2,x);plot(k,y1,k,y2,':');xlabel('k');ylabel('幅度');legend('y1''y2');图1_2_2结论：H1（z）的两个极点都在单位圆内，所以系统稳定，从图中可以看出响应曲线升高后有回落，系统最终趋向于0；H2（z）的一个极点在单位圆内，另一个在单位圆上，所以系统最终临界稳定，从图中可以看出响应曲线上升后没有回落，系统最终趋向于6.7左右。
数字信号处理实验指导书(带源程序)
实验一离散时间系统与MA TLAB一. 实验目的1. 进一步加深对离散时间系统的理解。
2. 学习在MATLAB中怎样表示离散时间信号。
3. 熟悉离散时间信号的作图。
二. 实验步骤1. 复习离散时间系统的有关容。
2. 复习MA TLAB的基本语法。
3. 按实验容熟悉stem。
4. 编写程序。
5. 输出结果，总结结论，按要求写出实验报告。
三. 实验容1．掌握stem函数STEM(Y) plots the data sequence Y as stems from the x axis terminated with circles for the data value.STEM(X,Y) plots the data sequence Y at the values specified in X.例：t=[0:0.1:2]; x=cos(pi*t+0.6); stem(t,x);xn=[4,2,2,3,6,7]; stem(xn);思考：STEM(Y)与STEM(X,Y)有什么不同？STEM与PLOT函数有什么不同？2．掌握subplot函数H = SUBPLOT(m,n,p), or SUBPLOT(mnp), breaks the Figure window into an m-by-n matrix of small axes, selects the p-th axes for the current plot, and returns the axis handle. The axes are counted along the top row of the Figure window, then the second row, etc.例：n1=0:3;x1=[1,1,1,1];subplot(221);stem(n1,x1);title('x1序列');n2=0:7;x2=[1,2,3,4,4,3,2,1];subplot(222);stem(n2,x2);title('x2序列');n3=0:7;x3=[4,3,2,1,1,2,3,4];subplot(223);stem(n3,x3);title('x3序列');n4=0:7;x41=cos((pi/4)*n4);subplot(224);stem(n4,x41);title('x4序列');思考：subplot是怎样分配各个作图分区的顺序号的？3．信号的运算]0,1.0,4.0,7.0,1[)(1=n x ，]9.0,7.0,5.0,3.0,1.0[)(2=n x ，请作出)()(21n x n x +，)()(21n x n x 的图形。
《数字信号处理》实验指导书(正文)
实验一 离散时间信号分析一、实验目的1．掌握各种常用的序列，理解其数学表达式和波形表示。
2．掌握在计算机中生成及绘制数字信号波形的方法。
3．掌握序列的相加、相乘、移位、反褶等基本运算及计算机实现与作用。
4．掌握线性卷积软件实现的方法。
5．掌握计算机的使用方法和常用系统软件及应用软件的使用。
6．通过编程，上机调试程序，进一步增强使用计算机解决问题的能力。
二、实验原理1．序列的基本概念离散时间信号在数学上可用时间序列来表示，其中代表序列的第n 个数字，n 代表时间的序列，n 的取值范围为∞<<∞-n 的整数，n 取其它值)(n x 没有意义。
离散时间信号可以是由模拟信号通过采样得到，例如对)(t x a 模拟信号进行等间隔采样，采样间隔为T ，得到一个{})(nT x a 有序的数字序列就是离散时间信号，简称序列。
2．常用序列常用序列有：单位脉冲序列（单位采样）)(n δ、单位阶跃序列)(n u 、矩形序列)(n R N 、实指数序列、复指数序列、正弦型序列等。
3．序列的基本运算序列的运算包括移位、反褶、和、积、标乘、累加、差分运算等。
4．序列的卷积运算∑∞∞-*=-=)()()()()(n h n x m n h m x n y上式的运算关系称为卷积运算，式中代表两个序列卷积运算。
两个序列的卷积是一个序列与另一个序列反褶后逐次移位乘积之和，故称为离散卷积，也称两序列的线性卷积。
其计算的过程包括以下4个步骤。
(1)反褶：先将)(n x 和)(n h 的变量n 换成m ，变成)(m x 和)(m h ，再将)(m h 以纵轴为对称轴反褶成)(m h -。
(2)移位：将)(m h -移位n ，得)(m n h -。
当n 为正数时，右移n 位；当n 为负数时，左移n 位。
(3)相乘：将)(m n h -和)(m x 的对应点值相乘。
(4)求和：将以上所有对应点的乘积累加起来，即得)(n y 。
数字信号处理实验指导书
《数字信号处理》实验指导书编写：刘梦亭审核：司玉娟阎维和适用专业：电子信息工程电子信息科学与技术通信工程等电子信息与工程系2009年9月目录实验一：离散时间信号分析 (1)实验二：离散时间系统分析 (3)实验三：离散系统的Z域分析 (6)实验四：FFT频谱分析及应用 (9)实验五：IIR数字滤波器的设计 (12)实验六：FIR数字滤波器的设计 (16)附录： MATLAB基本操作及常用命令 (20)实验一：离散时间信号分析实验学时：2学时 实验类型：验证 实验要求：必修 一、实验目的1） 掌握离散卷积计算方法； 2） 学会差分方程的迭代解法；3） 了解全响应、零输入响应、零状态响应和初始状态的物理意义和具体求法； 二、实验内容 1、信号的加数学描述 )()()(21n x n x n x += MATLAB 实现 21X X X +=设[ x10=[1 0.7 0.4 0.1 0]; x20=[0.1 0.3 0.5 0.7 0.9 1];]2、信号的乘数学描述 )()()(21n x n x n x *= MATLAB 实现 2.1X X X *=设[ x10=[1 0.7 0.4 0.1 0]; x20=[0.1 0.3 0.5 0.7 0.9 1];]3、计算卷积用MATLAB 计算序列{-2 0 1 –1 3}和序列{1 2 0 -1}的离散卷积。
首先用手工计算，然后用MATLAB 编程验证。
三、实验组织运行要求1、学生在进行实验前必须进行充分的预习，熟悉实验内容；2、学生根据实验要求，读懂并理解相应的程序；3、学生严格遵守实验室的各项规章制度，注意人身和设备安全，配合和服从实验室人员管理；4、教师在学生实验过程中予以必要的辅导，独立完成实验；5、采用集中授课形式。
四、实验条件1、具有WINDOWS 98/2000/NT/XP 操作系统的计算机一台； 2.、MATLAB 编程软件。
《数字信号处理》实验讲义(信息计算)
《数字信号处理》实验指导书实验一 常见离散信号的产生一、实验目的1. 加深对离散信号的理解。
2. 掌握典型离散信号的Matlab 产生和显示。
二、实验原理及方法在MATLAB 中，序列是用矩阵向量表示，但它没有包含采样信息，即序列位置信息，为此，要表示一个序列需要建立两个向量；一是时间序列n,或称位置序列，另一个为取值序列x ，表示如下： n=[…,-3,-2,-1,0,1,2,3,…]x=[…,6,3,5,2,1,7,9,…]一般程序都从0 位置起始，则x= [x(0), x(1), x(2),…]对于多维信号需要建立矩阵来表示，矩阵的每个列向量代表一维信号。
数字信号处理中常用的信号有指数信号、正弦信号、余弦信号、方波信号、锯齿波信号等，在MATLAB 语言中分别由exp, sin, cos, square, sawtooth 等函数来实现。
三、实验内容1. 用MATLAB 编制程序，分别产生长度为N(由输入确定)的序列：①单位冲击响应序列：()n δ可用MATLAB 中zeros 函数来实现； ②单位阶跃序列：u(n)可用MATLAB 中ones 函数来实现； ③正弦序列：()sin()x n n ω=； ④指数序列：(),nx n a n =-∞<<+∞⑤复指数序列：用exp 函数实现()0()a jb nx n K e+=，并给出该复指数序列的实部、虚部、幅值和相位的图形。
(其中00.2,0.5,4,40a b K N =-===.)参考流程图：四、实验报告要求1. 写出实验程序，绘出单位阶跃序列、单位阶跃序列、正弦序列、指数序列的图形以及绘 出复指数序列的实部、虚部、幅值和相位的图形。
2. 序列信号的实现方法。
3. 在计算机上实现正弦序列0()sin(2)x n A fn πϕ=+。
实验二 离散信号的运算一、实验目的1. 掌握离散信号的时域特性。
2. 用MATLAB 实现离散信号的各种运算。
数字信号处理实验指导书(wcx)
1．线性和非线性系统
例2－1设系统为
y[n]-0.4y[n-1]+0.75y[n-2]=2.2403x[n]+2.4908x[n-1]+2.2403x[n-2]
要求用MATLAB程序仿真系统，输入三个不同的输入序列x1(n),x2(n)和
x(n)=a.x1(n)+b.x2(n)，计算并求出相应的输出响应y1[n],y2[n]和y[n]。
数字信号处理应用的一个常见例子是从被加性噪声污染的信号中移除噪声。假定信号s[n]被噪声d[n]所污染，得到一个含有噪声的信号x[n]=s[n]+d[n]。我们需要对x[n]进行运算，产生一个合理的逼近s[n]，对时刻n的样本求平均，产生输出信号是一种简单有效的方法。如：三点滑动平均的信号。
程序1－3实现三点滑动平均的信号运算：
(2)程序1－2：正弦序列的产生和绘制
% Program P1_2
% Generation of a sinusoidal sequence
n = 0:40;
f = 0.1;
phase = 0;
A = 1.5;
arg = 2*pi*f*n - phase;
x = A*cos(arg);
clf;% Clear old graph
由固冇频率wn把模拟低通滤波器原型转换为低通高通带通带阻滤运用脉冲响应不变法或双线性变换法把模拟滤波器转换成数字滤波器matlab信号处理工具箱提供了儿个用于直接设计iir数字滤波器的函数这一直接设计iir数字滤波器1butterworth模拟和数字滤波器设计数字域
数字信号处理
实验指导书
王创新文卉
长沙理工大学电气与信息工程学院电子信息工程教研室
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
实验一参考程序：1 产生10点的单位抽样序列δ(n)；function unit_pulse(N)% unit_pulse.mN=10;x=zeros(1,N);x(1)=1;n=0:N-1;figure(1);stem(n,x);xlabel('单位抽样序列')axis([-1 20 0 1.1])2产生10点同时移位3位的单位抽样序列δ(n-3)；function shift_unit_pulse(N,k)% shift_unit_pulse.mN=10;k=3;x=zeros(1,N);x(k+1)=1;n=0:N-1;figure(2);stem(n,x);xlabel('移位3位的单位抽样序列')axis([-1 10 0 1.1])function [x, n]=i shift_unit_pulse (n0,ns,nf)n=[0:9]；x=[(n-3)==0]3产生任意序列f(n)=8δ(n)+7δ(n-1)+6δ(n-2) +5δ(n-3)+ 4δ(n-4)+7δ(n-5)；function arbitrary_pulse(N)% arbitrary_pulse.mN=10x=zeros(1,N);x(1)=8;x(2)=7;x(3)=6;x(4)=5;x(5)=4;x(6)=7;n=0:N-1;figure(3);stem(n,x);xlabel('任意序列f(n)')axis([-1 10 0 9])4产生N＝10点的单位阶跃序列function unit_step(N)% unit_step.mN=10;x=ones(1,N);n=0:N-1;figure(4);stem(n,x);xlabel('单位阶跃序列')axis([-1 15 0 1.1])5产生斜率为3，n0=4,点数为20点的斜坡序列g(n)=B(n-n0)function slope(N,k,B)% slope.mN=20;k=4;B=3;x=[zeros(1,k) ones(1,N-k)];for i=1:Nx(i)=B*x(i)*(i-k);endn=0:N-1;figure(5);stem(n,x);xlabel('斜坡序列')axis([-1 20 0 90])6产生幅度A=3，频率f＝100，初始相位 ＝1.2，点数为32点的正弦序列。
function sine(N,A,f,fai)% sine.mA=3;f=100;fai=1.2;N=32;n=0:N-1;x=A*sin(2*pi*f*(n/N)+fai);figure(6);stem(n,x);xlabel('正弦序列')axis([-1 32 -3.2 3.2])7产生幅度A=3，角频率ω＝314，点数为32点的复正弦序列。
function complex_sine(N,A,w)% complex_sine.mA=3; w=314;N=32n=0:N-1;x=A*exp(j*w*n);figure(7);stem(n,x);xlabel('复正弦序列')axis([-1 32 -3.2 3.2])8产生幅度A=3，a＝0.7，点数为32点的实指数序列。
function real_exponent(N,A,a)% real_exponent.mN=32;a=0.7;A=3;n=0:N-1;x=A*a.^n;figure(8);stem(n,x);xlabel('实指数序列');axis([-1 32 0 3.2])9 观察并分析采用不同频率时，对函数50()218.2sin(50)()t a x t e t u t ππ-=的频谱影响。
（a ）：以1000s f Hz =，对其进行采样得到x 1(n)。
（b ）：以300s f Hz =，对其进行采样得到x 2(n)（c ）：以200s f Hz =，对其进行采样得到x 3(n)a=218.2;b=50*pi;fs=1000; %采样频率1000hzT=1/fs;Tp=50*0.001; %观察时间50微秒M=Tp*fs; %采样点数n=0:M-1y=a*exp(-b*n*T).*sin(b*n*T ) %函数表达式subplot(3,2,1)stem(n,y)xlabel('n');ylabel('xa(nT)');title('fs=1000HZ');axis([0,M,-2,1.2*max(abs(y))])yk=T*fft(y,M) %M 点FFTK=0:M-1;fk=K/Tp;subplot(3,2,2)plot(fk,abs(yk))xlabel('f(Hz)');ylabel('幅度');title('T*FFT£¨xa(nt)),fs=1000HZ'); axis([0,fs,0,1.2*max(abs(yk))])fs=300; %采样频率300hzT=1/fs;Tp=50*0.001;M=Tp*fs;n=0:M-1y=a*exp(-b*n*T).*sin(b*n*T )subplot(3,2,3)stem(n,y)xlabel('n');ylabel('xa(nT)');title('fs=300HZ');axis([0,M,-2,1.2*max(abs(y))])yk=T*fft(y,M)K=0:M-1;fk=K/Tp;subplot(3,2,4)plot(fk,abs(yk))xlabel('f(Hz)');ylabel('幅度');title('T*FFT,fs=300HZ'); axis([0,fs,0,1.2*max(abs(yk))]);fs=200; %采样频率200hzT=1/fs;Tp=50*0.001;M=Tp*fs;n=0:M-1y=a*exp(-b*n*T).*sin(b*n*T )subplot(3,2,5)stem(n,y)xlabel('n');ylabel('xa(nT)');title('fs=200HZ');axis([0,M,-2,1.2*max(abs(y))]);yk=T*fft(y,M) %M点FFTK=0:M-1;fk=K/Tp;subplot(3,2,6)plot(fk,abs(yk))xlabel('f(Hz)');ylabel('幅度');title('T*FFT,fs=200HZ'); axis([0,fs,0,1.2*max(abs(yk))])实验二参考程序：%时域直接计算线形卷积clear all;XN=[1 2 3 4 5];HN=[1 2 1 2];YN=COVN(XN,HN);M=LENGTH(YN)-1;N=0:1:M;SUBPLOT(2,1,1)STEM(N,YN);XLABEL('N'); YLABEL('幅度'); TITLE('线性卷积')%用DFT计算卷积clear all;xn=[1 2 3 4 5];hn=[1 2 1 2];m1=length(xn);m2=length(hn);M=m1+m2-1;X=fft(xn,M);H=fft(hn,M);Y=X.*Hyn=ifft(Y,M)n=0:1:M-1;subplot(2,1,1)stem(n,yn);xlabel('n'); ylabel('幅度'); title('线性卷积')%计算循环卷积clear all;xn=[1 2 3 4 5];hn=[1 2 1 2];n=5n1=0:1:n-1;X=fft(xn,n);H=fft(hn,n)Y=H.*X;y=ifft(Y,n)stem(n1,y)xlabel('n'); ylabel('幅度'); title('5点循环卷积')实验三参考程序：（1）直接设计法设计BW（巴特沃斯）高通数字滤波器fp=40; %带通截止频率fs=30; %阻通截止频率ft=200; %采样频率rp=0.5;rs=40;wp=fp/(ft/2); %利用Nyquist频率进行归一化ws=fs/(ft/2);[n,wc]=buttord(wp,ws,rp,rs); %求数字滤波器的最小阶数和截止频率[b,a]=butter(n,wc, 'high'); %设计高通数字滤波器系数b，a[H,W]=freqz(b,a); %求系统频响特性，W为数字角频率，单位rad plot(W*ft/(2*pi),abs(H));grid; %绘出频率响应曲线xlabel('频率/Hz');ylabel('幅值');（2）脉冲响应不变法设计数字滤波器clear all;fs=100wp=0.5*pi;ws=0.8*pi;rp=1;rs=15;Wp=wp*fs; %由数字角频率转换为模拟角频率Ws=ws*fs;[n,wc]=buttord(Wp,Ws,rp,rs,'s'); %选择滤波器的最小阶数[z,p,k]=buttap(n);[Bp,Ap]=zp2tf(z,p,k); %将零极点增益转换为分子分母参数[b,a]=lp2lp(Bp,Ap,wc) % 将低通原型转换为模拟低通[bz,az]=impinvar(b,a,fs); %脉冲相应不变法变换为数字滤波器[H,W]= freqz(bz,az); %求解数字滤波器的频率响应plot(W*fs/(2*pi),abs(H));grid;xlabel('频率/hz');ylabel('幅值');（3）双线性变换法设计数字滤波器clear all;ft=1000;fpl=100;fph=250;wp1= fpl *2*pi; %临界频率采用模拟角频率表示wph= fph*2*pi; %临界频率采用模拟角频率表示wp=[ wp1,wph];wpb=wp/ ft; %求数字频率rp=3;rs=30;fsl=50;fsh=300;ws1= fsl *2*pi; %临界频率采用模拟角频率表示wsh= fsh *2*pi; %临界频率采用模拟角频率表示ws=[ ws1, wsh];wsb=ws/ ft; %求数字频率OmegaP=2* ft *tan(wpb/2);%频率预畸OmegaS=2* ft *tan(wsb/2);%频率预畸%选择滤波器的最小阶数[n,Wn]=buttord(OmegaP,OmegaS, rp, rs,'s'); %此处是代入经预畸变后获得的归一化模拟频率参数[bt,at]=butter(n,Wn,'s'); % 设计一个n阶的巴特沃思模拟滤波器[bz,az]=bilinear(bt,at, ft); %双线性变换为数字滤波器[H,W] = freqz(bz,az); %求解数字滤波器的频率响应plot(W*ft/(2*pi),abs(H));grid;xlabel('频率/Hz');ylabel('幅值');。