Source code for sknano.core.geometric_regions._3D_regions

# -*- coding: utf-8 -*-
"""
============================================================================
3D geometric regions (:mod:`sknano.core.geometric_regions._3D_regions`)
============================================================================

.. currentmodule:: sknano.core.geometric_regions._3D_regions

"""
from __future__ import absolute_import, division, print_function, \
    unicode_literals
__docformat__ = 'restructuredtext en'

from abc import ABCMeta, abstractmethod
import numbers
import numpy as np

from sknano.core.math import Point, Vector, vector as vec
from .base import GeometricRegion, GeometricTransformsMixin, ndim_errmsg

__all__ = ['Geometric3DRegion', 'Parallelepiped', 'Cuboid', 'Cube',
           'Ellipsoid', 'Sphere', 'Cylinder', 'Cone', 'geometric_3D_regions']


[docs]class Geometric3DRegion(GeometricRegion, GeometricTransformsMixin, metaclass=ABCMeta): """Abstract base class for representing 3D geometric regions.""" @property def ndim(self): """Return the dimensions.""" return 3 @property @abstractmethod def volume(self): """Volume of 3D geometric region.""" raise NotImplementedError @property def measure(self): """Alias for :attr:`~Geometric3DRegion.volume`, which is the measure \ of a 3D geometric region.""" return self.volume @property def bounding_box(self): """Bounding :class:`Cuboid`.""" return Cuboid(pmin=self.pmin.tolist(), pmax=self.pmax.tolist())
[docs]class Parallelepiped(Geometric3DRegion): """`Geometric3DRegion` for a parallelepiped. .. versionadded:: 0.3.0 Represents a parallelepiped with origin :math:`(o_x, o_y, o_z)` and direction vectors :math:`\\mathbf{u}=(u_x, u_y, u_z)`, :math:`\\mathbf{v}=(v_x, v_y, v_z)`, and :math:`\\mathbf{w}=(w_x, w_y, w_z)`. Parameters ---------- o : array_like, optional Parallelepiped origin. If `None`, it defaults to `o=[0, 0, 0]`. u, v, w : array_like, optional Parallelepiped direction vectors stemming from origin `o`. If `None`, then the default values are `u=[1, 0, 0]`, `v=[1, 1, 0]`, and `w=[0, 1, 1]`. Notes ----- :class:`Parallelepiped` represents the bounded region :math:`\\left\\{o+\\lambda_1\\mathbf{u}+\\lambda_2\\mathbf{v}+ \\lambda_3\\mathbf{w}\\in R^3|0\\le\\lambda_i\\le 1\\right\\}`, where :math:`\\mathbf{u}`, :math:`\\mathbf{v}`, and :math:`\\mathbf{w}` have to be linearly independent. Calling `Parallelepiped` with no parameters is equivalent to `Parallelepiped`\ `(o=[0, 0, 0], u=[1, 0, 0], v=[1, 1, 0], w=[0, 1, 1])`. """ def __init__(self, o=None, u=None, v=None, w=None): super().__init__() if o is None: o = [0, 0, 0] self._o = Point(o) if u is None: u = [1, 0, 0] self._u = Vector(u, p0=self.o) if v is None: v = [1, 1, 0] self._v = Vector(v, p0=self.o) if w is None: w = [0, 1, 1] self._w = Vector(w, p0=self.o) self.points.append(self.o) self.vectors.extend([self.u, self.v, self.w]) self.fmtstr = "o={o!r}, u={u!r}, v={v!r}, w={w!r}" @property def o(self): """3D point coordinates :math:`(o_x, o_y, o_z)` of origin. Returns ------- :class:`Point` 3D :class:`Point` coordinates :math:`(o_x, o_y, o_z)` of origin. """ return self._o @o.setter def o(self, value): if not isinstance(value, (tuple, list, np.ndarray)) or len(value) != 3: raise TypeError('Expected a 3-element array_like object') self.translate(Vector(p0=self._o, p=Point(value))) # self._o[:] = Point(value) @property def u(self): """3D direction vector :math:`\\mathbf{u}=(u_x, u_y, u_z)`, with \ origin :attr:`~Parallelepiped.o`""" return self._u @u.setter def u(self, value): if not isinstance(value, (tuple, list, np.ndarray)) or len(value) != 3: raise TypeError('Expected a 3-element array_like object') self._u[:] = Vector(value) self._u.p0 = self.o @property def v(self): """3D direction vector :math:`\\mathbf{v}=(v_x, v_y, v_z)`, with \ origin :attr:`~Parallelepiped.o`""" return self._v @v.setter def v(self, value): if not isinstance(value, (tuple, list, np.ndarray)) or len(value) != 3: raise TypeError('Expected a 3-element array_like object') self._v[:] = Vector(value) self._v.p0 = self.o @property def w(self): """3D direction vector :math:`\\mathbf{w}=(w_x, w_y, w_z)`, with \ origin :attr:`~Parallelepiped.o`""" return self._w @w.setter def w(self, value): if not isinstance(value, (tuple, list, np.ndarray)) or len(value) != 3: raise TypeError('Expected a 3-element array_like object') self._w[:] = Vector(value) self._w.p0 = self.o @property def lengths(self): """:class:`~python:tuple` of side lengths""" return self.a, self.b, self.c @property def abc(self): """Alias for :attr:`Parallelepiped.lengths`.""" return self.lengths @property def a(self): """Alias for :attr:`Parallelepiped.u.length`.""" return self.u.length @property def b(self): """Alias for :attr:`Parallelepiped.v.length`.""" return self.v.length @property def c(self): """Alias for :attr:`Parallelepiped.w.length`.""" return self.w.length @property def centroid(self): """Parallelepiped centroid, :math:`(c_x, c_y, c_z)`. Computed as the 3D point :math:`(c_x, c_y, c_z)` with coordinates: .. math:: c_x = o_x + \\frac{u_x + v_x + w_x}{2} c_y = o_y + \\frac{u_y + v_y + w_y}{2} c_z = o_z + \\frac{u_z + v_z + w_z}{2} where :math:`(o_x, o_y, o_z)`, :math:`(u_x, u_y, u_z)`, :math:`(v_x, v_y, v_z)`, and :math:`(w_x, w_y, w_z)` are the :math:`(x, y, z)` coordinates of the origin :math:`o` and :math:`(x, y, z)` components of the direction vectors :math:`\\mathbf{u}`, :math:`\\mathbf{v}`, and :math:`\\mathbf{w}`, respectively. Returns ------- :class:`Point` 3D :class:`Point` of centroid. """ ox, oy, oz = self.o ux, uy, uz = self.u vx, vy, vz = self.v wx, wy, wz = self.w cx = ox + (ux + vx + wx) / 2 cy = oy + (uy + vy + wy) / 2 cz = oz + (uz + vz + wz) / 2 return Point([cx, cy, cz]) @property def volume(self): """Parallelepiped volume, :math:`V`. Computed as: .. math:: V = |\\mathbf{u}\\cdot\\mathbf{v}\\times\\mathbf{w}| """ return np.abs(vec.scalar_triple_product(self.u, self.v, self.w))
[docs] def contains(self, point): """Test region membership of `point` in :class:`Parallelepiped`. Parameters ---------- point : array_like Returns ------- :class:`~python:bool` `True` if `point` is within :class:`Parallelepiped`, `False` otherwise Notes ----- A `point` :math:`(p_x, p_y, p_z)` is within the bounded region of a parallelepiped with origin :math:`(o_x, o_y, o_z)` and direction vectors :math:`\\mathbf{u}=(u_x, u_y, u_z)`, :math:`\\mathbf{v}=(v_x, v_y, v_z)`, and :math:`\\mathbf{w}=(w_x, w_y, w_z)` if the following is true: .. math:: 0\\le\\frac{ v_z (w_x (p_y - o_y) + w_y (o_x - p_x)) + w_z (v_x (o_y - p_y) + v_y (p_x - o_x)) + o_z (v_y w_x - v_x w_y) + p_z (v_x w_y - v_y w_x)}{ u_z (v_x w_y - v_y w_x) + u_y (v_z w_x - v_x w_z) + u_x (v_y w_z - v_z w_y)}\\le 1 \\land 0\\le\\frac{ u_z (w_x (p_y - o_y) + w_y (o_x - p_x)) + w_z (u_x (o_y - p_y) + u_y (p_x - o_x)) + o_z (u_y w_x - u_x w_y) + p_z (u_x w_y - u_y w_x)}{ u_z (v_y w_x - v_x w_y) + u_y (v_x w_z - v_z w_x) + u_x (v_z w_y - v_y w_z)}\\le 1 \\land 0\\le\\frac{ u_z (v_x (p_y - o_y) + v_y (o_x - p_x)) + v_z (u_x (o_y - p_y) + u_y (p_x - o_x)) + o_z (u_y v_x - u_x v_y) + p_z (u_x v_y - u_y v_x)}{ u_z (v_x w_y - v_y w_x) + u_y (v_z w_x - v_x w_z) + u_x (v_y w_z - v_z w_y)}\\le 1 """ px, py, pz = Point(point) ox, oy, oz = self.o ux, uy, uz = self.u vx, vy, vz = self.v wx, wy, wz = self.w q1 = (vz * (wx * (py - oy) + wy * (ox - px)) + wz * (vx * (oy - py) + vy * (px - ox)) + oz * (vy * wx - vx * wy) + pz * (vx * wy - vy * wx)) / \ (uz * (vx * wy - vy * wx) + uy * (vz * wx - vx * wz) + ux * (vy * wz - vz * wy)) q2 = (uz * (wx * (py - oy) + wy * (ox - px)) + wz * (ux * (oy - py) + uy * (px - ox)) + oz * (uy * wx - ux * wy) + pz * (ux * wy - uy * wx)) / \ (uz * (vy * wx - vx * wy) + uy * (vx * wz - vz * wx) + ux * (vz * wy - vy * wz)) q3 = (uz * (vx * (py - oy) + vy * (ox - px)) + vz * (ux * (oy - py) + uy * (px - ox)) + oz * (uy * vx - ux * vy) + pz * (ux * vy - uy * vx)) / \ (uz * (vx * wy - vy * wx) + uy * (vz * wx - vx * wz) + ux * (vy * wz - vz * wy)) return q1 >= 0 and q1 <= 1 and q2 >= 0 and q2 <= 1 and \ q3 >= 0 and q3 <= 1
[docs] def todict(self): """Returns a :class:`~python:dict` of the :class:`Parallelepiped` \ constructor parameters.""" return dict(o=self.o, u=self.u, v=self.v, w=self.w)
[docs]class Cuboid(Geometric3DRegion): """`Geometric3DRegion` for a cuboid. .. versionadded:: 0.3.0 Represents an axis-aligned cuboid with lower corner :math:`p_{\\mathrm{min}}= (x_{\\mathrm{min}},y_{\\mathrm{min}},z_{\\mathrm{min}})` and upper corner :math:`p_{\\mathrm{max}}= (x_{\\mathrm{max}},y_{\\mathrm{max}},z_{\\mathrm{max}})`. Parameters ---------- pmin, pmax : array_like, optional The minimum and maximum 3D point coordinates of the axis-aligned cuboid from `pmin=[xmin, ymin, zmin]` to `pmax=[xmax, ymax, zmax]`. xmin, ymin, zmin : float, optional The minimum :math:`(x, y, z)` coordinates of the axis-aligned cuboid. xmax, ymax, zmax : float, optional The maximum :math:`(x, y, z)` coordinates of the axis-aligned cuboid. Notes ----- :class:`Cuboid` represents the region :math:`\\left\\{\\{x, y, z\\}| x_{\\mathrm{min}}\\le x\\le x_{\\mathrm{max}}\\land y_{\\mathrm{min}}\\le y\\le y_{\\mathrm{max}}\\land z_{\\mathrm{min}}\\le z\\le z_{\\mathrm{max}}\\right\\}` Calling `Cuboid` with no parameters is equivalent to `Cuboid`\ `(pmin=[0, 0, 0], pmax=[1, 1, 1])`. """ def __init__(self, pmin=None, pmax=None, xmin=0, ymin=0, zmin=0, xmax=1, ymax=1, zmax=1): super().__init__() if pmin is None: pmin = [xmin, ymin, zmin] self._pmin = Point(pmin) if pmax is None: pmax = [xmax, ymax, zmax] self._pmax = Point(pmax) assert np.all(np.less_equal(self.pmin, self.pmax)) self.points.extend([self._pmin, self._pmax]) # assert self.pmin <= self.pmax self.fmtstr = "pmin={pmin!r}, pmax={pmax!r}" def _update_points(self): pmin = self.pmin.copy() pmax = self.pmax.copy() self._pmin = np.minimum(pmin, pmax) self._pmax = np.maximum(pmin, pmax) assert np.all(np.less_equal(self._pmin, self._pmax)) self.points.clear() self.points.extend([self.pmin, self.pmax]) @property def pmin(self): """:class:`Point` at minimum extent.""" return super().pmin @pmin.setter def pmin(self, value): if not isinstance(value, (tuple, list, np.ndarray)) or \ len(value) != self.ndim: raise TypeError(ndim_errmsg.format(self.ndim)) self._pmin[:] = Point(value) self._update_points() @property def pmax(self): """:class:`Point` at maximum extent.""" return super().pmax @pmax.setter def pmax(self, value): if not isinstance(value, (tuple, list, np.ndarray)) or \ len(value) != self.ndim: raise TypeError(ndim_errmsg.format(self.ndim)) self._pmax[:] = Point(value) self._update_points() @property def xmin(self): """:math:`x_{\\mathrm{min}}` coordinate.""" return self.pmin.x @xmin.setter def xmin(self, value): if not isinstance(value, numbers.Number): raise TypeError('Expected a number') self._pmin.x = float(value) self._update_points() @property def xmax(self): """:math:`x_{\\mathrm{max}}` coordinate.""" return self.pmax.x @xmax.setter def xmax(self, value): if not isinstance(value, numbers.Number): raise TypeError('Expected a number') self._pmax.x = float(value) self._update_points() @property def ymin(self): """:math:`y_{\\mathrm{min}}` coordinate.""" return self.pmin.y @ymin.setter def ymin(self, value): if not isinstance(value, numbers.Number): raise TypeError('Expected a number') self._pmin.y = float(value) self._update_points() @property def ymax(self): """:math:`y_{\\mathrm{max}}` coordinate.""" return self.pmax.y @ymax.setter def ymax(self, value): if not isinstance(value, numbers.Number): raise TypeError('Expected a number') self._pmax.y = float(value) self._update_points() @property def zmin(self): """:math:`z_{\\mathrm{min}}` coordinate.""" return self.pmin.z @zmin.setter def zmin(self, value): if not isinstance(value, numbers.Number): raise TypeError('Expected a number') self._pmin.z = float(value) self._update_points() @property def zmax(self): """:math:`z_{\\mathrm{max}}` coordinate.""" return self.pmax.z @zmax.setter def zmax(self, value): if not isinstance(value, numbers.Number): raise TypeError('Expected a number') self._pmax.z = float(value) self._update_points() @property def lengths(self): """:class:`~python:tuple` of side lengths""" return self.lx, self.ly, self.lz @property def lx(self): """Distance between :math:`x_{\\mathrm{max}}-x_{\\mathrm{min}}`.""" return self.xmax - self.xmin @property def ly(self): """Distance between :math:`y_{\\mathrm{max}}-y_{\\mathrm{min}}`.""" return self.ymax - self.ymin @property def lz(self): """Distance between :math:`z_{\\mathrm{max}}-z_{\\mathrm{min}}`.""" return self.zmax - self.zmin @property def abc(self): """Alias for :attr:`Cuboid.lengths`.""" return self.lengths @property def a(self): """Alias for :attr:`Cuboid.lx`.""" return self.lx @property def b(self): """Alias for :attr:`Cuboid.ly`.""" return self.ly @property def c(self): """Alias for :attr:`Cuboid.lz`.""" return self.lz @property def volume(self): """:class:`Cuboid` volume, :math:`V=\\ell_x\\ell_y\\ell_z`.""" return self.lx * self.ly * self.lz @property def centroid(self): """:class:`Cuboid` centroid, :math:`(c_x, c_y, c_z)`. Computed as the 3D :class:`Point` :math:`(c_x, c_y, c_z)` with coordinates: .. math:: c_x = \\frac{x_{\\mathrm{min}} + x_{\\mathrm{max}}}{2} c_y = \\frac{y_{\\mathrm{min}} + y_{\\mathrm{max}}}{2} c_z = \\frac{z_{\\mathrm{min}} + z_{\\mathrm{max}}}{2} Returns ------- :class:`Point` 3D :class:`Point` of centroid. """ cx = (self.xmin + self.xmax) / 2 cy = (self.ymin + self.ymax) / 2 cz = (self.zmin + self.zmax) / 2 return Point([cx, cy, cz])
[docs] def contains(self, point): """Test region membership of `point` in :class:`Cuboid`. Parameters ---------- point : array_like Returns ------- :class:`~python:bool` `True` if `point` is within :class:`Cuboid`, `False`, otherwise. Notes ----- A point :math:`(p_x, p_y, p_z)` is within the bounded region of a cuboid with lower corner at :math:`p_{\\mathrm{min}}= (x_{\\mathrm{min}}, y_{\\mathrm{min}}, z_{\\mathrm{min}})` and upper corner at :math:`p_{\\mathrm{max}}= (x_{\\mathrm{max}}, y_{\\mathrm{max}}, z_{\\mathrm{max}})` if the following is true: .. math:: x_{\mathrm{min}}\\le x\\le x_{\\mathrm{max}}\\land y_{\mathrm{min}}\\le y\\le y_{\\mathrm{max}}\\land z_{\mathrm{min}}\\le z\\le z_{\\mathrm{max}} """ px, py, pz = Point(point) return (px >= self.xmin) and (px <= self.xmax) and \ (py >= self.ymin) and (py <= self.ymax) and \ (pz >= self.zmin) and (pz <= self.zmax)
[docs] def todict(self): """Returns a :class:`~python:dict` of the :class:`Cuboid` \ constructor parameters.""" return dict(pmin=self.pmin, pmax=self.pmax)
[docs]class Cube(Geometric3DRegion): """`Geometric3DRegion` for a cube. .. versionadded:: 0.3.0 :class:`Cube` represents the region :math:`\\left\\{\\{c_i\\pm\\frac{a}{2}\\}|a>0\\forall i\\in\\{x,y,z\\}\\right\\}`. Parameters ---------- center : array_like, optional The :math:`(x,y,z)` coordinate of the axis-aligned cube center :math:`(c_x, c_y, c_z)`. a : float, optional Side length :math:`a` of axis-aligned cube. Notes ----- Calling :class:`Cube` with no parameters is equivalent to :class:`Cube`\ `(center=[0, 0, 0], a=1)`. """ def __init__(self, center=None, a=1): super().__init__() if center is None: center = [0, 0, 0] self._center = Point(center) self.a = a self.points.append(self.center) self.fmtstr = "center={center!r}, a={a:.3f}" @property def center(self): """Center point :math:`(c_x, c_y, c_z)` of axis-aligned cube.""" return self._center @center.setter def center(self, value): if not isinstance(value, (tuple, list, np.ndarray)) or len(value) != 3: raise TypeError('Expected a 3-element array_like object') self._center[:] = Point(value) @property def a(self): """Side length :math:`a` of the axis-aligned cube.""" return self._a @a.setter def a(self, value): self._a = float(value) @property def centroid(self): """Alias for :attr:`~Cube.center`.""" return self.center @property def volume(self): """:class:`Cube` volume, :math:`V=a^3`.""" return self.a ** 3
[docs] def contains(self, point): """Test region membership of `point` in :class:`Cube`. Parameters ---------- point : array_like Returns ------- :class:`~python:bool` `True` if `point` is within :class:`Cube`, `False` otherwise Notes ----- A `point` :math:`(p_x, p_y, p_z)` is within the bounded region of a cube with center :math:`(c_x, c_y, c_z)` and side length :math:`a` if the following is true: .. math:: c_i-\\frac{a}{2}\\le p_i\\le c_i+\\frac{a}{2}\\forall i\\in \\{x, y, z\\} """ px, py, pz = Point(point) cx, cy, cz = self.center a = self.a xmin = cx - a / 2 xmax = cx + a / 2 ymin = cy - a / 2 ymax = cy + a / 2 zmin = cz - a / 2 zmax = cz + a / 2 return (px >= xmin) and (px <= xmax) and \ (py >= ymin) and (py <= ymax) and \ (pz >= zmin) and (pz <= zmax)
[docs] def todict(self): """Returns a :class:`~python:dict` of the :class:`Cube` \ constructor parameters.""" return dict(center=self.center, a=self.a)
[docs]class Ellipsoid(Geometric3DRegion): """`Geometric3DRegion` for an ellipsoid. .. versionadded:: 0.3.0 Represents an axis-aligned ellipsoid centered at the point :math:`(c_x, c_y, c_z)` with semi-axes lengths :math:`r_x, r_y, r_z`. Parameters ---------- center : array_like or :class:`~sknano.core.Point` 3-tuple of floats or an instance of the :class:`~sknano.core.Point` class specifying the center point coordinates :math:`(x,y,z)` of the axis-aligned :class:`Ellipsoid` with semi-axes lengths `rx`, `ry`, `rz`. rx, ry, rz : float, optional Semi-axes lengths :math:`r_x, r_y, r_z` of axis-aligned :class:`Ellipsoid` centered at `center`. Notes ----- :class:`Ellipsoid` represents the axis-aligned ellipsoid: .. math:: \\left\\{\\{x, y, z\\}\\in R^3| \\left(\\frac{x-c_x}{r_x}\\right)^2 + \\left(\\frac{y-c_y}{r_y}\\right)^2 + \\left(\\frac{z-c_z}{r_z}\\right)^2\\le 1\\right\\} Calling :class:`Ellipsoid` with no parameters is equivalent to :class:`Ellipsoid`\ `(center=[0, 0, 0], rx=1, ry=1, rz=1)`. """ def __init__(self, center=None, rx=1, ry=1, rz=1): super().__init__() if center is None: center = [0, 0, 0] self._center = Point(center) self.rx = rx self.ry = ry self.rz = rz self.points.append(self.center) self.fmtstr = \ "center={center!r}, rx={rx:.3f}, ry={ry:.3f}, rz={rz:.3f}" @property def center(self): """:class:`Ellipsoid` center point :math:`(c_x, c_y, c_z)`.""" return self._center @center.setter def center(self, value): if not isinstance(value, (tuple, list, np.ndarray)) or len(value) != 3: raise TypeError('Expected a 3-element array_like object') self._center[:] = Point(value) @property def rx(self): """Length of semi-axis :math:`r_x`.""" return self._rx @rx.setter def rx(self, value): self._rx = float(value) @property def ry(self): """Length of semi-axis :math:`r_y`.""" return self._ry @ry.setter def ry(self, value): self._ry = float(value) @property def rz(self): """Length of semi-axis :math:`r_z`.""" return self._rz @rz.setter def rz(self, value): self._rz = float(value) @property def centroid(self): """Alias for :attr:`~Ellipsoid.center`.""" return self.center @property def volume(self): """:class:`Ellipsoid` volume, \ :math:`V=\\frac{4}{3}\\pi r_x r_y r_z`.""" return 4 / 3 * np.pi * self.rx * self.ry * self.rz
[docs] def contains(self, point): """Test region membership of `point` in :class:`Ellipsoid`. Parameters ---------- point : array_like Returns ------- :class:`~python:bool` `True` if `point` is within :class:`Ellipsoid`, `False` otherwise Notes ----- A `point` :math:`(p_x, p_y, p_z)` is within the bounded region of an ellipsoid with center :math:`(c_x, c_y, c_z)` and semi-axes lengths :math:`r_x, r_y, r_z` if the following is true: .. math:: \\left(\\frac{p_x - c_x}{r_x}\\right)^2 + \\left(\\frac{p_y - c_y}{r_y}\\right)^2 + \\left(\\frac{p_z - c_z}{r_z}\\right)^2\\le 1 """ px, py, pz = Point(point) cx, cy, cz = self.center rx, ry, rz = self.rx, self.ry, self.rz q1 = (px - cx) ** 2 / rx ** 2 q2 = (py - cy) ** 2 / ry ** 2 q3 = (pz - cz) ** 2 / rz ** 2 return q1 + q2 + q3 <= 1.0
[docs] def todict(self): """Returns a :class:`~python:dict` of the :class:`Ellipsoid` \ constructor parameters.""" return dict(center=self.center, rx=self.rx, ry=self.ry, rz=self.rz)
[docs]class Sphere(Geometric3DRegion): """`Geometric3DRegion` for a sphere. .. versionadded:: 0.3.0 Parameters ---------- center : array_like, optional Either a 3-tuple of floats or an instance of the :class:`~sknano.core.Point` class specifying the :math:`(x,y,z)` coordinates of the :class:`Sphere` center. r : float, optional Sphere radius :math:`r` """ def __init__(self, center=None, r=1): super().__init__() if center is None: center = [0, 0, 0] self._center = Point(center) self.r = r self.points.append(self.center) self.fmtstr = "center={center!r}, r={r:.3f}" @property def center(self): """:class:`Sphere` center point :math:`(h, k, l)`.""" return self._center @center.setter def center(self, value): if not isinstance(value, (tuple, list, np.ndarray)) or len(value) != 3: raise TypeError('Expected a 3-element array_like object') self._center[:] = Point(value) @property def centroid(self): """Alias for :attr:`~Sphere.center`.""" return self.center @property def r(self): """:class:`Sphere` radius, :math:`r`.""" return self._r @r.setter def r(self, value): self._r = float(value) @property def volume(self): """:class:`Sphere` volume, :math:`V=\\frac{4}{3}\\pi r^3`.""" return 4 / 3 * np.pi * self.r ** 3
[docs] def contains(self, point): """Test region membership of `point` in :class:`Sphere`. Parameters ---------- point : array_like Returns ------- :class:`~python:bool` `True` if `point` is within :class:`Sphere`, `False` otherwise. Notes ----- A `point` :math:`(p_x, p_y, p_z)` is within the bounded region of a sphere with center :math:`(h, k, l)` and radius :math:`r` if the following is true: .. math:: (p_x - h)^2 + (p_y - k)^2 + (p_z - l)^2 \\le r^2 """ x, y, z = Point(point) h, k, l = self.center r = self.r return (x - h) ** 2 + (y - k) ** 2 + (z - l) ** 2 <= r ** 2
[docs] def todict(self): """Returns a :class:`~python:dict` of the :class:`Sphere` \ constructor parameters.""" return dict(center=self.center, r=self.r)
[docs]class Cylinder(Geometric3DRegion): """`Geometric3DRegion` for a cylinder. .. versionadded:: 0.3.10 Represents a cylinder of radius :math:`r` around the line from :math:`(x_1, y_1, z_1)` to :math:`(x_2, y_2, z_2)`. Parameters ---------- p1, p2 : array_like, optional 3-tuples or :class:`~sknano.core.Point` class instances specifying the `Cylinder` axis from point :math:`p_1=(x_1,y_1,z_1)` to :math:`p_2=(x_2,y_2,z_2)`. r : float, optional `Cylinder` radius :math:`r` Notes ----- :class:`Cylinder` represents a cylinder region :math:`\\left\\{p_1+\\rho\\cos(\\theta)\\mathbf{v}_1 + \\rho\\sin(\\theta)\\mathbf{v}_2 + \\mathbf{v}_3 z| 0\\le\\theta\\le 2\\pi\\land 0\\le\\rho\\le 1\\land 0\\le z\\le 1\\right\\}` where :math:`\\mathbf{v}_3=p_2 - p_1` and the vectors :math:`\\{\\mathbf{v}_1,\\mathbf{v}_2,\\mathbf{v}_3\\}` are orthogonal with :math:`|\\mathbf{v}_1|=|\\mathbf{v}_2|=1`, and :math:`p_1=(x_1, y_1, z_1)` and :math:`p_2=(x_2,y_2,z_2)`. Calling :class:`Cylinder` with no parameters is equivalent to :class:`Cylinder`\ `(p1=[0, 0, -1], p2=[0, 0, 1], r=1)`. """ def __init__(self, p1=None, p2=None, r=1): super().__init__() if p1 is None: p1 = [0, 0, -1] self._p1 = Point(p1) if p2 is None: p2 = [0, 0, 1] self._p2 = Point(p2) self.r = r self.points.extend([self.p1, self.p2]) self.fmtstr = "p1={p1!r}, p2={p2!r}, r={r:.3f}" @property def r(self): """:class:`Cylinder` radius :math:`r`.""" return self._r @r.setter def r(self, value): self._r = float(value) @property def p1(self): """:class:`Cylinder` axis point :math:`p_1=(x_1, y_1, z_1)`.""" return self._p1 @p1.setter def p1(self, value): if not isinstance(value, (tuple, list, np.ndarray)) or len(value) != 3: raise TypeError('Expected a 3-element array_like object') self._p1[:] = Point(value) @property def p2(self): """:class:`Cylinder` axis point :math:`p_2=(x_2, y_2, z_2)`.""" return self._p2 @p2.setter def p2(self, value): if not isinstance(value, (tuple, list, np.ndarray)) or len(value) != 3: raise TypeError('Expected a 3-element array_like object') self._p2[:] = Point(value) @property def axis(self): """:class:`Cylinder` axis :class:`Vector` \ :math:`\\boldsymbol{\\ell}=p_2 - p_1`. Returns ------- :class:`Vector` 3D :class:`Vector` along :class:`Cylinder` axis from :class:`Point` :math:`p_1=(x_1, y_1, z_1)` to :class:`Point` :math:`p_2=(x_2, y_2, z_2)`. """ return Vector(p0=self.p1, p=self.p2) @property def centroid(self): """:class:`Cylinder` centroid, :math:`(c_x, c_y, c_z)`. Computed as: .. math:: c_x = \\frac{x_1 + x_2}{2} c_y = \\frac{y_1 + y_2}{2} c_z = \\frac{z_1 + z_2}{2} """ return (self.p1 + self.p2) / 2 @property def volume(self): """:class:`Cylinder` volume, :math:`V=\\pi r^2 \\ell`.""" return np.pi * self.r ** 2 * self.axis.length
[docs] def contains(self, point): """Test region membership of `point` in :class:`Cylinder`. Parameters ---------- point : array_like Returns ------- :class:`~python:bool` `True` if `point` is within :class:`Cylinder`, `False` otherwise. Notes ----- A `point` :math:`(p_x, p_y, p_z)` is within the bounded region of a cylinder with radius :math:`r` around the line from :math:`p_1=(x_1, y_1, z_1)` to :math:`p_2 = (x_2, y_2, z_2)` if the following is true: .. math:: 0\\le q\\le 1\\land (x_1 - p_x + (x_2 - x_1) q)^2 + (y_1 - p_y + (y_2 - y_1) q)^2 + (z_1 - p_z + (z_2 - z_1) q)^2 \\le r^2 where :math:`q` is: .. math:: q = \\frac{(p_x - x_1)(x_2 - x_1) + (p_y - y_1)(y_2 - y_1) + (p_z - z_1)(z_2 - z_1)}{(x_2 - x_1)^2 + (y_2 - y_1)^2 + (z_2 - z_1)^2} """ px, py, pz = Point(point) x1, y1, z1 = self.p1 x2, y2, z2 = self.p2 r = self.r q1 = ((px - x1) * (x2 - x1) + (py - y1) * (y2 - y1) + (pz - z1) * (z2 - z1)) / \ ((x2 - x1) ** 2 + (y2 - y1) ** 2 + (z2 - z1) ** 2) q2 = (x1 - px + (x2 - x1) * q1) ** 2 + \ (y1 - py + (y2 - y1) * q1) ** 2 + \ (z1 - pz + (z2 - z1) * q1) ** 2 return (not np.allclose(x1, x2) or not np.allclose(y1, y2) or not np.allclose(z1, z2)) and r > 0 and q1 >= 0 and q1 <= 1 \ and q2 <= r ** 2
[docs] def todict(self): """Returns a :class:`~python:dict` of the :class:`Cylinder` \ constructor parameters.""" return dict(p1=self.p1, p2=self.p2, r=self.r)
[docs]class Cone(Geometric3DRegion): """`Geometric3DRegion` for a cone. .. versionadded:: 0.3.10 Represents a cone with a base of radius :math:`r` centered at :math:`p_1=(x_1,y_1,z_1)` and a tip at :math:`p_2=(x_2, y_2, z_2)`. Parameters ---------- p1, p2 : array_like, optional 3-tuples or :class:`~sknano.core.Point` class instances for a :class:`Cone` with a base of radius `r` centered at :math:`p_1=(x_1,y_1,z_1)` and a tip at :math:`p_2=(x_2,y_2,z_2)`. r : float, optional Radius :math:`r` of :class:`Cone` base Notes ----- :class:`Cone` represents a bounded cone region :math:`\\left\\{p_1+\\rho(1-z)\\cos(\\theta)\\mathbf{v}_1 + \\rho(1-z)\\sin(\\theta)\\mathbf{v}_2+\\mathbf{v}_3 z| 0\\le\\theta\\le 2\pi\\land 0\\le\\rho\\le 1\\land 0\\le z\\le 1\\right\\}` where :math:`\\mathbf{v}_3=p_2-p_1` and vectors :math:`(\\mathbf{v}_1,\\mathbf{v}_2,\\mathbf{v}_3)` are orthogonal with :math:`|\\mathbf{v}_1|=|\\mathbf{v}_2|=1` and :math:`p_1=(x_1,y_1,z_1)` and :math:`p_2=(x_2, y_2, z_2)`. Calling :class:`Cone` with no parameters is equivalent to :class:`Cone`\ `(p1=[0, 0, 0], p2=[0, 0, 2], r=1)`. """ def __init__(self, p1=None, p2=None, r=1): super().__init__() if p1 is None: p1 = [0, 0, 0] self._p1 = Point(p1) if p2 is None: p2 = [0, 0, 2] self._p2 = Point(p2) self.r = r self.points.extend([self.p1, self.p2]) self.fmtstr = "p1={p1!r}, p2={p2!r}, r={r:.3f}" @property def r(self): """Radius :math:`r` of :class:`Cone` base.""" return self._r @r.setter def r(self, value): self._r = float(value) @property def p1(self): """Center point :math:`(x_1, y_1, z_1)` of :class:`Cone` base.""" return self._p1 @p1.setter def p1(self, value): if not isinstance(value, (tuple, list, np.ndarray)) or len(value) != 3: raise TypeError('Expected a 3-element array_like object') self._p1[:] = Point(value) @property def p2(self): """Point :math:`(x_2, y_2, z_2)` of :class:`Cone` tip.""" return self._p2 @p2.setter def p2(self, value): if not isinstance(value, (tuple, list, np.ndarray)) or len(value) != 3: raise TypeError('Expected a 3-element array_like object') self._p2[:] = Point(value) @property def axis(self): """:class:`Cone` axis :class:`Vector` from \ :math:`\\boldsymbol{\\ell}=p_2 - p_1` Returns ------- :class:`Vector` """ return Vector(p0=self.p1, p=self.p2) @property def centroid(self): """:class:`Cone` centroid, :math:`(c_x, c_y, c_z)`. Computed as: .. math:: c_x = \\frac{3 x_1 + x_2}{4} c_y = \\frac{3 y_1 + y_2}{4} c_z = \\frac{3 z_1 + z_2}{4} Returns ------- :class:`Point` 3D :class:`Point` of :class:`Cone` centroid. """ return (3 * self.p1 + self.p2) / 4 @property def volume(self): """:class:`Cone` volume, :math:`V=\\frac{1}{3}\\pi r^2 \\ell`.""" return np.pi * self.r ** 2 * self.axis.length / 3
[docs] def contains(self, point): """Test region membership of `point` in :class:`Cone`. Parameters ---------- point : array_like Returns ------- :class:`~python:bool` `True` if `point` is within :class:`Cone`, `False` otherwise. Notes ----- A `point` :math:`(p_x, p_y, p_z)` is within the bounded region of a cone with a base of radius :math:`r` centered at :math:`p_1=(x_1, y_1, z_1)` and tip at :math:`p_2 = (x_2, y_2, z_2)` if the following is true: .. math:: 0\\le q\\le 1\\land (x_1 - p_x + (x_2 - x_1) q)^2 + (y_1 - p_y + (y_2 - y_1) q)^2 + (z_1 - p_z + (z_2 - z_1) q)^2 \\le r^2 q^2 where :math:`q` is: .. math:: q = \\frac{(p_x - x_1)(x_2 - x_1) + (p_y - y_1)(y_2 - y_1) + (p_z - z_1)(z_2 - z_1)}{(x_2 - x_1)^2 + (y_2 - y_1)^2 + (z_2 - z_1)^2} """ px, py, pz = Point(point) x1, y1, z1 = self.p1 x2, y2, z2 = self.p2 r = self.r q1 = ((px - x1) * (x2 - x1) + (py - y1) * (y2 - y1) + (pz - z1) * (z2 - z1)) / \ ((x2 - x1) ** 2 + (y2 - y1) ** 2 + (z2 - z1) ** 2) q2 = (x1 - px + (x2 - x1) * q1) ** 2 + \ (y1 - py + (y2 - y1) * q1) ** 2 + \ (z1 - pz + (z2 - z1) * q1) ** 2 q3 = r ** 2 * q1 ** 2 return (not np.allclose(x1, x2) or not np.allclose(y1, y2) or not np.allclose(z1, z2)) and r > 0 and q1 >= 0 and q1 <= 1 \ and q2 <= q3
[docs] def todict(self): """Returns a :class:`~python:dict` of the :class:`Cone` \ constructor parameters.""" return dict(p1=self.p1, p2=self.p2, r=self.r)
geometric_3D_regions = \ [Parallelepiped, Cuboid, Cube, Ellipsoid, Sphere, Cylinder, Cone]