Source code for irfpy.util.triangle2d

'''Triangle in 2D.

.. codeauthor:: Yoshifumi Futaana

Triangle, defined by x-y plane, is implemented.

:class:`Triangle2D` provides the main functionality.
:meth:`create_triangle` will create a :class:`Triangle2D` object.

Triangle P0=(1,3), P1=(3,3), P2=(1,8) is created by

>>> tri = create_triangle([1, 3], [3, 3], [1, 8])
>>> print(tri.get_area())  # Calculate the area.
5.0
>>> print(tri.isinside([2, 4]))  # P=(2,4) is inside the triangle
True
>>> print(tri.isinside([1, 2]))  # P=(1,2) is on the triangle
True
>>> print(tri.isinside([1, 2], bound=False))  # P=(1,2) is on the triangle
False

'''

import numpy as np

[docs]def create_triangle(p0, p1, p2): return Triangle2D([p0[0], p1[0], p2[0]], [p0[1], p1[1], p2[1]])
[docs]class Triangle2D: def __init__(self, xlist, ylist): ''' Create 2D triangle. ''' xs = np.array(xlist, dtype=float) ys = np.array(ylist, dtype=float) if xs.shape != (3,): raise ValueError('xlist should have (3,) shape') if ys.shape != (3,): raise ValueError('ylist should have (3,) shape') self.p0 = np.array([xs[0], ys[0]]) ''' The first point ''' self.p1 = np.array([xs[1], ys[1]]) ''' The second point ''' self.p2 = np.array([xs[2], ys[2]]) ''' The third point ''' if self.get_area() == 0: raise ValueError('The given points do not form triangle.')
[docs] def get_area(self): ''' Return the area of the triangle. Triangle, P0 is origin, P1 is x=1, P2 is y=1. >>> tri = Triangle2D([1, 0, 0], [0, 1, 0]) >>> print(tri.get_area()) 0.5 ''' v1 = self.p1 - self.p0 v2 = self.p2 - self.p0 return np.abs(v1[0] * v2[1] - v1[1] * v2[0]) / 2.
[docs] def isinside(self, point, bound=True): ''' Return true if the given point is inside or on the bound of the triangle >>> tri = Triangle2D([0.1, 1, 0.1], [0.1, 0.1, 1]) >>> tri.isinside([0.1, 0.1]) # It is on the bound. True >>> tri.isinside([0.1, 0.1], bound=False) # It is on the bound. False >>> tri.isinside([0.3, 0.3]) True ''' l0 = self.p1 - self.p0 l1 = self.p2 - self.p1 l2 = self.p0 - self.p2 p = np.array(point) c0 = np.cross(self.p1 - self.p0, p - self.p0) c1 = np.cross(self.p2 - self.p1, p - self.p1) c2 = np.cross(self.p0 - self.p2, p - self.p2) cc = np.array([c0, c1, c2]) if (cc == 0).any(): return bound if (cc > 0).all(): return True if (cc < 0).all(): return True return False
def __str__(self): return "<Triangle2D: (P0) %g, %g (P1) %g, %g (P2) %g %g>" % ( self.p0[0], self.p0[1], self.p1[0], self.p1[1], self.p2[0], self.p2[1], ) def __repr__(self): return self.__str__()
[docs]class Polygon2D: ''' Polygon, extending :class:`Triangle2D` using triangulation. >>> poly = Polygon2D([-1, 0, 1, 0], [0, -1, 0, 1]) >>> print(poly.isinside([0, 0])) # Is (0,0) inside? Yes. True ''' def __init__(self, xs, ys): self.xs = np.array(xs) self.ys = np.array(ys) if self.xs.shape != self.ys.shape: raise ValueError('Shape in xs and ys differs')
[docs] def isinside(self, point): tris = self.triangulate() for tri in tris: if tri.isinside(point): return True return False
[docs] def triangulate(self): ''' Triangulation. :returns: A set of :class:`Triangle2D` Indeed, this is not yet fully debugged. c.f. http://sonson.jp/?p=205 ''' xs = self.xs.copy().tolist() ys = self.ys.copy().tolist() n_points = len(self.xs) triangles = [] n_tri = n_points - 2 while len(triangles) < n_tri: xn = np.array(xs) yn = np.array(ys) len2 = xn * xn + yn * yn pi1 = len2.argmax() # The furtherset point. pi0 = (pi1 - 1) % len(xs) pi2 = (pi1 + 1) % len(xs) # Triangle p0 = np.array([xn[pi0], yn[pi0]]) p1 = np.array([xn[pi1], yn[pi1]]) p2 = np.array([xn[pi2], yn[pi2]]) try: tri = create_triangle(p0, p1, p2) # Triangle2D object. except ValueError: # If p0, p1, p2 does not form triangle. n_tri -= 1 xs.pop(pi1) ys.pop(pi1) continue # Outer product masterouterproduct = np.cross(p1 - p0, p2 - p1) while True: outerproduct = np.cross(p1 - p0, p2 - p1) flag = True if outerproduct * masterouterproduct < 0: flag = False for x, y in zip(xn, yn): if tri.isinside([x, y], bound=False): flag = False if flag: triangles.append(tri) xs.pop(pi1) ys.pop(pi1) break else: pi0 = (pi0 - 1) % len(xs) pi1 = (pi1 - 1) % len(xs) pi2 = (pi2 - 1) % len(xs) p0 = np.array([xn[pi0], yn[pi0]]) p1 = np.array([xn[pi1], yn[pi1]]) p2 = np.array([xn[pi2], yn[pi2]]) tri = create_triangle(p0, p1, p2) # Triangle2D object. return triangles
import unittest import doctest
[docs]def doctests(): return unittest.TestSuite(( doctest.DocTestSuite(), ))
if __name__ == '__main__': unittest.main(defaultTest='doctests')