ar1psd.m

蒙托卡罗模拟奇异谱分析

function [P,c0]=ar1psd(g,a,N)

% AR1PSD - returns the spectral density of an AR(1) process.
%
% Syntax: [P,c0]=ar1psd(g,a,N); P=ar1psd(g,a,N);
%
% Input:   g - lag-one autocorrelation.
%          a - noise variance.
%          N - number of points desired in spectral density vector.
%
% Output:  P - power spectral density, N points from 0 to .5 in a column vector.
%         c0 - zero-lag autocovariance.
%
%  AR(1) process: x(t) = g * x(t-1) + a * z(t),
%        where z(t) is Gaussian noise.
%
% The spectrum is roughly normalized, so that sum(P)=c0*N;
%
% See e.g Jenkins and Watts 1968, pg. 228.
%
% Written by Eric Breitenberger.      Version 4/22/97
% Please send comments and suggestions to [e-mail]       
%

P=linspace(0,.5, N)';
P=a./(1+g^2-2*g*cos(2*pi*P));

c0=a^2/(1-g^2);

P=P*c0*N/sum(P);