Modifica di Asse Trasversale Cerniera (sezione)

==Trace Resolution and Point Fitting==

The approximate calculation presented in the previous section considers only the initial and final points of a movement that follows an arc of a circle, along with the segment distance <math>h</math> (or, to remain consistent with other calculations, three points on the circumference). Modern measuring devices can record the positions of mandibular movement at very short time intervals, but even in this case, not everything is perfect. The problem shifts to the sampling capabilities of the phenomenon and the recording channels (topics that will be addressed in due course). Therefore, not only the initial and final points are known, but also many other points along the trace can be measured. It is assumed that this additional information can improve the accuracy of determining the hinge axis. Under ideal conditions, all these points should lie on a segment of a circle around the hinge axis. By fitting a circle to these points, the hinge axis can be determined with precision. However, the imprecision of each measurement introduces uncertainty at each point (Fig. 7). This "noise" affects the detection of the hinge axis.

The best algorithms for reconstructing rotational movement are those that use a least-squares optimization approach to fit a circle. The most effective methods are the Taubin fitting<ref name="Taubin1991">Taubin G. Estimation of Planar Curves, Surfaces and Nonplanar Space Curves defined by Implicit Equations, with Applications to Edge and Range Image Segmentation. IEEE Trans 1991. PAMI;13:1115–1138.</ref> or Pratt fitting methods<ref name="Pratt1987">Pratt V. Direct least-squares fitting of algebraic surfaces. Computer Graphics 1987;21:145–152.</ref>, both available in MATLAB. In our case, we used the Levenberg-Marquardt method.

