Matlab绘制对数啁啾的对数

Question

我创建了一个与matlab帮助页面完全相同的对数啁啾.

t = 0:0.001:10;      % 10 seconds @ 1kHz sample rate
fo = 10; f1 = 400;   % Start at 10Hz, go up to 400Hz
X = chirp(t,fo,10,f1,'logarithmic');
figure(2);
spectrogram(X,256,200,256,1000,'yaxis');

然后我使用以下代码将其带到频域,该代码适用于我的其他应用程序.

fft_prep = fftshift(fft(X));
fft_mag = abs(fft_prep);
pos_fft = fft_mag(1:ceil(length(fft_mag)/2));
db_fft = 20*log10(pos_fft);
figure(1);
plot(db_fft);

我很惊讶地发现以下图表看起来令人兴奋1kHz-5kHz:

我对matlab中的啁啾功能并不熟悉,并且想知道是否有人看到了我遗漏的明显东西.欢迎任何其他指针.

Answer 1

chirp功能没有错......

您只需要根据频率值绘制db_fft,而不是矢量索引=).

plot(linspace(fo,f1,length(db_fft)), db_fft);

我还测试了使用其他FFT方法计算信号的FFT,它们也指示了0到400 Hz之间的范围.

更新:

IMO,我发现在视觉上更容易不以dB或功率(周期图)绘图.这是一个很好的例子和我计算时域信号FFT的goto方法:mathworks.se/help/matlab/ref/fft.html

响应:

经过一番思考,我同意我的答案不正确,但不是因为你说的原因.频域中的x轴不应该达到啁啾的实际长度(或一半,或dubbel或类似的东西).频域中的x轴应该达到信号采样率的一半(Fs/2),并且您有义务确保采样频率是您希望/希望的最大频率的两倍.解决.

换句话说,假设您的FFT与时域信号的长度相同/两倍/一半是不正确的,因为我们可以选择任意数量的频率区来表示FFT,最佳实践是长度= N ^ 2( 2)的功率,用于快速计算.想一想,为什么在计算FFT时甚至需要知道时间值？你不!您只需要采样频率(应设置为Fs = 1000 btw,而不是Fs = 0.001).

我上面的回答是不正确的,它应该是:

plot(linspace(0, Fs/2, length(db_fft)), db_fft)

而不是Fs/2,你有写长度(t)/(2*Tfinal).它(几乎)与Fs/2的值相同,但它不是正确的方式=).

这是我的goto FFT方法(值不是以dB为单位).

function [X,f] = myfft(x,Fs,norm)
    % usage: [X, f] = myfft(x,Fs,norm);
    %        figure(); plot(f,X)
    % norm: 'true' normalizes max(amplitude(fft))=1, default=false.
    if nargin==2
        norm=false;
    end
    L = length(x); NFFT = 2^nextpow2(L);
    f = Fs/2*linspace(0,1,NFFT/2+1);
    %f =0:(Fs/NFFT):Fs/2;
    X = fft(x,NFFT)/L; X = 2*abs(X(1:NFFT/2+1));
    if norm==true; X = X/max(abs(X)); end
end

这是[Xfft,f] = myfft(X,Fs)的结果图; 情节(F,Xfft); 注意,根据NyQuist定理,返回频率bin向量具有max(f)= Fs/2(任何比Fs/2更高的频率都不能被解析).