Skip to content

Sphere Fit

MAGYC - Benchmark Methods - Sphere Fit

This module contains sphere fit appraches for magnetometer calibration.

Functions:

Name Description
sphere_fit

Standard sphere fit method.

sphere_fit(magnetic_field)

The sphere fit method fits a sphere to a collection of data using a closed form for the solution. With this purpose, propose an optimization problem that seeks to minimize the sum:

\[\sum_i ((x_i-x_c)^2+(y_i-y_c)^2+(z_i-z_c)^2-r^2)^2\]

Where x, y, and z is the data; \(x_c\), \(y_c\), and \(z_c\) are the sphere center; and r is the radius.

The method assumes that points are not in a singular configuration and are real numbers to solve this problem. If you have coplanar data, use a circle fit with svd for determining the plane, recommended Circle Fit (Pratt method), by Nikolai Chernov

Inspired by Alan Jennings, University of Dayton, implementation (source)

Parameters:

Name Type Description Default
magnetic_field ndarray or list

Magnetic field measurements in a 3xN or Nx3 numpy array or list.

required

Returns:

Name Type Description
hard_iron ndarray

Hard iron bias.

calibrated_magnetic_field ndarray

Calibrated magnetic field measurements

Raises:

Type Description
TypeError

If the input is not a numpy array or a list.

ValueError

If the input is not a 3xN or Nx3 numpy array.

Source code in magyc/benchmark_methods/spherefit.py
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
def sphere_fit(magnetic_field: Union[np.ndarray, list]) -> Tuple[np.ndarray, np.ndarray]:
    """
    The sphere fit method fits a sphere to a collection of data using a closed
    form for the solution. With this purpose, propose an optimization problem that
    seeks to minimize the sum:

    $$\\sum_i ((x_i-x_c)^2+(y_i-y_c)^2+(z_i-z_c)^2-r^2)^2$$

    Where x, y, and z is the data; $x_c$, $y_c$, and $z_c$ are the sphere center;
    and r is the radius.

    The method assumes that points are not in a singular configuration and are
    real numbers to solve this problem. If you have coplanar data, use a circle
    fit with svd for determining the plane, recommended [Circle Fit (Pratt method),
    by Nikolai Chernov](http://www.mathworks.com/matlabcentral/fileexchange/22643)

    Inspired by Alan Jennings, University of Dayton, implementation ([source](
    https://www.mathworks.com/matlabcentral/fileexchange/34129-sphere-fit-least-squared))

    Args:
        magnetic_field (numpy.ndarray or list): Magnetic field measurements in a
            3xN or Nx3 numpy array or list.

    Returns:
        hard_iron (numpy.ndarray): Hard iron bias.
        calibrated_magnetic_field (numpy.ndarray): Calibrated magnetic field measurements

    Raises:
        TypeError: If the input is not a numpy array or a list.
        ValueError: If the input is not a 3xN or Nx3 numpy array.
    """
    # Check if the input is a list and convert it to a numpy array
    if isinstance(magnetic_field, list):
        magnetic_field = np.array(magnetic_field)

    # Check if the input is a numpy array
    if not isinstance(magnetic_field, np.ndarray):
        raise TypeError("The input must be a numpy array or a list.")

    # Check if the input is a 3xN or Nx3 numpy array
    if magnetic_field.ndim != 2 or (magnetic_field.shape[0] != 3 and magnetic_field.shape[1] != 3):
        raise ValueError("The input must be a 3xN or Nx3 numpy array.")

    # Force the array to be a Nx3 numpy array
    if magnetic_field.shape[0] == 3:
        magnetic_field = magnetic_field.T

    # Compute magnetic field calibration
    mf = magnetic_field
    a_matrix = np.array(
        [
            [
                np.mean(mf[:, 0] * (mf[:, 0] - np.mean(mf[:, 0]))),
                2 * np.mean(mf[:, 0] * (mf[:, 1] - np.mean(mf[:, 1]))),
                2 * np.mean(mf[:, 0] * (mf[:, 2] - np.mean(mf[:, 2]))),
            ],
            [
                0,
                np.mean(mf[:, 1] * (mf[:, 1] - np.mean(mf[:, 1]))),
                2 * np.mean(mf[:, 1] * (mf[:, 2] - np.mean(mf[:, 2]))),
            ],
            [0, 0, np.mean(mf[:, 2] * (mf[:, 2] - np.mean(mf[:, 2])))],
        ]
    )

    a_matrix = a_matrix + a_matrix.T
    b_matrix = np.array(
        [
            [np.mean((mf[:, 0] ** 2 + mf[:, 1] ** 2 + mf[:, 2] ** 2) * (mf[:, 0] - np.mean(mf[:, 0])))],
            [np.mean((mf[:, 0] ** 2 + mf[:, 1] ** 2 + mf[:, 2] ** 2) * (mf[:, 1] - np.mean(mf[:, 1])))],
            [np.mean((mf[:, 0] ** 2 + mf[:, 1] ** 2 + mf[:, 2] ** 2) * (mf[:, 2] - np.mean(mf[:, 2])))],
        ]
    )

    hard_iron = np.array(np.linalg.lstsq(a_matrix, b_matrix, rcond=None)[0])

    # Calibrate magnetic field
    calibrated_magnetic_field = magnetic_field - hard_iron.flatten()

    return hard_iron.flatten(), calibrated_magnetic_field