function Par = RiemannSWFLa(XY) %-------------------------------------------------------------------------- % % Riemann circle fit, the SWFL version % A. Strandlie, J. Wroldsen, R. Fruhwirth, and B. Lillekjendlie, % "Particle tracks fitted on the Riemann sphere", % Computer Physics Commun., Vol. 131, pp. 95-108, (2000) % % Input: XY(n,2) is array of coordinates of n points % % Output: Par = [a b R] is the fitting circle: % center (a,b) and radius R % % This version is modified from the original one; it uses algebraic parameters % %-------------------------------------------------------------------------- n = size(XY,1); % number of data points centroid = mean(XY); X = XY(:,1) - centroid(1); Y = XY(:,2) - centroid(2); Z = X.*X + Y.*Y; factor = 2*sqrt(sum(Z)/n); X = X/factor; Y = Y/factor; Z = Z/factor/factor; Zp = Z+1; Zm = Z-1; Zpp = Zp'*Zp; Zpm = Zp'*Zm; ZpX = Zp'*X; ZpY = Zp'*Y; A1 = (Zm*Zpp-Zp*Zpm)/2; A2 = X*Zpp - Zp*ZpX; A3 = Y*Zpp - Zp*ZpY; AAA = [A1 A2 A3]; [U,S,V]=svd(AAA,0); P = V(:,3); Q = -2*(P(1)*Zpm/2 + P(2)*ZpX + P(3)*ZpY)/Zpp; A = (P(1) + Q)/2; D = (Q - P(1))/2; Par = [-P(2)/A/2 -P(3)/A/2 sqrt(P(2:3)'*P(2:3)-4*A*D)/2/abs(A)]; Par = Par*factor; Par(1:2) = Par(1:2) + centroid; end % RiemannSWFLa