Add test suite for triangulation code, fix edge cases

[hannobraun]

The triangulation algorithm first does a Delaunay triangulation, then filters out triangles that are not part of the face, by performing a series of point-in-polygon tests. It does this, by shooting a ray from the tested point to another point outside of the polygon, and counting how often that ray hits a polygon edge.

I believe that the algorithm does not address all edges cases. I can think of the following ones that are currently missing.

Ray hits edge that is parallel to it

1. The ray would be tested against the edge before (from the ray's perspective) the parallel edge, hitting it in the vertex those two edges share. The count would be 1.
2. The ray would be tested against the parallel edge, hitting it in the same vertex. This would be recognized as a duplicate hit and not be counted again.
3. The ray would be tested against the edge after the parallel edge, hitting it in the vertex those two edges share. The count would be 2.

A count of 2 would indicate that the ray has not crossed the polygon boundary. This would only be correct, if the polygon was concave in that location. For a convex polygon, or one convex in that location at least, this would be an incorrect count.

Ray hits vertex at a concave location

1. The ray is tested against one of the edges that share the vertex. The count would be 1.
2. The ray is tested against the other edge that shares the vertex. This would be recognized as a duplicate hit and not counted.

A count of 1 indicates that the ray has crossed the polygon boundary, which is not correct in this case. It just touched the polygon, but staid inside.

---

The triangulation algorithm was somewhat hacked together by me, and there is no test suite. I think it's high time to move the algorithm to a dedicated module and add an extensive test suite. Once covering the existing functionality, the test suite should be extended to verify and fix the edge cases I have presented here.

[hannobraun]

I've started working on this.

[alexoro412]

I came across what might be another edge case in the triangulation algorithm so I thought I'd add it here. I was playing around with the star model, and I modified it to make every 3rd point slightly longer:

let vertex_iter = (0..num_vertices).map(|i| {
        let angle = 2. * PI / num_vertices as f64 * i as f64;
        let radius = if i % 2 == 0 {
            r1
        } else {
            r2 + (if i % 3 == 0 { 0.2 } else { 0.0 })
        };
        (angle, radius)
    });

But instead of making every 3rd point slightly longer, every 3rd point just disappears (I'm using num_points=12 here):

I commented out the filtering step in the triangulation algorithm, and you can see below the mesh it generates for my example.

I think the issue here is that the Delaunay step is producing some triangles that are partially inside and partially outside the shape. When this happens, an edge that was in the original shape does not appear in the triangulation. Then, these triangles either get filtered out, which removes a chunk of the shape, or they get left in, which adds a new piece to the shape.

One way to approach this might be to have a step after Delaunay that cuts all of the triangles along the edges of the shape, although I'm not sure how efficient/robust that would be.

[hannobraun]

Thank you for reporting this, @alexoro412!

I've been working towards a better point-in-polygon test recently, using an external library (geo), which will hopefully fix a lot of issues that currently originate in my own crummy code. This will take a bit more time though, as some cleanups have been, and still are, necessary, to integrate geo, and the computer that has all my in-progress work is currently out for repair.

I hope that this will take care of a lot of these issues. Once geo is integrated, I'll revisit all the known triangulation bugs, including yours.

I think the issue here is that the Delaunay step is producing some triangles that are partially inside and partially outside the shape.

Oh, that is a very interesting idea! I didn't consider it before, but I suppose it could happen.

There's a feature in spade (the library we're using for Delaunay triangulation) which lets you define a bunch of edges that have to be triangle edges in the final triangulation. I think we have to pass the face edges there to prevent that case.

Thank you, this is very interesting! I don't know if your specific issue is indeed caused by this (I'd have to take a closer look at the triangles to be sure), but I think this is definitely a real bug.

[hannobraun]

@alexoro412 I've decided to open a separate issue (#430) to track the corner case you reported, as the corner cases tracked in this issue are related to the point-in-polygon test. Thanks again!