SystemEvolution_findk(inter,inter2,intra,intra2,N,Ns,Hsep)

function [ters,km,g1,g2,ind,ind2] = SystemEvolution_findk(inter,inter2,intra,intra2,N,Ns,Hsep)

%(2)  collecting partition and merge energies gained
ind = [inter(:,3)'; max(inter(:,1:2)')];
g1 = inter(:,4)';
A1 = find(inter2(:,4)>0);
if length(A1) > 0
   B1= find(g1>0);
   A2 = setdiff(A1-1,B1);
   B2 = inter2(:,4)';
   g1(A2) = B2(A2+1);
   B3 = [inter2(:,3)'; max(inter2(:,1:2)')];
   ind(:,A2) = B3(:,A2+1);
end
ind2 = [intra(:,3)'; max(intra(:,1:2)')];
g2 = intra(:,4)';

B2 = intra2(:,4)';
A1 = find(B2>0);
B1= find(g2>0);
if length(A1) > 0 && length(B1) > 0
   B3 = B1(end);
   A3 = setdiff(A1-1,B1);
   B3 = find(A3<B3);
   if length(B3) > 0
      A2 = A3(B3);
      g2(A2) = B2(A2+1);
      B3 = [intra2(:,3)'; max(intra2(:,1:2)')];
      ind2(:,A2) = B3(:,A2+1);
   end
end
ters = ind(1:2,:);

%(3) System Evolution via Partition Energy > Merge Energy
ko = 1;                         % initial number of clusters
kn = 1;                         % initial number of clusters
A3 = 3;                        % falk conditions allowed A3 times
A2 = g2>0;
%Hsep = 0;                 % threshold for seperation
for i = Ns
  if g1(i) > Hsep
     ko = i;
     A3 = 3;
  elseif g2(i) > Hsep
     kn = i;
     A3 = 3;
     ters(:,i)=ind2(:,i);
  else
     A3 = A3-1;
  end
   if A3 < 1
       B1 = min([i+3 Ns(end)]);
       A1 = A2(i:B1);
       if sum(A1) >= 3
          A3 = A3+1;
          ters(:,i-2:i)=ind2(:,i-2:i);
       else
          break; 
       end
    end
end

%(4) plotting results
if kn > ko
   km = [kn ko];
else
   km = ko;
end
g1 = Ns(1);                                   % starting point in plotting
g2=min([km(1)+10 N]);                % ending point in plotting
B3 = NaN*ones(2,g2);
B3(:,Ns) = ters(:,Ns);
ters = B3;