alge_odr.m

椭圆拟合的相关介绍与数学运算方法

function [z, a, b, alpha, step] = ...
        alge_odr (X, z, a, b, alpha, show);
%ALGE_ODR       
%       
%       [z, a, b, alpha, step] = ...
%       alge_odr (X, z, a, b, alpha, show);
%       
%       computes the best fit ellipse in parameterform
%       x = z(1) + a  cos(phi-alpha), y = z(2) + b sin(phi-alpha)
%       using the odr algorithm.  %
%       X: given points <X(i,1), X(i,2)>
%       z, a, b, alpha: starting values
%       show: if (show), test output
%
%       z, a, b, alpha: ellipse found
%       step: nof iterations
%
%       See also ODR, ALGE_ODR_HOOK

%       This implicit problem is solved by successive
%       solution of weighted explicit problems.
%       For further discussion see
%       Gill/Murray/Wright
%       "Practical Optimization"
%       Academic Press, New York, 1981.

  OPTI_ONEW = 9;
  OPTI_ONED = 10;

  wate = sqrt(eps);
  awate = sqrt(sqrt(wate))+1e-3;

  c = cos(alpha); s = sin(alpha);
  Q = [c -s; s c];
  A = Q*[1/(a^2) 0; 0 1/(b^2)]*Q';
  beta = [z(1); z(2); A(1,1); A(1,2); A(2,2)];
  delta = zeros(size(X));

  step = 0;
  while (awate >= wate),
    OPTIONS = [[OPTI_ONED, awate]];
    [beta, delta, err, astep] = odr ('alge_odr_hook', ...
      X, zeros(size(X,1),1), beta, OPTIONS, [], [], delta);
    awate = awate^2;
    step = step + astep;
  end

  [zno, a, b, alpha, err] = ellipse_params (...
    [beta(3); 2*beta(4); beta(5); 0; 0; -1]);
  z = [beta(1); beta(2)];

end % alge_odr