function Par = Landau(XY,ParIni) %-------------------------------------------------------------------------- % % Geometric circle fit (minimizing orthogonal distances) by Landau algorithm % M. Landau, "Estimation of a circular arc center and its radius", % Computer Vision, Graphics and Image Processing, Vol. 38, pages 317-326, (1987) % % Input: XY(n,2) is the array of coordinates of n points x(i)=XY(i,1), y(i)=XY(i,2) % ParIni = [a b R] is the initial guess (supplied by user) % % Output: Par = [a b R] is the fitting circle: % center (a,b) and radius R % %-------------------------------------------------------------------------- epsilon=0.00001; % tolerance (small threshold) IterMAX = 800; % maximal number of iterations; with a bad initial guess it may take >100 iterations centroid = mean(XY); % the centroid of the data set X = XY(:,1) - centroid(1); % centering data Y = XY(:,2) - centroid(2); % centering data ParNew = ParIni - [centroid , 0]; % centering the initial guess for iter=1:IterMAX % main loop, each run is one iteration ParOld = ParNew; Dx = X - ParOld(1); Dy = Y - ParOld(2); D = sqrt(Dx.*Dx + Dy.*Dy); ParNew = [-mean(Dx./D) , -mean(Dy./D) , 1]*mean(D); progress = (norm(ParNew-ParOld))/(norm(ParOld)+epsilon); if progress < epsilon, break, end % stopping rule end % the end of the main loop (over iterations) Par = ParOld + [centroid , 0]; % stopping rule end % Landau