close all
clear all

%% start
addpath('../infrastructure/tools');

%% read model
load('../infrastructure/model/model.mat');

%% load data
d = dbload('data/history.csv');

%% set dates
sfilt = qq(1990,1);
ehist = qq(2017,4);
efilt = qq(2022,4);
rpt_rng = sfilt:efilt;

%% Load the observed data
%% Prices
dbobs.obs_lpettax = d.lpettax;
dbobs.obs_lbfp    = d.lbfp   ;
dbobs.obs_w_bfp   = d.w_bfp  ;

% Consumer price components 
dbobs.obs_lcpi_food      = d.lcpi_food    ;
dbobs.obs_lcpi_petr      = d.lcpi_petr    ;
dbobs.obs_lcpi_elec      = d.lcpi_elec    ;
dbobs.obs_lcpi_serv      = d.lcpi_serv    ;
dbobs.obs_lcpi_goodsx    = d.lcpi_goodsx  ;
dbobs.obs_lcpi_core    = d.lcpi_core ;

dbobs.obs_target         = d.target       ;

% Consumer price weights
dbobs.obs_w_food         = d.w_food       ;
dbobs.obs_w_elec         = d.w_elec       ;
dbobs.obs_w_petr         = d.w_petr       ;
dbobs.obs_w_serv         = d.w_serv       ;
dbobs.obs_w_goodsx       = d.w_goodsx     ;

% Wages
dbobs.obs_dot_w          = d.dot_w        ;

%% Relative Prices

dbobs.obs_lrp_food_eq    = hpf(tseries(qq(2000,1):qq(2016,2),d.lrp_food(qq(2000,1):qq(2016,2)))) ;

dbobs.obs_lrp_petr_eq    = hpf(tseries(qq(2000,1):qq(2016,2),d.lrp_petr(qq(2000,1):qq(2016,2)))) ;

dbobs.obs_lrp_serv_eq    = hpf(tseries(qq(2000,1):qq(2016,3),d.lrp_serv(qq(2000,1):qq(2016,3)))) ;

dbobs.obs_lrp_goodsx_eq  = hpf(tseries(qq(2000,1):qq(2017,1),d.lrp_goodsx(qq(2000,1):qq(2017,1)))) ;

%% Rest of World

% Consumer prices
dbobs.obs_lx_cpi         = d.lx_cpi       ;

% Real economy
dbobs.obs_lx_gdp_gap     = d.lx_gdp_gap   ;
dbobs.obs_lx_gdp         = d.lx_gdp       ;
dbobs.obs_g_x            = d.g_x          ; 

% Financial
dbobs.obs_x_rn           = d.x_rn         ;

theTunex = tseries; 
theTunex(qq(2017,4))     = -0.2; %-0.53; % Weighted G3 based on US 0.04, Euro area -0.5, and Japan -1
x_rr_eq_temp            = tseries(qq(2000,1):qq(2017,4),d.x_rr(qq(2000,1):qq(2017,4)));
dbobs.obs_x_rr_eq       = hpf(x_rr_eq_temp,qq(2000,1):qq(2012,4),'level=',theTunex, 'lambda=', 10000);


%% Commodity Prices

dbobs.obs_lfood          = d.lfood        ;
dbobs.obs_lforcmd        = d.lforcmd      ;
dbobs.obs_loil           = d.loil         ;

lqforcmd_eq             = hpf(d.lqforcmd,qq(2000,1):qq(2018,4),'lambda=',835); %see d1_OptimiseCommodities for method used
dbobs.obs_lqforcmd_eq   = lqforcmd_eq{qq(2000,1):qq(2018,4)};

%% Real Economy

dbobs.obs_lemp          = d.lemp        ;
dbobs.obs_lgdp_eq       = d.lgdp_bar    ; 
dbobs.obs_lgdp          = d.lgdp        ;
dbobs.obs_gg            = d.dot_gdp_pot ;

%% Financial 
% Interest rates
dbobs.obs_rn             = d.rn          ;
dbobs.obs_rn8            = d.rn8         ;
dbobs.obs_rn20           = d.rn20        ;
dbobs.obs_rn40           = d.rn40        ;

theTune = tseries; 
theTune(qq(2000,1)) = 5.0;
theTune(qq(2009,2)) = 1.8;
theTune(qq(2012,4)) = 1.5;
theTune(qq(2017,4)) = 1.79;
rr_eq_temp          = tseries(qq(2000,1):qq(2017,4),d.rr(qq(2000,1):qq(2017,4)));
dbobs.obs_rr_eq     = hpf(rr_eq_temp,qq(2001,1):qq(2009,1),'level=',theTune);

% Country Risk Premium
theTune = tseries;
theTune(qq(2000,1)) = 1.8277;
theTune(qq(2006,1)) = 1.5;
theTune(qq(2014,1)) = 1.8;
theTune(qq(2017,2)) = 2.0;
dbobs.obs_prem      = (d.prem5+d.prem6+d.prem7+d.prem8+d.prem9+d.prem10)/6;
dbobs.obs_prem_eq   = hpf(dbobs.obs_prem,qq(2000,1):qq(2017,4),'level=',theTune); 
dbobs.obs_prem_alt    = d.prem5; 

% Exchange rate
dbobs.obs_ls             = d.ls           ;

