Fix internal trapezoid triangulation & triangle angle calculations

- Fix Triangle.getAngles
- Fix Triangle.getMinAngle
- Fix optimal bisector triangle computations
- Fix trapezoid computations

Triangulation is (mostly) working now. Still need to handle the edge
case where the base of a trapezoid corresponds to multiple upper bound
vectors.

src/pslg.ts

@@ -1,7 +1,19 @@
-import {Matrix} from "./linear";
-
 export type Angle = number
 
+/**
+ * A better rounding function than .toFixed(...).
+ * @param value
+ * @param precision
+ */
+export function safeRound(value: number, precision = 12) {
+    return parseFloat(
+        Math.round(
+            // @ts-ignore
+            value.toFixed(precision + 1) + 'e' + precision
+        ) + 'e-' + precision
+    );
+}
+
 export function deg2rad(degrees: number): number {
     return degrees * (Math.PI / 180)
 }
@@ -30,8 +42,8 @@ export class Point {
         public readonly name?: string,
     ) {
         this.coordinate = {
-            x: parseFloat(coordinate.x.toFixed(12)),
+            x: safeRound(coordinate.x),
-            y: parseFloat(coordinate.y.toFixed(12)),
+            y: safeRound(coordinate.y),
         }
     }
 
@@ -179,10 +191,7 @@ export class Segment {
             )
         }
 
