close all

AddPostFix=[];
if DO_TREND
    AddPostFix=[AddPostFix '_hp'];
end
if StartYear>0
    AddPostFix=[AddPostFix '_' num2str(StartYear)];
end
if EndYear<inf
    AddPostFix=[AddPostFix '_' num2str(StartYear)];
end

DescriptiveTable=cell(Params.p,1);
for pp=1:length(VarNames)
%     
    if pp<=Params.p
        DescriptiveTable{pp}=[VarLabels{pp} ' & ' num2str(nanmean(reshape(Y(pp,:,:),1,[])), '%1.2f ') ' ' num2str(nanstd(reshape(Y(pp,:,:),1,[])), ' (%1.2f) \\\\') ];
        YearFE=mean(squeeze(Y(pp,:,:))');
        PlotData=squeeze(Y(pp,:,:))'-ones(Params.N,1)*YearFE;
    else
        DescriptiveTable{pp}=[VarLabels{pp} ' & ' num2str(nanmean(reshape(Z(pp-Params.p,:,:),1,[])), '%1.2f ') ' ' num2str(nanstd(reshape(Z(pp-Params.p,:,:),1,[])), ' (%1.2f) \\\\') ];
        YearFE=mean(squeeze(Z(pp-Params.p,:,:))');
        PlotData=squeeze(Z(pp-Params.p,:,:))'-ones(Params.N,1)*YearFE;
    end
    PlotData=PlotData-mean(PlotData,2)*ones(1,T-StartObs+1);
%     hndl=plot(vYear(StartObs:T),prctile(PlotData,[5 25 50 75 95])'); 
%     hold all
%     plot(vYear(StartObs:T),YearFE,'-k');
%     
%     close all
    if pp<=Params.p
        if VarScale(pp)==1 
            DisplayUnit='p.p.';
        else
            DisplayUnit='%';
        end
        hndl=my_make_IRF(vYear(StartObs:T),prctile(PlotData,[5 25 50 75 95])',YearFE'-mean(YearFE), 'DoAreaPlot', true, 'YLabel',DisplayUnit,'XLabel','');
        legend(hndl([3 4]), {'MSAs net of fixed effects'; 'US fixed effects'}, 'Location', 'Best');
    else
        DisplayUnit='% / p.p.';
        hndl=my_make_IRF(vYear(StartObs:T),prctile(PlotData,[5 25 50 75 95])', 'DoAreaPlot', true, 'YLabel',DisplayUnit,'XLabel','');
        legend(hndl([3 ]), {'MSAs net of fixed effects'}, 'Location', 'Best');
    end
    axis([vYear(StartObs) vYear(T) ylim]);
    legend('boxoff');
    title(VarLabels{pp});
    
    print('-depsc', ['Graphs/Descriptive_' VarNames{pp} AddPostFix '.eps']);
%     hold all;
%     plot(); 
end
