Re: MATLAB code for Power Spectral Density
- From: "Glennaebad" <fitlike@xxxxxxxxxxxxxxxxx>
- Date: Sun, 11 Dec 2005 11:45:46 +1300
"Mike Yarwood" <mpyarwood@xxxxxxxxxxxxxxx> wrote in message
news:dnet0d$d1k$1@xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
>
> "john" <johns@xxxxxxxxxx> wrote in message
> news:1134210735.813800.60750@xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
> >
> > ngeva0 wrote:
> >> Can anyone please check if my code is correct? SINE512 is just a file
> >> name.
> >> This file has 512 data points and it's a sine wave. Thank you!!
> >>
> >> %Sampling frequency
> >> Fs=50000;
> >>
> >> %# of samples in the data
> >> datasize=size(SINE512);
> >> numsample=datasize(1);
> >>
> >> %Windowing
> >> H=hann(numsample);
> >> W=H.*(SINE512(:,2));
> >>
> >> %Fourier Transform
> >> FFTX=fft(W,numsample);
> >>
> >> %Power: magnitude^2
> >> X=FFTX(1:floor(numsample/2)).*conj(FFTX(1:floor(numsample/2)));
> >>
> >> %Bandwidth
> >> BW=1.5*Fs*numsample;
> >>
> >> %PSD=magnitude^2/bandwidth
> >> PSD=X/BW;
> >>
> >> %Computing the corresponding frequency values
> >> Omega=Fs*(0:numsample-1)/numsample;
> >> Omega=Omega(1:floor(numsample/2));
> >>
> >> %Plot PSD
> >> H=plot(Omega,PSD);
> >> set(H,'Color','BLACK');
> >
> > I am happy to report that your code behaves exactly as programmed.
> > Congratulations.
> >
> But you may want to look at :
> http://www.mathworks.com/access/helpdesk/help/techdoc/ref/fft.html#998360
>
> where you can see that they normalise power in the bin by dividing by the
> number of points in the DFT and get powers of about 2*pi for their unit
> peak amplitude sinewaves.
>
> A 1Volt peak amplitude sinewave across a 1 ohm resistor dissipates 1/2
watt
> and your PSD display may either register a 1 volt (peak) input sinewave as
0
> dBV (w.r.t. 1V pk), -3.01 dBV (w.r.t. 1 V r.m.s) or -6.02 dBV (w.r.t. 1
> Vpk-pk but this is very unlikely). So put a known amplitude sinewave into
> your instrument (with as little extra noise as possible) and see what it
> displays.
>
> None of the above? has it scaled by a further 10*log10(512/50000) ~ -19.9
> dB ?
>
> Now you know what the instrument is doing you can decide what your matlab
> routine should do.
>
>
> Best of Luck - Mike
>
>
Remember its the AREA under the PSD that gives power - not the peak values.
For a sine wave on a bin this means you need to multiply by the width of the
bin.
Glen
.
- References:
- MATLAB code for Power Spectral Density
- From: ngeva0
- Re: MATLAB code for Power Spectral Density
- From: john
- Re: MATLAB code for Power Spectral Density
- From: Mike Yarwood
- MATLAB code for Power Spectral Density
- Prev by Date: Re: Extracting a signal from noise?
- Next by Date: Re: Audio Compressor/AGC in C or C++
- Previous by thread: Re: MATLAB code for Power Spectral Density
- Next by thread: Re: MATLAB code for Power Spectral Density
- Index(es):
Relevant Pages
|
