信号与系统频域研讨

1、信号与系统课程研究性学习手册姓名 纪良源 学号 12212052 同组成员 贾皓新 曹世杰 杨欣 高昀峰 指导教师 时间 2014.12.25 信号的频域分析专题研讨【目的】(1) 建立工程应用中有效带宽的概念，了解有限次谐波合成信号及吉伯斯现象。(2) 掌握带限信号，带通信号、未知信号等不同特性的连续时间信号的抽样，以及抽样过程中的参数选择与确定。认识混叠误差，以及减小混叠误差的措施。(3) 加深对信号频域分析基本原理和方法的理解。(4) 锻炼学生综合利用所学理论和技术，分析与解决实际问题的能力。【研讨内容】基础题题目1：吉伯斯现象(1)以定义信号的有效带宽，试确定下图所示信号的有效带宽，取

2、A=1，T=2。(2)画出有效带宽内有限项谐波合成的近似波形，并对结果加以讨论和比较。(3)增加谐波的项数，观察其合成的近似波形，并对结果加以讨论和比较。(a) 周期矩形信号 (b) 周期三角波信号【知识点】连续周期信号的频域分析，有效带宽，吉伯斯现象【信号频谱及有效带宽计算】（a） 周期矩形信号：由于此周期矩形信号的周期T=2，所以W0=2/T=周期信号在区间-1,1的表达式为由于x(t)是奇对称信号，因此有 C0=0根据傅里叶系数的计算 公式，有计算可得周期矩形信号的频谱为 由功率计算公式：得 P=1/4根据有效带宽的定义： 得 （b） 周期三角波信号的频谱及有效带宽的计算：由于此周期矩形

