%% 1. PRELIMINARIES
% =======================================================================
clear all; clear session; close all; clc
warning off

addpath './VAR_Toolbox/VAR'
addpath './VAR_Toolbox/Utils'
addpath './VAR_Toolbox/Auxiliary'
addpath './VAR_Toolbox/ExportFig'
addpath './VAR_Toolbox/Figure'
addpath './VAR_Toolbox/Stats'

country_list = ["Australia", "Austria", "Belgium", "Brazil", "Canada", "Chile", "China", ...
    "Colombia", "Czech Republic", "Denmark", "Estonia", "Finland", ...
    "France", "Germany", "Greece", "Hungary", "Iceland", ...
    "Ireland", "Israel", "Italy", "Japan", "Korea", ...
    "Latvia", "Lithuania", "Luxembourg", "Mexico", "Netherlands", ...
    "New Zealand", "Norway", "Poland", "Portugal", "Slovak Republic", ...
    "Slovenia", "Spain", "Sweden", "Switzerland", "Turkey", ...
    "United Kingdom","United States", ];


country_list_unc =["Australia", "Austria", "Belgium", "Brazil", "Canada", "Chile", "China", ...
    "Colombia", "Czech", "Denmark", "Estonia", "Finland", ...
    "France", "Germany", "Greece", "Hungary", "Iceland", ...
    "Ireland", "Israel", "Italy", "Japan", "Korea", ...
    "Latvia", "Lithuania", "Luxembourg", "Mexico", "Netherlands", ...
    "New Zealand", "Norway", "Poland", "Portugal", "Slovak Republic", ...
    "Slovenia", "Spain", "Sweden", "Switzerland", "Turkey", ...
    "US", "UK"];

%load US and control variables
[xlsdata,text] = xlsread('./Data/main_data.xlsx','United States');

ip_us = (xlsdata(:,2));
reu_us = (xlsdata(:,1));
ff  = (xlsdata(:,3));
emp  = (xlsdata(:,4));
sp500  = (xlsdata(:,5));
cpi  = (xlsdata(:,6));
vxo = (xlsdata(:,7));

report_out = [];
for jj = 1:size(country_list,2)
    disp(jj);
    
    % Load
    [xlsdata,text] = xlsread('./Data/main_data.xlsx',country_list{jj});
    [xlsdata_unc,text_unc] = xlsread('./Data/uncertainty_index_wps.xlsx',country_list_unc{jj});
        %dates = 1973+5/12 : 1/12 : 2018+12/12;
%          if jj == 4|8|12
%      REU = [xlsdata(1:end-4,1);NaN;NaN;NaN;NaN];
%          else
%      
%             REU = xlsdata(1:end,1);
%          end
    
    dates_1 = text(2,1);
    temp = char(dates_1(1));
    year = str2double(temp(1:4));
    month = str2double(temp(6:end));
    
    dates_first = year + month/12;
    
    dates_f = 1960+7/12;
    
    index = floor((dates_first - dates_f)*12);
    
    
    %  REU = xlsdata(1:end,1);
    REU = [xlsdata_unc(1:end,1);NaN(0,1)];
    ip = xlsdata(:,2);
    
    
    REU_index = [NaN(index,1);REU];
    
    
    ip_index = [NaN(index,1);ip];

    ind1 =isnan(REU_index) ==0;
    ind2 =isnan(ip_index) ==0;
    ind3 = isnan(vxo) == 0;
    ind = ind1&ind2&ind3;
    
    REU = REU_index(ind,1);
    ip =  ip_index(ind,1);
    vxo_index = vxo(ind);
    
    
    T = length(ip);
    index = [1,2];
    
    if jj == 17
    X = [vxo_index,REU,-100*log(ip)];
    
    else
     X = [vxo_index,REU,100*log(ip)];   
    end
    
    %FigName = ['./Figure_2020/IRF_bi_validation_',country_list{jj}];
    vnames = text(1,2:3);
    % Define number of variables and of observations
    [nobs, nvar] = size(X);
    
    
    %% VAR ESTIMATION
    % =======================================================================
    % Set the case for the VARout (0, 1, or 2)
    det = 2;
    % Set number of nlags
    %nlags = 12;
    nlags = 6;
    % Estimate
    [VAR, VARopt] = VARmodel(X,nlags,det);
    % Print at screen and create table
    VARopt.vnames = ['A','B',country_list(jj)];
    [beta, tstat, TABLE] = VARprint_notable(VAR,VARopt);
    
    
    %% IMPULSE RESPONSE
    % =======================================================================
    % Set options some options for IRF calculation
    VARopt.nsteps = 12*5;
    VARopt.ident = 'oir';
    VARopt.quality = 0;
    VARopt.pick = 2;
    VARopt.pctg = 68;
    % Compute IRF
    [IRF, VAR] = VARir(VAR,VARopt);
    % Compute error bands
    [IRFINF,IRFSUP,IRFMED] = VARirband(VAR,VARopt);
    % Plot
    
    
    if jj <= 20
        subplot(4,5,jj)
        VARirplot_val(IRFMED,VARopt,IRFINF,IRFSUP,'test');
    else
        %figure
        subplot(4,5,jj-20)
        VARirplot_val(IRFMED,VARopt,IRFINF,IRFSUP,'test');
    end
    
    if jj == 20
        
        dim = [8,6];
        set(gcf,'papersize',dim,'paperposition',[0,0,dim]);
        figure(gcf)
        print(gcf,'-dpdf','./figure/IRF_validation_first20') %Change Figure Name Here
        close all
    elseif jj == 39
        
        dim = [8,6];
        set(gcf,'papersize',dim,'paperposition',[0,0,dim]);
        figure(gcf)
        print(gcf,'-dpdf','./figure/IRF_validation_next19') %Change Figure Name Here
        close all
        
    end
    %
    report_temp= [IRFINF(:,3,2),IRFMED(:,3,2),IRFSUP(:,3,2)];
    report_out = [report_out, report_temp];
   % pause
    
end

