function [parent] = findkloop(parent,varargin)
if size(varargin,2)==2;
disp('Arguments mutually exclusive, please try again.');
return
end
%Set Up Variables
    sparsity=0.5;
    kmax=15;
    repeats=5;
    score1=1;
%varargin arguments
if isempty(varargin)
    data=parent.score;
    disp('System: Default to scores matrix.')
    J=size(data,1);
    L=round(J*(sparsity/100));
    score2=3;
else
    if strcmp(varargin{1},'norm')==1;
    data=parent.normscore;
    disp('User: normalised scores.')
    J=size(data,1);
    L=round(J*(sparsity/100));
    score2=3;
    else 
        if strcmp(varargin{1},'data')==1;
        data=parent.data;
        disp('User: parent data, brute-force in N-dimensions.')
        J=size(data,1);
        L=round(J*(sparsity/100));
        score2=size(data,2);
        else
            if strcmp(varargin(1),'data')==0;
            data=parent.score;
            disp('System: Options invalid. Default to scores matrix.')
            J=size(data,1);
            L=round(J*(sparsity/100));
            score2=3;
            end
        end
    end
end
SAMPS=5;
disp(sprintf('Begin %0.0f-dimensional separation, %0.0f loops over %0.0f samplings of %0.2f percent; %0.0f/%0.0f:',score2-score1+1,SAMPS*repeats,SAMPS,sparsity,L,J));
mf=zeros(kmax,repeats,SAMPS);
pf=zeros(kmax,repeats,SAMPS);
tf=zeros(kmax,repeats,SAMPS);
idxf=zeros(kmax,repeats,L,SAMPS);
scf=zeros(kmax,SAMPS);
for a=1:SAMPS
    fprintf('> %1.0f of %0.0f:',a,SAMPS);
    [mf(:,:,a),pf(:,:,a),tf(:,:,a),idxf(:,:,:,a),scf(:,a),entb(:,:,:,a),silnt(:,:,a)] = findkint(data,score1,score2,kmax,repeats,sparsity);
end
parent.findk.name=sparsity;
parent.findk.means=reshape(mf,kmax,repeats*SAMPS);
parent.findk.quotients=reshape(pf,kmax,repeats*SAMPS);
parent.findk.processtimes=reshape(tf,kmax,repeats*SAMPS);
parent.findk.indices=reshape(idxf,kmax,repeats*SAMPS,L);
parent.findk.entr=scf;
parent.findk.ent_long=entb;
parent.findk.sent=reshape(silnt,kmax,repeats*SAMPS);
%disp('Rendering...');
%allplotsfindk(parent);
disp(sprintf('Complete. Total time %0.3f seconds',sum(sum(sum(tf)))));
end
function [means,props,timerecs,idx,ent,entb,sent] = findkint(data,a,b,kmax,repeats,sparsity)
means=ones(kmax,repeats);
% means(1,:)=NaN;
props=ones(kmax,repeats);
% props(1,:)=NaN;
timerecs=zeros(kmax,repeats);
sent=zeros(kmax,repeats);
% timerecs(1,:)=NaN;
data=data(:,a:b);
J=size(data,1);
L=round(J*(sparsity/100));
idx=zeros(kmax,repeats,L);
datasamp=data(randsample(J,L),:);
for n=1:repeats
   for k=2:kmax
   tic;
   idx(k,n,:)=kmeans(datasamp,k,'Maxiter',200,'emptyaction','drop','replicates',5,'Distance','cityblock');
   sil=silhouette(datasamp,reshape(idx(k,n,:),1,[]),'cityblock');
   sent(k,n)=entropy(sil);
   means(k,n)=mean(sil);
   props(k,n)=((size(sil(sil>0),1)/size(sil,1)));
   timerecs(k,n)=toc;
   end
   fprintf('|');
end
% for k=1:kmax
% a(k)=(mean(mean(abs(corr(reshape(idx(k,:,:),repeats,[]))))));
% end
fprintf(' ');
ent=zeros(kmax,1);
entb=zeros(kmax,repeats,repeats);
for k=1:kmax
I=reshape(idx(k,:,:),L,repeats);
for m=1:repeats
    for n=m:repeats
        [~,~,B] = unique(I(:,[m n]),'rows','stable');
        [xp] = plotgram(B,L);
        p=xp(:,2);
        p(p==0) = [];
        p = p ./ numel(B);
        entb(k,m,n) = -sum(p.*log10(p));
        entb(k,n,m) = entb(k,m,n);
    end
end
ent(k)=mean(mean(entb(k,:,:)));
end
ent(isnan(ent))=0;
ent=ent./max(ent);
fprintf(' ');
close all
end