demo.m

无参数功率谱估计

clc
disp(['Power Spectrum Estimation Via Thresholded Cepstrum Demonstration'])
disp(['Code by Steve Lutes, concept from "Smoothed Nonparametric'])
disp(['Spectral Estimation Via Cepstrum Thresholding" by Perte Stocia'])
disp(['and Niclas Sandgren, IEEE Signal Processing Magazine, November 2006 vol 23, number 6'])

disp([''])

disp(['There are 4 sets of data in this demo.']);
disp([''])
disp(['Press space bar to continue.'])
pause

clc
disp(['Demo 1-Power Spectrum of White Noise'])


Input = randn(1,200);

plot(Input)
title(['Time Plot Of Gaussian White Noise-Demo 1'])

Power1 = abs(fft(Input).^2);
[Power2, Message, StatusFlag] = CepSpec(Input, 'wide', 'raw');
MidPoint = length(Power2);

figure
plot(Power1(1:1:MidPoint))
title(['Power Spectrum Of Gaussian White Noise By Wiener-Kinchine Method-Demo 1'])

figure
plot(Power2)
title(['Power Spectrum Of Gaussian White Noise By Thresholded Cepstrum Method Method-Demo 1'])

disp([''])
disp(['Press space bar to continue.'])
pause
close all
clc

disp(['Demo 2-Power Spectrum of Multisine With Added White Noise'])

TimeValues = [0:1:299];
Noise = randn(size(TimeValues));

Input = sin(2*pi*(TimeValues/50)) + (1.0*Noise);
MidPoint = fix(length(TimeValues)/2);

plot(Input)
title(['Time Plot Of Single Sine w/ Added Noise-Demo 2'])

Power1 = abs(fft(Input-mean(Input))).^2;
[Power2, Message, StatusFlag] = CepSpec(Input-mean(Input), 'narrow', 'raw');

figure
plot(Power1(1:1:MidPoint))
title(['Power Spectrum Of Single Sine w/ Added Noise By Wiener-Kinchine Method-Demo 2'])

figure
plot(Power2)
title(['Power Spectrum Of Single Sine w/ Added Noise By Thresholded Cepstrum Method Method-Demo 2'])


disp([''])
disp(['Press space bar to continue.'])
pause
close all
clc

disp(['Demo 4-Power Spectrum of US$-Japanese Yen Exchange Rate i.e. number of Yen/US$'])

Input = dlmread('CanadianLynxPop.txt','\t');
MidPoint = fix(length(Input)/2);

plot(Input)
title(['Time Plot Of Candian Lnyx Population Data'])

Power1 = abs(fft(Input-mean(Input))).^2;
[Power2, Message, StatusFlag] = CepSpec(Input-mean(Input), 'narrow', 'raw');

figure
plot(Power1(1:1:MidPoint))
title(['Power Spectrum Of Candian Lnyx Population Data By Wiener-Kinchine Method-Demo 3'])

figure
plot(Power2)
title(['Power Spectrum Of Candian Lnyx Population Data By Thresholded Cepstrum Method Method-Demo 3'])


disp([''])
disp(['Press space bar to continue.'])
pause
close all
clc

disp(['Demo 4-Power Spectrum of US$-Japanese Yen Exchange Rate Data'])
disp(['Exchange rate is Number of Yen per US dollar'])

Input = dlmread('USDYenExchange.txt','\t');
MidPoint = fix(length(Input)/2);

plot(Input)
title(['Time Plot Of US$-Japanese Yen Exchange Rate Data'])

Power1 = abs(fft(Input)).^2;
[Power2, Message, StatusFlag] = CepSpec(Input, 'narrow', 'raw');

figure
loglog(Power1(1:1:MidPoint))
title(['Power Spectrum Of US$-Japanese Yen Exchange Rate By Wiener-Kinchine Method-Demo 3'])

figure
loglog(Power2)
title(['Power Spectrum Of US$-Japanese Yen Exchange Rate By Thresholded Cepstrum Method Method-Demo 3'])


disp([''])
disp(['Press space bar to continue.'])

pause
close all
clc

disp(['All demos finished.'])