现代信号分析课程大作业

课程名称：现代信号分析考试形式：□专题研究报告□论文大作业□综合考试学生姓名：学号：序号分项类别得分1题目12题目23题目34题目45题目56题目67题目78题目89题目910题目10总分评阅人：时间：年月目录一、FFT算法的特点。.........................................-3-二、平稳随机过程.............................................-8-三、功率谱估计..............................................-10-四、自适应滤波器（LMS和RLS）...............................-16-五、维纳滤波器..............................................-23-六、FIR维纳滤波器..........................................-25-七、AR/MA滤波器............................................-28-八、卡尔曼滤波器............................................-31-九、小波分析方法及应用......................................-37-十、其他现代信号分析方法....................................-43-现代信号分析-3-一、FFT算法的特点。题目：请举例说明FFT算法的特点计算离散傅里叶变换的一种快速算法，简称FFT。快速傅里叶变换是1965年由J.W.库利和T.W.图基提出的。采用这种算法能使计算机计算离散傅里叶变换所需要的乘法次数大为减少，特别是被变换的抽样点数N越多，FFT算法计算量的节省就越显著。根据对序列分解与选取方法的不同而产生了FFT的多种算法，基本算法是基2DIT和基2DIF。FFT在离散傅里叶反变换、线性卷积和线性相关等方面也有重要应用。快速傅氏变换（FFT），是离散傅氏变换的快速算法，它是根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的。有限列长为N的序列x（n）的DFT变换公式在MATLAB中的表达式为11)()(NnnkNWnxkX1,...,2,1Nk其逆变换为NknkNWkXNnx1)(1)(Nn,...,2,1FFT主要就是利用了nkNW以下两个特性使长序列的DFT分解为更小点数的DFT来实现的。（1）利用nkNW的对称性使DFT运算中有些项合并*)()(knNknNnNkN（2）利用nkNW的周期性和对称性使长序列的DFT分解为更小点数的DFTnNkNNnkNnkN)()(FFT算法在应用时能够大大减少计算量，但是也存在着和DFT同样的一些问题，主要有：频谱混叠，泄露问题以及栅栏效应的问题，下面通过运用MATLAB的FFT相关函数对信号进行变换来进行实际分析。（1）频谱混叠奈奎斯特定理已被众所周知了，所以几乎所有人的都知道为了不让频谱混叠，理论上采样频谱大于等于信号的最高频率。采样周期的倒数是频谱分辨率，最高频率的倒数是采样周期。设定采样点数为N，采样频率fs，最高频率fh，故频谱分辨率f=fs/N,而fs=2fh，所以可以看出最高频率与频谱分辨率是相互矛盾的，提高频谱分辨率f的同时，在N确定的情况下必定会导致最高频率fh的减小；同样的，提高最高频率fh的同时必会引起f的增大，即分辨率变大,下面举例说明。设输入信号为tts300sin5240cos3，即它的主要频率为120Hz和150Hz。分别取采样频率为400Hz，300Hz，200Hz，100Hz，程序如下：L=1024;%采样点数NFFT=2^nextpow2(L);%Nextpowerof2fromlengthofyFs1=400;%采样频率-4-t1=(0:L-1)/Fs1;s1=3*cos(2*pi*120*t1)+5*sin(2*pi*150*t1);Y1=fft(s1,NFFT)/L;f1=Fs1/2*linspace(0,1,NFFT/2+1);Fs2=300;t2=(0:L-1)/Fs2;s2=3*cos(2*pi*120*t2)+5*sin(2*pi*150*t2);Y2=fft(s2,NFFT)/L;f2=Fs2/2*linspace(0,1,NFFT/2+1);Fs3=200;t3=(0:L-1)/Fs3;s3=3*cos(2*pi*120*t3)+5*sin(2*pi*150*t3);Y3=fft(s3,NFFT)/L;f3=Fs3/2*linspace(0,1,NFFT/2+1);Fs4=100;t4=(0:L-1)/Fs4;s4=3*cos(2*pi*120*t4)+5*sin(2*pi*150*t4);Y3=fft(s4,NFFT)/L;f4=Fs4/2*linspace(0,1,NFFT/2+1);subplot(2,2,1);plot(f1,2*abs(Y1(1:NFFT/2+1)));title('采样频率为400Hz');gridon;xlabel('频率(Hz)');subplot(2,2,2);plot(f2,2*abs(Y2(1:NFFT/2+1)));title('采样频率为300Hz');gridon;xlabel('频率(Hz)');subplot(2,2,3);plot(f3,2*abs(Y3(1:NFFT/2+1)));title('采样频率为200Hz');gridon;xlabel('频率(Hz)');subplot(2,2,4);plot(f4,2*abs(Y4(1:NFFT/2+1)));title('采样频率为100Hz');gridon;xlabel('频率(Hz)');gridon;结果如图1-1所示，可见当采样频率为300Hz时，没有发生频率混叠，但当频率小于两倍的信号频率时，经过采样就会发生混叠现象，输入信号中的某些频率就无法辨识出。现代信号分析-5-图1-1FFT的混淆现象（2）泄漏现象实际中的信号序列往往很长，甚至是无限长序列。为了方便，我们往往用截短的序列来近似它们。这样可以使用较短的DFT来对信号进行频谱分析。这种截短等价于给原信号序列乘以一个窗函数。而窗函数的频谱不是有限带宽的，从而它和原信号的频谱进行卷积以后会扩展原信号的频谱。值得一提的是，泄漏是不能和混叠完全分离开的，因为泄露导致频谱的扩展，从而造成混叠。为了减小泄漏的影响，可以选择适当的窗函数使频谱的扩散减到最小。下面举例说明。设输入信号为tts300sin5240cos3，即它的主要频率为120Hz和150Hz。，用MATLAB自带的信号处理工具箱（通过命令sptool调出）分别对信号加矩形窗（boxcar），海明窗（Hamming），汉宁窗（Hanning）和布莱克曼窗（Blackman），分析不同的窗函数对减小泄露现象的影响。具体操作如下：先生成信号s，然后在sptool工具箱的file菜单中import信号s，Fs为400Hz，如图1-2所示。图1-2sptool工具箱界面然后点create按钮，在patameters栏选择welch方法，Nfft点数1024，Nwind为256，分别加布-6-莱克曼窗,矩形窗，海明窗，和三角窗，得如下PSD图。(a)blackman(b)boxcar(c)hamming(d)triang图1-3FFT的泄露现象由以上四图可见，布莱克曼窗主瓣宽，旁瓣小，频率识别精度最低，但幅值识别精度最高。矩形窗主瓣窄，旁瓣大，频率识别精度最高，幅值识别精度最低；海明窗是余弦窗的一种，又称改进的升余弦窗，海明窗与汉宁窗都是余弦窗，只是加权系数不同。海明窗加权的系数能使旁瓣达到更小，海明窗的频谱是由3个矩形时窗的频谱合成，但其旁瓣衰减速度为20dB／（10oct）。三角窗与矩形窗比较,主瓣宽约等于矩形窗的两倍,但旁瓣小,而且无负旁瓣。所以对于窗函数的选择，应考虑被分析信号的性质与处理要求。如果仅要求精确读出主瓣频率，而不考虑幅值精度，则可选用主瓣宽度比较窄而便于分辨的矩形窗，如果分析窄带信号，且有较强的干扰噪声，则应选用旁瓣幅度小的窗函数，如汉宁窗。（3）栅栏效应栅栏效应是因为用计算频谱只限制为基频的整数倍而不可能将频谱视为一连续函数而产生的。设信号的最高频率为hf，要求的频率分辨率为F（即在频率轴上所能得到的最小频率间隔），则采样点数N需要满足以下条件：FfNh2若不满足，就会产生栅栏效应。减小栅栏效应通常有两个方法：一个方法是在原纪录末端添加一些零值点来变动时间周期内的点数，并保持记录不变。从而在保持原有频谱连续形式不变的情况下，变更了谱线的位置。这样，原来看不到的频谱分量就能移到可见的位置上了。但是，补零并没有增加序列的有效长度，所以并不能提高分辨率。但补零同时可以使数据长度N为2的整数幂，以便于应用FFT。另一个方法是增加采样点数，从而减小频率间隔，使原来被挡住的一些频谱的峰点或谷点显露出来。注意，这时候每根谱线多对应的频率和原来的已经不相同了。下面通过实际的函数来分析两种方法的效果。现代信号分析-7-程序如下：L=90;%采样点数NFFT=2^nextpow2(L);%Nextpowerof2fromlengthofyFs=400;%²采样频率t=(0:L-1)/Fs;s=cos(200*pi*t)+2*cos(204*pi*t);Y1=fft(s,NFFT)/L;f1=Fs/2*linspace(0,1,NFFT/2+1);plot(f1,2*abs(Y1(1:NFFT/2+1)));xlabel('频率');%----原序列末端补零-------------------------------------------L3=1024;NFFT=2^nextpow2(L3);s3=[szeros(1,L3-L)];Y3=fft(s3,NFFT)/L3;f3=Fs/2*linspace(0,1,NFFT/2+1);figure;plot(f3,2*abs(Y3(1:NFFT/2+1)));xlabel('频率');('频率');%----增加采样点数---------------------------------------------L2=1024;NFFT=2^nextpow2(L2);t=(0:L2-1)/Fs;s2=cos(200*pi*t)+2*cos(204*pi*t);Y2=fft(s2,NFFT)/L2;f2=Fs/2*linspace(0,1,NFFT/2+1);figure;plot(f2,2*abs(Y2(1:NFFT/2+1)));xlabel('频率');设输入信号为tts204cos2200cos，采样频率为400Hz。可见，信号含有2个频率分量，分别为100Hz和102Hz，则要求的频率分辨率为2Hz，可算得N为202，即至少采样202个点，才不会发生栅栏现象。首先取采样点数为90，对其进行FFT，得到的图形如图1-4所示。图1-4FFT的栅栏效应-8-可见，虽然采样频率满足香农采样定理，但由于采样点数较少，频率为20Hz和20.5Hz的信号无法分辨。下面采用补零的方法，将N补满1024，得到的FFT图形如图1-5所示。图1-5信号末端补零改善栅栏效应由上图可以看出，补零对分辨率没有影响，但是对频谱起到了平滑的作用。同时，由于增加了大量零值，使得信号的平均幅频响应大大减小。然后采用增加采样点的方法，将采样点数增加到1024，FFT后的图形如图1-6所示。图1-6增加采样点数改善栅栏效应由上图可以看出，增加采样点数后就能分辨出100Hz和102Hz的信号，有效地改善了栅栏效应。从时域看，这种方法相当于对信号进行整周期采样，实际中常