lz_eq_tune              = tseries;
lz_eq_tune(qq(2008,3))  = 37; 
lz_eq_tune(qq(2018,1))  = 51.3;
dbobs.obs_lz_eq         = hpf(d.lz,qq(2000,1):qq(2017,4),...
                            'lambda=',3200,'level=',lz_eq_tune);

dbsave(dbobs,'data\history_transformed.csv');

% return; 

%% filtration
p.std_e_vathike    = 0;
p.std_e_vat_food   = 0;
p.std_e_vat_serv   = 0;
p.std_e_vat_goodsx = 0;

p.std_e_lpetrol    = 0;
p.std_e_dot_bfp    = 3;
p.std_e_dot_pettax = 3;
p.std_e_w_bfp      = 0.1;

p.std_e_lgdp_gap           = 0.596; %0.2443; %0.6;
p.std_e_lcpi               = 0;
p.std_e_lcpi_core          = 1;
p.std_e_dot_cpi_goodsx     = 1.2420; %0.6453; %1;
p.std_e_dot_cpi_serv       = 0.9410; %1;
p.std_e_dot_cpi_food       = 1.6045; %1; %3; %2;
p.std_e_dot_cpi_petr       = 4.7012; %5; %7.6175; %2;
p.std_e_dot_cpi_elec       = 3.4096; %2;
p.std_e_dot_pm_serv        = 0.5;
p.std_e_dot_pm_goodsx      = 0.5;
p.std_e_dot_w              = 1.7407; %2;
% p.std_e_lrw_gap            = 0.2; %0;
p.std_e_lemp_gap           = 0.5; %1;
p.std_e_dot_emp_eq         = 0.3; %1
p.std_e_ls                 = 2.5670; %2;
p.std_e_rn                 = 0.2473; %1;
p.std_e_target             = 0.4;
p.std_e_lx_gdp_gap         = 0.2879; %1;
p.std_e_x_rn               = 0.2380; %1;
p.std_e_dot_x_cpi          = 0.5466; %1;
p.std_e_x_rr_eq            = 0.3; %1; %0.6 %0.3
p.std_e_rr_eq              = 0.6; %1; %0.6 %0.3
p.std_e_dot_z_eq           = 0.6; %1.5
p.std_e_prem               = 0.1860; %0.3;
p.std_e_prem_eq            = 0.01;
p.std_e_dot_gdp_eq         = 0.491; %0.6;
p.std_e_dot_rw_eq          = 0.0; %0; %1;
p.std_e_dot_rp_food_eq     = 1; %1
p.std_e_dot_rp_food_eq_ss  = 1;
p.std_e_dot_rp_petr_eq     = 2; %0; %1
p.std_e_dot_rp_petr_eq_ss  = 2;
p.std_e_dot_rp_elec_eq     = 1;
p.std_e_dot_rp_goodsx_eq   = 1;
p.std_e_dot_rp_goodsx_eq_ss= 1;

p.std_e_dot_rp_serv_eq     = 1;
p.std_e_dot_qfood_eq       = 1;
p.std_e_lqfood_gap         = 3;
p.std_e_dot_qforcmd_eq     = 0.6; %0.6
p.std_e_lqforcmd_gap       = 3;
p.std_e_dot_qoil_eq        = 1;
p.std_e_lqoil_gap          = 3;
% p.std_e_term               = 1;

p.std_e_term8               = 1;
p.std_e_term20              = 1;
p.std_e_term40              = 1;

p.std_mes_lgdp_eq          = 0; %3;
p.std_mes_lz_eq            = 0; %3;
p.std_mes_rr_eq            = 0.2; %0.3; %0.5
p.std_mes_x_rr_eq          = 0.1; %0.5;
p.std_mes_prem_eq          = 0.0;

p.std_mes_g                = 0.0; %0.5;
p.std_mes_g_x              = 0.0; %0.5;

p.std_e_lgdp_eq            = 0.1; %0.2; %0.3; %0.3; %0.2; %0.3;
p.std_e_g                  = 0.05; %0.1; %0.15; %0.3 %0.7;

p.std_e_lx_gdp_eq          = 0.15; %0.5; %0.3;
p.std_e_g_x                = 0.1; %0.2; %0.7;

p.std_e_w_food             = 0.001; %0.7;
p.std_e_w_elec             = 0.001; %0.7;
p.std_e_w_petr             = 0.001; %0.7;
p.std_e_w_serv             = 0.001; %0.7;
p.std_e_w_goodsx           = 0.001; %0.7;

m = assign(m,p);

[mfilt,filtdb,v_factor,delta,pe,smcov] = filter(m,dbobs,sfilt:efilt,'returnCont=',false);
f = filtdb.mean;

save('data\filter_db.mat','filtdb');

fprintf('\nv_factor = %6.4f\n\n',v_factor);

% shocks' contributions to log-likelihood
eList = get(m,'eList');
disp('Contributions of shocks to the overall minus log likelihood');
for i=1:length(eList)
  sh = filtdb.mean.(eList{i});
  ssh(i) = sum((sh/mfilt.(['std_' eList{i}])).^2); %#ok<SAGROW>
  fprintf('%30s: %8.4f\n',eList{i},ssh(i));
end

%% Graphing and report 
Plotrng = qq(2002,1):ehist;
d = filtdb.mean;
d.zeroline = tseries(Plotrng,0);

% make report
commandwindow;
graphsFilter

%% finish
rmpath('../infrastructure/tools');