matlab实现DFT

基于DFTMATLAB的实现一、实验目的1 .掌握DFT函数的使用方法。2 .基于DFT的信号检测和频谱分析。3 .了解信号截取长度对频谱分析的影响。二、实验内容1 .使用DFT计算信号的功率谱。实验步骤：t=0:0.001:0.6；x=sin (2* pi * 50 * t ) sin (2* pi * 120 * t ) rann (1，长度(t )；Y=dft(x，512 )P=Y.*conj(Y)/512；f=1000*(0愚人节55)/512；plot(f，P(1:256 ) )2 .进行信号检测。 分析对应于信号频谱的频率轴的数字频率与频率的关系。 对模拟信号进行采样，求出n点DFT的振幅频谱。实验步骤：subplot (2，2，1 )N=45； n=0:N-1； t=0.01*n；q=n*2*pi/N；x=2*sin(4*pi*t) 5*cos(8*pi*t )y=dft(x，n )plot(q，abs(y ) )； title(DFT N=45 )subplot (2，2，2 )N=50； n=0:N-1； t=0.01*n； q=n*2*pi/N；x=2*sin(4*pi*t) 5*cos(8*pi*t )y=dft(x，n )plot(q，abs(y ) )； title(DFT N=50 )subplot (2，2，3 )N=55； n=0:N-1； t=0.01*n；q=n*2*pi/N；x=2*sin(4*pi*t) 5*cos(8*pi*t )y=dft(x，n )plot(q，abs(y ) )； title(DFT N=55 )subplot (2，2，4 )N=60； n=0:N-1； t=0.01*n；q=n*2*pi/N；x=2*sin(4*pi*t) 5*cos(8*pi*t )y=dft(x，n )plot(q，abs(y ) )； title(DFT N=60 )针对3.2，将限幅长度和DFT点数进一步增加，例如将n增加至256，观察信号频谱的变化并分析这一变化的原因。 在剪辑长度不变化的条件下改变采样频率，观察信号频谱的变化，分析这种变化的原因。将n扩展到256的步骤：N=256； n=0:N-1； t=0.01*n；q=n*2*pi/N；x=2*sin(4*pi*t) 5*cos(8*pi*t )y=dft(x，n )plot(q，abs(y ) )； title(DFT N=256 )分析原因：由于如果T=0.01s，第一序列的周期为100，第二序列的周期为50，因此在采样点数小于100的情况下，频率分辨率不足，无法区分两个信号。 当采样点足够多(256 )时，频率分辨率增加，从而能够区分两个频率的信号。在将采样间隔设为T=0.1s的情况下，n保持为45的程序：N=45； n=0:N-1； t=0.1*n；q=n*2*pi/N；x=2*sin(4*pi*t) 5*cos(8*pi*t )y=dft(x，n )plot(q，abs(y ) )； title(DFT N=45 )分析原因：如果T=0. 1s，则第一序列的周期为10，第二序列的周期为5，因此在采样点数为45的情况下能够区分两个信号。参数是相同的，n取64，将噪声w(t )与信号相加。figure(2)subplot (2，1，1 )N=64； n=0:N-1； t=0.01*n；q=n*2*pi/N；x=2* sin (4* pi * t )5* cos (8* pi * t ) y=DFT (x，n )plot(q，abs(y ) )； title(DFT N=64 )subplot (2，1，2 )N=64； n=0:N-1； t=0.01*n；q=n*2*pi/N；x=2* sin (4* pi * t )5* cos (8* pi * t ) 0.8 * rann (1，n) y=dft(x，n )plot(q，abs(y ) )； title(DFT N=64(with noise ) )可以看出这个噪声不影响信号检测。对于4.3，将噪声增大到2*randn(1，n )和8*randn(1，n )，并且以不同噪声绘制时域波形和频谱进行比较。subplot (2，1，1 )N=64； n=0:N-1； t=0.01*n；q=n*2*pi/N；x=2* sin (4* pi * t )5* cos (8* pi * t )2* rann (1，n) y=dft(x，n )plot(q，abs(y ) )； title(DFT N=64(with noise2) )subplot (2，1，2 )N=64； n=0:N-1； t=0.01*n；q=n*2*pi/N；x=2* sin (4* pi * t )5* cos (8* pi * t )8* rann (1，n) y=dft(x，n )plot(q，abs(y ) )； title(DFT N=64(with noise8) )subplot (3，2，1 )N=64； n=0:N-1； t=0.01*n；q=n*2*pi/N；x=2* sin (4* pi * t )5* cos (8* pi * t )0.8* rann (1，n )plot(x) title (噪声为0.8*w的信号)y=dft(x，n )subplot (3，2，2 )plot(q，abs(y ) )； title (噪声为0.8*w时的频谱)subplot (3，2，3 )N=64； n=0:N-1； t=0.01*n；q=n*2*pi/N；x=2* sin (4* pi * t )5* cos (8* pi * t )2* rann (1，n )plot(x) title (噪声为2*w时的信号)y=dft(x，n )subplot (3，2，4 )plot(q，abs(y ) )； title (噪声为2*w时的频谱)subplot (3，2，5 )N=64； n=0:N-1； t=0.01*n；q=n*2*pi/N；x=2* sin (4* pi * t )5* cos (8* pi * t )8* rann (1，n )plot(x) title (噪声为8*w时的信号)y=dft(x，n )subplot (3，2，6 )plot(q，abs(y ) )； title (噪声为8*w时的频谱)实验分析：噪声较小时，不影响信号检测，但噪声较大时，原信号的频率成分消失，继续加大噪声，频谱看起来杂乱。5 .使用一个n点DFT计算两个长度为n的实序列n点离散傅立叶变换，并将结果与直接使用两个n点DFT的结果进行比较。x=1 2 3 4 5 6；y=6 5 4 3 2 1；a，b=dft_2(x，y )a=Columns 1 through 321.0000-3.00005.1962 I-3.00001.7321 Icolumns 4路径6- 3.0000-3.0000-1.7321 I-3.0000-5.1962 Ib=Columns 1 through 321.00003.0000-5.1962 I3. 0000-1.7321 Icolumns 4路径63.00003.00001.7321 I3. 00005.1962 I函数文件如下：函数 y 1，y2=dft_2(a，b )N=length(a )x=zeros(1，n )x=a j*b；X=dft(x，n )X0=conj(fliplr(X ) )X0=X0(N) X0(1:N-1)；y1=(X X0)./2；y2=(X-X0)./2./j；要执行直接计算：dft(x )ans=Columns 1 through 321.0000-3.00005.1962 I-3.00001.7321 Icolumns 4路径6- 3.0000-3.0000-1.7321 I-3.0000-5.1962 Idft(y )ans=Columns 1 through 321.00003.0000-5.1962 I3. 0000-1.7321 Icolumns 4路径63.00003.00001.7321 I3. 00005.1962 I6 .将DFT与DFT的计算时间进行比较。 (定时函数tic，toc )n分别为256、512、1024、2048、4096、方案包括：N=256；N=4096；x=randn(1，n )ticy=dft(x，n )tocticz=dft(x )tocn=256:elapsedtimeis 0.172000 seconds。Elapsed time is 0.015000 seconds。n=512:elapsedtimeis 0.687000第二代。Elapsed time is 0.000000 seconds。n=1024:elapsedtimeis 3.031000 seconds。Elapsed time is 0.047000 seconds。n=2048:elapsedtimeis 13.3753000 seconds。Elapsed time is 0.063000 seconds。n=4096:elapsedtimeis 59.250000 seconds。Elapsed time is 0.125000 seconds。7 .对给定声音信号进行频谱分析，写出采样频率，描绘声音信号的波形和频谱，分析声音信号的频率分布特征。(1)绘制时域波形，对整个声音序列施加DFTx，fs=wav读取(c :ai1. wav )subplot (2，1，1 )N=length(x )n=0:N-1；plot(n，x )xlabel(n )ylabel(x )title (时域波形)subplot (2，1，2 )N=length(x )n=0:N-1；t=0.01*n；q=n*2*pi/N；y=dft(x，n )plot(q，abs(y ) )；xlabel(n )ylabel(ai1)title(DFT )(2)分别求出与k=300，3500对应的信号频率(Hz )。x，fs=wav读取(c :ai1. wav )N=length(x )n=0:N-1；t=n*(1/fs )q=n*2*pi/N；n1=300；q1=n1*2*pi/N；f1=q1*fs/(2*pi )n2=3500；q2=n2*2*pi/N；f2=q2*fs/(2*pi )fs=16000(3)从声音中提取一部分声音(256点)，进行DFT，得到频谱，描绘时域的波形和频谱。 分别求出与k=5、60情况对应的信号频率(Hz )方案包括：x，fs=wav读取(c :ai1. wav )subplot (2，1，1 )N=256；n=0:N-1；x=x(1:256 )plot(n，x )xlabel(n )ylabel(x )title(256点时域波形)subplot (2，1，2 )N=256；n=0:N-1；t=0.01*n；q=n*2*pi/N；x=x(1:256 )y=dft(x，n )plot(q，abs(y )