% est_pchirp2weighteddirect...                         Steve Mann, 1992-1993
%
% estimate the 8 projective parameters from a pair of images, simplistic weighting
%
% sign convention: m (or X) down, n (or Y) to right
% properly handles situation when one or both images contain some NaN points
%
% See also: est_affine.m,  est_pseudo.m,  etc...

function [p1,p2,p3,p4,p5,p6,p7,p8]=f(E1,E2);  % two images in, 0,1,or 8 outs

%%%%mask = isnan(E1).*isnan(E2);  % NaN mask: regions where either is undefined
% doesnt work because nanmask makes it still nan

[El,Em,En] = derivatives(E1,E2);
[M,N]= size(El);                    % size of any of the 3 derivative matrices
[n m]=meshdom(1:N,M:-1:1);          % variables ``x'' and ``y''
%disp('est_pseudo: converting from pixels to units on [0,1); doesn''t work???')
%m=m/M;    disp('doesnt work?')
%n=n/N;    disp('doesnt work?')

DERIVATIVES = ...
[...
  [sum2nan(m.*Em                     .* m.*Em) sum2nan(m.*Em                     .*  n.*Em) sum2nan(m.*Em                    .* Em                       ) sum2nan(m.*Em                     .*                     m.*En)...
   sum2nan(m.*Em                     .* n.*En) sum2nan(m.*Em                     .*  En   ) sum2nan(m.*Em                    .* (m.*El-m.*m.*Em-m.*n.*En)) sum2nan(m.*Em                     .* (n.*El-m.*n.*Em-n.*n.*En))];...

  [sum2nan(n.*Em                     .* m.*Em) sum2nan(n.*Em                     .*  n.*Em) sum2nan(n.*Em                    .* Em                       ) sum2nan(n.*Em                     .*                     m.*En)...
   sum2nan(n.*Em                     .* n.*En) sum2nan(n.*Em                     .*  En   ) sum2nan(n.*Em                    .* (m.*El-m.*m.*Em-m.*n.*En)) sum2nan(n.*Em                     .* (n.*El-m.*n.*Em-n.*n.*En))];...

  [sum2nan(Em                        .* m.*Em) sum2nan(Em                        .*  n.*Em) sum2nan(Em                       .* Em                       ) sum2nan(Em                        .*                     m.*En)...
   sum2nan(Em                        .* n.*En) sum2nan(Em                        .*  En   ) sum2nan(Em                       .* (m.*El-m.*m.*Em-m.*n.*En)) sum2nan(Em                        .* (n.*El-m.*n.*Em-n.*n.*En))];...

  [sum2nan(m.*En                     .* m.*Em) sum2nan(m.*En                     .*  n.*Em) sum2nan(m.*En                    .* Em                       ) sum2nan(m.*En                     .*                     m.*En)...
   sum2nan(m.*En                     .* n.*En) sum2nan(m.*En                     .*  En   ) sum2nan(m.*En                    .* (m.*El-m.*m.*Em-m.*n.*En)) sum2nan(m.*En                     .* (n.*El-m.*n.*Em-n.*n.*En))];...

  [sum2nan(n.*En                     .* m.*Em) sum2nan(n.*En                     .* n.*Em) sum2nan(n.*En                     .* Em                       ) sum2nan(n.*En                     .*                     m.*En)...
   sum2nan(n.*En                     .* n.*En) sum2nan(n.*En                     .* En   ) sum2nan(n.*En                     .* (m.*El-m.*m.*Em-m.*n.*En)) sum2nan(n.*En                     .* (n.*El-m.*n.*Em-n.*n.*En))];...

  [sum2nan(En                        .* m.*Em) sum2nan(En                        .* n.*Em) sum2nan(En                        .* Em                       ) sum2nan(En                        .*                     m.*En)...
   sum2nan(En                        .* n.*En) sum2nan(En                        .* En   ) sum2nan(En                        .* (m.*El-m.*m.*Em-m.*n.*En)) sum2nan(En                        .* (n.*El-m.*n.*Em-n.*n.*En))];...

  [sum2nan((m.*El-m.*m.*Em-m.*n.*En) .* m.*Em) sum2nan((m.*El-m.*m.*Em-m.*n.*En) .* n.*Em) sum2nan((m.*El-m.*m.*Em-m.*n.*En) .* Em                       ) sum2nan((m.*El-m.*m.*Em-m.*n.*En) .*                     m.*En)...
   sum2nan((m.*El-m.*m.*Em-m.*n.*En) .* n.*En) sum2nan((m.*El-m.*m.*Em-m.*n.*En) .* En   ) sum2nan((m.*El-m.*m.*Em-m.*n.*En) .* (m.*El-m.*m.*Em-m.*n.*En)) sum2nan((m.*El-m.*m.*Em-m.*n.*En) .* (n.*El-m.*n.*Em-n.*n.*En))];...

  [sum2nan((n.*El-m.*n.*Em-n.*n.*En) .* m.*Em) sum2nan((n.*El-m.*n.*Em-n.*n.*En) .* n.*Em) sum2nan((n.*El-m.*n.*Em-n.*n.*En) .* Em                       ) sum2nan((n.*El-m.*n.*Em-n.*n.*En) .*                     m.*En)...
   sum2nan((n.*El-m.*n.*Em-n.*n.*En) .* n.*En) sum2nan((n.*El-m.*n.*Em-n.*n.*En) .* En   ) sum2nan((n.*El-m.*n.*Em-n.*n.*En) .* (m.*El-m.*m.*Em-m.*n.*En)) sum2nan((n.*El-m.*n.*Em-n.*n.*En) .* (n.*El-m.*n.*Em-n.*n.*En))]...
 ];

Derivatives = ...
  [sum2nan((m.*Em+n.*En-El) .* m.*Em); ...
   sum2nan((m.*Em+n.*En-El) .* n.*Em); ...
   sum2nan((m.*Em+n.*En-El) .*    Em); ...
   sum2nan((m.*Em+n.*En-El) .* m.*En); ...
   sum2nan((m.*Em+n.*En-El) .* n.*En); ...
   sum2nan((m.*Em+n.*En-El) .*    En); ...
   sum2nan((m.*Em+n.*En-El) .* (m.*El-m.*m.*Em-m.*n.*En)); ...
   sum2nan((m.*Em+n.*En-El) .* (n.*El-m.*n.*Em-n.*n.*En)) ... 
  ];                         % vector of derivatives

parameters = DERIVATIVES\Derivatives;

if (nargout == 8) | (nargout == 0)
  % p1 = parameters(1);     % if called with no semicolon, do not want to see any ans=[], etc, so make no reference to p1
  p2 = parameters(2);
  p3 = parameters(3);
  p4 = parameters(4);
  p5 = parameters(5);
  p6 = parameters(6);
  p7 = parameters(7);
  p8 = parameters(8);
end%if
if nargout == 8
  p1 = parameters(1);
end%if

if nargout == 1
  p1 = parameters([1 2 3 4 5 6 7 8]); % no longer permuted
end%if

if nargout == 0
  disp('est_pseudo: 0 output args given; displaying result to screen:')
  disp(sprintf('%g %g %g %g %g %g %g %g',parameters(1),p2,p3,p4,p5,p6,p7,p8))
end%if

if (nargout ~=0) & (nargout ~= 1) & (nargout ~= 8)
  error('est_pseudo: number of output arguments must be 0, 1, or 8')
end%if