【问题标题】:Matlab, FFT frequency range differences or are they the same?Matlab,FFT频率范围差异还是它们相同?
【发布时间】:2016-03-01 23:14:44
【问题描述】:

我试图了解 matlab 中的 FFT 是如何工作的,特别是如何定义频率范围来绘制它。碰巧我已经阅读了 matlab 帮助链接和这里的其他讨论,我认为(猜想)我对此感到困惑。 在 matlab 链接中: http://es.mathworks.com/help/matlab/math/fast-fourier-transform-fft.html 他们将这样的频率范围定义为:

f = (0:n-1)*(fs/n)

nfs 为:

n = 2^nextpow2(L); % Next power of 2 from length of signal x
fs = N/T; % N number of samples in x and T the total time of the recorded signal

但是,另一方面,在之前的帖子中Understanding Matlab FFT example (基于以前版本的matlab),得到的频率范围定义为:

f = fs/2*linspace(0,1,NFFT/2+1);

NFFT 与上述n 相同(信号长度x 的下一个2 的幂)。 那么,基于此,这些不同的向量(方程 1 和最终方程)如何相同? 如果你能看到,向量是不同的,因为前者有n 点,而后者有NFFT/2 点!事实上,(fs/n) 的因子与fs/2 不同。

【问题讨论】:

    标签: matlab fft dft


    【解决方案1】:

    那么,基于此,这些不同的向量(方程 1 和最终方程)如何相同?

    example in the documentation from Mathworks 绘制 FFT 的整个 n 点输出。这涵盖了从 0 到接近 fs(确切地说是 (n-1)/n * fs)的频率。然后他们进行以下观察(对 FFT 的实际输入有效):

    频率范围的前半部分(从 0 到奈奎斯特频率fs/2)足以识别数据中的分量频率,因为后半部分只是前半部分的反映。

    other post you refer to 只是选择不显示多余的后半部分。然后它使用一半的点数,也覆盖一半的频率范围。

    其实因子(fs/n)和fs/2是不一样的。

    也许理解它的最简单方法是比较两个表达式的输出,以获得一些小的值n,如n=8 和设置fs=1(因为fs 将两个表达式相乘)。一方面,第一个表达式[0:n-1]*(fs/n) 的输出将是:

    0.000  0.125  0.250  0.500  0.625  0.750  0.875
    

    fs/2*linspace(0,1,n/2+1) 的输出将是:

    0.000  0.125  0.250  0.500
    

    如您所见,频率集与奈奎斯特频率 fs/2 完全相同。

    【讨论】:

    • 平面和简单。感谢您的解释 SleuthEye
    【解决方案2】:

    混淆可能是由于您引用的两个示例以不同方式绘制 fft 的结果。请参阅下面的代码以获取此说明中的引用。

    在第一个示例中,绘图是频率范围内的功率谱(周期图)。请注意,在第一个图中,周期图的中心不是 0,这意味着频率范围似乎是奈奎斯特采样频率的两倍。如 mathworks 链接中所述,通常的做法是将周期图居中于 0 以避免这种混淆(图 2)。

    对于第二个示例,采用相同的参数,原始图是傅立叶谱的幅度,其归一化与第一个示例不同(图 3)。使用 Matlab 的全频率排序的语法(如代码中所述),将这个看似不同的 fft 结果转换为示例 1 的结果是微不足道的;图 4 中复制了以 0 为中心的周期图的相同结果。

    因此,为了具体回答您的问题,两种情况下的频率范围相同,最大频率等于奈奎斯特采样频率,如下所示:

    f = fs/2*linspace(0,1,NFFT/2+1);
    

    理解 dfft 如何工作(也在 Matlab 中)的关键是要了解您只是将离散数据集投影到傅立叶空间中,其中 matlab 中 fft() 函数返回的是每个频率分量的展开和系数的顺序由以下方式给出(在 Matlab 中如示例 2 所示):

    f = [f(1:end-1) -fliplr(f(1,2:end))];
    

    有关更多详细信息,请参阅 DFT 上的 Wikipedia 页面: https://en.wikipedia.org/wiki/Discrete_Fourier_transform

    将 fft 省略为 2 的幂参数作为

    y = fft(x).
    

    在这种情况下,您只会看到 y 中的几个非零分量对应于输入信号的确切系数。 mathworks 页面声称以下是使用或不使用此长度的动机:

    “使用 2 的幂作为变换长度可以优化 FFT 算法,尽管在实践中,使用 n = m 的执行时间通常差别不大。”

    %% First example:
    % http://www.mathworks.com/help/matlab/math/fast-fourier-transform-fft.html
    
    fs = 10;                                % Sample frequency (Hz)
    t = 0:1/fs:10-1/fs;                      % 10 sec sample
    x = (1.3)*sin(2*pi*15*t) ...             % 15 Hz component
      + (1.7)*sin(2*pi*40*(t-2));            % 40 Hz component 
    % Removed the noise
    
    m = length(x);          % Window length
    n = pow2(nextpow2(m));  % Transform length
    y = fft(x,n);           % DFT
    f = (0:n-1)*(fs/n);     % Frequency range
    power = y.*conj(y)/n;   % Power of the DFT
    
    subplot(2,2,1)
    plot(f,power,'-o')
    xlabel('Frequency (Hz)')
    ylabel('Power')
    title('{\bf Periodogram}')
    
    y0 = fftshift(y);          % Rearrange y values
    f0 = (-n/2:n/2-1)*(fs/n);  % 0-centered frequency range
    power0 = y0.*conj(y0)/n;   % 0-centered power
    
    subplot(2,2,2)
    plot(f0,power0,'-o')
    % plot(f0,sqrt_power0,'-o')
    xlabel('Frequency (Hz)')
    ylabel('Power')
    title('{\bf 0-Centered Periodogram} Ex. 1')
    
    %% Second example:
    % http://stackoverflow.com/questions/10758315/understanding-matlab-fft-example
    
    % Let's redefine the parameters for consistency between the two examples
    
    Fs = fs;                      % Sampling frequency
    % T = 1/Fs;                   % Sample time (not required)
    L = m;                        % Length of signal
    % t = (0:L-1)*T;              % Time vector (as above)
    % % Sum of a 3 Hz sinusoid and a 2 Hz sinusoid
    % x = 0.7*sin(2*pi*3*t) + sin(2*pi*2*t); %(as above)
    
    NFFT = 2^nextpow2(L); % Next power of 2 from length of y
                          % NFFT == n (from above)
    Y = fft(x,NFFT)/L;
    f = Fs/2*linspace(0,1,NFFT/2+1);
    
    % Plot single-sided amplitude spectrum.
    subplot(2,2,3)
    plot(f,2*abs(Y(1:NFFT/2+1)),'-o') 
    title('Single-Sided Amplitude Spectrum of y(t)')
    xlabel('Frequency (Hz)')
    ylabel('|Y(f)|')
    
    
    % Get the 0-Centered Periodogram using the parameters of the second example
    f = [f(1:end-1) -fliplr(f(1,2:end))]; % This is the frequency ordering used
                                          % by the full fft in Matlab
    
    power = (Y*L).*conj(Y*L)/NFFT;
    
    % Rearrange for nicer plot
    ToPlot = [f; power]; [~,ind] = sort(f); 
    ToPlot = ToPlot(:,ind);
    
    subplot(2,2,4)
    plot(ToPlot(1,:),ToPlot(2,:),'-o')
    xlabel('Frequency (Hz)')
    ylabel('Power')
    title('{\bf 0-Centered Periodogram} Ex. 2')
    

    【讨论】:

    • 亲爱的 PZwan,谢谢。实际上,这是一个完整的答案和非常有用的代码,可以帮助您理解 DFT 和 matlab 的 fft 的解释和一般内容。
    • 亲爱的 Gohann,这是我的荣幸。很高兴它对您有所帮助。
    猜你喜欢
    • 1970-01-01
    • 1970-01-01
    • 2014-05-04
    • 1970-01-01
    • 2022-07-23
    • 2011-12-16
    • 2013-12-04
    • 1970-01-01
    • 2011-05-04
    相关资源
    最近更新 更多