1、信号处理MATLAB的应用 机电工程学院工程测试与信号处理 学生姓名： 学 号： 专业班级： 任课教师：大学硕士研究生课程试卷（结课作业）题号一二三四五六七八九总分成绩每个同学根据所具备的实验测量条件，采集一个多频率成分的信号，要求：1介绍试验测量方案，并利用MATLAB对采集的信号进行预处理（掐头去尾等）。根据试验测量方案、现场照片和采集的信号，评价实验设计和实验过程。（共10分）2对信号进行你认为必要的时域分析，并写出分析结论。（共10分）3对信号做快速傅里叶变换（FFT）。写出MATLAB程序（5分）；画出信号的FFT结果图形（5分）；对FFT结果进行分析，并写出分析结论（10分）。（共
2、20分）4对信号做窗口傅里叶变换。写出MATLAB程序（5分）；画出分析结果图形（5分）；对变换结果进行分析，并写出分析结论（10分）。（共20分）5对信号做小波变换。写出MATLAB程序（5分）；画出分析结果图形（5分）；对变换结果进行分析，并写出分析结论（10分）。（共20分）6设计能滤除信号中你认为不需要的频率成分的数字滤波器，并求信号通过该滤波器后的输出。写出MATLAB程序（5分）；绘出输出的时域结果和频域结果（5分）；对结果进行分析，并写出分析结论（10分)。（共20分）1. 信号的采集：模仿老师课堂上所录“河南工业大学”六字音频的例子，我在下面录制“羊年快乐”并保存为 WAV 格
3、式音频文件。此实验信号采集为实验室内操作，考虑到白天噪声干扰，将采集时间选在夜里。现场照相如图1。对录制音频进行测试，得fs=44100，bits=16，符合实验对设备要求。图12. 利用MATLAB在时域对音频文件进行掐头去尾处理。 1）在时域对音频文件进行分析，播放音频文件，并绘制图形如图2，使用程序如下：x,fs,bits=wavread(jinfeng.wav);sound(x(:,1),fs,16);figure(1)N,M=size(x);t=(0:N-1)/fs;subplot(2,1,1);plot(t,x);grid ontitle(原声第一声道，fs=44100 Hz);x
4、label(time/s);ylabel(Amplitude);图2由图形可以看出，音频信号集中于四处，分别为“羊年快乐”四字发音，中间停顿处幅值接近于零，说明所使用音频文件干扰信号微弱，此后进行的一系列分析可以以此信号为基础。 2）对音频文件进行掐头去尾处理，绘制时域图，如图3。程序如下：y=x(30000:N-39000,1);t=t(30000:N-39000);%掐头去尾N,M=size(y);save(wenyan.mat,y);%保存为MATLAB文件%sound(y,fs);%回放t=0:N-1;t=t*1/fs;%计算时间序列subplot(2,1,2);plot(t,y);g
5、rid on %绘制波形图title(第一声道fs=44100 Hz);xlabel(time/s);ylabel(Amplitude);图33. 利用MATLAB 对信号做快速傅里叶变换（FFT），绘制频域分析图形，为便于分析，将功率谱图一同绘出，如图4。程序如下：figure(2)Y=fft(y);magy=abs(Y);f=(0:N-1)*fs/N;subplot(2,1,1);plot(f(1:N),2*magy(1:N)/N);grid onaxis(0 4000 0 0.04)title(原声第一声道，fs=44100 Hz);xlabel(frequency/Hz);ylabel
6、(magnitude);subplot(2,1,2);plot(f(1:N),(2*magy(1:N)/N).2);grid on %功率谱title(第一声道的功率谱);xlabel(frequency/Hz);ylabel(magnitude);axis(0 4000 0 0.0015)图4由于4000Hz以上幅值接近于零，为便于观察数据，将频率轴截取（0,4000）。功率谱图中只保留了较强信号在6000Hz处有微弱信号，但功率谱图中也并没有显示。两图对比可知，功率谱图中幅值与上图并不一致，但增加了对比度，从定性的角度分析频率成分更加直观。由图可以看出，音频文件频率主要分布于200Hz15
7、00Hz，有四个峰值，分别为250Hz、500Hz、800Hz、1400Hz。4. 利用MATLAB 对音频文件做窗口傅里叶变换。分析结果图如图5。程序如下：figure(3)nfft=512;%nfft是FFT的点数wn=nfft;%窗函数的宽度应小于或等于nfftB,f,t=specgram(y,nfft,fs,kaiser(wn,5),wn/2);%反变换的STFTmesh(t,f,2*abs(B)/wn);grid ontitle(第一声道反变换的STFT);xlabel(time/s);ylabel(frequency/Hz);zlabel(magnitude);axis(0 3.5
8、 0 7000 0 0.3)图5 由图5可以看出，“羊年快乐”四个字频率集中分布在02000Hz，但每个字，每个时间点的频率成分并不一样。6000Hz左右频率成分集中分布在“羊”字时段，且此时段低频成分所占比重较大。因另外三个字低频0200Hz低频成分所占比重不大，且高频成分几乎可以忽略，初步假设，使此音频信号通过2001600带通滤波器，会对“羊”字时段影响较大。5. 对信号做小波变换，并对结果图形和变换结果进行分析。clear all;clc;load wsy;whosx=x(1:261507); len = length(x);cw1 = cwt(x,1:32,db1,plot); ti
9、tle(Continuous Transform, absolute coefficients.) ylabel(Scale)cw1,sc = cwt(x,1:32,db1,scal);title(Scalogram) 小波变换结果图 通过分析小波变换的图形结果，随着时间的增加，采集声音信号的能量同样能够直观显示，可知声音频率成分与被放大的小波相似度越高，在声音发生的时间段内，尺度系数越大，颜色越亮白，可知声音频率成分与被放大的小波相似度越高，即以低频成分为主。6. 设计能滤除信号中你认为不需要的频率成分的数字滤波器，并求信号通过该滤波器后的输出。 依据上步假设，使音频信号通过带通滤波器，通带
10、频率为2001600Hz，程序如下：figure(4);Fs=4500;wp=200,1600/(Fs/2);ws=100,1700/(Fs/2);Rp=0.35;Rs=40;N,Wn=cheb1ord(wp,ws,Rp,Rs);b,a=cheby1(N,Rp,Wn);H,W=freqz(b,a);plot(W*Fs/(2*pi),abs(H);y2 = filter(b,a,y);figure(5);plot(y2);grid ontitle(信号经带通滤波器（时域）);sound(200*y2,fs);%回放figure(6)derta_fs = fs/length(y);plot( -f
11、s/2:derta_fs:fs/2-derta_fs,abs(fftshift(fft(y2);grid on title(信号经过带通滤波器（频域）);axis(-10000 10000 0 1200)图6图7图8图9图6为滤波器的幅频、相频特性曲线，图7为信号通过滤波器的时域分析图，图8为信号通过滤波器的频域分析图，图9为信号通过滤波器的窗口傅里叶变换分析图。综合上图知，阻带内频率成分已被滤掉。由图9可以看出，“羊”字时段通过滤波器后信号变得非常微弱，而另外三个字信号相对较强。听声音回放可以明显分辨出后三个字，所受影响不大，而“羊”字虽可听到声音，但已分辨不清具体读音。由此可以验证上步假设的合理性。通过上述分析可知，采集到一组信号后，可对此进行时域、频域各方面进行分析，并作出相应处理。通过时域截断，可将无用信号剔除，压缩数据量；通过滤波，可滤除干扰信号（或任意不需要信号），使结果更加清晰。
