...
 
Commits (30)
......@@ -83,37 +83,6 @@ test Windows:
- shared-win
deploy docs to staging:
stage: deploy
script:
# Note: the group/username directory must already exist on the server before calling this command
- rsync -avz -e'ssh -v' --numeric-ids --delete doc/build/html/* staging_docs_rsync:WEISS/SoftwareRepositories/scikit-surgerycore 2>&1
tags:
- docs-staging
environment:
name: staging
url: http://weisslab-lin.cs.ucl.ac.uk/staging/WEISS/SoftwareRepositories/scikit-surgerycore
only:
- master
dependencies:
- build docs
deploy docs to production:
stage: deploy
script:
# Note: the group/username directory must already exist on the server before calling this command
- rsync -avz -e'ssh -v' --numeric-ids --delete doc/build/html/* production_docs_rsync:WEISS/SoftwareRepositories/scikit-surgerycore 2>&1
tags:
- docs-production
environment:
name: production
only:
- public
dependencies:
- build docs
deploy pip to production:
stage: deploy
......
......@@ -6,7 +6,7 @@ scikit-surgerycore
:width: 128px
:target: https://weisslab.cs.ucl.ac.uk/WEISS/SoftwareRepositories/SNAPPY/scikit-surgerycore
.. image:: https://weisslab.cs.ucl.ac.uk/WEISS/SoftwareRepositories/SNAPPY/scikit-surgerycore/badges/master/build.svg
.. image:: https://weisslab.cs.ucl.ac.uk/WEISS/SoftwareRepositories/SNAPPY/scikit-surgerycore/badges/master/pipeline.svg
:target: https://weisslab.cs.ucl.ac.uk/WEISS/SoftwareRepositories/SNAPPY/scikit-surgerycore/pipelines
:alt: GitLab-CI test status
......
# -*- coding: utf-8 -*-
"""Functions for point based registration using Orthogonal Procrustes."""
""" Functions for point based registration using Orthogonal Procrustes. """
import numpy as np
import sksurgerycore.algorithms.vector_math as vm
def validate_procrustes_inputs(fixed, moving):
......@@ -65,3 +66,59 @@ def compute_fre(fixed, moving, rotation, translation):
mean_squared_error = sum_squared_error / fixed.shape[0]
fre = np.sqrt(mean_squared_error)
return fre
def compute_tre_from_fle(fiducials, mean_fle_squared, target_point):
"""
Computes an estimation of TRE from FLE and a list of fiducial locations.
See:
`Fitzpatrick (1998), equation 46 <http://dx.doi.org/10.1109/42.736021>`_.
:param fiducials: Nx3 ndarray of fiducial points
:param mean_fle_squared: expected (mean) FLE squared
:param target_point: a point for which to compute TRE.
:return: mean TRE squared
"""
# pylint: disable=literal-comparison
if not isinstance(fiducials, np.ndarray):
raise TypeError("fiducials is not a numpy array'")
if not fiducials.shape[1] == 3:
raise ValueError("fiducials should have 3 columns")
if fiducials.shape[0] < 3:
raise ValueError("fiducials should have at least 3 rows")
if not isinstance(target_point, np.ndarray):
raise TypeError("target_point is not a numpy array'")
if not target_point.shape[1] == 3:
raise ValueError("target_point should have 3 columns")
if not target_point.shape[0] == 1:
raise ValueError("target_point should have 1 row")
number_of_fiducials = fiducials.shape[0]
centroid = np.mean(fiducials, axis=0)
covariance = np.cov(fiducials.T)
assert covariance.shape[0] == 3
assert covariance.shape[1] == 3
_, eigen_vectors_matrix = np.linalg.eig(covariance)
f_array = np.zeros(3)
for axis_index in range(3):
sum_f_k_squared = 0
for fiducial_index in range(fiducials.shape[0]):
f_k = vm.distance_from_line(centroid,
eigen_vectors_matrix[axis_index],
fiducials[fiducial_index])
sum_f_k_squared = sum_f_k_squared + f_k * f_k
f_k_rms = np.sqrt(sum_f_k_squared / number_of_fiducials)
f_array[axis_index] = f_k_rms
inner_sum = 0
for axis_index in range(3):
d_k = vm.distance_from_line(centroid,
eigen_vectors_matrix[axis_index],
target_point)
inner_sum = inner_sum + (d_k * d_k / f_array[axis_index])
mean_tre_squared = (mean_fle_squared / fiducials.shape[0]) * \
(1 + (1./3.) * inner_sum)
return mean_tre_squared
# -*- coding: utf-8 -*-
"""Functions for pivot calibration."""
""" Functions for pivot calibration. """
import random
import numpy as np
# pylint: disable=literal-comparison
def pivot_calibration(matrices4x4):
def pivot_calibration(tracking_matrices):
"""
Performs Pivot Calibration and returns Residual Error.
Performs Pivot Calibration, using Algebraic One Step method,
and returns Residual Error.
1. matricesNx4x4 and moving must be numpy array
2. matricesNx4x4 and moving should have 4 columns
3. matricesNx4x4 and moving should have at least 4 rows
See `Yaniv 2015 <https://dx.doi.org/10.1117/12.2081348>`_.
:param matrices4x4: point set, N x 4 x 4 ndarray
:returns: residual_error
:param tracking_matrices: N x 4 x 4 ndarray, of tracking matrices.
:returns: pointer offset, pivot point and RMS Error about centroid of pivot.
:raises: TypeError, ValueError
"""
if not isinstance(matrices4x4, np.ndarray):
raise TypeError("matrices4x4 is not a numpy array'")
if not isinstance(tracking_matrices, np.ndarray):
raise TypeError("tracking_matrices is not a numpy array'")
if not tracking_matrices.shape[1] == 4:
raise ValueError("tracking_matrices should have 4 rows per matrix")
if not tracking_matrices.shape[2] == 4:
raise ValueError("tracking_matrices should have 4 columns per matrix")
if not matrices4x4.shape[1] == 4: # pylint: disable=literal-comparison
raise ValueError("matrices4x4 should have 4 columns per matrix")
number_of_matrices = tracking_matrices.shape[0]
number_of_matrices = len(matrices4x4)
# See equation in section 2.1.2 of Yaniv 2015.
# Ax = b.
size_a = 3 * number_of_matrices, 6
a_values = np.zeros(size_a, dtype=np.float64)
size_x = 6, 1
x_values = np.zeros(size_x, dtype=np.float64)
size_b = 3 * number_of_matrices, 1
b_values = np.zeros(size_b, dtype=np.float64)
for i in range(number_of_matrices):
b_values[i * 3 + 0, 0] = -1 * matrices4x4[i, 0, 3]
b_values[i * 3 + 1, 0] = -1 * matrices4x4[i, 1, 3]
b_values[i * 3 + 2, 0] = -1 * matrices4x4[i, 2, 3]
a_values[i * 3 + 0, 0] = matrices4x4[i, 0, 0]
a_values[i * 3 + 1, 0] = matrices4x4[i, 1, 0]
a_values[i * 3 + 2, 0] = matrices4x4[i, 2, 0]
a_values[i * 3 + 0, 1] = matrices4x4[i, 0, 1]
a_values[i * 3 + 1, 1] = matrices4x4[i, 1, 1]
a_values[i * 3 + 2, 1] = matrices4x4[i, 2, 1]
a_values[i * 3 + 0, 2] = matrices4x4[i, 0, 2]
a_values[i * 3 + 1, 2] = matrices4x4[i, 1, 2]
a_values[i * 3 + 2, 2] = matrices4x4[i, 2, 2]
# Column vector containing -1 * translation from each tracking matrix.
b_values[i * 3 + 0, 0] = -1 * tracking_matrices[i, 0, 3]
b_values[i * 3 + 1, 0] = -1 * tracking_matrices[i, 1, 3]
b_values[i * 3 + 2, 0] = -1 * tracking_matrices[i, 2, 3]
# A contains rotation matrix from each tracking matrix.
a_values[i * 3 + 0, 0] = tracking_matrices[i, 0, 0]
a_values[i * 3 + 1, 0] = tracking_matrices[i, 1, 0]
a_values[i * 3 + 2, 0] = tracking_matrices[i, 2, 0]
a_values[i * 3 + 0, 1] = tracking_matrices[i, 0, 1]
a_values[i * 3 + 1, 1] = tracking_matrices[i, 1, 1]
a_values[i * 3 + 2, 1] = tracking_matrices[i, 2, 1]
a_values[i * 3 + 0, 2] = tracking_matrices[i, 0, 2]
a_values[i * 3 + 1, 2] = tracking_matrices[i, 1, 2]
a_values[i * 3 + 2, 2] = tracking_matrices[i, 2, 2]
# A contains -I for each tracking matrix.
a_values[i * 3 + 0, 3] = -1
a_values[i * 3 + 1, 3] = 0
a_values[i * 3 + 2, 3] = 0
......@@ -63,13 +73,10 @@ def pivot_calibration(matrices4x4):
u_values, s_values, v_values = np.linalg.svd(a_values, full_matrices=False)
c_values = np.dot(u_values.T, b_values)
w_values = np.linalg.solve(np.diag(s_values), c_values)
w_values = np.dot(np.diag(1 / s_values), c_values)
x_values = np.dot(v_values.T, w_values)
# Back substitution
# Calculating the rank
# Calculating the rank, and removing close to zero singular values.
rank = 0
for i, item in enumerate(s_values):
if item < 0.01:
......@@ -77,34 +84,97 @@ def pivot_calibration(matrices4x4):
if item != 0:
rank += 1
if rank < 6: # pylint: disable=literal-comparison
if rank < 6:
raise ValueError("PivotCalibration: Failed. Rank < 6")
# Residual Matrix
# Compute RMS error.
residual_matrix = (np.dot(a_values, x_values) - b_values)
residual_error = 0.0
for i in range(number_of_matrices * 3):
residual_error = residual_error + \
np.dot(residual_matrix[i, 0],
residual_matrix[i, 0])
residual_error = np.sum(residual_matrix * residual_matrix)
residual_error = residual_error / float(number_of_matrices * 3)
residual_error = np.sqrt(residual_error)
# Output
# MakeIdentity matrix
# output_matrix = np.identity(4)
#
# output_matrix[0, 3] = x_values[0, 0]
# output_matrix[1, 3] = x_values[1, 0]
# output_matrix[2, 3] = x_values[2, 0]
# print("pivotCalibration=(", x_values[3, 0], ","
# , x_values[4, 0], ",", x_values[5, 0],
# "),residual=", residual_error)
return x_values, residual_error
# print(x_values)
return np.matrix(x_values), residual_error
def pivot_calibration_with_ransac(tracking_matrices,
number_iterations,
error_threshold,
concensus_threshold,
early_exit=False
):
"""
Written as an exercise for implementing RANSAC.
:param tracking_matrices: N x 4 x 4 ndarray, of tracking matrices.
:param number_iterations: the number of iterations to attempt.
:param error_threshold: distance in millimetres from pointer position
:param concensus_threshold: the minimum percentage of inliers to finish
:param early_exit: If True, returns model as soon as thresholds are met
:returns: pointer offset, pivot point and RMS Error about centroid of pivot.
:raises: TypeError, ValueError
"""
if number_iterations < 1:
raise ValueError("The number of iterations must be > 1")
if error_threshold < 0:
raise ValueError("The error threshold must be a positive distance.")
if concensus_threshold < 0 or concensus_threshold > 1:
raise ValueError("The concensus threshold must be [0-1] as percentage")
if not isinstance(tracking_matrices, np.ndarray):
raise TypeError("tracking_matrices is not a numpy array'")
number_of_matrices = tracking_matrices.shape[0]
population_of_indices = range(number_of_matrices)
minimum_matrices_required = 3
highest_number_of_inliers = -1
best_model = None
best_rms = -1
for iter_counter in range(number_iterations):
indexes = random.sample(population_of_indices,
minimum_matrices_required)
sample = tracking_matrices[indexes]
try:
model, _ = pivot_calibration(sample)
except ValueError:
print("RANSAC, iteration " + str(iter_counter) + ", failed.")
continue
# Need to evaluate the number of inliers.
# Slow, but it's written as a teaching exercise.
world_point = model[3:6]
number_of_inliers = 0
inlier_indices = []
for matrix_counter in range(number_of_matrices):
offset = np.vstack((model[0:3], 1))
transformed_point = tracking_matrices[matrix_counter] @ offset
diff = world_point - transformed_point[0:3]
norm = np.linalg.norm(diff)
if norm < error_threshold:
number_of_inliers = number_of_inliers + 1
inlier_indices.append(matrix_counter)
percentage_inliers = number_of_inliers / number_of_matrices
# Keep the best model so far, based on the highest number of inliers.
if percentage_inliers > concensus_threshold \
and number_of_inliers > highest_number_of_inliers:
highest_number_of_inliers = number_of_inliers
inlier_matrices = tracking_matrices[inlier_indices]
best_model, best_rms = pivot_calibration(inlier_matrices)
# Early exit condition, as soon as we find model with enough fit.
if percentage_inliers > concensus_threshold and early_exit:
return best_model, best_rms
if best_model is None:
raise ValueError("Failed to find a model using RANSAC.")
print("RANSAC Pivot, from " + str(number_of_matrices)
+ " matrices, used " + str(highest_number_of_inliers)
+ " matrices, with error threshold = " + str(error_threshold)
+ " and consensus threshold = " + str(concensus_threshold)
)
return best_model, best_rms
......@@ -6,7 +6,6 @@ import numpy as np
from sksurgerycore.algorithms.errors \
import validate_procrustes_inputs, compute_fre
# pylint: disable=invalid-name, line-too-long
......@@ -53,8 +52,9 @@ def orthogonal_procrustes(fixed, moving):
# Note: numpy factors h = u * np.diag(s) * v
svd = np.linalg.svd(H)
# Arun equation 13
X = np.matmul(svd[2].transpose(), svd[0].transpose())
# Replace Arun Equation 13 with Fitzpatrick, chapter 8, page 470,
# to avoid reflections, see issue #19
X = _fitzpatricks_X(svd)
# Arun step 5, after equation 13.
det_X = np.linalg.det(X)
......@@ -67,7 +67,6 @@ def orthogonal_procrustes(fixed, moving):
" and no singular values are close enough to zero")
if det_X < 0 and np.any(np.isclose(svd[1], np.zeros((3, 1)))):
# Implement 2a in section VI in Arun paper.
v_prime = svd[2].transpose()
v_prime[0][2] *= -1
......@@ -85,3 +84,18 @@ def orthogonal_procrustes(fixed, moving):
fre = compute_fre(fixed, moving, R, T)
return R, T, fre
def _fitzpatricks_X(svd):
"""This is from Fitzpatrick, chapter 8, page 470.
it's used in preference to Arun's equation 13,
X = np.matmul(svd[2].transpose(), svd[0].transpose())
to avoid reflections.
"""
VU = np.matmul(svd[2].transpose(), svd[0])
detVU = np.linalg.det(VU)
diag = np.eye(3, 3)
diag[2][2] = detVU
X = np.matmul(svd[2].transpose(), np.matmul(diag, svd[0].transpose()))
return X
# -*- coding: utf-8 -*-
""" Various small maths utilities. """
import numpy as np
def distance_from_line(p_1, p_2, p_3):
"""
Computes distance of a point p_3, from a line defined by p_1 and p_2.
See `here <https://en.wikipedia.org/wiki/Distance_from_a_point_to_a_line>`_.
:return: euclidean distance
"""
# I want notation same as wikipedia page, so disabling warning.
# pylint: disable=invalid-name
n = p_2 - p_1
n = n / np.linalg.norm(n)
a_minus_p = p_1 - p_3
vector_to_line = a_minus_p - (np.dot(a_minus_p, n) * n)
distance = np.linalg.norm(vector_to_line)
return distance
......@@ -13,8 +13,10 @@ Design principles:
| its up to the consumer to know where to find the data.
"""
import os
import json
import copy
import sksurgerycore.utilities.file_utilities as fu
import sksurgerycore.utilities.validate_file as f
......@@ -30,18 +32,37 @@ class ConfigurationManager:
def __init__(self, file_name,
write_on_setter=False
):
""" Constructor. """
f.validate_is_file(file_name)
abs_file = fu.get_absolute_path_of_file(file_name)
f.validate_is_file(abs_file)
if write_on_setter:
f.validate_is_writable_file(file_name)
f.validate_is_writable_file(abs_file)
with open(file_name, "r") as read_file:
with open(abs_file, "r") as read_file:
self.config_data = json.load(read_file)
self.file_name = file_name
self.file_name = abs_file
self.write_on_setter = write_on_setter
def get_file_name(self):
"""
Returns the absolute filename that was used when
the ConfigurationManager was created.
:return: str absolute file name
"""
return self.file_name
def get_dir_name(self):
"""
Returns the directory name of the file that was used when
creating the ConfigurationManager.
:return: str dir name
"""
return os.path.dirname(self.file_name)
def get_copy(self):
""" Returns a copy of the data read from file.
......
# -*- coding: utf-8 -*-
"""Functions to load MITK's mps point set file. """
import xml.etree.ElementTree as ET
import numpy as np
def load_mps(file_name):
"""
Load a pointset from a .mps file. For now, just loads points,
without geometry information.
:param file_name: string representing file path.
:return: ids (length N), points (Nx3)
"""
tree = ET.parse(file_name)
root = tree.getroot()
points = root[1]
time_series = points[0]
tmp = {}
for point in time_series.findall('point'):
i = point.find('id').text
x_c = point.find('x').text
y_c = point.find('y').text
z_c = point.find('z').text
tmp[int(i)] = (float(x_c), float(y_c), float(z_c))
ids = np.zeros(len(tmp.keys()), dtype=int)
result = np.zeros((len(tmp.keys()), 3))
counter = 0
for k in tmp:
ids[counter] = k
result[counter][0] = tmp[k][0]
result[counter][1] = tmp[k][1]
result[counter][2] = tmp[k][2]
counter = counter + 1
return ids, result
......@@ -63,7 +63,7 @@ def construct_rz_matrix(angle, is_in_radians=True):
return rot_z
# pylint: disable=too-many-branches
def construct_rotm_from_euler(
angle_a, angle_b, angle_c,
sequence, is_in_radians=True):
......@@ -98,6 +98,18 @@ def construct_rotm_from_euler(
rot_c = construct_rz_matrix(angle_c, is_in_radians)
break
if sequence == 'zxy':
rot_a = construct_rz_matrix(angle_a, is_in_radians)
rot_b = construct_rx_matrix(angle_b, is_in_radians)
rot_c = construct_ry_matrix(angle_c, is_in_radians)
break
if sequence == 'zyx':
rot_a = construct_rz_matrix(angle_a, is_in_radians)
rot_b = construct_ry_matrix(angle_b, is_in_radians)
rot_c = construct_rx_matrix(angle_c, is_in_radians)
break
if sequence == 'xyx':
rot_a = construct_rx_matrix(angle_a, is_in_radians)
rot_b = construct_ry_matrix(angle_b, is_in_radians)
......@@ -110,6 +122,18 @@ def construct_rotm_from_euler(
rot_c = construct_rx_matrix(angle_c, is_in_radians)
break
if sequence == 'xyz':
rot_a = construct_rx_matrix(angle_a, is_in_radians)
rot_b = construct_ry_matrix(angle_b, is_in_radians)
rot_c = construct_rz_matrix(angle_c, is_in_radians)
break
if sequence == 'xzy':
rot_a = construct_rx_matrix(angle_a, is_in_radians)
rot_b = construct_rz_matrix(angle_b, is_in_radians)
rot_c = construct_ry_matrix(angle_c, is_in_radians)
break
if sequence == 'yxy':
rot_a = construct_ry_matrix(angle_a, is_in_radians)
rot_b = construct_rx_matrix(angle_b, is_in_radians)
......@@ -122,15 +146,15 @@ def construct_rotm_from_euler(
rot_c = construct_ry_matrix(angle_c, is_in_radians)
break
if sequence == 'xyz':
rot_a = construct_rx_matrix(angle_a, is_in_radians)
rot_b = construct_ry_matrix(angle_b, is_in_radians)
if sequence == 'yxz':
rot_a = construct_ry_matrix(angle_a, is_in_radians)
rot_b = construct_rx_matrix(angle_b, is_in_radians)
rot_c = construct_rz_matrix(angle_c, is_in_radians)
break
if sequence == 'zyx':
rot_a = construct_rz_matrix(angle_a, is_in_radians)
rot_b = construct_ry_matrix(angle_b, is_in_radians)
if sequence == 'yzx':
rot_a = construct_ry_matrix(angle_a, is_in_radians)
rot_b = construct_rz_matrix(angle_b, is_in_radians)
rot_c = construct_rx_matrix(angle_c, is_in_radians)
break
......
......@@ -236,12 +236,14 @@ class TransformManager:
return
for candidate in candidates:
if candidate in list_of_nodes:
continue
else:
list_of_nodes.append(candidate)
self.__get_list(candidate, after, list_of_nodes)
if list_of_nodes[-1] == after:
break
else:
list_of_nodes.pop()
list_of_nodes.append(candidate)
self.__get_list(candidate, after, list_of_nodes)
if list_of_nodes[-1] == after:
break
list_of_nodes.pop()
# -*- coding: utf-8 -*-
""" File processing utils. """
import os
def get_absolute_path_of_file(file_name, dir_name=None):
"""
Filenames in our .json config could be absolute or relative to the
current working dir. This method tries to find the valid, full file path.
:param file_name:
:param dir_name: prefix, for example, the dirname of our .json file.
:return: absolute path name of file if found, otherwise None.
"""
if not file_name:
raise ValueError("Empty file_name.")
file = None
if os.path.isabs(file_name):
if os.path.isfile(file_name):
file = file_name
elif dir_name is not None:
joined = os.path.join(dir_name, file_name)
if os.path.isfile(joined):
file = joined
elif os.getcwd() is not None:
joined = os.path.join(os.getcwd(), file_name)
if os.path.isfile(joined):
file = joined
return file
......@@ -7,7 +7,7 @@ import sksurgerycore.algorithms.pivot as p
from glob import glob
def test_empty_matrices4x4():
def test_empty_matrices():
with pytest.raises(TypeError):
p.pivot_calibration(None)
......@@ -19,9 +19,9 @@ def test_rank_lt_six():
file_names = glob('tests/data/PivotCalibration/1378476416922755200.txt')
arrays = [np.loadtxt(f) for f in file_names]
matrices = np.concatenate(arrays)
numberOf4x4Matrices = int(matrices.size/16)
matrices4x4 = matrices.reshape(numberOf4x4Matrices, 4, 4)
p.pivot_calibration(matrices4x4)
number_of_matrices = int(matrices.size/16)
matrices = matrices.reshape(number_of_matrices, 4, 4)
p.pivot_calibration(matrices)
def test_four_columns_matrices4x4():
......@@ -41,9 +41,9 @@ def test_return_value():
file_names = glob('tests/data/PivotCalibration/*')
arrays = [np.loadtxt(f) for f in file_names]
matrices = np.concatenate(arrays)
numberOf4x4Matrices = int(matrices.size/16)
matrices4x4 = matrices.reshape(numberOf4x4Matrices, 4, 4)
x_values, residual_error =p.pivot_calibration(matrices4x4)
number_of_matrices = int(matrices.size/16)
matrices = matrices.reshape(number_of_matrices, 4, 4)
x_values, residual_error =p.pivot_calibration(matrices)
assert 1.838 == round(residual_error, 3)
assert -14.476 == round(x_values[0, 0], 3)
assert 395.143 == round(x_values[1, 0], 3)
......@@ -62,8 +62,26 @@ def test_rank_if_condition():
file_names = glob('tests/data/test_case_data.txt')
arrays = [np.loadtxt(f) for f in file_names]
matrices = np.concatenate(arrays)
numberOf4x4Matrices = int(matrices.size/16)
matrices4x4 = matrices.reshape(numberOf4x4Matrices, 4, 4)
p.pivot_calibration(matrices4x4)
number_of_matrices = int(matrices.size/16)
matrices = matrices.reshape(number_of_matrices, 4, 4)
p.pivot_calibration(matrices)
def test_pivot_with_ransac():
file_names = glob('tests/data/PivotCalibration/*')
arrays = [np.loadtxt(f) for f in file_names]
matrices = np.concatenate(arrays)
number_of_matrices = int(matrices.size/16)
matrices = matrices.reshape(number_of_matrices, 4, 4)
model_1, residual_1 = p.pivot_calibration(matrices)
print("Without RANSAC:" + str(model_1) + ", RMS=" + str(residual_1))
model_2, residual_2 = p.pivot_calibration_with_ransac(matrices, 10, 4, 0.25)
print("With RANSAC:" + str(model_2) + ", RMS=" + str(residual_2))
assert residual_2 < residual_1
model_3, residual_3 = p.pivot_calibration_with_ransac(matrices, 10, 4, 0.25, early_exit=True)
print("With Early Exit RANSAC:" + str(model_3) + ", RMS=" + str(residual_3))
......@@ -65,7 +65,10 @@ def test_identity_result():
def test_reflection_data():
"""This seems to be testing that a rotation
is prefered to a reflection. To transform these
points you can either reflect through yz, or
rotate 180 about y"""
fixed = np.zeros((4, 3))
fixed[0][1] = 1
fixed[2][0] = 2
......@@ -80,11 +83,60 @@ def test_reflection_data():
expected_rotation = np.eye(3)
expected_rotation[0][0] = -1
expected_rotation[2][2] = -1
print (np.linalg.det(rotation))
print (rotation)
assert np.allclose(rotation, expected_rotation, 0.0000001)
assert np.allclose(translation, np.zeros((3, 1)), 0.0000001)
assert error < 0.0000001
def test_reflection_BARD_data():
"""This is a reflection test using data taken from using
BARD at the 2019 summer school. This data set will give
you reflections, unless you replace equation 13 from Arun
with Fitzpatrick's modification. This was done at issue #19"""
ct_fids = np.zeros((4,3))
world_fids = np.zeros((4,3))
ct_fids[0][0] = 37.82
ct_fids[0][1] = 52.99
ct_fids[0][2] = -191.00
ct_fids[1][0] = 310.90
ct_fids[1][1] = 70.36
ct_fids[1][2] = -196.50
ct_fids[2][0] = 171.98
ct_fids[2][1] = 57.59
ct_fids[2][2] = -376.96
ct_fids[3][0] = 176.56
ct_fids[3][1] = 71.30
ct_fids[3][2] = -7.42
world_fids[0][0] = -7.060625292517708829e+01
world_fids[0][1] = -4.340731346419998715e+01
world_fids[0][2] = -7.678145283688323275e+01
world_fids[1][0] = -8.290290497289761618e+01
world_fids[1][1] = 2.365685716505672360e+02
world_fids[1][2] = -6.248063455571644909e+01
world_fids[2][0] = 1.061984137321565527e+02
world_fids[2][1] = 1.054365783777915624e+02
world_fids[2][2] = -6.347752061138085367e+01
world_fids[3][0] = -2.658595567434794020e+02
world_fids[3][1] = 9.556619277927123335e+01
world_fids[3][2] = -7.335614476329364209e+01
rotation, translation, error = p.orthogonal_procrustes(world_fids, ct_fids)
expected_rotation = np.eye(3)
expected_rotation[0][0] = -0.04933494
expected_rotation[0][1] = -0.01438556
expected_rotation[0][2] = -0.99867869
expected_rotation[1][0] = 0.99212873
expected_rotation[1][1] = 0.11451658
expected_rotation[1][2] = -0.05066094
expected_rotation[2][0] = 0.11509405
expected_rotation[2][1] = -0.99331717
expected_rotation[2][2] = 0.00862266
assert np.allclose(rotation, expected_rotation, 0.0000001)
def test_translation_result():
......
# -*- coding: utf-8 -*-
""" Tests concerning computing TRE from FLE. """
import pytest
import numpy as np
import sksurgerycore.algorithms.errors as err
def measure_tre_1(mean_fle_squared, target):
fiducials = np.zeros((4, 3))
fiducials[0][0] = 1
fiducials[0][1] = 1
fiducials[1][0] = -0.5
fiducials[1][1] = 0.5
fiducials[2][0] = -1
fiducials[2][1] = -1
fiducials[3][0] = 0.5
fiducials[3][1] = -0.5
error = err.compute_tre_from_fle(fiducials, mean_fle_squared, target)
return error
def test_tre_origin_zero_fle():
target = np.zeros((1, 3))
error = measure_tre_1(0, target)
assert np.isclose(error, 0)
def test_tre_origin_1mm_fle():
target = np.zeros((1, 3))
error = measure_tre_1(1, target)
assert np.isclose(error, 0.25)
def test_tre_origin_1mm_fle_non_zero_target():
target = np.zeros((1, 3))
target[0][0] = 2
error = measure_tre_1(1, target)
assert np.isclose(error, 1.048142396999972)
def test_invalid_because_fiducials_wrong_type():
with pytest.raises(TypeError):
err.compute_tre_from_fle("not an arrray", 1, np.ones((1, 3)))
def test_invalid_because_fiducials_wrong_columns():
with pytest.raises(ValueError):
err.compute_tre_from_fle(np.ones((1,4)), 1, np.ones((1, 3)))
def test_invalid_because_fiducials_wrong_rows():
with pytest.raises(ValueError):
err.compute_tre_from_fle(np.ones((2, 3)), 1, np.ones((1, 3)))
def test_invalid_because_target_wrong_type():
with pytest.raises(TypeError):
err.compute_tre_from_fle(np.ones((3, 3)), 1, "not an array")
def test_invalid_because_target_wrong_columns():
with pytest.raises(ValueError):
err.compute_tre_from_fle(np.ones((3, 3)), 1, np.ones((1, 4)))
def test_invalid_because_fiducials_wrong_rows():
with pytest.raises(ValueError):
err.compute_tre_from_fle(np.ones((3, 3)), 1, np.ones((2, 3)))
# -*- coding: utf-8 -*-
import os
import pytest
import sksurgerycore.configuration.configuration_manager as cm
......@@ -20,6 +21,8 @@ def test_constructor_with_valid_file():
manager = cm.ConfigurationManager("tests/data/FordPrefect.json")
assert manager is not None
assert manager.get_file_name() is not None
assert manager.get_dir_name() == os.path.dirname(manager.get_file_name())
def test_setter_getter_loop():
......
<?xml version="1.0" encoding="UTF-8" ?>
<point_set_file>
<file_version>0.1</file_version>
<point_set>
<time_series>
<time_series_id>0</time_series_id>
<Geometry3D ImageGeometry="false" FrameOfReferenceID="0">
<IndexToWorld type="Matrix3x3" m_0_0="1" m_0_1="0" m_0_2="0" m_1_0="0" m_1_1="1" m_1_2="0" m_2_0="0" m_2_1="0" m_2_2="1" />
<Offset type="Vector3D" x="0" y="0" z="0" />
<Bounds>
<Min type="Vector3D" x="-134.87921472311" y="-87.0495708179474" z="1570.031022310257" />
<Max type="Vector3D" x="111.4532385277748" y="-19.0247973632812" z="1773.7020261287689" />
</Bounds>
</Geometry3D>
<point>
<id>0</id>
<specification>0</specification>
<x>10.2352026367</x>
<y>-87.0495708179</y>
<z>1773.70202613</z>
</point>
<point>
<id>1</id>
<specification>0</specification>
<x>-134.879214723</x>
<y>-19.0247973633</y>
<z>1574.53710377</z>
</point>
<point>
<id>2</id>
<specification>0</specification>
<x>111.453238528</x>
<y>-22.3497973633</y>
<z>1570.03102231</z>
</point>
</time_series>
</point_set>
</point_set_file>
# coding=utf-8
"""Tests for load_mps"""
import six
import pytest
import sksurgerycore.io.load_mps as lmps
def test_load_mps():
ids, points = lmps.load_mps('tests/data/pointset.mps')
assert points.shape[0] == 3
assert points.shape[1] == 3
assert ids.shape[0] == 3
# -*- coding: utf-8 -*-
""" File processing utils. """
import os
import pytest
import sksurgerycore.utilities.file_utilities as fu
def test_get_absolute_file():
cwd = os.getcwd()
this_file = "tests/utilities/test_file_utilities.py"
fp = fu.get_absolute_path_of_file(this_file)
assert fp == os.path.join(cwd, fp)
fp2 = fu.get_absolute_path_of_file(fp)
assert fp2 == fp
fp3 = fu.get_absolute_path_of_file(this_file, cwd)
assert fp3 == fp
# content of: tox.ini , put in same dir as setup.py
[tox]
envlist = py27,py36,lint
envlist = py36,lint
skipsdist = True
[travis]
python =
2.7: py27
3.6: py36, docs, lint
[testenv]
......@@ -31,13 +30,6 @@ commands = sphinx-build -M html . build
basepython=python3.6
commands=pyinstaller --onefile sksurgerycore.py --noconfirm --windowed
[testenv:pip2]
basepython=python2.7
changedir=pip_test
skip_install=True
commands = pip install {posargs}
sksurgerycore --help
[testenv:pip3]
basepython=python3.6
changedir=pip_test
......