北航数字信号处理上机实验一实验报告.docx
- 文档编号:7012635
- 上传时间:2023-01-16
- 格式:DOCX
- 页数:20
- 大小:283.98KB
北航数字信号处理上机实验一实验报告.docx
《北航数字信号处理上机实验一实验报告.docx》由会员分享,可在线阅读,更多相关《北航数字信号处理上机实验一实验报告.docx(20页珍藏版)》请在冰豆网上搜索。
北航数字信号处理上机实验一实验报告
离散时间信号处理实验报告
实验一信号的采样与重构
班级
学号
姓名
同组者
日期
实验介绍
连续时间信号采样是获得离散时间信号的一种重要方式,但是时域上的离散化会带来信号在频域上发生相应的变化。
在本实验中,我们将分别看到低通信号和带通信号在不同的采样率下得到的离散信号波形与连续信号波形在时域和频域上的对应关系。
同时,离散信号的二次采样在实际的应用中可能是必须的,有时甚至是非常重要的。
在实验的最后,我们也会看到离散信号的抽取和内插所带来的频谱变化。
由于matlab语言无法表达连续信号,实验中我们采用足够密的采样点来模拟连续信号(远大于奈奎斯特采样的要求),即:
t=0:
Ts:
T(Ts=1/fs<<奈奎斯特采样频率)
实验中,为了分析离散信号与连续信号之间的频谱关系,加深对采样定理的理解,了解模拟频谱、数字频谱、以及离散信号被加窗后各自的频谱,从而直观的理解采样频率对频谱的影响和加窗后对频谱的影响。
由此可以掌握数字处理方法对模拟信号进行频谱分析的基本原则,即:
如何选择合适的信号长度、采样周期以使得对模拟信号的频谱分析的误差达到分析的要求。
在该实验中,用到的Matlab函数有:
plot(x,y),其作用是在坐标中以x为横坐标、y为纵坐标的曲线,注意x和y都是长度相同的离散向量;
xlabel(‘xxx’),其作用是对x轴加上坐标轴说明“xxx”;
ylabel(‘yyy’),其作用是对y轴加上坐标轴说明“yyy”;
title(‘ttt’),其作用是对坐标系加上坐标轴说明“ttt”;
subplot(m,n,w),其作用是当需要在同一显示面板中显示多个不同的坐标系时,m、n分别指明每行和每列的坐标系个数,w为当前显示坐标系的流水号(1到m*n之间)。
在实验中我们需要画出信号的频谱,对于连续信号频谱的逼近需要你自己编写,原理如下:
连续时间非周期信号
的傅里叶变换对为:
用DFT方法对该变换逼近的方法如下:
1、将
在t轴上等间隔(宽度为T)分段,每一段用一个矩形脉冲代替,脉冲的幅度为其起始点的抽样值
,然后把所有矩形脉冲的面积相加。
该方法实际为平顶处理,利用采样和零阶保持器就可以完成,则有:
2、将序列
截断成从t=0开始长度为
的有限长序列,包含有N个采样,即时域加矩形窗,则上式又可以进一步近似为:
由于时域采样,采样频率为
,则频域产生以
为周期的周期延拓,如果是限带信号,则有可能不产生混叠,成为连续周期频谱序列,频域周期为
。
3、由于数值计算的限制,在频域上也只能计算离散点(频域抽样)上的数值。
我们将频域的一个周期
中也分成N段,即
。
每个频域采样点间的间隔为
。
则上式可以进一步化简为:
由此我们看到了DFT与CTFT之间的近似关系。
如果将T和
取得尽量小,则我们可以得到模拟信号的较精确的时频特性。
实验解答
一、低通采样
考虑以下模拟信号xa(t),xa(t)=cos(2*pi*f0*t)+2*cos(2*pi*10*f0*t),令f0=100Hz。
显然该信号的最高频率1kHz。
1.采样频率fs大于信号最高频率10f0的2倍是信号的恢复。
1)对信号以fs1=2.2kHz的频率进行采样,得到抽样信号x(n)=xa(t)|t=nT;
2)利用内插公式xr(t)=x(n)*[sin(pi*(t-n*Ts)/Ts)/(pi*(t-n*Ts)/Ts)]恢复出原始信号x^a(t);
3)绘出Δx(t)=xa(t)-x^a(t)的时域波形;
4)利用前面介绍的画频谱的方法画出xa(t)和x(n)的频谱,对应3)中的误差曲线和频谱图说明是否有新的频率分量产生,是何原因?
2.采样频率fs低于信号最高频率10f0的2倍是信号的恢复。
分别取fs=6*f0和fs=14*f0,重复1中的操作,注意观察时域波形和频谱的对应变化。
(一)、MATLAB源代码
clear;
f0=100;
Fs=3000;
Ts=1/Fs;
ts=0:
Ts:
0.1;
X=xa(ts);
t=0:
Ts/10:
0.1;
L=length(t);
Xa=xa(t);
Xr=zeros(1,L);
forn=1:
L
Xr(n)=xr(X,Ts,t(n));
end
subplot(6,1,1);
plot(t,Xa);
title('原始信号');
subplot(6,1,2);
stem(ts,X);
title('采样信号');
subplot(6,1,3);
plot(t,Xr);
title('恢复信号');
subplot(6,1,4);
plot(t,Xa-Xr);
title('原信号与恢复信号之差');
fs=20000;
y1=fft(Xa,fs);
mag=abs(y1);
f=(0:
length(y1)-1)'*(10*Fs)/length(y1);
subplot(6,1,5)
plot(f(1:
fs/2),mag(1:
fs/2));
axis([0120003000]);
fs=3000;
y2=fft(X,fs);
title('原信号频谱');
mag=abs(y2);
f=(0:
length(y2)-1)'*(Fs)/length(y2);
subplot(6,1,6)
plot(f(1:
fs/2),mag(1:
fs/2));
title('抽样信号频谱');
functionf=xa(t)
f0=100;
f=cos(2*pi*f0*t)+2*cos(2*pi*10*f0*t);
end
functionx=xr(X,Ts,t)
x=0;
forn=1:
length(X)
x=x+X(n)*sinc((t-(n-1)*Ts)/Ts);
end
(二)、运行结果
1、
fs=2200Hz
2、fs=600Hz
3、fs=1400Hz
(三)、分析及结论
1(4):
当fs=2200Hz时,基本可认为没有新的频率分量产生,当fs=600Hz或fs=1400Hz时,有新的频率分量产生。
可能原因之一是当不满足fs>2f0时,采样信号不能完全恢复原始信号,因此在频域上会产生新的频率分量。
2、当不满足fs>2f0时,可以看到时域恢复信号出现失真,且fs越小,恢复信号与原信号差距越大。
在频域上表现为产生了新的频率分量(200/400Hz)。
综上所述,欲使采样信号能无失真地恢复原信号,要求是fs>=2f0,若不满足该条件,会出现恢复波形失真的情况。
二、带通采样
通过该题目进一步加深对带通信号采样定理的认识,即不需要象低通采样那样需要信号最高频率的2倍才能恢复出原始信号。
带通信号的采样定理可描述如下:
如果模拟信号f(t)为带通信号,其角频率限制在fL和fm之间,则必须的最低采样频率fs>=2(fm-fL)。
且当采样率满足fs=4f0/(2n+1)(f0=(fm+fL)/2为信号的中心频率,n=0,1,2,……),则可以无失真的从采样信号中恢复原始信号。
注意当fs的取值大于2fm时,与低通采样定理意义相同。
考虑模拟信号xa(t),xa(t)=cos(2*pi*f0*t)+2*cos(2*pi*f1*t),令f0=1.6kHz,f1=2kHz。
该信号带宽0.4k< 1.请你根据带通采样定理,选取适当的采样速率以得到抽样信号x(n),然后同样用内插公式恢复成模拟信号,并思考如何由该模拟信号得到原来的带通信号。 2.请绘出带通信号的频谱及x(n)的频谱,观察并比较它们。 (一)、MATLAB源代码 1、fc=400;%信号带宽400HZ fmin=2*fc; fc0=(1600+2000)/2; fc1=1600;%信号一的频率 fc2=2000;%信号二的频率 fs1=400;%欠采样的采样频率 fs2=800;%临界采样采样频率 n=0; fs3=4*fc0/(2*n+1);%n取1时的采样速率 %fs3=5000;%n取1时的采样速率 f0=80000;%用以模拟连续信号的离散信号采样速率 %(大nyquist率) t=0: 1/f0: 1; N=1*f0; xt0=cos(2*pi*fc1*t)+2*cos(2*pi*fc2*t); figure (1); subplot(3,1,1); plot(t,xt0); axis([00.03-44]); title('待采样信号波形'); xlabel('x');ylabel('xt'); yjw0=fft(xt0,N); absy0=abs(yjw0); f=(0: N-1)*f0/N; %subplot(2,1,2); %plot(f,absy0); %axis([08000-410000]); %title('正常采样信号频谱') %xlabel('数字频率');ylabel('幅度谱'); t=0: 1/fs3: 1; N=1*fs3; xts=cos(2*pi*fc1*t)+2*cos(2*pi*fc2*t); figure (1); subplot(3,1,2); plot(t,xts); axis([00.03-44]); title('采样信号波形'); xlabel('x');ylabel('xt'); yjw0=fft(xts,N); absy0=abs(yjw0); f=(0: N-1)*f0/N; fs=fs3;%采样频率 Ts=1/fs3;%采样周期 T0=0.01; n=0: (3*T0)/Ts; t1=0: Ts: 3*T0; %disp('t1='+t1); %x1=sin(2*pi*f0*n/fs)+1/3*sin(6*pi*f0*n/fs);%采样信号 x1=cos(2*pi*fc1*n/fs)+2*cos(2*pi*fc2*n/fs); T_N=ones(length(n),1)*t1-n'*Ts*ones(1,length(t1));%表示t-nT xa=x1*sinc(fs*T_N); subplot(3,1,3); plot(t1,xa); axis([00.03-44]); title('恢复信号波形'); 2、fc=400;%信号带宽400HZ fmin=2*fc; fc0=(1600+2000)/2; fc1=1600;%信号一的频率 fc2=2000;%信号二的频率 fs1=400;%欠采样的采样频率 fs2=800;%临界采样采样频率 n=0; fs3=4*fc0/(2*n+1);%n取1时的采样速率 f0=80000;%用以模拟连续信号的离散信号采样速率 t=0: 1/f0: 1; N=1*f0; xt0=cos(2*pi*fc1*t)+2*cos(2*pi*fc2*t); figure (1); yjw0=fft(xt0,N); absy0=abs(yjw0); f=(0: N-1)*f0/N; subplot(2,1,1); plot(f,absy0); %axis([08000-410000]); title('待采样信号频谱') xlabel('数字频率');ylabel('幅度谱'); %(大nyquist率) t=0: 1/fs3: 1; N=1*fs3; xt0=cos(2*pi*fc1*t)+2*cos(2*pi*fc2*t); figure (1); yjw0=fft(xt0,N); absy0=abs(yjw0); f=(0: N-1)*f0/N; subplot(2,1,2); plot(f,absy0); axis([080000-410000]); title('正常采样信号频谱') xlabel('数字频率');ylabel('幅度谱'); (二)、运行结果 1、 2、 (三)、分析及结论 1、对采样信号及频谱,这里n取0;对恢复信号: 通过内插公式由采样信号得到恢复信号,可以看到,采样信号比较完整的恢复了原始信号。 2、可以看到,采样信号的频谱发生了一定搬移,但包含和了原始信号的所有频谱信息。 总的来说,当采样率满足fs=4f0/(2n+1)(f0=(fm+fL)/2为信号的中心频率,n=0,1,2,……)时,即正常采样下的信号可以恢复出原信号一个完整周期的频谱。 可能是由于采样频率和信号长度的问题,恢复信号有一定的失真,但这已经能够证明带通采样定理的正确性,即采样后的信号没有丢失原信号的信息。 三、加窗后信号幅频特性的变化 考虑一个以指数率衰减的信号xa(t)=e-atcos(2*pi*f0*t),为简便起见,衣复述形式表示为xa(t)=e-ate2*pi*f0。 现在以fs=1/T抽样,则得到抽样信号x(n)=xa(nT),n=…,-3,-2,-1,0,1,2,…。 这样的无限长序列计算机是无法存储的,通常的做法是令xL(nT)=x(nT),n=0,1,.…,L-1,它是长度为L的加窗信号。 模拟、采样、加窗后信号的频谱幅度如下: 我们有极限 , 其中: α=0.2s-1,f0=0.5Hz,fs=1Hz。 1、请在同一图中绘出模拟信号频谱|Xa(f)|2及采样信号频谱T|Xa(f)|2,在另外一张图中绘出三个频谱|Xa(f)|2,T|Xa(f)|2,T|XL(f)|2(频率范围: 0 2、改变xL(nT)=x(nT)的长度L并重复1,观察其中的变化。 (一)、MATLAB源代码 f=[0: 0.01: 3]; a=0.2; f0=0.5; fs=1; T=1/fs; L=25; xa=1./(a.^2+(2*pi*(f-f0)).^2); xf=1./(1-2*(exp(-a*T))*cos(2*pi*(f-f0))+exp(-2*a*T)); xl=(1-2*(exp(-a*T*L))*cos(2*pi*(f-f0)*L)+exp(-2*a*T*L))./(1-2*(exp(-a*T))*cos(2*pi*(f-f0))+exp(-2*a*T)); subplot(3,1,1); plot(f,xa); title('模拟信号频谱幅度'); holdon subplot(3,1,2); plot(f,xf); title('采样信号频谱幅度'); subplot(3,1,3) plot(f,xl); title('长度为100的加窗信号') holdoff (二)、运行结果 取α=0.2s-1,f0=0.5Hz,fs=1Hz 1、模拟信号频谱,采样后频谱和加窗长度为2时处理的信号频谱 2、模拟信号频谱,采样后频谱和加窗长度为25时处理的信号频谱 3、模拟信号频谱,采样后频谱和加窗长度为100时处理的信号频谱 (三)、分析及结论 对模拟信号进行数字处理时,首先要对模拟信号进行采样,采样频率由奈奎斯特采样定理决定。 对采样而来的数字信号进行DTFT处理得到其频谱。 由DTFT的计算公式可知,DTFT的计算需要用到信号的所有采样点,当信号为无限长或者是相当长时,这样的计算不可行也没有实际意义。 因此会把信号分成许多一定长度的数据段,然后分段处理。 如果把数据进行分段,相当于对信号进行了加矩形窗的处理。 当窗信号长度较低时,这时候频谱会有较大的失真,这是因为此时由于窗函数是一个频带无限的函数,所以即使原信号x(t)是限带宽信号,而在截断以后也必然成为无限带宽的函数,即信号在频域的能量与分布被扩展了。 又从采样定理可知,无论采样频率多高,只要信号一经截断,就不可避免地引起混叠,因此信号截断必然导致一些误差。 那么加长截断长度T,即矩形窗口加宽,则窗谱W(ω)将被压缩变窄(π/T减小)。 虽然理论上讲,其频谱范围仍为无限宽,但实际上中心频率以外的频率分量衰减较快,因而泄漏误差将减小。 当窗口宽度T趋于无穷大时,则谱窗W(ω)将变为δ(ω)函数,而δ(ω)与X(ω)的卷积仍为X(ω),这说明,如果窗口无限宽,即不截断,就不存在泄漏误差。 这些在实验结果中得到了体现。 四、信号的抽取和内插 设离散时间信号为汉宁窗信号(升余弦): 当N=32时: 1、给出 的幅频特性曲线; 2、给出x(n)经两倍抽取之后的幅频特性; 3、给出x(n)经两倍内插之后的幅频特性; 4、当N=16时,重复1、2、3。 (一)、MATLAB源代码 n=1: 32; x=0.5*(1-cos(2*pi*n/31)); X=fft(x); fs=2000;%采样频率 f=n*fs/32/1000; subplot(3,2,1); plot(f,abs(X)); xlabel('相位(单位pi)','fontsize',13); ylabel('幅值','fontsize',13); title('(N=32)幅频特性曲线1','fontsize',13); x2=dyaddown(x); n2=1: 16; X2=fft(x2); f2=n2*fs/16/1000; subplot(3,2,3); plot(f2,abs(X2)); xlabel('相位(单位pi)','fontsize',13); ylabel('幅值','fontsize',13); title('(N=32)幅频特性曲线2','fontsize',13); x3=dyadup(x); n3=1: 65; X3=fft(x3); f3=n3*fs/65/1000; subplot(3,2,5); plot(f3,abs(X3)); xlabel('相位(单位pi)','fontsize',13); ylabel('幅值','fontsize',13); title('(N=32)幅频特性曲线3','fontsize',13); N=500; n=1: 16; y=0.5*(1-cos(2*pi*n/15)); Y=fft(y); f=n*fs/16/1000; subplot(3,2,2); plot(f,abs(Y)); xlabel('相位(单位pi)','fontsize',13); ylabel('幅值','fontsize',13); title('(N=16)幅频特性曲线1','fontsize',13); y2=dyaddown(y); n2=1: 8; Y2=fft(y2); f2=n2*fs/8/1000; subplot(3,2,4); plot(f2,abs(Y2)); xlabel('相位(单位pi)','fontsize',13); ylabel('幅值','fontsize',13); title('(N=16)幅频特性曲线2','fontsize',13); y3=dyadup(y); n3=1: 33; Y3=fft(y3); f3=n3*fs/33/1000; subplot(3,2,6); plot(f3,abs(Y3)); xlabel('相位(单位pi)','fontsize',13); ylabel('幅值','fontsize',13); title('(N=16)幅频特性曲线3','fontsize',13); (二)、运行结果及分析 1、N=32时的仿真结果 (1)(N=32)幅频特性曲线1为x(n)的幅频特性曲线。 可以看出,横轴零处为直流分量的频谱,0-2pi内两个上升沿为余弦信号的频谱体现。 (2)(N=32)幅频特性曲线2为x(n)经两倍抽取之后的幅频特性。 可以看出,经过抽取之后,频谱图形变宽,上升沿和直流分量图形更加明显。 (3)(N=32)幅频特性曲线3为x(n)经两倍内插之后的幅频特性。 可以看出,经过内插,频域被挤压缩小,频谱图形被压缩,1pi处分别是直流分量和余弦信号经过内插后的频域图形表示。 2、N=16时的仿真结果 (1)(N=16)幅频特性曲线1为x(n)的幅频特性曲线。 可以看出,横轴零处为直流分量的频谱,0-2pi内两个上升沿为余弦信号的频谱体现。 (2)(N=16)幅频特性曲线2为x(n)经两倍抽取之后的幅频特性。 可以看出,经过抽取之后,频谱图形变宽,上升沿和直流分量图形更加明显。 (3)(N=16)幅频特性曲线3为x(n)经两倍内插之后的幅频特性。 可以看出,经过内插,频域被挤压缩小,频谱图形被压缩,1pi处分别是直流分量和余弦信号经过内插后的频域图形表示。
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 北航 数字信号 处理 上机 实验 报告