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
41 changes: 19 additions & 22 deletions pygem/affine.py
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
"""
Utilities for the affine transformations of the bounding box of the Free Form Deformation.
Utilities for the affine transformations of the bounding box of the Free Form
Deformation.
"""
import math
import sys
Expand All @@ -9,22 +10,23 @@

def angles2matrix(rot_z=0, rot_y=0, rot_x=0):
"""
This method returns the rotation matrix for given rotations around z, y and x axes.
The output rotation matrix is equal to the composition of the individual rotations.
Rotations are counter-clockwise. The default value of the three rotations is zero.
This method returns the rotation matrix for given rotations around z, y and
x axes. The output rotation matrix is equal to the composition of the
individual rotations. Rotations are counter-clockwise. The default value of
the three rotations is zero.

:param float rot_z: rotation angle (in radians) around z-axis.
:param float rot_y: rotation angle (in radians) around y-axis.
:param float rot_x: rotation angle (in radians) around x-axis.

:return: rot_matrix: rotation matrix for the given angles. The matrix shape is always (3, 3).
:return: rot_matrix: rotation matrix for the given angles. The matrix shape
is always (3, 3).
:rtype: numpy.ndarray

:Example:

>>> import pygem.affine as at
>>> import numpy as np

>>> # Example of a rotation around x, y, z axis
>>> rotz = 10*np.pi/180
>>> roty = 20*np.pi/180
Expand All @@ -34,8 +36,8 @@ def angles2matrix(rot_z=0, rot_y=0, rot_x=0):
.. note::

- The direction of rotation is given by the right-hand rule.
- When applying the rotation to a vector, the vector should be column vector
to the right of the rotation matrix.
- When applying the rotation to a vector, the vector should be column
vector to the right of the rotation matrix.
"""
rot_matrix = []
if rot_z:
Expand Down Expand Up @@ -63,9 +65,10 @@ def angles2matrix(rot_z=0, rot_y=0, rot_x=0):

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.
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.

Expand All @@ -75,7 +78,6 @@ def to_reduced_row_echelon_form(matrix):
:Example:

>>> import pygem.affine as at

>>> matrix = [[1., 1., 1.], [1., 1., 1.], [1., 1., 1.]]
>>> rref_matrix = at.to_reduced_row_echelon_form(matrix)

Expand Down Expand Up @@ -118,9 +120,9 @@ def affine_points_fit(points_start, points_end):
:param numpy.ndarray points_start: set of starting points.
:param numpy.ndarray points_end: set of ending points.

:return: transform_vector: function that transforms a vector according to the
affine map. It takes a source vector and return a vector transformed
by the reduced row echelon form of the map.
:return: transform_vector: function that transforms a vector according to
the affine map. It takes a source vector and return a vector transformed
by the reduced row echelon form of the map.
:rtype: function

:Example:
Expand Down Expand Up @@ -164,7 +166,7 @@ def affine_points_fit(points_start, points_end):

