Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
79 changes: 41 additions & 38 deletions pygem/affine_trans.py
Original file line number Diff line number Diff line change
Expand Up @@ -44,6 +44,45 @@ def angles2matrix(rot_z=0, rot_y=0, rot_x=0):
return np.eye(3)


def to_reduced_row_echelon_form(matrix):
"""
This method computes the reduced row echelon form (a.k.a. row canonical form) of a matrix.
The code is taken from https://rosettacode.org/wiki/Reduced_row_echelon_form#Python and
edited with minor changes.

:param matrix matrix: matrix to be reduced.

:return matrix: the reduced matrix.
:rtype: matrix

.. note::
`matrix` will change after calling this function.
"""
lead = 0
row_count = len(matrix)
column_count = len(matrix[0])
for r in range(row_count):
if lead >= column_count:
return matrix
i = r
while matrix[i][lead] == 0:
i += 1
if i == row_count:
i = r
lead += 1
if column_count == lead:
return matrix
matrix[i], matrix[r] = matrix[r], matrix[i]
lv = matrix[r][lead]
matrix[r] = [mrx / float(lv) for mrx in matrix[r]]
for i in range(row_count):
if i != r:
lv = matrix[i][lead]
matrix[i] = [iv - lv*rv for rv, iv in zip(matrix[r], matrix[i])]
lead += 1
return matrix


def affine_points_fit(points_start, points_end):
"""
Fit an affine transformation from starting points to ending points through a
Expand All @@ -60,7 +99,7 @@ def affine_points_fit(points_start, points_end):
:Example:

>>> import pygem.affine_trans as at

>>> # Example of a rotation
>>> p_start = np.array([[1,0,0], [0,1,0], [0,0,1], [0,0,0]])
>>> p_end = np.array([[0,1,0], [-1,0,0], [0,0,1], [0,0,0]])
Expand Down Expand Up @@ -92,43 +131,7 @@ def affine_points_fit(points_start, points_end):
Q[i][j] += qt[i] * qt[j]


def to_reduced_row_echelon_form(matrix):
"""
This method computes the reduced row echelon form (a.k.a. row canonical form) of a matrix.
The code is taken from https://rosettacode.org/wiki/Reduced_row_echelon_form#Python and
edited with minor changes.

:param matrix matrix: matrix to be reduced.

:return matrix: the reduced matrix.
:rtype: float
"""
lead = 0
row_count = len(matrix)
column_count = len(matrix[0])
for r in range(row_count):
if lead >= column_count:
return matrix
i = r
while matrix[i][lead] == 0:
i += 1
if i == row_count:
i = r
lead += 1
if column_count == lead:
return matrix
matrix[i], matrix[r] = matrix[r], matrix[i]
lv = matrix[r][lead]
matrix[r] = [mrx / float(lv) for mrx in matrix[r]]
for i in range(row_count):
if i != r:
lv = matrix[i][lead]
matrix[i] = [iv - lv*rv for rv, iv in zip(matrix[r], matrix[i])]
lead += 1
return matrix


# Augement Q with c and solve Q * a' = c by Gauss-Jordan
# Augement Q with c and get the reduced row echelon form of the result
M = [Q[i] + c[i] for i in range(dim+1)]

if np.linalg.cond(M) < 1/sys.float_info.epsilon:
Expand Down
51 changes: 50 additions & 1 deletion tests/test_affine_trans.py
Original file line number Diff line number Diff line change
Expand Up @@ -53,7 +53,56 @@ def test_angles2matrix_rot_xyz(self):

mat_test = at.angles2matrix(rotz, roty, rotx)
np.testing.assert_array_almost_equal(mat_exact, mat_test)



def test_to_reduced_row_echelon_form_1(self):
matrix = [[1., 1., 1.], [1., 1., 1.], [1., 1., 1.]]
rref_matrix = at.to_reduced_row_echelon_form(matrix)
rref_matrix_exact = [[1.0, 1.0, 1.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]]
assert rref_matrix == rref_matrix_exact


def test_to_reduced_row_echelon_form_2(self):
matrix = [[1., -1., 1.], [-1., 1., -1.], [3., 4., 5.]]
rref_matrix = at.to_reduced_row_echelon_form(matrix)
rref_matrix_exact = [[1.0, 0.0, 1.2857142857142856], [0.0, 1.0, 0.2857142857142857], [0.0, 0.0, 0.0]]
assert rref_matrix == rref_matrix_exact


def test_to_reduced_row_echelon_form_3(self):
matrix = [[0., 0., 0.], [0., 0., 0.], [0., 0., 0.]]
rref_matrix = at.to_reduced_row_echelon_form(matrix)
rref_matrix_exact = [[0., 0., 0.], [0., 0., 0.], [0., 0., 0.]]
assert rref_matrix == rref_matrix_exact


def test_to_reduced_row_echelon_form_4(self):
matrix = [[0., 0., 0.], [-3., 6., -3.], [0., 0., 0.]]
rref_matrix = at.to_reduced_row_echelon_form(matrix)
rref_matrix_exact = [[1.0, -2.0, 1.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]]
assert rref_matrix == rref_matrix_exact


def test_to_reduced_row_echelon_form_5(self):
matrix = [[0., 0., 0.], [0., 0., 0.], [2., 4., 6.]]
rref_matrix = at.to_reduced_row_echelon_form(matrix)
rref_matrix_exact = [[1.0, 2.0, 3.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]]
assert rref_matrix == rref_matrix_exact


def test_to_reduced_row_echelon_form_6(self):
matrix = [[0., 0., 0.], [0., 0., 0.], [2., 4., 6.], [1., 0., 0.]]
rref_matrix = at.to_reduced_row_echelon_form(matrix)
rref_matrix_exact = [[1.0, 0.0, 0.0], [-0.0, 1.0, 1.5], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]]
assert rref_matrix == rref_matrix_exact


def test_to_reduced_row_echelon_form_7(self):
matrix = [[0., 0., 0., 4], [0., 0., -3., 2], [2., 4., 6., 0]]
rref_matrix = at.to_reduced_row_echelon_form(matrix)
rref_matrix_exact = [[1.0, 2.0, 0.0, 0.0], [0.0, 0.0, 1.0, 0.0], [0.0, 0.0, 0.0, 1.0]]
assert rref_matrix == rref_matrix_exact


def test_affine_points_fit_identity_1(self):
p_start = np.array([[1,0,0], [0,1,0], [0,0,1], [0,0,0]])
Expand Down