Thread Subject: PSD and sptool

Hi,

I am working with sptool for analyze the PSD (Power Spectral
Density) of signals.

Now, I need to implement the code to do the analysis in
automatic manner.

This is the code (x contains the time serie):

y = fft(y, 2048);
Pyy = y * conj(y);

Now, I have searched how to calculare normalization factor
(to render plot of my Pyy the same of sptool), but I haven't
found nothing.

Could you help me?

On Aug 4, 9:39 am, "Vincenzo Manzoni" <ervin...@tiscali.it> wrote:
> Hi,
>
> I am working with sptool for analyze the PSD (Power Spectral
> Density) of signals.
>
> Now, I need to implement the code to do the analysis in
> automatic manner.
>
> This is the code (x contains the time serie):
>
> y = fft(y, 2048);
> Pyy = y * conj(y);
>
> Now, I have searched how to calculare normalization factor
> (to render plot of my Pyy the same of sptool), but I haven't
> found nothing.
>
> Could you help me?

Google Groups

greg-heath power-spectrum-scaling

Hope this helps.

Greg

> greg-heath power-spectrum-scaling

In this link
(http://www.mathworks.com/support/tech-notes/1700/1702.html)
the magnitude of fft is scaled by NFFT, while you, in your
post, scale the fft by sqrt(NFFT).

What is the difference? Thanks you for your help...

Uhm Greg, another question: if I calculate the power of this
signal

Fs = 1024; % Sample frequency: 1 kHz
t=0:1/Fs:10;
x = sin(2*pi*10*t); % Sin with frequency = 10 Hz

with fft following the article which I was cited the post
before, I obtain an armonic at 10 Hz with power 0.5... that
is theoretically corret.

If I calculate the power of signal with sptool, I obtain an
armonic at 10 Hz with power 1.5625... which change if I
change NFFT.

I am extremely confused...

On Aug 5, 7:35 am, "Vincenzo Manzoni" <ervin...@tiscali.it> wrote:
> >greg-heathpower-spectrum-scaling
>
> In this link
> (http://www.mathworks.com/support/tech-notes/1700/1702.html)
> the magnitude of fft is scaled by NFFT, while you, in your
> post, scale the fft by sqrt(NFFT).
>
> What is the difference? Thanks you for your help...

The difference is my convention yields the relations

Pxav = sum(PSD)/N = sumsqr(x)/N


Hope this helps.

Greg

On Aug 5, 8:32 am, "Vincenzo Manzoni" <ervin...@tiscali.it> wrote:
> UhmGreg, another question: if I calculate thepowerof this
> signal
>
> Fs = 1024; % Sample frequency: 1 kHz

UMM ... 1KHz = 1000Hz

> t=0:1/Fs:10;

dt = 1/Fs = 0.0009765625

T = 10 + dt = 10.0009765625

N = T/dt = 10*Fs+1 = 10241

> x = sin(2*pi*10*t); % Sin with frequency = 10 Hz

f0 = 10
T0 = 1/f0 = 0.1

T/T0 = 100.009765625

> with fft following the article which I was cited the post
> before, I obtain an armonic

correction: a harmonic

> at 10 Hz with power 0.5... that is theoretically corret.

Show how you obtained that value.

> If I calculate the power of signal with sptool, I obtain an
> armonic (sic) at 10 Hz with power 1.5625... which change if I
> change NFFT.
>
> I am extremely confused...

There might be a difference in scaling conventions. However,
unless we know exactly what you did and exactly what sptool does,
we can't do anything but guess.

Hope this helps.

Greg

>> at 10 Hz with power 0.5... that is theoretically corret.
 
> Show how you obtained that value.

I know the area below the power density spectrum is the
signal variance. Then:

Fs = 1024;
t = 0:1/Fs:1;
x = sin(2*pi*t*10);

var(x) -> 0.5

If I calculate the PSD with normalization of the article
which I have cited,

sum(mx(:)) -> 0.5

My assumption is correct?

> However,
> unless we know exactly what you did and exactly what
sptool does,
> we can't do anything but guess.

In sptool, I import the same signal and I create the
spectrum with fft as choose (with the same number of sample
I use in when calculate the PSD without sptool).

If your purpose is to obtain an estimate of PSD using periodogram, then I'll
suggest to use pwelch or even better, the spectrum object instead.

For example, according to your post, the PSD can be simply calculated as

Hspec = spectrum.periodogram; % define psd estimator
Hopts = psdopts(Hspec, y); % construct psd
estimator options
set(Hopts,'NFFT',2048); % set NFFT
Hpsd = psd(Hspec,y); % calculate psd
plot(Hpsd); % plot psd


The tech note 1702 is a very simple example of calculating power spectrum
using FFT, it is not a power spectrum DENSITY. If you do want a power
spectrum, you can use spectrum object too. All you have to do is to replace
psd function with msspectrum function, i.e.,

Hpsd = msspectrum(Hspec,y);

There are two demos that you may be interested in, one of them is a brief
introduction of using spectrum object and the other is about how to measure
the power of a given signal.

http://www.mathworks.com/products/signal/demos.html?file=/products/demos/shipping/signal/spectralanalysisobjsdemo.html
http://www.mathworks.com/products/signal/demos.html?file=/products/demos/shipping/signal/deterministicsignalpower.html

Hope this helps,

Honglei

"Vincenzo Manzoni" <ervincio@tiscali.it> wrote in message
news:f91vj2$2td$1@fred.mathworks.com...
> Hi,
>
> I am working with sptool for analyze the PSD (Power Spectral
> Density) of signals.
>
> Now, I need to implement the code to do the analysis in
> automatic manner.
>
> This is the code (x contains the time serie):
>
> y = fft(y, 2048);
> Pyy = y * conj(y);
>
> Now, I have searched how to calculare normalization factor
> (to render plot of my Pyy the same of sptool), but I haven't
> found nothing.
>
> Could you help me?
>
>