''' Cone is a cone, a geometry.
.. codeauthor:: Yoshifumi Futaana
'''
import numpy as _np
[docs]def get_surface_vectors(axis_vector, cone_angle_deg, ndiv=180):
''' Make vectors of cone surfaces.
:param axis_vector: A 3-D vector.
:type axis_vector: ``np.array``
:param cone_angle_deg: Cone angle in degrees
:type cone_angle_deg: ``float``
:returns: Vectors in (3, 180) shaped ndarray. Normalized.
:rtype: ``np.array``
The cone vectors, with its center acis vextor of (1, 1, -1)
and the angle 60 degrees can be calculated as follows.
>>> vectors = get_surface_vectors((1., 1., -1.), 60.0)
>>> print(vectors.shape)
(3, 180)
>>> import numpy as np
>>> cos_angles = np.array((1., 1., -1)).dot(vectors) / np.sqrt(3)
>>> np.allclose(cos_angles, 0.5) # All the values are 0.5, i.e. 60 deg
True
'''
axis = _np.array(axis_vector)
naxis = _np.linalg.norm(axis)
nvaxis = axis / naxis
# Original vectors (normalized)
phi = _np.linspace(0, 2 * _np.pi, ndiv)
cone_angle = _np.radians(cone_angle_deg)
vec = _np.array([_np.sin(cone_angle) * _np.cos(phi),
_np.sin(cone_angle) * _np.sin(phi),
_np.zeros_like(phi) + _np.cos(cone_angle) ])
# zeros... is for broadcast.
lat = _np.arccos(nvaxis[2])
lon = _np.arctan2(nvaxis[1], nvaxis[0])
mat1 = _np.array([[_np.cos(lon), _np.sin(lon), 0],
[-_np.sin(lon), _np.cos(lon), 0],
[0, 0, 1]]).T
mat2 = _np.array([[_np.cos(lat), 0, -_np.sin(lat)],
[0, 1, 0],
[_np.sin(lat), 0, _np.cos(lat)]]).T
return mat1.dot(mat2).dot(vec)
[docs]class Cone:
""" A cone. Single direction.
"""
def __init__(self, apex, axis, cone_angle):
"""
:param apex: The location of the apex. (3,) shaped numpy array
:param axis: The direction. (3,) shaped numpy array
:param cone_angle: Cone angle in degrees
"""
self.apex = _np.array(apex)
"""The coordinate of the apex"""
self.axis = _np.array(axis) / _np.linalg.norm(axis)
""""Normal vector of the cone central axis"""
self.cone_angle = cone_angle
"""Cone angle in degrees"""
[docs] def get_surface_vectors(self, length=1.0, ndiv=180):
"""
:param ndiv:
:return: (ndiv, 3) shaped vector
.. note::
The returned dimension is different (transposed) from :func:`get_surface_vectors`!!
It is because the (N, 3) shape is more intuitive for specifying multiple 3-D vectors.
>>> cone = Cone([-1, 0, 0], [1, 0, 0], 45)
>>> surf = cone.get_surface_vectors(ndiv=4)
>>> print(surf)
[[-2.92893219e-01 0.00000000e+00 -7.07106781e-01]
[-2.92893219e-01 6.12372436e-01 3.53553391e-01]
[-2.92893219e-01 -6.12372436e-01 3.53553391e-01]
[-2.92893219e-01 -1.73191211e-16 -7.07106781e-01]]
"""
vecs = get_surface_vectors(self.axis, self.cone_angle, ndiv=ndiv) * length
vecs = vecs + self.apex[:, _np.newaxis]
return _np.copy(vecs.T)
[docs] def get_surface_vectors_from_apex(self, length=1.0, ndiv=180):
vecs = get_surface_vectors(self.axis, self.cone_angle, ndiv=ndiv) * length
return vecs.T
def __contains__(self, coord):
""" Returns if a single point coordinate is inside the cone or not.
:param coord: (3,) shaped numpy array.
:return: Boolean. True if coord in the cone.
>>> c = Cone([0, 0, 0], [1, 0, 0], 45)
>>> [3, 0, 0] in c # On the axis
True
>>> [-3, 0, 0] in c # On the axis but different direction
False
>>> [0, 3, 0] in c # Off the axis
False
>>> [3, 3, 0] in c # On the surface
True
>>> [3, 3.1, 0] in c # Beyond the surface
False
>>> [0, 0, 0] in c # At apex
True
"""
rvec = _np.array(coord) - self.apex
lrvec = _np.linalg.norm(rvec)
if lrvec == 0:
return True # The vertex
rvec = rvec / lrvec
inner = _np.inner(rvec, self.axis) # self.axis has been normalized
return inner >= _np.cos(_np.deg2rad(self.cone_angle))
[docs] def within(self, points):
""" Return a list of in-out.
Multipoint-version of in-out determinator.
For a single point, you can use
>>> cone = Cone([0, 0, 0], [1, 0, 0], 45) # Apex at origin, axis in a-direction, 45 deg cone-angle
>>> (3, 0, 0) in cone # Return if (3, 0, 0) is in the cone.
True
This :meth:`within` method get multiple point.
>>> points = [[0, 0, 0], [1, 0, 0], [2, 0, 0], [-1, 0, 0]]
>>> cone.within(points)
array([ True, True, True, False])
It is equivalent to
>>> [p in cone for p in points]
[True, True, True, False]
:param points: (N, 3) shaped array representing 3-D points
:return: (N,) shpaed boolean array
"""
pts = _np.array(points) # (N, 3)
rvec = pts - self.apex[_np.newaxis, :] # (N, 3)
lrvec = _np.linalg.norm(rvec, axis=1) # (N,)
rvec = rvec / lrvec[:, _np.newaxis] # (N, 3)
inner = (rvec * self.axis[_np.newaxis, :]).sum(axis=1) # Inner product
return (inner >= _np.cos(_np.deg2rad(self.cone_angle))) | (lrvec == 0)