3、信号的周期T=2，所以W0=2/T=周期信号在区间0,2的表达式为由于这个周期三角波信号的偶对称信号，所以根据傅里叶系数的计算 公式，在时有计算可得周期矩形信号的频谱为 由功率计算公式 得 P=1/3 根据有效带宽的定义： 得 【仿真程序】周期方波的频谱t=-2:0.001:2;N=input(N=);c0=0;xN=c0*ones(1,length(t);for n=1:2:NxN=xN+sin(n*pi*t)/(n*pi)-cos(n*pi)*sin(n*pi*t)/(n*pi);endplot(t,xN);grid on;三角波的频谱t=-2:0.001:2;N=input(N=);c0

4、=0.5;xN=c0*ones(1,length(t);for n=1:2:NxN=xN+2*cos(n*pi)*(1-cos(n*pi)*cos(n*pi*t)/(n*n*pi*pi);endplot(t,xN);grid on【仿真结果】周期方波的吉普斯现象N=10N=20N=50N=100三角波的吉普斯现象N=4N=8N=20【结果分析】提示：应从以下几方面对结果进行分析：(1) 图(a) 和图(b)信号有效带宽内有限项谐波合成波形与原波形的近似度比较。答：随着N值的不断增加，图(a) 和图(b)信号有效带宽内有限项谐波合成波形与原波形越来越接近。到一定的程度，两个图形基本达到一致。(2

5、) 分析图(a) 和图(b)信号的时域特性与有效带宽内谐波次数的关系。答：a) 和图(b)信号的时域特性与有效带宽内谐波次数的关系。谐波会随着其次数的增加，而减弱其对信号时域的影响 。(3) 谐波次数增加，图(a) 和图(b)信号合成波形分别有什么变化，从中能得出什么结论？答：谐波次数增加，（a）在不连续点附近部分和 x(t)所呈现的起伏,这个起伏的峰值大小似乎不随N 增大而下降，（b）图中也是一样也就是说，一个不连续信号x(t)的傅里叶级数的截断近似xN(t),一般来说,在接近不连续点处将呈现高频起伏和超量,而且,若在实际情况下利用这样一个近似式的话,就应该选择足够大的N ,以保证这些起伏拥

6、有的总能量可以忽略.当然,在极限情况下,近似误差的能量是零,而且一个不连续的信号(如方波)的傅里叶级数表示是收敛的。这就是吉普斯现象。【自主学习内容】周期方波以及周期三角波的matlab表示，吉普斯现象与傅里叶级数的关系，。【阅读文献】陈后金，胡健，薛健.信号与系统（第二版）M.北京：清华大学出版社，北京交通大学出版社，2005.【发现问题】不管是什么样的信号，其频谱图若用傅里叶级数来模拟信号的话，就会出现吉普斯现象。【问题探究】 【研讨内容】中等题题目2：分析音阶的频谱(1) 录制你所喜欢乐器（如钢琴、小提琴等）演奏的音阶，并存为wav格式。(2) 画出各音阶的时域波形，并进行比较。(3)

7、对所采集的音阶信号进行频谱分析，比较各音阶的频谱。【知识点】连续时间信号的频域分析【温馨提示】利用MATLAB提供的函数fft计算频谱。【题目分析】我选的是钢琴的其中几个音阶，利用fft函数画出他们的各自的频谱，来进行比较。【仿真程序】（1）我由于没有钢琴，所以我就在网上下了几个古筝的音阶。分别是d调的6和e调的6.即d6.ogg和e6.ogg首先用格式转换器将其转换成wav格式音频。再使用wavread函数读取音频。程序如下：fs=44000;nbits=16;y1,fs,nbits=wavread(D:d6.wav);%sound(y1,fs);plot(y1);xlabel(t/s);y

8、label(x(t);title(d6);N=length(y1);n1 = 0:(N-1)/2;n2 = (N-1)/2:-1:0;f=n*fs/N;f =-n2*fs/N n1*fs/N;Xabs=abs(fft(y1,N)/N;Xabs=Xabs/max(Xabs);plot(f,fftshift(Xabs);xlabel(f/Hz);ylabel(X(jw);title(d6);axis(-5000 5000 0 1.5);fs=44000;nbits=16;y1,fs,nbits=wavread(D:e6.wav);%sound(y1,fs);plot(y1);xlabel(t/s)

9、;ylabel(x(t);title(e6);N=length(y1);n1 = 0:(N-1)/2;n2 = (N-1)/2:-1:0;f=n*fs/N;f =-n2*fs/N n1*fs/N;Xabs=abs(fft(y1,N)/N;Xabs=Xabs/max(Xabs);plot(f,fftshift(Xabs);xlabel(f/Hz);ylabel(X(jw);title(e6);axis(-5000 5000 0 1.5);复制过来注释变成了乱码，我也不知道为什么，下面是截图：【仿真结果】 【结果分析】提示：应从以下几方面对结果进行分析：(1) 你所选择乐器演奏的音阶，其时域波形的

10、包络有何特点？答：发现振幅呈现整体衰减的趋势，并有小突起，或者衰减为零，两者图像都差不多，看不出来很大区别。(2) 你所选择乐器演奏的音阶，其频谱有何特点？基波是多少？谐波是多少？答：从频谱来看，每一个音阶的频率在03000Hz左右，而当频率为1500Hz左右时，其幅度可以瞬间升到1左右，频率最高为3500Hz左右。如果仔细看的话，d6达到峰值的频率在1100Hz左右，e6达到峰值的频率在1400Hz左右。【自主学习内容】fft函数的用法【阅读文献】陈后金，胡健，薛健.信号与系统（第二版）M.北京：清华大学出版社，北京交通大学出版社，2005.【发现问题】(1) 改变音阶的包络，相应音阶听起来

11、会有什么变化？答：音色，响度发生变化。(2) 音阶频谱中的谐波分量有什么作用？答：影响音色(3) 你所分析的乐器各音阶对应的频率是多少，之间存在什么关系？答：钢琴各级音阶频率有所不同，通一1、2、3、4、5、6、7在不同调次达到峰值对应频率与调数相关。调越高，达到峰值对应频率越高。【问题探究】【研讨内容】拓展题题目3：连续时间信号的抽样(1) 对带限信号（如,等），确定合适的抽样间隔T，分析的频谱和抽样所得到离散信号的频谱X(ejW)，并将两者进行比较。 (2) 将正弦信号按抽样频率fs=8kHz进行1 秒钟抽样，得离散正弦序列xk为比较f0=2kHz, 2.2 kHz, 2.4 kHz, 2

12、.6 kHz和 f0=7.2 kHz, 7.4 kHz, 7.6 kHz, 7.8 kHz 两组信号抽样所得离散序列的声音，解释所出现的现象。(3) 对于许多具有带通特性的信号，举例验证可否不需要满足?【知识点】连续非周期信号的频谱，离散非周期信号的频谱，时域抽样，频域抽样【温馨提示】 (1) 利用MATLAB提供的函数fft计算抽样所得序列xk的频谱。(2) 利用MATLAB函数 sound(x, fs)播放正弦信号和声音信号。(3) 可以利用仪器或仿真软件产生具有带通特性的信号。【题目分析】(1)一拿到题目我们就遇到了一个问题，正弦信号为什么可以转变为声音，后来经过小组内讨论，大家一致认为

13、正弦信号是一种机械波的模拟，是由振动产生的，因此，可以通过sound函数播放出声音。(2)要想做到抽样1秒钟，就必须注意k的取值范围，即k的取值范围为08000。(3)用stem函数画出抽样后离散的点。【仿真程序】（1）以带限信号Sa(t)为研究对象，x(t)= Sa(t)调用形式为xt=sinc(t);X（）=：w=-5 -pi -pi pi pi 5;xw=0 0 1 1 0 0;plot(w,xw)xlabel(w);ylabel(X(jw)axis(-5 5 -0.5 1.5)【仿真结果】将x(t)按间隔进行T=0.5s,得离散信号xk为xk= Sa(2k)【仿真程序】k=-100:1

14、:100;f=2;x=sinc(pi/f*k);stem(k,x)axis(-*pi 5*pi -1 1)xlabel(k);ylabel( xk)【仿真结果】离散信号xk的频谱X(ejW)=1/T(w-n)=0.1-n 【仿真程序】t=-100:0.0001:100;T=4*pi;w0=2*pi/T;ft=0.5+0.5*square(w0*(t+pi),50);plot(t,ft)axis(-20,20,-0.5,1.5)xlabel(w);ylabel( X(ej)【仿真结果】比较：将x(t)在时域离散化得到xk,则在频域内X(ejW) 是x(jw)周期化（2）【仿真程序】k=0:1:8

15、000;f0=input(f0=);f=8000;x=sinc(pi*f0/f*k);sound(x,f);subplot(211)stem(k,x)axis(0 20 -1 1)subplot(212)plot(k,x)axis(0 20 -1 1)【仿真结果】（a）2KHZ2.2KHZ2.4KHZ2.6KHZ现象：随着声音频率的变大，声音在逐渐变得尖细且高。 (b)7.2KHZ7.4KHZ7.6KHZ7.8KHZ现象：随着声音频率的变大，声音逐渐变得粗犷且低沉。（3）以（1）中的信号x(t)= Sa(t)为例，若，则得到频谱X(ejW)如下：此时，频谱恰好无法区分，如果频谱将会混叠，所以必

16、须满足。【结果分析】频率同在增大，产生两个不同的现象，是因为离散时域信号的频域是以2为周期的，在抽样频率较小时，在频域上会出现频域的非零值叠加(即混叠)，此时随着抽样频率的增加，叠加部分减少，并且使频谱向的奇数倍移动，高频分量增加，声音变得尖一些；而在抽样频率接近信号的频率时，由于频谱的周期性，此时频谱是会随着抽样频率的增加向的偶数倍方向移动，高频分量减少，声音变得低沉。只有符合时域抽样定理的抽样，才不会发生混叠，即抽样的最小频率为2fm。【自主学习内容】（1） Input函数、sound函数、stem函数的用法。（2） 奈奎斯特抽样定理在MATLAB仿真中的应用。（3） MATLAB中对信号

17、进行傅里叶变换的函数fft的用法。【阅读文献】陈后金，胡健，薛健.信号与系统（第二版）M.北京：清华大学出版社，北京交通大学出版社，2005.【发现问题】 若连续时间信号的最高频率未知，该如何确定对信号进行抽样的最大间隔？ 【问题探究】带通信号抽样频率确定的理论分析。【研讨内容】扩展题题目4：男生女生声音信号的转换(1) 采集wav格式的男女生语音信号。(2) 对所采集的男女生信号进行频谱分析，比较男女生频谱的特点。(3) 实现男生女生声音信号的转换。【知识点】连续非周期信号的频谱，离散非周期信号的频谱，时域抽样，频域抽样【温馨提示】可以根据傅里叶变换公式，利用数值积分计算；也可以利用MATL

18、AB提供的函数fft计算。【男生女生声音信号转换的方法】（1） 更改基频通过资料的查询和整理， 可以知道不同人的基频不同 统计如下：正常成年男声：0200Hz；正常成年女声：200450Hz；小孩声音的基频要比女声的高，老年人的基频要比男声的低。经过整理统计可知女声基频=男声基频*1.5。本程序使用的是通过抽样与插值的方式来达到基频的改变。以女变男为例：用整数D对语音信号X（n）进行抽取Xd=X（Dn）；然后将X（n）的抽样频率提高到I（整数）倍 ，即为对X（n）的插值。D/I=3/2；（2） 时长规整通过抽样插值来改变基频 也使播放速度，播放时间发生改变， 因此通过时长规整的方式来使播放速度

19、 和时间恢复到原来。本程序使用的是用重叠叠加算法来达到时长规整。重叠叠加算法原理：它分为两个阶段分解和合成将原始信号以帧长N，帧间距sa进行分解，然后以帧间距ss进行合成。sa与ss的的比值决定了时长规整因子F=sa/ss。为保证重叠区域幅度不变【仿真程序】（2）对所采集的男女生信号进行频谱分析，比较男女生频谱的特点fs=44000;nbits=16;y1,fs,nbits=wavread(D:.wav);%sound(y1,fs);% plot(y1);xlabel(t/s);ylabel(x(t);title();N=length(y1);n=0:N-1;n1 = 0:(N-1)/2;n2

20、 = (N-1)/2:-1:0;f=n*fs/N;f =-n2*fs/N n1*fs/N;Xabs=abs(fft(y1,N)/N;Xabs=Xabs/max(Xabs);plot(f,fftshift(Xabs);xlabel(f/Hz);ylabel(X(jw);title();axis(-5000 5000 0 1.5); fs=44000;nbits=16;y1,fs,nbits=wavread(D:.wav);%sound(y1,fs);% plot(y1);xlabel(t/s);ylabel(x(t);title();N=length(y1);n=0:N-1;n1 = 0:(N-

21、1)/2;n2 = (N-1)/2:-1:0;f=n*fs/N;f =-n2*fs/N n1*fs/N;Xabs=abs(fft(y1,N)/N;Xabs=Xabs/max(Xabs);plot(f,fftshift(Xabs);xlabel(f/Hz);ylabel(X(jw);title();axis(-5000 5000 0 1.5); （3）function Y=voice1(x)%更改采样率使基频改变d=resample(x,3,2);%时长整合使语音文件恢复原来时长W=400;Wov=W/2;Kmax=W*2;Wsim=Wov;xdecim=8;kdecim=2;X=d;F=1.5

22、;Ss =W-Wov;xpts = size(X,2);ypts = round(xpts / F);Y = zeros(1, ypts);xfwin = (1:Wov)/(Wov+1);ovix = (1-Wov):0;newix = 1:(W-Wov);simix = (1:xdecim:Wsim) - Wsim;padX = zeros(1, Wsim), X, zeros(1,Kmax+W-Wov);Y(1:Wsim) = X(1:Wsim);xabs = 0;lastxpos = 0;km = 0;for ypos = Wsim:Ss:(ypts-W);xpos = F * ypos

23、;kmpred = km + (xpos - lastxpos);lastxpos = xpos;if (kmpred = Kmax)km = kmpred;elseysim = Y(ypos + simix);rxy = zeros(1, Kmax+1);rxx = zeros(1, Kmax+1);Kmin = 0;for k = Kmin:kdecim:Kmaxxsim = padX(Wsim + xpos + k + simix);rxx(k+1) = norm(xsim);rxy(k+1) = (ysim * xsim);endRxy = (rxx = 0).*rxy./(rxx+(

24、rxx=0);km = min(find(Rxy = max(Rxy)-1);endxabs = xpos+km;Y(ypos+ovix) = (1-xfwin).*Y(ypos+ovix) + (xfwin.*padX(Wsim+xabs+ovix);Y(ypos+newix) = padX(Wsim+xabs+newix);endend【仿真结果】【结果分析】通过女声与男声的频谱比较，可以看出，女声频谱高频分量较多，而男声高频分量较少，这也就是男声听起来比较浑厚，女声听起来比较尖锐的原因。【自主学习内容】声音的采集及频谱分析【阅读文献】陈后金，胡健，薛健.信号与系统（第二版）M.北京：清华

25、大学出版社，北京交通大学出版社，2005.【发现问题】函数fft计算出的离散频谱Xm和原连续信号的频谱X(jw)之间有何对应关系？答： 利用fft计算出的离散频谱Xm的包络线就是X(jw)的波形一致。【问题探究】系统的频域分析专题研讨【目的】(1) 加深对系统频域分析基本原理和方法的理解。(2) 加深对信号幅度调制与解调基本原理和方法的理解。(3) 锻炼学生综合利用所学理论和技术，分析与解决工程实际问题的能力。【研讨内容】题目1幅度调制和连续信号的Fourier变换本题研究莫尔斯码的幅度调制与解调。本题中信号的形式为其中信号x(t)由文件ctftmod.mat定义，可用命令Load ctftm

26、od 将文件ctftmod.mat定义的变量装入系统内存。运行命令Load ctftmod后，装入系统的变量有af bf dash dot f1 f2 t x其中bf af： 定义了一个连续系统H(s)的分子多项式和分母多项式。可利用freqs(bf,af,w)求出该系统的频率响应，也可用sys=tf(bf,af)得到系统的模型，从而用lsim求出信号通过该系统的响应。dash dot： 给出了莫尔斯码中的基本信号dash和dot的波形f1 f2： 载波频率t: 信号x(t)的抽样点x: 信号x(t)的在抽样点上的值信号x(t)含有一段简单的消息。Agend 007的最后一句话是The fut

27、ure of technology lies in 还未说出最后一个字，Agend 007就昏倒了。你(Agend 008)目前的任务就是要破解Agend 007的最后一个字。该字的信息包含在信号x(t)中。信号x(t)具有式0的形式。式中的调制频率分别由变量f1和f2给出，信号m1(t)，m2(t)和m3(t)对应于字母表中的单个字母，这个字母表已用国际莫尔斯码进行编码，如下表所示：A - H O - - -V -B - I P -W -C -J - - -Q - -X -D -K -R -Y -E L -S Z -F -M - -T -G - -N -U -(1) 字母B可用莫尔斯码表示为

28、b=dash dot dot dot，画出字母B莫尔斯码波形；(2) 用freqs(bf,af,w)画出系统的幅度响应；(3) 利用lsim求出信号dash通过由sys=tf(bf,af)定义的系统响应，解释你所获得的结果；(4)用解析法推导出下列信号的Fourier变换(5)利用(4)中的结果，设计一个从x(t)中提取信号m1(t)的方案，画出m1(t)的波形并确定其所代表的字母；(6)对信号m2(t)和m3(t)重复(5)。请问Agent 008The future of technology lies in 【题目分析】【仿真程序】(1）w=linspace(-5,5,200);RC=0

29、.5;b=1;a=RC 1;H=freqs(b,a,w);subplot(4,1,1);plot(w,abs(H);xlabel=(omega);ylabel=(|H(jomega)|);RC=1;a=RC 1;H=freqs(b,a,w);subplot(4,1,2);plot(w,abs(H);xlabel=(omega);ylabel=(|H(jomega)|);RC=2;a=RC 1;H=freqs(b,a,w);subplot(4,1,3);plot(w,abs(H);xlabel=(omega);ylabel=(|H(jomega)|);RC=4;a=RC 1;H=freqs(b,

30、a,w);subplot(4,1,4);plot(w,abs(H);xlabel=(omega);ylabel=(|H(jomega)|);（2）t=linspace(0,0.2,1024);RC=0.4w1=100;w2=3000;H1=1/(j*w1*RC+1);H2=1/(j*w2*RC+1);x=cos(w1*t)+cos(w2*t);y=abs(H1)*cos(w1*t) +abs(H2)*cos(w2*t);subplot(2,1,1);plot(t,x);ylabel(x(t);xlabel(时间（秒）);subplot(2,1,2)plot(t,y);ylabel(y(t);x

31、label(时间（秒）);(3)【1】t=linspace(0,0.08,1024);RC=input(RC=);H=1/(j*100*pi*RC+1);x=10*abs(sin(100*pi*t);y=abs(H)*10*abs(sin(100*pi*t);plot(t,y);ylabel(y(t);xlabel(时间（秒）);【2】t=0:1/1800:0.1;xt=10*abs(sin(100*pi*t);RC=1;w1=100;t=linspace(0,0.1,1024);H1=1/(j*w1*RC+1);y=abs(H1)*abs(sin(100*pi*t);zx=sum(xt);z

32、y=sum(y);subplot(2,1,1);plot(t,zx);legend(x的直流分量);subplot(2,1,2);plot(t,zy);legend(y的直流分量);【仿真结果】（1）（2）（3）【1】RC=0.4RC=1RC=2RC=5【2】【结果分析】(1)随着RC取值的增加，幅度响应曲线越集中。 （2）起初的波形很混乱，经过低频滤波后变成规则的图形【自主学习内容】利用sum函数进行直流分量的计算，用freqs计算连续系统的频率响应。 如何滤除高频信号。【阅读文献】信号与系统 Matlab原理及简单应用【发现问题】 (专题研讨或相关知识点学习中发现的问题)：在滤波的过程中，

33、随着RC的变化，波形呈现规律性变化。【问题探究】在进行正弦信号的频谱分析时.发现正弦信号的频谱的冲击会随着t=linspace(a,b,f)中的f的值而发生改变，如下列程序：t=linspace(0,0.2,2000)w1=100;x=cos(w1*t)n=length(x);X=fft(x,n);subplot(2,1,1);plot(t,x);仿真结果为：这是当f=2000时的仿真结果。下面是f=4000时的仿真结果： 后来我在网上搜了一下关于正弦信号的频谱分析问题：得到下面的结果：程序：t=-1:0.01:1;x=5*sin(2*pi*10*t);N=length(x);fx=fft(x

34、);df=100/N;n=0:N/2;f=n*df;subplot(211);plot(t,x);grid on; subplot(212);plot(f,abs(fx(n+1)*2/N); grid;仿真结果：这是对f=1/T=w/2而言的,得到如下程序：t=-1:0.01:1;x=5*sin(2*pi*10*t);N=length(x);fx=fft(x);df=100/N;n=0:N/2;w=n*df*2*pi;subplot(211);plot(t,x);grid on; subplot(212);plot(w,abs(fx(n+1)*2/N); grid;仿真结果为：【研讨内容】题目2分析实际物理系统的频率响应。【数学模型】接入端。输出端。通过分析：发现（1）和（2）的系统描述还有点相似，所以把它弄他们弄到一块讨论，（3）和（4）同点是（3）是对称的，而（4）中xk-2没有，且两边的数是相反数，所以也把他们弄到一块讨论。此题考查的是离散系统的频率响应的求解，可根据描述方程的下列频率响应：（1）(2) (3) (4) 再照着课本求离散系统频率响应。【仿真程序】（1） 和（2）a1=1;a2=1;b1=1.69 1.05 0.45 0.45 1.05 1.69;b2=1.69 1.05 0.45 -0.45 -1.05 -1.69;h1,w=