if np.linalg.cond(affine_matrix) < 1 / sys.float_info.epsilon:
rref_aff_matrix = to_reduced_row_echelon_form(affine_matrix)
rref_aff_matrix = np.array(rref_aff_matrix)
rref_aff_matrix = np.array(rref_aff_matrix)[:, 4:]
else:
raise RuntimeError(
"Error: singular matrix. Points are probably coplanar."
Expand All @@ -179,12 +181,7 @@ def transform_vector(source):
:return destination: numpy.ndarray representing the transformed vector.
:rtype: numpy.ndarray
"""
destination = np.zeros(dim)
for i in range(dim):
for j in range(dim):
destination[j] += source[i] * rref_aff_matrix[i][j + dim + 1]
# Add the last line of the rref
destination[i] += rref_aff_matrix[dim][i + dim + 1]
destination = source.dot(rref_aff_matrix[:-1]) + rref_aff_matrix[-1]
return destination

return transform_vector
88 changes: 42 additions & 46 deletions pygem/freeform.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,65 +3,67 @@

:Theoretical Insight:

Free Form Deformation is a technique for the efficient, smooth and accurate geometrical
parametrization. It has been proposed the first time in *Sederberg, Thomas W., and Scott
R. Parry. "Free-form deformation of solid geometric models." ACM SIGGRAPH computer
graphics 20.4 (1986): 151-160*. It consists in three different step:
Free Form Deformation is a technique for the efficient, smooth and accurate
geometrical parametrization. It has been proposed the first time in
*Sederberg, Thomas W., and Scott R. Parry. "Free-form deformation of solid
geometric models." ACM SIGGRAPH computer graphics 20.4 (1986): 151-160*. It
consists in three different step:

- Mapping the physical domain to the reference one with map :math:`\\boldsymbol{\psi}`.
In the code it is named *transformation*.
- Mapping the physical domain to the reference one with map
:math:`\\boldsymbol{\psi}`. In the code it is named *transformation*.

- Moving some control points to deform the lattice with :math:`\\hat{T}`.
The movement of the control points is basically the weight (or displacement)
:math:`\\boldsymbol{\mu}` we set in the *parameters file*.
The movement of the control points is basically the weight (or displacement)
:math:`\\boldsymbol{\mu}` we set in the *parameters file*.

- Mapping back to the physical domain with map :math:`\\boldsymbol{\psi}^-1`.
In the code it is named *inverse_transformation*.
- Mapping back to the physical domain with map
:math:`\\boldsymbol{\psi}^-1`. In the code it is named
*inverse_transformation*.

FFD map (:math:`T`) is the composition of the three maps, that is

.. math::
T(\\cdot, \\boldsymbol{\\mu}) = (\\Psi^{-1} \\circ \\hat{T} \\circ \\Psi)
(\\cdot, \\boldsymbol{\\mu})
.. math:: T(\\cdot, \\boldsymbol{\\mu}) = (\\Psi^{-1} \\circ \\hat{T} \\circ
\\Psi) (\\cdot, \\boldsymbol{\\mu})

In this way, every point inside the FFD box changes it position according to

.. math::
\\boldsymbol{P} = \\boldsymbol{\psi}^-1 \\left( \\sum_{l=0} ^L \\sum_{m=0} ^M
\\sum_{n=0} ^N \\mathsf{b}_{lmn}(\\boldsymbol{\\psi}(\\boldsymbol{P}_0))
\\boldsymbol{\\mu}_{lmn} \\right)
.. math:: \\boldsymbol{P} = \\boldsymbol{\psi}^-1 \\left( \\sum_{l=0} ^L
\\sum_{m=0} ^M \\sum_{n=0} ^N
\\mathsf{b}_{lmn}(\\boldsymbol{\\psi}(\\boldsymbol{P}_0))
\\boldsymbol{\\mu}_{lmn} \\right)

where :math:`\\mathsf{b}_{lmn}` are Bernstein polynomials.
We improve the traditional version by allowing a rotation of the FFD lattice in order
to give more flexibility to the tool.
where :math:`\\mathsf{b}_{lmn}` are Bernstein polynomials. We improve the
traditional version by allowing a rotation of the FFD lattice in order to
give more flexibility to the tool.

You can try to add more shapes to the lattice to allow more and more involved transformations.
You can try to add more shapes to the lattice to allow more and more
involved transformations.

"""
import numpy as np
from scipy import special
import pygem.affine as at


class FFD(object):
"""
Class that handles the Free Form Deformation on the mesh points.

:param FFDParameters ffd_parameters: parameters of the Free Form Deformation.
:param numpy.ndarray original_mesh_points: coordinates of the original points of the mesh.
:param FFDParameters ffd_parameters: parameters of the Free Form
Deformation.
:param numpy.ndarray original_mesh_points: coordinates of the original
points of the mesh.

:cvar FFDParameters parameters: parameters of the Free Form Deformation.
:cvar numpy.ndarray original_mesh_points: coordinates of the original points of the mesh.
The shape is `n_points`-by-3.
:cvar numpy.ndarray modified_mesh_points: coordinates of the points of the deformed mesh.
The shape is `n_points`-by-3.
:cvar numpy.ndarray original_mesh_points: coordinates of the original points
of the mesh. The shape is `n_points`-by-3.
:cvar numpy.ndarray modified_mesh_points: coordinates of the points of the
deformed mesh. The shape is `n_points`-by-3.

:Example:

>>> import pygem.freeform as ffd
>>> import pygem.params as ffdp
>>> import numpy as np

>>> ffd_parameters = ffdp.FFDParameters()
>>> ffd_parameters.read_parameters('tests/test_datasets/parameters_test_ffd_sphere.prm')
>>> original_mesh_points = np.load('tests/test_datasets/meshpoints_sphere_orig.npy')
Expand All @@ -77,11 +79,8 @@ def __init__(self, ffd_parameters, original_mesh_points):

def perform(self):
"""
This method performs the deformation on the mesh points. After the execution
it sets `self.modified_mesh_points`.

.. todo::

This method performs the deformation on the mesh points. After the
execution it sets `self.modified_mesh_points`.
"""
# translation and then affine transformation
translation = self.parameters.origin_box
Expand Down Expand Up @@ -113,7 +112,8 @@ def perform(self):
) & (reference_frame_mesh_points[:, 2] <= 1.)]
(n_rows_mesh, n_cols_mesh) = mesh_points.shape

# Initialization. In order to exploit the contiguity in memory the following are transposed
# Initialization. In order to exploit the contiguity in memory the
# following are transposed
(dim_n_mu, dim_m_mu, dim_t_mu) = self.parameters.array_mu_x.shape
bernstein_x = np.zeros((dim_n_mu, n_rows_mesh))
bernstein_y = np.zeros((dim_m_mu, n_rows_mesh))
Expand Down Expand Up @@ -172,19 +172,15 @@ def perform(self):
@staticmethod
def _transform_points(original_points, transformation):
"""
This private static method transforms the points according to the affine transformation taken from
affine_points_fit method.
This private static method transforms the points according to the affine
transformation taken from affine_points_fit method.

:param numpy.ndarray original_points: coordinates of the original points.
:param function transformation: affine transformation taken from affine_points_fit method.
:param numpy.ndarray original_points: coordinates of the original
points.
:param function transformation: affine transformation taken from
affine_points_fit method.

:return: modified_points: coordinates of the modified points.
:rtype: numpy.ndarray
"""
n_rows_mesh, n_cols_mesh = original_points.shape
modified_points = np.zeros((n_rows_mesh, n_cols_mesh))

for i in range(0, n_rows_mesh):
modified_points[i, :] = transformation(original_points[i])

return modified_points
return transformation(original_points)