'''Levenberg-Marquardt Method:''' This method is extremely versatile and can be applied to a wide range of non-linear fitting problems, not only for fitting circles but also other curves or surfaces. It is known for its ability to converge even when the initial estimate is not particularly close to the optimal solution. This method is a compromise between the gradient method and the Gauss-Newton method, balancing efficiency and stability. The Levenberg-Marquardt method is well-supported by many numerical optimization libraries, such as SciPy in Python, making it a practical and accessible choice for most users. Its implementation is relatively simple and well-documented, making it suitable for anyone wishing to perform curve or surface fitting quickly. Since the problem of fitting a circle is non-linear (due to the quadratic dependence of the coordinates on the circle's center), the Levenberg-Marquardt method is well-equipped to handle these complexities.

'''Taubin Method:''' The Taubin method, on the other hand, is particularly recommended when a circle fitting that is less sensitive to noise and measurement errors is desired. This method minimizes the orthogonal error, which often better represents real errors in experimentally collected data. If the priority were to robustly minimize the error, regardless of the direction of fitting, the Taubin method might be more appropriate.

The choice between Levenberg-Marquardt and Taubin depends on the specific needs of the problem and the characteristics of the data. For general non-linear fitting problems with good initial estimates and the need for a quick and practical solution, Levenberg-Marquardt is a solid choice. However, for more specialized applications, such as precise circle or conic fitting, the Taubin method may offer advantages in terms of accuracy and robustness.

In conclusion, the simulation showed that the fitting method used was able to estimate both the center and radius of the circle with good accuracy, despite the presence of noise in the data. The overall error is low, indicating that the effect of noise on the final result was minimal. An error of <math>\cong 0.63 </math> <small><math display="inline">mm</math></small> is rather small, which implies that the fitting was accurate despite the noise introduced into the data. This error represents the deviation of the calculated circle's center from the true center.

Figure 7 illustrates the dependency of the hinge axis determination on the distance of the measurement position from the center of rotation. It also highlights that a measurement close to the TMJ (with a small radius) does not offer greater accuracy with the circle fitting method compared to a measurement taken further from the TMJ or closer to the jaw (e.g., at the incisal point). On the contrary, radii smaller than <math>20</math> <small><math display="inline">mm</math></small> significantly increase the imprecision in determining the hinge axis. This result is positive in a clinical context, as it indicates that, even in the presence of noise (as can happen in real measurements), it is possible to accurately estimate the <math>_tHA</math> and the associated radius, which are crucial for patient diagnosis and rehabilitative treatment.

<Center>
<gallery widths="350" heights="282" perrow="2" mode="slideshow">
File:Official HA 1.jpg| Figura 7: 
</gallery>
</Center>

===Mathematical Formalism: Fitting Error=== 
To illustrate the concept of error due to noise in determining the points of the curve, the preliminary mathematical steps are as follows (refer to Figure 7, represented in Geogebra).

The script starts by defining the original circle center with the coordinates:

<math>O_{HA} = \begin{pmatrix} 11.55 \ 61.14\end{pmatrix}</math>

This center represents the central point around which a circle will be constructed.

The radius of the circle is set to <math>r = 26.38</math> mm. This value represents the distance from the center to the points on the edge of the circle. Using the center and radius, 10 points are uniformly distributed along the circumference of the circle. The angles <math>\theta</math> used to position the points are distributed from <math>0</math> to <math>2\pi</math>. For each angle, the coordinates of the points are calculated using the parametric equations of the circle:

<math>x = x_c + r \cdot \cos(\theta)</math>

<math>y = y_c + r \cdot \sin(\theta)</math>

A function called add_noise_to_points is used to add noise to the generated points. The noise is generated from a normal (or Gaussian) distribution with a mean of zero and a standard deviation of <math>0.1</math> units. This step simulates the measurement errors that may occur in practice. The points disturbed by noise are stored in the noisy_points variable.

The function residuals is defined to calculate the residuals, that is, the difference between the distance of each noisy point from the estimated center and the estimated radius. This function is essential for optimization since the circle fitting aims to minimize these residuals:

<math>R_i = \sqrt{(x_i - x_c)^2 + (y_i - y_c)^2} - r</math>

where <math>R_i</math> is the residual for the <math>i</math>-th point.

An initial estimate of the noisy circle's center is made by calculating the mean of the noisy points' coordinates:

<math>x_m = \frac{1}{n} \sum_{i=1}^{n} x_i</math>, <math>y_m = \frac{1}{n} \sum_{i=1}^{n} y_i</math>

This average is used as a starting point for the fitting.

'''Initial Radius'''

The initial radius is estimated as the mean distance between each noisy point and the estimated center:

<math>r_m = \frac{1}{n} \sum_{i=1}^{n} \sqrt{(x_i - x_m)^2 + (y_i - y_m)^2}</math>

This approximated value is used as the initial radius for the optimization.

'''Circle Fitting with Levenberg-Marquardt Method'''

*Fitting Optimization:

The least_squares function from the SciPy library is used to perform the circle fitting. This method minimizes the residuals calculated earlier, aiming to find the optimal values for the center coordinates <math>(x_c, y_c)</math> and the radius <math>r</math>. The optimization results are returned as the calculated center <math>(x_c, y_c)</math> and the radius <math>r</math>.

'''Error Calculation'''

*Error between Centers:

The distance between the calculated center and the original center is computed. This distance represents the localization error, indicating how far the calculated center deviates from the true center due to noise:

<math>\text{error} = \sqrt{(x_c - x_{\text{original}})^2 + (y_c - y_{\text{original}})^2}</math>

The error is then converted into millimeters if the original units were in centimeters, for greater interpretative precision.

===='''Script Python: Errore di fitting'''====
<syntaxhighlight lang="python">
import numpy as np
from scipy.optimize import least_squares
import matplotlib.pyplot as plt

# Original data provided
original_center = np.array([11.56, 61.21])
radius = 26.38
angles = np.linspace(0, 2 * np.pi, 10)
points = np.array([
    original_center + radius * np.array([np.cos(angle), np.sin(angle)])
    for angle in angles
])

# Function to add noise to the points
def add_noise_to_points(points, noise_level):
    noisy_points = points + np.random.normal(0, noise_level, points.shape)
    return noisy_points

# Add noise to the points
noise_level = 0.1  # Reduced noise level
noisy_points = add_noise_to_points(points, noise_level)

# Residual function for non-linear fitting
def residuals(c, points):
    xc, yc, r = c
    Ri = np.sqrt((points[:, 0] - xc)**2 + (points[:, 1] - yc)**2)
    return Ri - r

# Calculate the center and radius using noisy points
x_m = np.mean(noisy_points[:, 0])
y_m = np.mean(noisy_points[:, 1])
r_m = np.mean(np.sqrt((noisy_points[:, 0] - x_m)**2 + (noisy_points[:, 1] - y_m)**2))
initial_guess = [x_m, y_m, r_m]

# Fitting the circle using the Levenberg-Marquardt method
result = least_squares(residuals, initial_guess, args=(noisy_points,))
xc, yc, r = result.x

calculated_center = np.array([xc, yc])

# Calculate the error between the centers
error = np.linalg.norm(calculated_center - original_center)
error_mm = error * 10  # Convert to millimeters (assuming the points might be in centimeters)

# Verify the calculated center
print(f'Calculated Center: {calculated_center}')
print(f'Radius: {r:.4f}')
print(f'Error: {error_mm:.4f} mm')

# Visualization of the fitted circle
theta = np.linspace(0, 2 * np.pi, 100)
x_fit = calculated_center[0] + r * np.cos(theta)
y_fit = calculated_center[1] + r * np.sin(theta)

plt.figure(figsize=(8, 8))
plt.scatter(points[:, 0], points[:, 1], color='b', label='Original Points')
plt.scatter(noisy_points[:, 0], noisy_points[:, 1], color='r', label='Noisy Points')
plt.scatter(original_center[0], original_center[1], color='r', marker='x', s=100, label='Original Center')
plt.plot(x_fit, y_fit, 'g--', label='Fitted Circle')
plt.scatter(calculated_center[0], calculated_center[1], color='g', label='Calculated Center')
plt.xlabel('X')
plt.ylabel('Y')
plt.legend()
plt.title('Circle Fitting with Noisy Points')
plt.grid(True)
plt.show()

</syntaxhighlight>