function Out = CutLevelSets(z,numsplit)
%Take level sets of a function with more samples in the middle
%Level set height sampling done with equal area under sine function
%
%INPUTS
%  z = 1000 x 1000 matrix, input data set from Continuum Model
%  h = Number of level sets to be taken (heights to sample from)
%
%OUTPUTS
%   Out = A z x 2 array which details the Area Fraction and the binary 
%         variable PathExist which indicates the presence of percolation
%         for each level set taken in melt
%
%By: Anthony Cheng (Mar 2015)
%===================================================================

%-------------------------------------------------------------------
%Take level sets of a function with equal area under sine function
%-------------------------------------------------------------------
%1 to -1 as integral of sin(x) is -cos(x)
equalarea = linspace(1,-1,numsplit); 
inv = acos(equalarea);

maxval = max(max(z));
minval = min(min(z));
aveval = (maxval+minval)/2;

%Normalize "inv" to a distribution centered at 0 rather than pi/2
norm = (inv-mean(inv))/std(inv);
%"Translate" norm to a distribution centered at the mean value of the data
%set rather than 0
dev = (minval-aveval)/norm(1,1);
spread = dev*norm+aveval;

%-------------------------------------------------------------------
%Take a binary sampling of the data set to simulate points of melt 
%ponds (1) on sea ice (0) and concatenate in variable "melt"
%-------------------------------------------------------------------
melt = [];
for n = 1:numsplit
    k = spread(1,n);
    binary = k >= z;
    binary = +binary;
    melt = cat(3,melt,binary);
end

%-------------------------------------------------------------------
%Determine if percolation occurs at each sampled level set, and the
%Area Fraction at each level set
%-------------------------------------------------------------------
Out = CriticalArea(melt);

end