%% data generation process:

% T = 500;
% 
% A = [   0.99    0.0074;
%    -0.0136    0.99];
% 
% C = [ 1 1 ; -1 +1];
% 
% x(:,1) = [3;1];
% 
% Sigma_w = %% you will find out by  running EM
% Sigma_v = %% you will find out by  running EM
% 
% w = randn(2,T); w = sqrtm(Sigma_w)*w;
% v = randn(2,T); v = sqrtm(Sigma_v)*v;
% 
% for t=1:T-1
% 	x(:,t+1) = A * x(:,t) + sqrtm(Sigma_w)*w(:,t);
% 	y(:,t) = C*x(:,t) + sqrtm(Sigma_v) * v(:,t);
% end
% y(:,T) = C*x(:,T) + sqrtm(Sigma_v) * v(:,T);
% 
% clear Sigma_w Sigma_v w v
%
% 
%
% save kf_data.mat

load kf_data.mat

%% note: x is the true state; y are the measurements; in a real
%% application, we would typically not have x available

figure; plot(x');
figure; plot(y');


P_0 = diag([.01 .01]);
x0 = x(:,1);

% we don't know Sigma_v, Sigma_w, let's use the following estimates:
Sigma_w_hat = diag([1 1]);
Sigma_v_hat = diag([.5 .5]);


%%% your code here %%%