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
6 changes: 4 additions & 2 deletions docs/source/code.rst
Original file line number Diff line number Diff line change
Expand Up @@ -5,9 +5,11 @@ Code Documentation
.. toctree::
:maxdepth: 2

file_handler
ffd_parameters
affine_trans
ffd
ffd_parameters
file_handler


.. automodule:: pygem
:members:
8 changes: 8 additions & 0 deletions docs/source/ffd.rst
Original file line number Diff line number Diff line change
@@ -0,0 +1,8 @@
FFD
=================

.. automodule:: pygem.ffd

.. autoclass:: FFD
:members:
:private-members:
1 change: 1 addition & 0 deletions pygem/affine_trans.py
Original file line number Diff line number Diff line change
Expand Up @@ -60,6 +60,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
68 changes: 35 additions & 33 deletions pygem/ffd.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,6 @@
"""
import numpy as np
from scipy import special
import pygem.ffd_parameters as ffdp
import pygem.affine_trans as at


Expand All @@ -23,60 +22,63 @@ def perform(self):
"""
(n_rows_mesh, n_cols_mesh) = self.original_mesh_points.shape

## translation and then affine transformation
translation = np.array([self.parameters.origin_box_x, self.parameters.origin_box_y, self.parameters.origin_box_z])
# translation and then affine transformation
translation = self.parameters.origin_box

fisical_frame = np.array([self.parameters.position_vertex_1 - translation, \
self.parameters.position_vertex_2 - translation, \
self.parameters.position_vertex_3 - translation, \
[0,0,0]])
reference_frame = np.array([[1,0,0], [0,1,0], [0,0,1], [0,0,0]])
[0, 0, 0]])
reference_frame = np.array([[1, 0, 0], [0, 1, 0], [0, 0, 1], [0, 0, 0]])

transformation = at.affine_points_fit(fisical_frame, reference_frame)
inverse_transformation = at.affine_points_fit(reference_frame, fisical_frame)
##

# initialization
(dim_n_mu, dim_m_mu, dim_t_mu) = self.parameters.array_mu_x.shape
bernstein_x = np.zeros((n_rows_mesh, dim_n_mu))
bernstein_y = np.zeros((n_rows_mesh, dim_m_mu))
bernstein_z = np.zeros((n_rows_mesh, dim_t_mu))
shift_mesh_points = np.zeros((n_rows_mesh, n_cols_mesh))

reference_frame_mesh_points = self._transform_points(self.original_mesh_points - translation, transformation)
reference_frame_mesh_points = self._transform_points(self.original_mesh_points - translation, \
transformation)

for i in range (0, dim_n_mu):
aux1 = np.power((1-reference_frame_mesh_points[:,0]),dim_n_mu-1-i)
aux2 = np.power(reference_frame_mesh_points[:,0],i)
bernstein_x[:,i] = special.binom(dim_n_mu-1,i) * np.multiply(aux1, aux2)
for i in range(0, dim_n_mu):
aux1 = np.power((1-reference_frame_mesh_points[:, 0]), dim_n_mu-1-i)
aux2 = np.power(reference_frame_mesh_points[:, 0], i)
bernstein_x[:, i] = special.binom(dim_n_mu-1, i) * np.multiply(aux1, aux2)

for i in range (0, dim_m_mu):
aux1 = np.power((1-reference_frame_mesh_points[:,1]),dim_m_mu-1-i)
aux2 = np.power(reference_frame_mesh_points[:,1],i)
bernstein_y[:,i] = special.binom(dim_m_mu-1,i) * np.multiply(aux1, aux2)
for i in range(0, dim_m_mu):
aux1 = np.power((1-reference_frame_mesh_points[:, 1]), dim_m_mu-1-i)
aux2 = np.power(reference_frame_mesh_points[:, 1], i)
bernstein_y[:, i] = special.binom(dim_m_mu-1, i) * np.multiply(aux1, aux2)

for i in range (0, dim_t_mu):
aux1 = np.power((1-reference_frame_mesh_points[:,2]),dim_t_mu-1-i)
aux2 = np.power(reference_frame_mesh_points[:,2],i)
bernstein_z[:,i] = special.binom(dim_t_mu-1,i) * np.multiply(aux1, aux2)
for i in range(0, dim_t_mu):
aux1 = np.power((1-reference_frame_mesh_points[:, 2]), dim_t_mu-1-i)
aux2 = np.power(reference_frame_mesh_points[:, 2], i)
bernstein_z[:, i] = special.binom(dim_t_mu-1, i) * np.multiply(aux1, aux2)

for i in range (0, dim_n_mu):
for j in range (0, dim_m_mu):
for k in range (0, dim_t_mu):
aux = np.multiply(bernstein_x[:,i], np.multiply(bernstein_y[:,j], bernstein_z[:,k]))
for i in range(0, dim_n_mu):
for j in range(0, dim_m_mu):
for k in range(0, dim_t_mu):
aux = np.multiply(bernstein_x[:, i], np.multiply(bernstein_y[:, j], bernstein_z[:, k]))
aux_x = aux * self.parameters.array_mu_x[i, j, k]
aux_y = aux * self.parameters.array_mu_y[i, j, k]
aux_z = aux * self.parameters.array_mu_z[i, j, k]
shift_mesh_points[:,0] += aux_x
shift_mesh_points[:,1] += aux_y
shift_mesh_points[:,2] += aux_z
shift_mesh_points[:, 0] += aux_x
shift_mesh_points[:, 1] += aux_y
shift_mesh_points[:, 2] += aux_z

# Splitting points inside and outside the lattice: TODO not very efficient
for i in range (0, n_rows_mesh):
if (reference_frame_mesh_points[i,0] < 0) or (reference_frame_mesh_points[i,1] < 0) or (reference_frame_mesh_points[i,2] < 0) or \
(reference_frame_mesh_points[i,0] > 1) or (reference_frame_mesh_points[i,1] > 1) or (reference_frame_mesh_points[i,2] > 1):
shift_mesh_points[i,:] = 0
for i in range(0, n_rows_mesh):
if (reference_frame_mesh_points[i, 0] < 0) or (reference_frame_mesh_points[i, 1] < 0) \
or (reference_frame_mesh_points[i, 2] < 0) or (reference_frame_mesh_points[i, 0] > 1) \
or (reference_frame_mesh_points[i, 1] > 1) or (reference_frame_mesh_points[i, 2] > 1):
shift_mesh_points[i, :] = 0

return self._transform_points(shift_mesh_points + reference_frame_mesh_points, inverse_transformation) + translation
return self._transform_points(shift_mesh_points + reference_frame_mesh_points, \
inverse_transformation) + translation


def _transform_points(self, original_points, transformation):
Expand All @@ -88,6 +90,6 @@ def _transform_points(self, original_points, transformation):

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

return modified_points
Loading