-        return (
+        return Point.from(this.getXAtY(x.y), this.getYAtX(x.x)).is(x)
-            this.getYAtX(x.x) === x.y
-            && this.getXAtY(x.y) === x.x
-        )
     }
 
     yValueIsWithinRange(y: number, inclusive = true) {
@@ -395,15 +404,15 @@ export class Trapezoid {
  */
 export class Triangle {
     get a(): Point {
-        return this.sides[0].from
+        return this.getPoints()[0]
     }
 
     get b(): Point {
-        return this.sides[1].from
+        return this.getPoints()[1]
     }
 
     get c(): Point {
-        return this.sides[2].from
+        return this.getPoints()[2]
     }
 
     get orderedSides(): [Segment, Segment, Segment] {
@@ -492,12 +501,22 @@ export class Triangle {
 
     /** Get the points of the triangle a, b, c, respectively. */
     getPoints(): [Point, Point, Point] {
-        return [this.a, this.b, this.c]
+        let points: Point[] = []
+        this.sides.some(side => {
+            if ( !points.some(point => point.is(side.from)) ) points.push(side.from)
+            if ( !points.some(point => point.is(side.to)) ) points.push(side.to)
+        })
+
+        points = points.sort((a, b) => {
+            if ( a.x === b.x ) return a.y - b.y
+            return a.x - b.x
+        })
+
+        return [points[0], points[1], points[2]]
     }
 
     getCircumcenter(): Point {
         const [pointA, pointB, pointC] = this.getPoints()
-        const [angleA, angleB, angleC] = this.getAngles()
 
         const [sin2A, sin2B, sin2C] = this.getAngles().map(x => Math.sin(2 * x))
 
@@ -524,7 +543,6 @@ export class Triangle {
         const numerator = p2.y * (p1.x - p3.x) + p1.y * (p3.x - p2.x) + p3.y * (p2.x - p1.x)
         const denominator = (p2.x - p1.x) * (p1.x - p3.x) + (p2.y - p1.y) * (p1.y - p3.y)
         const radio = numerator / denominator
-
         return Math.atan(radio)
     }
 }
@@ -735,9 +753,15 @@ export class Graph {
         if ( existing ) return existing
 
         this.segments.push(x)
+        this.findExistingPointOrAdd(x.from)
+        this.findExistingPointOrAdd(x.to)
         return x
     }
 
+    hasExistingSegment(x: Segment) {
+        return this.segments.some(segment => segment.is(x))
+    }
+
     findExistingTriangleOrAdd(x: Triangle): Triangle {
         const existing = this.triangles.find(triangle => triangle.is(x))
         if ( existing ) return existing
@@ -767,6 +791,7 @@ export class Graph {
             const newFrom = newPoints.find(point => point.is(segment.from))
             const newTo = newPoints.find(point => point.is(segment.to))
             if ( !newFrom || !newTo ) {
+                console.log({from: segment.from.coordinate, to: segment.to.coordinate})
                 throw new Error('Tried to clone segment, but could not match all points')
             }
 

src/trapezoidTriangulation.ts

@@ -9,7 +9,7 @@ import {
     SegmentWithIntersection, Triangle
 } from "./pslg";
 
-export function getFirstIntersectingSegmentInDirection(raySegment: Segment, boundary: GraphBoundary, graph: Graph, direction: GraphDirection): [Segment, Point] {
+export function getFirstIntersectingSegmentInDirection(raySegment: Segment, boundary: GraphBoundary, graph: Graph, direction: GraphDirection, inclusive = false): [Segment, Point] {
     const intersectingSegment = boundary.getBoundary(direction)
     const intersectingPoint = raySegment.getIntersectionWith(intersectingSegment)
     if ( !intersectingPoint ) {
@@ -29,7 +29,7 @@ export function getFirstIntersectingSegmentInDirection(raySegment: Segment, boun
 
             return {
                 segment,
-                intersect: segment.getIntersectionWithin(raySegment)
+                intersect: segment[inclusive ? 'getIntersectionWith' : 'getIntersectionWithin'](raySegment)
             }
         })
         .filter(x => x && x.intersect) as SegmentWithIntersection[])
@@ -45,12 +45,13 @@ export function getFirstIntersectingSegmentInDirection(raySegment: Segment, boun
 
 export function triangulate(originalGraph: Graph): Graph {
     const graph = originalGraph.clone()
-    const boundary = addBoundingSquareTo(graph)
+    const trapezoidSegments: Segment[] = [] //graph.segments.filter(segment => segment.isHorizontal())
 
+    const boundary = addBoundingSquareTo(graph)
     const leftBound = boundary.getLeftBoundary()
     const rightBound = boundary.getRightBoundary()
 
-    const trapezoidSegments: Segment[] = []
+    // trapezoidSegments.push(boundary.getLowerBoundary())
 
     // For each vertex in the original graph, create a horizontal line that
     // extends in both directions until it intersects with either (1) the boundary
@@ -123,15 +124,21 @@ export function triangulate(originalGraph: Graph): Graph {
         }
     }
 
+    // Any horizontal segments present in the original graph will also be used for form
+    // trapezoids, so push them onto the list of trapezoid base segments.
+    originalGraph.segments
+        .filter(segment => segment.isHorizontal())
+        .forEach(segment => trapezoidSegments.push(segment))
+
     // Now, go through and identify trapezoids for all the horizontal segments we just added
     for ( const segment of trapezoidSegments ) {
-        // First, find the trapezoid formed with the segment as the bottom
+        // Find the trapezoid formed with the segment as the bottom
         // Create a vertical segment from the midpoint of the segment to the top boundary
         const horizontalMidpoint = segment.getMidpoint()
         let upperBoundaryPoint = Point.from(horizontalMidpoint.x, boundary.ymax)
         let upperBoundaryVerticalSegment = new Segment(horizontalMidpoint, upperBoundaryPoint)
 
-        const [upperIntersectSegment, upperIntersectPoint] = getFirstIntersectingSegmentInDirection(
+        let [upperIntersectSegment, upperIntersectPoint] = getFirstIntersectingSegmentInDirection(
             upperBoundaryVerticalSegment,
             boundary,
             graph,
@@ -151,10 +158,10 @@ export function triangulate(originalGraph: Graph): Graph {
             leftBoundaryHorizontalSegment,
             boundary,
             graph,
-            GraphDirection.LEFT
+            GraphDirection.LEFT,
+            true
         )
 
-        console.log('got leftIntersectSegment', leftBoundaryHorizontalSegment.toQuickDisplay(), leftIntersectSegment.toQuickDisplay(), leftIntersectPoint.coordinate)
         leftBoundaryHorizontalSegment = new Segment(verticalMidpoint, leftIntersectPoint)
 
         // Repeat to get the right boundary
@@ -166,10 +173,38 @@ export function triangulate(originalGraph: Graph): Graph {
             boundary,
             graph,
             GraphDirection.RIGHT,
+            true
         )
 
         rightBoundaryHorizontalSegment = new Segment(verticalMidpoint, rightIntersectPoint)
 
+        // Check if the upper boundary segment extends beyond the x-range of the left- and right-boundary segments
+        // If so, we need to split it to fit within the bounds of the current trapezoid, starting with the right side
+        if ( upperIntersectSegment.xmax > rightIntersectSegment.xmax ) {
+            let [upperIntersectSplit1, upperIntersectSplit2] = upperIntersectSegment.splitAt(
+                Point.from(rightIntersectSegment.xmax, upperIntersectSegment.ymax)
+            )
+
+            graph.removeSegment(upperIntersectSegment)
+            upperIntersectSplit1 = graph.findExistingSegmentOrAdd(upperIntersectSplit1)
+            upperIntersectSplit2 = graph.findExistingSegmentOrAdd(upperIntersectSplit2)
+
+            upperIntersectSegment = upperIntersectSplit1.xmax === rightIntersectSegment.xmax ? upperIntersectSplit1 : upperIntersectSplit2
+        }
+
+        // Repeat for the left side
+        if ( upperIntersectSegment.xmin < leftIntersectSegment.xmin ) {
+            let [upperIntersectSplit1, upperIntersectSplit2] = upperIntersectSegment.splitAt(
+                Point.from(leftIntersectSegment.xmin, upperIntersectSegment.ymax)
+            )
+
+            graph.removeSegment(upperIntersectSegment)
+            upperIntersectSplit1 = graph.findExistingSegmentOrAdd(upperIntersectSplit1)
+            upperIntersectSplit2 = graph.findExistingSegmentOrAdd(upperIntersectSplit2)
+
+            upperIntersectSegment = upperIntersectSplit1.xmin === leftIntersectSegment.xmin ? upperIntersectSplit1 : upperIntersectSplit2
+        }
+
         // Now, check if we actually have a 4-bound trapezoid, or if we have a triangle
         const points = Point.distinct([
             segment.from,
@@ -178,21 +213,6 @@ export function triangulate(originalGraph: Graph): Graph {
             upperIntersectSegment.to,
         ])
 
-        if ( points.length === 3 ) {
-            // We found a triangle! Less work.
-            // Create the triangle and push it onto the graph
-            const [p1, p2, p3] = points.map(x => graph.findExistingPointOrAdd(x))
-            const s12 = graph.findExistingSegmentOrAdd(new Segment(p1, p2))
-            const s23 = graph.findExistingSegmentOrAdd(new Segment(p2, p3))
-            const s31 = graph.findExistingSegmentOrAdd(new Segment(p3, p1))
-            graph.findExistingTriangleOrAdd(new Triangle([s12, s23, s31]))
-            continue  // FIXME - remove to handle below-segment case
-        }
-
-        if ( points.length !== 4 ) {
-            throw new RangeError('Found shape with invalid number of distinct points!')
-        }
-
         // Now, we have the 4 bounding segments of the trapezoid.
         // Let's find the segments that make up the trapezoid
         // We will do this by re-creating segments for the four sides of the trapezoid
@@ -205,33 +225,76 @@ export function triangulate(originalGraph: Graph): Graph {
 
         // TODO Account for the case where we don't need to split the segment.
         let trapezoidLeftBoundSegment = leftIntersectSegment
-        let leftSegment1: Segment | undefined
+        // let leftSegment1: Segment | undefined
-        let leftSegment2: Segment | undefined
+        // let leftSegment2: Segment | undefined
         if ( !leftIntersectSegment.hasPoint(leftSegmentIntersectPoint) ) {
-            let [localLeftSegment1, localLeftSegment2] = leftIntersectSegment.splitAt(leftSegmentIntersectPoint)
+            let [leftSegment1, leftSegment2] = leftIntersectSegment.splitAt(leftSegmentIntersectPoint)
             graph.removeSegment(leftIntersectSegment)
-            leftSegment1 = graph.findExistingSegmentOrAdd(localLeftSegment1)
+            leftSegment1 = graph.findExistingSegmentOrAdd(leftSegment1)
-            leftSegment2 = graph.findExistingSegmentOrAdd(localLeftSegment2)
+            leftSegment2 = graph.findExistingSegmentOrAdd(leftSegment2)
 
             // We care about the upper-segment from the split, as that is the bound of our trapezoid
             trapezoidLeftBoundSegment = leftSegment1.ymin === leftSegmentIntersectPoint.y ? leftSegment1 : leftSegment2
+
+            // Now, we need to consider the case where the upper segment we split extends beyond the upper bound of
+            // the trapezoid we are working with now. If so, split the upper segment again.
+            if ( trapezoidLeftBoundSegment.ymax > upperIntersectPoint.y && points.length > 3 ) {
+                // The left bound extends beyond the top of this trapezoid. So, split it.
+                let localLeftUpperSplitPoint = Point.from(trapezoidLeftBoundSegment.xmin, upperIntersectPoint.y)
+                let [leftUpperSegment1, leftUpperSegment2] = trapezoidLeftBoundSegment.splitAt(localLeftUpperSplitPoint)
+                graph.removeSegment(trapezoidLeftBoundSegment)
+                leftUpperSegment1 = graph.findExistingSegmentOrAdd(leftUpperSegment1)
+                leftUpperSegment2 = graph.findExistingSegmentOrAdd(leftUpperSegment2)
+
+                trapezoidLeftBoundSegment = leftUpperSegment1.ymax === upperIntersectPoint.y ? leftUpperSegment1 : leftUpperSegment2
             }
+        }
+
+        graph.findExistingSegmentOrAdd(trapezoidLeftBoundSegment)
 
         // Repeat this process for the right-side segment
         const rightSegmentIntersectPoint = rightIntersectSegment.getIntersectionWith(segment)
         if ( !rightSegmentIntersectPoint ) throw new Error('Unable to find trapezoid segment intersection')
 
         let trapezoidRightBoundSegment = rightIntersectSegment
-        let rightSegment1: Segment | undefined
-        let rightSegment2: Segment | undefined
         if ( !rightIntersectSegment.hasPoint(rightSegmentIntersectPoint) ) {
-            let [localRightSegment1, localRightSegment2] = rightIntersectSegment.splitAt(rightSegmentIntersectPoint)
+            let [rightSegment1, rightSegment2] = rightIntersectSegment.splitAt(rightSegmentIntersectPoint)
             graph.removeSegment(rightIntersectSegment)
-            rightSegment1 = graph.findExistingSegmentOrAdd(localRightSegment1)
+            rightSegment1 = graph.findExistingSegmentOrAdd(rightSegment1)
-            rightSegment2 = graph.findExistingSegmentOrAdd(localRightSegment2)
+            rightSegment2 = graph.findExistingSegmentOrAdd(rightSegment2)
 
             // We care about the upper-segment from the split, as that is the bound of our trapezoid
             trapezoidRightBoundSegment = rightSegment1.ymin === rightSegmentIntersectPoint.y ? rightSegment1 : rightSegment2
+
+            // Now, we need to consider the case where the upper segment we split extends beyond the upper bound of
+            // the trapezoid we are working with now. If so, split the upper segment again.
+            if ( trapezoidRightBoundSegment.ymax > upperIntersectPoint.y && points.length > 3 ) {
+                // The left bound extends beyond the top of this trapezoid. So, split it.
+                let localRightUpperSplitPoint = Point.from(trapezoidRightBoundSegment.xmin, upperIntersectPoint.y)
+                let [rightUpperSegment1, rightUpperSegment2] = trapezoidRightBoundSegment.splitAt(localRightUpperSplitPoint)
+                graph.removeSegment(trapezoidRightBoundSegment)
+                rightUpperSegment1 = graph.findExistingSegmentOrAdd(rightUpperSegment1)
+                rightUpperSegment2 = graph.findExistingSegmentOrAdd(rightUpperSegment2)
+
+                trapezoidRightBoundSegment = rightUpperSegment1.ymax === upperIntersectPoint.y ? rightUpperSegment1 : rightUpperSegment2
+            }
+        }
+
+        // break;
+
+        if ( points.length === 3 ) {
+            // We found a triangle! Less work.
+            // Create the triangle and push it onto the graph
+            // const [p1, p2, p3] = points.map(x => graph.findExistingPointOrAdd(x))
+            // const s12 = graph.findExistingSegmentOrAdd(new Segment(p1, p2))
+            // const s23 = graph.findExistingSegmentOrAdd(new Segment(p2, p3))
+            // const s31 = graph.findExistingSegmentOrAdd(new Segment(p3, p1))
+            graph.findExistingTriangleOrAdd(new Triangle([trapezoidLeftBoundSegment, trapezoidRightBoundSegment, segment]))
+            continue  // FIXME - remove to handle below-segment case
+        }
+
+        if ( points.length !== 4 ) {
+            throw new RangeError('Found shape with invalid number of distinct points!')
         }
 
         // Now we have all 4 bounding segments. We find the bisector that creates
@@ -243,29 +306,27 @@ export function triangulate(originalGraph: Graph): Graph {
 
         const bottomLeftBisectorSegment = new Segment(lowerLeftPoint, upperRightPoint)
         const bottomLeftBisectorUpperTriangle = new Triangle([bottomLeftBisectorSegment, upperIntersectSegment, trapezoidLeftBoundSegment])
-        // const bottomLeftBisectorLowerTriangle = new Triangle([bottomLeftBisectorSegment, segment, trapezoidRightBoundSegment])
+        const bottomLeftBisectorLowerTriangle = new Triangle([bottomLeftBisectorSegment, segment, trapezoidRightBoundSegment])
-        // const bottomLeftBisectorMinAngle = Math.min(bottomLeftBisectorUpperTriangle.getMinimumAngle(), bottomLeftBisectorLowerTriangle.getMinimumAngle())
+        const bottomLeftBisectorMinAngle = Math.min(bottomLeftBisectorUpperTriangle.getMinimumAngle(), bottomLeftBisectorLowerTriangle.getMinimumAngle())
 
-        // const upperLeftPoint = graph.findExistingPointOrAdd(Point.from(upperIntersectSegment.xmin, upperIntersectSegment.ymax))
+        const upperLeftPoint = graph.findExistingPointOrAdd(Point.from(upperIntersectSegment.xmin, upperIntersectSegment.ymax))
-        // const lowerRightPoint = graph.findExistingPointOrAdd(Point.from(segment.xmax, segment.ymin))
+        const lowerRightPoint = graph.findExistingPointOrAdd(Point.from(segment.xmax, segment.ymin))
-        //
+
-        // const topRightBisectorSegment = new Segment(upperLeftPoint, lowerRightPoint)
+        const topRightBisectorSegment = new Segment(upperLeftPoint, lowerRightPoint)
-        // const upperRightBisectorUpperTriangle = new Triangle([topRightBisectorSegment, upperIntersectSegment, trapezoidRightBoundSegment])
+        const upperRightBisectorUpperTriangle = new Triangle([topRightBisectorSegment, upperIntersectSegment, trapezoidRightBoundSegment])
-        // const upperRightBisectorLowerTriangle = new Triangle([topRightBisectorSegment, trapezoidLeftBoundSegment, segment])
+        const upperRightBisectorLowerTriangle = new Triangle([topRightBisectorSegment, trapezoidLeftBoundSegment, segment])
-        // const upperRightBisectorMinAngle = Math.min(upperRightBisectorUpperTriangle.getMinimumAngle(), upperRightBisectorLowerTriangle.getMinimumAngle())
+        const upperRightBisectorMinAngle = Math.min(upperRightBisectorUpperTriangle.getMinimumAngle(), upperRightBisectorLowerTriangle.getMinimumAngle())
-        //
+
-        // const optimalBisectorUpperTriangle = upperRightBisectorMinAngle > bottomLeftBisectorMinAngle ? upperRightBisectorUpperTriangle : bottomLeftBisectorUpperTriangle
+        const optimalBisectorUpperTriangle = upperRightBisectorMinAngle > bottomLeftBisectorMinAngle ? upperRightBisectorUpperTriangle : bottomLeftBisectorUpperTriangle
-        // const optimalBisectorLowerTriangle = upperRightBisectorMinAngle > bottomLeftBisectorMinAngle ? upperRightBisectorLowerTriangle : bottomLeftBisectorLowerTriangle
+        const optimalBisectorLowerTriangle = upperRightBisectorMinAngle > bottomLeftBisectorMinAngle ? upperRightBisectorLowerTriangle : bottomLeftBisectorLowerTriangle
-        //
+
-        // // Add the triangles to the graph
+        // Add the triangles to the graph
-        // const upperTriangleSegments = optimalBisectorUpperTriangle.sides.map(side => graph.findExistingSegmentOrAdd(side))
+        const upperTriangleSegments = optimalBisectorUpperTriangle.sides.map(side => graph.findExistingSegmentOrAdd(side))
-        // graph.findExistingTriangleOrAdd(new Triangle(upperTriangleSegments as [Segment, Segment, Segment]))
+        graph.findExistingTriangleOrAdd(new Triangle(upperTriangleSegments as [Segment, Segment, Segment]))
-        //
+
-        // const lowerTriangleSegments = optimalBisectorLowerTriangle.sides.map(side => graph.findExistingSegmentOrAdd(side))
+        const lowerTriangleSegments = optimalBisectorLowerTriangle.sides.map(side => graph.findExistingSegmentOrAdd(side))
-        // graph.findExistingTriangleOrAdd(new Triangle(lowerTriangleSegments as [Segment, Segment, Segment]))
+        graph.findExistingTriangleOrAdd(new Triangle(lowerTriangleSegments as [Segment, Segment, Segment]))
     }
 
-    // FIXME handle the lower-trapezoid case
-
     return graph
 }