Basics of Paleomagnetism and Rock Magnetism

Circle fitting: Classical method of Kasa (1976)

An observation and a circle with radius R. Consider \(n\) observation points \((x_i,y_i)\) (\(i=1,\cdots,n\)) which are fitted by a circle of radius \(R\) at the center \((a,b)\). The "geometric" distance of an observation point from the arc of the circle is denoted as, \[ d_i = r_i -R, \] where \[ r_i = \sqrt{ (x_i - a)^2 + (y_i - b)^2 }. \] The geometric method seeks \(a\), \(b\), and \(R\) which minimize the sum of the squares of the distance \(d_i\) as, \[ \sum_{i=1}^n d_i^2 = \sum_{i=1}^n \left(\sqrt{(x_i-a)^2 + (y_i-b)^2} - R\right)^2. \] However, this is a nonlinear problem which is quite difficult to solve. Here, less accurate but simpler, and widely used, "algebraic" methods are considered. This page introduces the classical method by Kasa (1976) and in the next page the improved one by Taubin (1991) is explained.

In the method of Kasa (1976), the algebraic distance is defined as, \[ f_i = r_i^2 - R^2. \] The quantity to be minimized is, \begin{eqnarray*} {\cal F_K} & = & \sum_{i=1}^n f_i^2, \\ & = & \sum_{i=1}^n (x_i^2 + y_i^2 - 2ax_i - 2by_i + a^2 + b^2 - R^2)^2 \end{eqnarray*} Here, we introduce the following new parameters, \begin{eqnarray} B & = & -2a, \nonumber \\ C & = & -2b, \label{eq01} \\ D & = & a^2 + b^2 - R^2. \nonumber \end{eqnarray} Notation avoiding "\(A\)" is to conform with another theory in the next page. Then, \[ {\cal F_K} = \sum_{i=1}^n (z_i + Bx_i + Cy_i + D)^2, \] where \[ z_i = x_i^2 + y_i^2. \] Setting \( \partial {\cal F_K}/\partial B = \partial {\cal F_K}/\partial C = \partial {\cal F_K}/\partial D =0 \), we obtain \begin{equation} \begin{array}{ccccccc} M_{xx}B & + & M_{xy}C & + & M_xD & = & -M_{xz} \\ M_{xy}B & + & M_{yy}C & + & M_yD & = & -M_{yz} \\ \label{eq02} M_xB & + & M_yC & + & nD & = & -M_z \end{array} \end{equation} where \(M_x = \sum_{i=1}^n x_i, M_{xx} = \sum_{i=1}^n x_i^2, M_{xy} = \sum_{i=1}^n x_i y_i\), etc. After solving \eqref{eq02} by a numerical method such as LU decomposition, \(a\), \(b\), and \(R\) are given by \eqref{eq01}.

Unfortunately this method sometimes gives an answer biased toward a small \(R\). According to Al-Sharadqah and Chernov (2009), \[ f_i = (r_i - R)(r_i + R) = d_i(2R + d_i) \approx 2Rd_i. \] Hence, Kasa method minimizes, \[ {\cal F_K} \approx 4 R^2 \sum_{i=1}^n d_i^2, \] and it often prefers to minimize \(R^2\) rather than \(\sum_{i=1}^n d_i^2\). To improve this, \[ {\cal F} = \frac{1}{4R^2} {\cal F_K} \] may be minimized as shown in the next page.

References:

↑ Top of this page   → Next page