Source code for irfpy.util.cone

''' 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)