format long;
set(0,'RecursionLimit',10000);
close all;

binnum = 25; %Number of bins
areabins = (1/binnum:1/binnum:1).';
allperccounts = zeros(binnum,1);
alltotalcounts = zeros(binnum,1);

directoryName = '../Data/pde_0.3_0.2/';
files = dir(directoryName);
fileIndex = find(~[files.isdir]);

%Create a cell structure for binning the results
areadata = horzcat(areabins, zeros(binnum,2));

%depth marker
depth = 20;

for i = 1:numel(fileIndex)
%Load the Image
fileName = files(fileIndex(i)).name;
load(strcat(directoryName,fileName))
threshdepth = 0.07;

%Analyze the level sets 
    for li = 1:size(hout,1)
        h = squeeze(hout(li,:,:));
        %h = rot90(h);
        Out = BinaryThresh(h,threshdepth); %Out = [AreaFraction,PathExist];
                                           
        %Extract Area Fraction data and the presence of percolation
        AreaFraction = Out(1,1); 
        PathExist = Out(1,2); 
        
        %Determine which bin to place the data
        bin = ceil(AreaFraction*binnum);
        
        if bin == 0
            bin = 1;
        end
        
        %Add data from level set to resultant bin
        areadata(bin,2) = areadata(bin,2) + PathExist;
        areadata(bin,3) = areadata(bin,3) + 1;

        perccounts = areadata(:,2);
        totalcounts = areadata(:,3);

        allperccounts = allperccounts + perccounts;
        alltotalcounts = alltotalcounts + totalcounts;
    end
end             

allavecounts = allperccounts./alltotalcounts;
x = areabins; y = allavecounts;

[beta,R,J,CovB,MSE,ErrorModelInfo] = ... 
nlinfit(x,y,@(c,x)((1/pi).*atan(c(1).*(x-c(2)))+c(3)),[40,0.5,0.5]);
ci = nlparci(beta,R,'jacobian',J); %for the 95% confidence intervals

hold on; plot(x,y,'k-');
line([beta(2) beta(2)],get(gca,'YLim'),'Color',[1 0 0]);
ylabel('Percent Percolation','fontsize',22);
xlabel(strcat('Area Fraction (Critical: ',num2str(beta(2)),')'),'fontsize',22);
title(strcat('PDE 0.2 0.2 Rotated Depth:',num2str(depth)),'fontsize',32);
disp(ci(2,1));
disp(ci(2,2));

hold off;

%filename = strcat('PDE_0.2_0.2.mat');
%save(filename,'areabins','allavecounts')
