tfar_fit_mse.m

来自「用于模拟时变非平稳的ARMA过程」· M 代码 · 共 53 行

function [AA, BB]= tfar_fit_mse(x, Mmax, Lmax, Psi, R)
% function [AA, BB]= tfar_fit_aic(x, Mmax, Lmax, Psi, R)
%   This file is part of the TFPM toolbox v1.0 (c)
%   [e-mail] and underlies the GPL.
% 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if(0)% TEST TEST TEST TEST TEST TEST TEST TEST TEST TEST
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear;tfpm
load('~/matlab/data/Motordaten/bmw1000.mat');
MM= 100;
N0= 1468;
N= 256;
offset= 64;
alpha= 1/2;

x= c1(offset:2:offset+2*N-1, 1)-mean(c1(offset:2:offset+2*N-1, 1));
plot(x)
Mmax= 4;
Lmax= 4;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
end;% TEST TEST TEST TEST TEST TEST TEST TEST TEST TEST
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

alpha= 1/2;
N= length(x);

%Psi= tf_multiwin(N, 3*Mmax, 4*Lmax, 0, 2, 1);
Psi= Psi(:, N/2-Mmax+1:N/2+Mmax+1);

Qalpha= corr_est(x, x, Mmax, alpha);
Aalpha= fft(Qalpha);
Aalpha= [ Aalpha(N/2+1:N, :); Aalpha(1:N/2, :)].*conj(Psi);
E= zeros(Mmax, Lmax);

for LAR= 1:Lmax
   Ayy= Aalpha(N/2-3*LAR+1:N/2+3*LAR+1, :);
   [AA, BB]= tfar_est_tfywu(Ayy, N);
   for MAR= 1:Mmax
      Rx= tfarma_wvsp(AA(:, 1:MAR+1, MAR), BB(:, :, MAR), N, alpha);
      E(MAR, LAR)= norm(Rx-Rengine)/norm(Rengine)
   end;
end;
[MAR, LAR]= find(E==min(min(E)));

Qalpha= corr_est(x, x, MAR, alpha);
Aalpha= fft(Qalpha);
Aalpha= [ Aalpha(N/2+1:N, :); Aalpha(1:N/2, :)];
Ayy= Aalpha(N/2-3*LAR+1:N/2+3*LAR+1, :);
[AA, BB]= tfar_est_tfywu(Ayy, N);