Autoregressive power spectral density estimate — modified covariance method
the power spectral density estimate,
pxx = pmcov(
a discrete-time signal,
x, found using the modified
covariance method. When
x is a vector, it is
treated as a single channel. When
x is a matrix,
the PSD is computed independently for each column and stored in the
corresponding column of
the distribution of power per unit frequency. The frequency is expressed
in units of rad/sample.
order is the order of
the autoregressive (AR) model used to produce the PSD estimate.
pxx = pmcov(
in the discrete Fourier transform (DFT). For real
even, and (
odd. For complex–valued
nfft. If you omit
or specify it as empty, then
pmcov uses a default
DFT length of 256.
returns a frequency vector,
f] = pmcov(___,
f, in cycles per unit time. The sampling
fs, is the number of samples per unit time. If the unit of time
is seconds, then
f is in cycles/second (Hz). For real-valued signals,
f spans the interval [0,
nfft is even and [0,
nfft is odd. For complex-valued signals,
f spans the
returns the two-sided AR PSD estimates at the frequencies specified in the vector,
f] = pmcov(
f. The vector
f must contain at least two elements,
because otherwise the function interprets it as
nfft. The frequencies in
f are in cycles per unit time. The sampling frequency,
fs, is the number of samples per unit time. If the unit of time is
f is in cycles/second (Hz).
pmcov(___) with no output arguments
plots the AR PSD estimate in dB per unit frequency in the current
Create a realization of an AR(4) wide-sense stationary random process. Estimate the PSD using the modified covariance 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
pmcov 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] = pmcov(y,4,1024,1); hold on plot(F,10*log10(Pxx)) legend('True Power Spectral Density','pmcov PSD Estimate')
Create a multichannel signal consisting of three sinusoids in additive white Gaussian noise. The sinusoids' frequencies are 100 Hz, 200 Hz, and 300 Hz. The sampling frequency is 1 kHz, and the signal has a duration of 1 s.
Fs = 1000; t = 0:1/Fs:1-1/Fs; f = [100;200;300]; x = cos(2*pi*f*t)'+randn(length(t),3);
Estimate the PSD of the signal using the modified covariance method with a 12th-order autoregressive model. Use the default DFT length. Plot the estimate.
morder = 12; pmcov(x,morder,,Fs)
Usage notes and limitations:
nfft argument is variable-size at compile time, then it must
not become a scalar or an empty array at runtime.