Re: FFT spectrum analyzer question

[Macy]

arrrggg, forgot to add the frequency!
make that
sig=sin(2*pi()*14e3*t); 
then you'll get 1/2

obviously just part of a 1Hz waveform won't have much energy. should have 
caught it with the non-zero DC value.

for your batlett's example you need to set the frequency to something integer 
of 32
So use 32kHz, or 3.2kHz


rewritten:
Ss=192000;
dt=1/Ss;
fc=3200;
for N = [192000 96000 48000 24000 12000 6000]
  max_sum = 0;
  sections = 0;
  for offset = [1:N:Ss]
    t=[offset:offset+N-1]*dt;

    sig=sin(2*pi*fc*t);

    sig_f=fft(sig)/N;
    max_sum = max_sum + max(abs(sig_f));
    sections = sections + 1;
    max_avg = max_sum / sections;
    printf('[%6d %6d] max: %f avg: %f\n', N, i, max(abs(sig_f)), max_avg)
  endfor
  printf('\n')
endfor

no idea what you're doing there. 


--- address@hidden wrote:

From: Bolin Hsu <address@hidden>
To: address@hidden
Subject: Re: FFT spectrum analyzer question
Date: Wed, 14 May 2014 10:33:27 -0700

Macy,

Thanks for your reply. Strangely I get 0.15578 instead of 1/2 from the
following:
Ss=192000;
dt=1/Ss;
N=9600;
t=[0:N-1]*dt;
sig=sin(2*pi()*t);
sig_f=fft(sig)/N;
max(abs(sig_f))

Anyway, I was studying Bartlett's method for power spectrum estimation.
That's why I fixed a run length, split the signal into multiple sections,
and apply fft on these sections. If I apply the same idea on your sample
code, I get the following:

Ss=192000;
dt=1/Ss;
for N = [192000 96000 48000 24000 12000 6000]
  max_sum = 0;
  sections = 0;
  for offset = [1:N:Ss]
    t=[offset:offset+N-1]*dt;
    sig=sin(2*pi*t);
    sig_f=fft(sig)/N;
    max_sum = max_sum + max(abs(sig_f));
    sections = sections + 1;
    max_avg = max_sum / sections;
    printf('[%6d %6d] max: %f avg: %f\n', N, i, max(abs(sig_f)), max_avg)
  endfor
  printf('\n')
endfor

The test shows that max(abs(sig_f)) is 1/2 only when N is equal to Ss. And
the average of the max over all sections doesn't converge to 1/2, either.
Is this inherent in Bartlett's method or caused by bugs in my
implementation?

Thanks,
Bolin


On Wed, May 14, 2014 at 8:31 AM, Macy <address@hidden> wrote:

> lot's of 'errors'
>
> given signal length as len [I prefer labeling run length as N]
> sig_f=fft(sig)/len;
> you do that you'll get closer to constant.
>
> next be sure to make the length 'correct' AND match sample lengths, must
> be ineteger related else you will have to use 'windowing' as partial values
> get shoved into adjacent bins
>
> your frequency, sampling rate, and packet lengths 'fight' each other.
>
> also, for time use one less than the length you set
> to understand try this:
> Ss=192000;
> dt=1/Ss;
> N=9600;
> t=[0:N-1]*dt;
> sig=sin(2*pi()*t);
> sig_f=fft(sig)/N;
> max(abs(sig_f))
> will always come up 1/2
> since the energy gets split during your fft.
> if you plot
> plot(20*log10(abs(sig_f)));
> you'll notice two halves to your spectrum which can get a bit confusing.
> just remember that your time waveform is ALWAYS real so the two havlesof
> the resulting spectrum will ALWAYS be related as conjugates, BUT the right
> half is inverted, coming down. see what i mean about confusing?
>
> now do again using
> t=[0:N]; and the plots will suck.
>
>
> I wrote a way to bounce back and forth between time and frequency domains:
> two custom functions: time2freq and freq2time
> [sig_f,ford]=time2freq(sig);
> so when you plot
> plot(ford/dt,20*log10(abs(sig_f)));
> you get an ACTUAL spectrum plot.
>
> if you want log-log
> plot(log10(ford(2:end)/dt/fc),20*log10(abs(sig_f(2:end))));
> why not index 1? because that is DC and it plays havoc with your log
> functions.
>
> The advanatage of this form is that you can do all kinds of spectral
> manipulations like run your signal through a filter by multiplying the two
> spectrums, then convert back to time and see the effect on the time
> waveform.
>
> sig_t=freq2time(sig_f.*filterspectrum);
> and get waveform back.
>
> ps: always use even length runs, gets really confusing if you don't.
>
>
>
>
>
>
> --- address@hidden wrote:
>
> From: Bolin Hsu <address@hidden>
> To: address@hidden
> Subject: FFT spectrum analyzer question
> Date: Tue, 13 May 2014 17:09:59 -0700
>
> I saw the thread about FFT spectrum analyzer at
> http://octave.1599824.n4.nabble.com/FFT-Spectrum-Analyzer-td4630647.html,
> and gave it a try. From http://en.wikipedia.org/wiki/Bartlett%27s_method,
> I
> got the idea that the energy is abs(fft(x))^2/M. So I ended up with the
> program below. The signal was sine wave with peak amplitude 1, so I thought
> I would get the RMS amplitude 1/sqrt(2). But the result I got was much
> smaller. And the result varies with different FFT size. What am I missing
> in this program?
>
> Fs=44100;
> w=14000;
> t=0:1/Fs:1;
> x=sin(2*pi*w*t);
> len=length(x);
> for N = [128 256 512 1024 2048 4096 len]
>   for i = [1:N:len]
>     if (i+N-1 <= len)
>       xdft=fft(x(i:i+N-1));
>       xdft=xdft(1:N/2+1);
>       psdx=(1/(Fs*N)).*abs(xdft).^2;
>       psdx(2:end-1)=2*psdx(2:end-1);
>       freq=0:Fs/N:Fs/2;
>       printf('[%5d, %5d] %f\n', N, i, max(psdx))
>     endif
>   endfor
> endfor
>
>
> Thanks,
> Bolin
>
>
> _______________________________________________
> Help-octave mailing list
> address@hidden
> https://lists.gnu.org/mailman/listinfo/help-octave
>
>
>


_______________________________________________
Help-octave mailing list
address@hidden
https://lists.gnu.org/mailman/listinfo/help-octave