Implement a function to calculate fitting RMSD.

cylinder_fitting/__init__.py

@@ -1,5 +1,6 @@
 from .fitting import fit
 
 from .visualize import show_fit
-
 from .visualize import show_G_distribution
+
+from .analysis import fitting_rmsd

cylinder_fitting/analysis.py

@@ -0,0 +1,9 @@
+import numpy as np
+
+from . import geometry
+
+
+def fitting_rmsd(w_fit, C_fit, r_fit, Xs):
+    '''Calculate the RMSD of fitting.'''
+    return np.sqrt(sum((geometry.point_line_distance(p, C_fit, w_fit) - r_fit) ** 2
+                    for p in Xs) / len(Xs)) 

cylinder_fitting/geometry.py

@@ -1,14 +1,29 @@
 import numpy as np
 
-def rotation_matrix_from_axis_and_angle(u, theta):
-  '''Calculate a rotation matrix from an axis and an angle.'''
 
-  x = u[0]
-  y = u[1]
-  z = u[2]
-  s = np.sin(theta)
-  c = np.cos(theta)
+def normalize(v):
+    '''Normalize a vector based on its 2 norm.'''
+    if 0 == np.linalg.norm(v):
+        return v
+    return v / np.linalg.norm(v)
+
+def rotation_matrix_from_axis_and_angle(u, theta):
+    '''Calculate a rotation matrix from an axis and an angle.'''
+
+    x = u[0]
+    y = u[1]
+    z = u[2]
+    s = np.sin(theta)
+    c = np.cos(theta)
+
+    return np.array([[c + x**2 * (1 - c), x * y * (1 - c) - z * s, x * z * (1 - c) + y * s],
+                     [y * x * (1 - c) + z * s, c + y**2 * (1 - c), y * z * (1 - c) - x * s ],
+                     [z * x * (1 - c) - y * s, z * y * (1 - c) + x * s, c + z**2 * (1 - c) ]])
 
-  return np.array([[c + x**2 * (1 - c), x * y * (1 - c) - z * s, x * z * (1 - c) + y * s],
-                   [y * x * (1 - c) + z * s, c + y**2 * (1 - c), y * z * (1 - c) - x * s ],
-                   [z * x * (1 - c) - y * s, z * y * (1 - c) + x * s, c + z**2 * (1 - c) ]])
+def point_line_distance(p, l_p, l_v):
+    '''Calculate the distance between a point and a line defined
+    by a point and a direction vector.
+    '''
+    l_v = normalize(l_v)
+    u = p - l_p
+    return np.linalg.norm(u - np.dot(u, l_v) * l_v)

setup.py

@@ -4,11 +4,11 @@
 
 setup(
     name='cylinder_fitting',
-    version='1.1.3',
+    version='1.1.4',
     author='Xingjie Pan',
     author_email='[email protected]',
     url='https://github.com/xingjiepan/cylinder_fitting',
-    download_url='https://github.com/xingjiepan/cylinder_fitting/archive/1.1.3.tar.gz',
+    download_url='https://github.com/xingjiepan/cylinder_fitting/archive/1.1.4.tar.gz',
     keywords = ['geometry', 'fitting-algorithm'],
     packages=[
         'cylinder_fitting',

tests/test_fitting.py

@@ -8,6 +8,7 @@
 from cylinder_fitting import show_fit
 from cylinder_fitting import show_G_distribution
 from cylinder_fitting import geometry
+from cylinder_fitting import fitting_rmsd
 
 
 def make_points_on_a_cylinder(theta, phi, C, r, N):
@@ -37,6 +38,8 @@ def test_fit():
 
     #show_G_distribution(data)
 
+    print("Fitting RMSD =", fitting_rmsd(w_fit, C_fit, r_fit, data))
+
     show_fit(w_fit, C_fit, r_fit, data)
 
     assert(np.absolute(1 - np.absolute(np.dot(w, w_fit))) < 1e-4)