%% Burg PSD Estimate of an AR(4) Process
% Create a realization of an AR(4) wide-sense stationary random process. Estimate the PSD using Burg's method. Compare the PSD estimate based on a single realization to the true PSD of the random process.
% Create an AR(4) system function. Obtain the frequency response and plot the PSD of the system.

A = [1 -2.7607 3.8106 -2.6535 0.9238];
[H,F] = freqz(1,A,[],1);
plot(F,20*log10(abs(H)))

xlabel('Frequency (Hz)')
ylabel('PSD (dB/Hz)')

%
% Create a realization of the AR(4) random process. Set the random number generator to the default settings for reproducible results. The realization is 1000 samples in length. Assume a sampling frequency of 1 Hz. Use pburg to estimate the PSD for a 4th-order process. Compare the PSD estimate with the true PSD.

rng default

x = randn(1000,1);
y = filter(1,A,x);
[Pxx,F] = pburg(y,4,1024,1);

hold on
plot(F,10*log10(Pxx))
legend('True Power Spectral Density','pburg PSD Estimate')

% Copyright 2015 The MathWorks, Inc.