1 个投票

Hi everybody! I'm trying to use dB scale in FFT. I've found the example http://www.dsplib.ru/content/winex/winex.html
I'm trying to get the same results in MATLAB. But I can't. I have different amplitude levels. What's wrong? So, I do:
Fs=1000;
T=1/Fs;
N=1024;
t=(0:N-1)*T;
f=(0:Fs/(N-1):Fs);
A0=1;
A1=0.003;
A2=10e-4;
f0=220;
f1=240;
f2=230;
s=A0*cos(2*pi*f0*t)+A1*cos(2*pi*f1*t)+A2*cos(2*pi*f2*t);
plot(f,mag2db(abs(fft(s))));

请先登录，再回答此问题。

 采纳的回答

David Goodmanson
David Goodmanson 2017-1-24

1 个投票

hi Alexander I assume you are looking for a couple of peaks at -6 dB. Otherwise this answer will not be of help. But to get -6 dB you need to divide the fft result by N before taking mag2dB. Before getting into that, I believe your frequency grid with the factor of N-1 is not quite right. It should have a factor of N and should not end at Fs:
f = (0:N-1)*Fs/N;
That change only scales the frequency axis though, and puts the peaks in the right spot. It does not affect the peak value. You will still not get -6 dB (right now it's more like -7) because f0 = 220 does not represent an exact periodic frequency to the fft. To get a sharp peak at -6 dB, the frequency must be a multiple of Fs/N = 1000/1024. You could replace 220 by 225*1000/1024 which is pretty close, but I think the best way is to just use N = 1000. The fft is blazingly fast anyway, so unless you are in some production situation doing a ten million of these, N = 2^n is not really necessary.

7 个评论
显示 5更早的评论 隐藏 5更早的评论

Hi, David Goodmanson! Thanks for your answer. About frequency grid - maybe you're right. But I'm looking for a couple of peaks at 0 dB(!), like in the example http://www.dsplib.ru/content/winex/winex.html Now I have 53 dB. Try to run my code and you'll see it.
David Goodmanson
David Goodmanson 2017-1-24
编辑：David Goodmanson 2017-1-24
Hi Alexandr (got your name right this time), I did run your code yesterday. Thank you for providing runnable code by the way, a lot of people posting questions here don't. I got +53 dB like you mentioned.
Dividing by N is necessary because of the way Matlab scales their fft. I would not be surprised that if they had it to do all over again they might include that factor in the fft (in which case ifft would also change). After dividing by N, the result is 53 - 20*log10(1024) = -7 dB. After going to a frequency that is a multiple of Fs/N, I got a value at the peak equal to mag2db(1/2) good to many decimal places.
As to the -6 dB, let w = (2pi f) and let's say you transform a cosine function with a phase angle g. Then
cos(wt+g) = (1/2) e^(iwt) e^(ig) + (1/2) e^(-iwt) e^(-ig)
If you report out all frequencies you get two complex peaks of magnitude 1/2, i.e. -6 dB. That's what your fft is trying to do. You are never going to get there, though, if you use a frequency that is not a multiple of Fs/N.
For the transform of a real function both peaks are the same size and contain the same information. Sometimes people choose to report out only the positive frequencies, throw away the phase information and double the size of the peak. So you get 0 dB. That's what they did in the link you provided.
Positive frequencies? Or, perhaps, you mean frequencies from 0 to Fs/2?
Fs=1000;
T=1/Fs;
N=1024;
t=(0:N-1)*T;
%f=(0:Fs/(N-1):Fs);
f = (0:N-1)*Fs/N;
A0=1;
A1=0.003;
A2=10e-4;
f0=220;
f1=240;
f2=230;
s=A0*cos(2*pi*f0*t)+A1*cos(2*pi*f1*t)+A2*cos(2*pi*f2*t);
val=2*abs(fft(s)/N);
f=f(1:length(f)/2);
val=val(1:length(f));
plot(f,mag2db(val));
Is it correct?
David Goodmanson
David Goodmanson 2017-2-4
编辑：David Goodmanson 2017-2-4
If you do the same thing as above and just change to N = 1000 rather than N = 1024 you will see a very large difference, three sharply defined peaks right where they should be, with the correct dB values.
Yes, by positive frequencies I meant those from 0 to Fs/2 (including 0 with the strictly positive ones). I think of frequencies from Fs/2 to Fs as negative since, due to aliasing or however you prefer to look at it, that range can be shifted down by Fs to produce the frequency range -Fs/2 to 0, with no change in the results.
Yes, the difference is great! Why? Could you explain, please! But not all values are correct, the third must be -80 dB but it is only -60 dB.
David Goodmanson
David Goodmanson 2017-2-7
For the fft to produce a sharp spike, there must be an exact integral number of oscillations in the time domain function. Otherwise what is repeating is not (say) a cosine wave but rather a cosine wave whose last oscillation is truncated. So (in your notation) f*t must equal an integer m. But for an N-point fft, t = NT and Fs = 1/T, so f = m Fs/N. and f must be a multiple of Fs/N. That's not the case for Fs = 1000, N = 1024, f = 220 but it is when you change N to 1000.
Also, you have A2 = 10e-4 not 1e-4, so -60 dB.
"Also, you have A2 = 10e-4 not 1e-4, so -60 dB."
Yes, it's my mistake. Now it's all right! Thank you very much!

请先登录，再进行评论。

更多回答（0 个）

类别

在 帮助中心 和 File Exchange 中查找有关 Parametric Spectral Estimation 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by