% Read in initial light results

clear all;
small = 1.0e-10;

% -- File Directories   
figdir = 'fig/';
outdir = 'out/';
matdir = 'mat/';

n = 250 ;         % Number of obs per draw

mod_str = ['mc_unitline_model_0_' num2str(n)];
load([matdir mod_str '_rejfreq_scpc']);
%load([matdir mod_str '_rejfreq_cscpc']);
load([matdir mod_str '_rejfreq_hr']);
rej_scpc_0 = packr(rejfreq_scpc);
%rej_cscpc_0 = packr(rejfreq_cscpc);
rej_hr_0 = packr(rejfreq_hr);
load([matdir mod_str '_rejfreq_conley']);
% Find Conley closest to 5%
tmp = mean(rejfreq_conley)';  % Average over x's
atmp = abs(tmp-0.05);
[~,ii]=min(atmp);
rej_conley_0 = rejfreq_conley(:,ii);


mod_str = ['mc_unitline_model_1_' num2str(n)];
load([matdir mod_str '_rejfreq_scpc']);
load([matdir mod_str '_rejfreq_cscpc']);
load([matdir mod_str '_rejfreq_hr']);
rej_scpc_1 = packr(rejfreq_scpc);
rej_cscpc_1 = packr(rejfreq_cscpc);
rej_hr_1 = packr(rejfreq_hr);
load([matdir mod_str '_rejfreq_conley']);
% Find Conley closest to 5%
tmp = mean(rejfreq_conley)';  % Average over x's
atmp = abs(tmp-0.05);
[~,ii]=min(atmp);
rej_conley_1 = rejfreq_conley(:,ii);

mod_str = ['mc_unitline_model_1_03_' num2str(n)];
load([matdir mod_str '_rejfreq_cscpc']);
rej_cscpc_1_03 = packr(rejfreq_cscpc);

mod_str = ['mc_unitline_model_2_dm_' num2str(n)];
load([matdir mod_str '_rejfreq_scpc']);
load([matdir mod_str '_rejfreq_cscpc']);
load([matdir mod_str '_rejfreq_hr']);
rej_scpc_2 = packr(rejfreq_scpc);
rej_cscpc_2 = packr(rejfreq_cscpc);
rej_hr_2 = packr(rejfreq_hr);
load([matdir mod_str '_rejfreq_conley']);
% Find Conley closest to 5%
tmp = mean(rejfreq_conley)';  % Average over x's
atmp = abs(tmp-0.05);
[~,ii]=min(atmp);
rej_conley_2 = rejfreq_conley(:,ii);

mod_str = ['mc_unitline_model_3_dm_' num2str(n)];
load([matdir mod_str '_rejfreq_scpc']);
load([matdir mod_str '_rejfreq_cscpc']);
load([matdir mod_str '_rejfreq_hr']);
rej_scpc_3 = packr(rejfreq_scpc);
rej_cscpc_3 = packr(rejfreq_cscpc);
rej_hr_3 = packr(rejfreq_hr);
load([matdir mod_str '_rejfreq_conley']);
% Find Conley closest to 5%
tmp = mean(rejfreq_conley)';  % Average over x's
atmp = abs(tmp-0.05);
[~,ii]=min(atmp);
rej_conley_3 = rejfreq_conley(:,ii);

% -- TABLE 1: -- %
outfile_name    = [outdir, 'Table_1.csv'];
fileID          = fopen(outfile_name,'w');
fprintf(fileID,'Model, HR, Kernel, SCPC \n');
fprintf(fileID,'Mod 1, %5.2f, %5.2f, %5.2f  \n',[mean(rej_hr_0) mean(rej_conley_0) mean(rej_scpc_0)]);
fprintf(fileID,'Mod 2, %5.2f, %5.2f, %5.2f  \n',[mean(rej_hr_1) mean(rej_conley_1) mean(rej_scpc_1)]);
fprintf(fileID,'Mod 3, %5.2f, %5.2f, %5.2f  \n',[mean(rej_hr_3) mean(rej_conley_3) mean(rej_scpc_3)]);
fprintf(fileID,'Mod 4, %5.2f, %5.2f, %5.2f  \n',[mean(rej_hr_2) mean(rej_conley_2) mean(rej_scpc_2)]);


% -- TABLE 2: -- %
pcvec = [0.05 0.25 0.50 0.75 0.95]';
q_scpc_1 = pctile(rej_scpc_1,pcvec);
q_scpc_2 = pctile(rej_scpc_2,pcvec);
q_scpc_3 = pctile(rej_scpc_3,pcvec);
outfile_name    = [outdir, 'Table_2.csv'];
fileID          = fopen(outfile_name,'w');
fprintf(fileID,'Model, 0.05, 0.25, 0.50, 0.75, 0.95 \n');
fprintf(fileID,'Mod 1, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f  \n',[mean(rej_scpc_0) mean(rej_scpc_0) mean(rej_scpc_0) mean(rej_scpc_0) mean(rej_scpc_0)]);
fprintf(fileID,'Mod 2, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f  \n', q_scpc_1');
fprintf(fileID,'Mod 3, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f  \n',q_scpc_3');
fprintf(fileID,'Mod 4, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f  \n',q_scpc_2');

% -- TABLE 3: -- %
pcvec = [0.05 0.25 0.50 0.75 0.95]';
q_cscpc_1 = pctile(rej_cscpc_1,pcvec);
q_cscpc_1_03 = pctile(rej_cscpc_1_03,pcvec);
q_cscpc_2 = pctile(rej_cscpc_2,pcvec);
q_cscpc_3 = pctile(rej_cscpc_3,pcvec);
outfile_name    = [outdir, 'Table_3.csv'];
fileID          = fopen(outfile_name,'w');
fprintf(fileID,'Model, Avg, 0.05, 0.25, 0.50, 0.75, 0.95 \n');
fprintf(fileID,'Mod 1, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f  \n',[mean(rej_scpc_0) mean(rej_scpc_0) mean(rej_scpc_0) mean(rej_scpc_0) mean(rej_scpc_0) mean(rej_scpc_0)]);
fprintf(fileID,'Mod 2a, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f  \n',[mean(rej_cscpc_1) q_cscpc_1']);
fprintf(fileID,'Mod 2b, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f  \n',[mean(rej_cscpc_1_03) q_cscpc_1_03']);
fprintf(fileID,'Mod 3, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f  \n',[mean(rej_cscpc_3) q_cscpc_3']);
fprintf(fileID,'Mod 4, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f, %5.2f  \n',[mean(rej_cscpc_2) q_cscpc_2']);



