Frames, Points and Point Clouds#

PointCloud is a numba-jittable class, dedicated to large sets of points in Euclidean space. It is a subclass of the Vector class of the Neumann library, thus being an estension of NumPy’s ndarray class equipped with a transformation mechanism and other goodies.

[1]:
from sigmaepsilon.mesh.space import PointCloud
from sigmaepsilon.mesh.triang import triangulate

coords, *_ = triangulate(size=(800, 600), shape=(3, 3))
points = PointCloud(coords)

A remarkable feature is that sliced objects retain their indices in the original pointcloud:

[2]:
points[3:6].inds

[2]:
array([3, 4, 5])

You can get the index of the closest point (according to the standard Euclidean metric)

[3]:
points.index_of_closest(points.center())

[3]:
4

or several indices of points being closest to an array of targets

[4]:
points.index_of_closest(points[:3])

[4]:
array([0, 1, 2], dtype=int64)

Likewise, for the indices of the points being the furthest from one or more targets:

[5]:
points.index_of_furthest(points[:3])

[5]:
array([8, 6, 6], dtype=int64)

Functions like closest and furthest return either a Point or a PointCloud, depending on the number of targets:

[6]:
points.closest(points.center())

[6]:
Array([400., 300.,   0.])
[7]:
type(points.closest(points.center()))

[7]:
sigmaepsilon.mesh.space.point.Point
[8]:
points.closest(points.center()).id

[8]:
4
[9]:
points.closest(points[:3])

[9]:
PointCloud([[  0.   0.   0.]
 [400.   0.   0.]
 [800.   0.   0.]])
[10]:
points.furthest(points[:3])

[10]:
PointCloud([[800. 600.   0.]
 [  0. 600.   0.]
 [  0. 600.   0.]])

You can get the bounds calling the bounds method on an instance:

[11]:
points.bounds()

[11]:
array([[  0., 800.],
       [  0., 600.],
       [  0.,   0.]])

Transformations#

The PointCloud class is an instance of the Vector class defined in sigmaepsilon.math and therefore is equipped with all the mechanisms that the library provides. For instance, to apply a 90 degree rotation about the Z axis and then to move the points along the X axis:

[12]:
import numpy as np

(
    points
    .rotate("Space", [0, 0, np.pi / 2], "XYZ")
    .move(np.array([10.0, 0., 0.]))
    .bounds()
)

[12]:
array([[-590.,   10.],
       [   0.,  800.],
       [   0.,    0.]])

Numba JIT compilation#

The data and the indices are both accessible inside numba-jitted functions, even in nopython mode. See code block below shows the available attributes:

[13]:
from numba import jit
import numpy as np


@jit(nopython=True)
def numba_nopython(arr):
    return arr[0], arr[0, 0], arr.x, arr.data, arr.inds


c = np.array([[0, 0, 0], [0, 0, 1.0], [0, 0, 0]])
pc = PointCloud(c, inds=np.array([0, 1, 2]))
numba_nopython(pc[1:])

[13]:
(array([0., 0., 1.]),
 0.0,
 array([0., 0.]),
 Array([[0., 0., 1.],
        [0., 0., 0.]]),
 array([1, 2]))

Reference frames#

sigmaepsilon.mesh relies on the mechanism introduced in sigmaepsilon.math to account for different reference frames and it introduces a new CartesianFrame class to be used for geometrical applications. You can then use these frames when defining a pointcloud. If you define a pointcloud without explicity specifying a frame, it is assumed that the points are embedded in the ambient frame.

[14]:
points = PointCloud(coords)
points.bounds()

[14]:
array([[  0., 800.],
       [  0., 600.],
       [  0.,   0.]])
[15]:
from sigmaepsilon.mesh.space import CartesianFrame

ambient_frame = CartesianFrame(dim=3)
definition_frame = ambient_frame.rotate(
    "Space", [0, 0, np.pi / 2], "XYZ", inplace=False)

points = PointCloud(coords, frame=definition_frame)
points.bounds()

[15]:
array([[-6.0000000e+02,  4.8985872e-14],
       [ 0.0000000e+00,  8.0000000e+02],
       [ 0.0000000e+00,  0.0000000e+00]])

When calling the bounds method on a pointcloud, we can specify a target reference frame. Of yourse, the bounds of the pointcloud in the definition frame is the same as before.

[16]:
points.bounds(definition_frame)

[16]:
array([[  0., 800.],
       [  0., 600.],
       [  0.,   0.]])

You can access the frame and the coordinates of a pointcloud through properties frame and array:

[17]:
points.frame

[17]:
array([[ 6.123234e-17,  1.000000e+00,  0.000000e+00],
       [-1.000000e+00,  6.123234e-17,  0.000000e+00],
       [ 0.000000e+00,  0.000000e+00,  1.000000e+00]])
[18]:
points.array

[18]:
Array([[  0.,   0.,   0.],
       [400.,   0.,   0.],
       [800.,   0.,   0.],
       [  0., 300.,   0.],
       [400., 300.,   0.],
       [800., 300.,   0.],
       [  0., 600.,   0.],
       [400., 600.,   0.],
       [800., 600.,   0.]])

If you want to get the coorinates of a pointcloud in a target frame, use the show method:

[19]:
target_frame = ambient_frame.rotate(
    "Space", [0, 0, np.pi], "XYZ", inplace=False)
points.show(target_frame)

[19]:
Array([[ 0.0000000e+00,  0.0000000e+00,  0.0000000e+00],
       [ 2.4492936e-14, -4.0000000e+02,  0.0000000e+00],
       [ 4.8985872e-14, -8.0000000e+02,  0.0000000e+00],
       [ 3.0000000e+02,  1.8369702e-14,  0.0000000e+00],
       [ 3.0000000e+02, -4.0000000e+02,  0.0000000e+00],
       [ 3.0000000e+02, -8.0000000e+02,  0.0000000e+00],
       [ 6.0000000e+02,  3.6739404e-14,  0.0000000e+00],
       [ 6.0000000e+02, -4.0000000e+02,  0.0000000e+00],
       [ 6.0000000e+02, -8.0000000e+02,  0.0000000e+00]])

An important difference between the CartesianFrame class in sigmaepsilon.mesh and Neumann is that the former also supports the concept of an origo:

[20]:
ambient_frame = CartesianFrame(dim=3)

origo = np.array([1.0, 0.0, 0.0])
definition_frame = CartesianFrame(dim=3, origo=origo)

coords = np.array([
    [1.0, 0.0, 0.0],
    [0.0, 1.0, 0.0],
    [0.0, 0.0, 1.0]
])
points = PointCloud(coords, frame=definition_frame)
points.show(ambient_frame)

[20]:
Array([[2., 0., 0.],
       [1., 1., 0.],
       [1., 0., 1.]])