% takes about 5 minutes on gigahertz p3
%
% function Finv=cls(J, kvalue)
%
% run least squares solution for 256 points of
% finverse directly from comparagram J and
% known 'k' value, 'kvalue'.
%
% sums up 65536 equations in 256 unknowns into groups of 128 equations, for
% total of 512 equations in 256 unknowns, in an attempt to get better rcond.
%
% see usage example by running 'runme.m'

function Finv=cls(J, kvalue);

   if (nargin < 2)
      error("use: try something like Finv=cls(J03,2)\n");
   end%if
   if (nargout < 1)
      error("use: try something like Finv=cls(J03,2)\n");
   end%if

  % J is the comparagram
  % b is rhs
  % n is the size of J

  % load comparagram
  disp("reading image"); fflush(1);
  %load("Jsum03.mat")

  % parameter initialization
  J=J.';
  % n=size(J,1); % assumes comparagram J is square
  [m,n]=size(J);
  if m~=n
    error('cls_fast: input comparagram must be square\n');
  end%if
  A=zeros(n*n,n);
  Finv=zeros(n,1); % zeros(n,1)*NaN; % just to declare as 0xDEADBEEF
  k=kvalue;
  K=log(k)/log(2);
  b=zeros(n*n,1); % zeros(n*n,1)*NaN; % as 0xDEADBEEF
  
  disp("creating system of equations"); fflush(1);
  for i=1:n
    for j=1:n
      A(i+n*(j-1), i) = A(i+n*(j-1), i) + J(i,j);
      A(i+n*(j-1), j) = A(i+n*(j-1), j) - J(i,j);
      b(i+n*(j-1))    = b(i+n*(j-1))    + K*J(i,j);
    end
    fprintf("%d\r",i); fflush(1);
  end
  Anew=NaN*zeros(512,256);
  bnew=NaN*zeros(512,1);
  for l=1:512
  fprintf("about to combine of 512 output rows, row %d\r",l); fflush(1);
    Anew(l,:)=sum(A(l*128-127:l*128,:));
    bnew(l,:)=sum(b(l*128-127:l*128));
  end%for
  fprintf("\n");
  % now A*Finv=b
  %solve it
  disp("solving..."); fflush(1);
  Finv=Anew\bnew;
  % finv = e.^(Finv);
  % finv = 2.^(Finv);
end%function