convert periodogram PSD estimation to signal amplitude
in [Matlab]
|
From: Dean on 31 Jul 2010 21:31 Hi guys, I've a signal consists of a frequency 3/32 Hz, sampled at 1 time per second (1 Hz). Since the signal amplitude is 1 volt, I am expecting the signal power is 0.5 v^2. Although I can obtain the frequency estimation perfectly, I don't know how to convert Periodogram PSD to the signal amplitude. If anyone knows the answer, please kindly let me know. Thank you all. Dean ____ Code as follows ________ clear all; close all; clc; N = 10000; Nfft = 2^nextpow2(N); A1 = 1; f1 = 3/32; fs = 1; t = (0:N-1)/fs; x = A1*cos(2*pi*f1*t); win = rectwin(N)'; [Pxx,f] = periodogram(x,win,Nfft,fs); figure; plot(f,Pxx,'k-o'); xlabel('Frequency, Hz'); ylabel('PSD, v^2/Hz'); grid on; xlim([0.05,0.15]);
From: Godzilla on 31 Jul 2010 23:37 "Dean " <wj.antispam(a)gmail.com> wrote in message <i32iop$m5r$1(a)fred.mathworks.com>... > Hi guys, > > I've a signal consists of a frequency 3/32 Hz, sampled at 1 time per second (1 Hz). Since the signal amplitude is 1 volt, I am expecting the signal power is 0.5 v^2. > > Although I can obtain the frequency estimation perfectly, I don't know how to convert Periodogram PSD to the signal amplitude. > > If anyone knows the answer, please kindly let me know. > > Thank you all. > > Dean > > ____ Code as follows ________ > clear all; close all; clc; > N = 10000; > Nfft = 2^nextpow2(N); > A1 = 1; > f1 = 3/32; > fs = 1; t = (0:N-1)/fs; > x = A1*cos(2*pi*f1*t); > win = rectwin(N)'; > [Pxx,f] = periodogram(x,win,Nfft,fs); > figure; plot(f,Pxx,'k-o'); > xlabel('Frequency, Hz'); ylabel('PSD, v^2/Hz'); > grid on; xlim([0.05,0.15]); My approach would be to 'calibrate' my code by using a signal with known amplitude and frequency. The output from the periodogram will then give you the scaling factor to convert periodgram units to real amplitude.
From: Wayne King on 1 Aug 2010 06:30 "Dean " <wj.antispam(a)gmail.com> wrote in message <i32iop$m5r$1(a)fred.mathworks.com>... > Hi guys, > > I've a signal consists of a frequency 3/32 Hz, sampled at 1 time per second (1 Hz). Since the signal amplitude is 1 volt, I am expecting the signal power is 0.5 v^2. > > Although I can obtain the frequency estimation perfectly, I don't know how to convert Periodogram PSD to the signal amplitude. > > If anyone knows the answer, please kindly let me know. > > Thank you all. > > Dean > > ____ Code as follows ________ > clear all; close all; clc; > N = 10000; > Nfft = 2^nextpow2(N); > A1 = 1; > f1 = 3/32; > fs = 1; t = (0:N-1)/fs; > x = A1*cos(2*pi*f1*t); > win = rectwin(N)'; > [Pxx,f] = periodogram(x,win,Nfft,fs); > figure; plot(f,Pxx,'k-o'); > xlabel('Frequency, Hz'); ylabel('PSD, v^2/Hz'); > grid on; xlim([0.05,0.15]); Hi Dean, You need to scale the periodogram estimate correctly. N = 10000; A1 = 1; f1 = 3/32; fs = 1; t = (0:N-1)/fs; x = A1*cos(2*pi*f1*t); h = spectrum.periodogram('Hamming'); msest = msspectrum(h,x,'Fs',1); plot(msest.Frequencies,msest.Data); Hope that helps, Wayne
From: Dean on 1 Aug 2010 08:54 "Wayne King" <wmkingty(a)gmail.com> wrote in message <i33ibt$7th$1(a)fred.mathworks.com>... > "Dean " <wj.antispam(a)gmail.com> wrote in message <i32iop$m5r$1(a)fred.mathworks.com>... > > Hi guys, > > > > I've a signal consists of a frequency 3/32 Hz, sampled at 1 time per second (1 Hz). Since the signal amplitude is 1 volt, I am expecting the signal power is 0.5 v^2. > > > > Although I can obtain the frequency estimation perfectly, I don't know how to convert Periodogram PSD to the signal amplitude. > > > > If anyone knows the answer, please kindly let me know. > > > > Thank you all. > > > > Dean > > > > ____ Code as follows ________ > > clear all; close all; clc; > > N = 10000; > > Nfft = 2^nextpow2(N); > > A1 = 1; > > f1 = 3/32; > > fs = 1; t = (0:N-1)/fs; > > x = A1*cos(2*pi*f1*t); > > win = rectwin(N)'; > > [Pxx,f] = periodogram(x,win,Nfft,fs); > > figure; plot(f,Pxx,'k-o'); > > xlabel('Frequency, Hz'); ylabel('PSD, v^2/Hz'); > > grid on; xlim([0.05,0.15]); > > Hi Dean, You need to scale the periodogram estimate correctly. > > N = 10000; > A1 = 1; > f1 = 3/32; > fs = 1; t = (0:N-1)/fs; > x = A1*cos(2*pi*f1*t); > h = spectrum.periodogram('Hamming'); > msest = msspectrum(h,x,'Fs',1); > plot(msest.Frequencies,msest.Data); > > Hope that helps, > Wayne Hi Wayne, Thank you for your code, it works perfectly. unfortunately, I've generated many graphs using power spectral density rather than mean-square spectrum. The questions is how can I convert PSD to the signal power/amplitude? Or is it possible at all? Cheers, Dean
From: Wayne King on 1 Aug 2010 14:18
"Dean " <jiangwei0515(a)gmail.com> wrote in message <i33qpd$mmp$1(a)fred.mathworks.com>... > "Wayne King" <wmkingty(a)gmail.com> wrote in message <i33ibt$7th$1(a)fred.mathworks.com>... > > "Dean " <wj.antispam(a)gmail.com> wrote in message <i32iop$m5r$1(a)fred.mathworks.com>... > > > Hi guys, > > > > > > I've a signal consists of a frequency 3/32 Hz, sampled at 1 time per second (1 Hz). Since the signal amplitude is 1 volt, I am expecting the signal power is 0.5 v^2. > > > > > > Although I can obtain the frequency estimation perfectly, I don't know how to convert Periodogram PSD to the signal amplitude. > > > > > > If anyone knows the answer, please kindly let me know. > > > > > > Thank you all. > > > > > > Dean > > > > > > ____ Code as follows ________ > > > clear all; close all; clc; > > > N = 10000; > > > Nfft = 2^nextpow2(N); > > > A1 = 1; > > > f1 = 3/32; > > > fs = 1; t = (0:N-1)/fs; > > > x = A1*cos(2*pi*f1*t); > > > win = rectwin(N)'; > > > [Pxx,f] = periodogram(x,win,Nfft,fs); > > > figure; plot(f,Pxx,'k-o'); > > > xlabel('Frequency, Hz'); ylabel('PSD, v^2/Hz'); > > > grid on; xlim([0.05,0.15]); > > > > Hi Dean, You need to scale the periodogram estimate correctly. > > > > N = 10000; > > A1 = 1; > > f1 = 3/32; > > fs = 1; t = (0:N-1)/fs; > > x = A1*cos(2*pi*f1*t); > > h = spectrum.periodogram('Hamming'); > > msest = msspectrum(h,x,'Fs',1); > > plot(msest.Frequencies,msest.Data); > > > > Hope that helps, > > Wayne > > Hi Wayne, > > Thank you for your code, it works perfectly. > > unfortunately, I've generated many graphs using power spectral density rather than mean-square spectrum. The questions is how can I convert PSD to the signal power/amplitude? Or is it possible at all? > > Cheers, > > Dean Hi Dean, in this case you can do: N = 10000; A1 = 1; f1 = 3/32; fs = 1; t = (0:N-1)/fs; x = A1*cos(2*pi*f1*t); [Pxx,F] = periodogram(x,[],[],1,'twosided'); plot(F,Pxx./length(x)); The power in the two complex exponentials is 1/4. If you want just the one-sided. [Pxx,F] = periodogram(x,[],[],1); plot(F,Pxx./length(x)); In general, if you have discrete spectra like your simple example here, then the msspectrum method is probably the best way. If you use the psd() method, there is an avgpower method for integrating under the PSD. h = spectrum.periodogram('Hamming'); psdest = psd(h,x,'Fs',1); psdest = psd(h,x,'Fs',1); avgpower(psdest) % returns 0.5 Hope that helps, Wayne |