arfilter.m

来自「AR滤波器的标准程序」· M 代码 · 共 37 行

% Title: AR filter in accordance of a power spectrum
% Author: Lilei
% Data: 2008.6.2
% Description: Input: the power spectrum (Spectrum),
%                     the variance of the white noise( Var )
%                     the order of the AR filter(Order)
%              Output: the AR filter coefficients(a)  
% Change: 
% Addition: 

function [ aa ] = ARFilter( Spectrum , Var , Order )
R = real( ifft( Spectrum ) );

figure;
plot(R);
xlim([1 200]);

%%%%%% compute the AR model coefficient %%%%%%%%%%
a_former = zeros( 1 , Order );
a_latter = zeros( 1 , Order );
var_former = 0.0;
var_latter = 0.0;

a_latter( 1 , 1 ) = - R( 1 , 2 ) / R( 1 , 1 );
var_latter = ( 1 - ( abs( a_latter( 1 , 1 ) ))^2) * R( 1 , 1 );

for( m = 2 : Order )
    var_former = var_latter;        
    a_former = a_latter;
    k = 1 : m - 1;
    a_latter( 1 , m ) = - ( sum ( a_former( 1 , k ) .* R( m - k + 1 ) ) + R( m + 1 ) ) / var_former;
    a_latter( 1 , k ) = a_former( 1 , k ) + a_latter( 1 , m ) .* a_former( 1 , m - k );
    var_latter = var_former * ( 1 - (abs(a_latter( 1 , m )))^2);
    plot( m , var_latter , 'r*' );
    hold on;
end
aa = a_latter;