Search code examples
algorithmgeometrypolygontraversal

Is there an algorithm to separate polygons that share an edge?


I have a list of vertices that define a polygon, a list of perimeter edges that connect those vertices and define the outline of the polygon, and a list of inner edges connecting vertices, effectively splitting the polygon. None of the edges intersect each other (they only meet at the start and end points).

I want to partition the bigger polygon into its smaller components by splitting it at the inner edges. Basically I need to know which sets of vertices are part of a polygon that has no intersecting edges.

Basically, this is the information I have:

enter image description here

The vertices [0, 1, 3, 5, 7, 9, 11, 13, 15] define the outside perimeter of my polygon and the blue edge (5, 13) is an inner edge splitting that perimeter into two polygons. (Disregard the horizontal purple lines, they are the result of the trapezoidalization of the polygon to find the edge (5, 13). They have no further meaning)

This is what I want to get to:

enter image description here

One polygon is defined by the vertices [0, 1, 3, 5, 13, 15] and the other is defined by [5, 7, 9, 11, 13].

I need a solution that works for any number of splitting inner edges. In the end, I would like to be able to partition a polygon like the following:

enter image description here

The sub-polygons are not necessarily convex. That might be of importance for any algorithm depending on it. The red sub-polygon has a concave vertex (13) for example.

My idea initially was to traverse the inside of each sub-polygon in a clockwise or counter-clockwise direction, keeping track of the vertices I encounter and their order. However I am having trouble finding a starting edge/vertex and guaranteeing that the next cw or ccw point is in fact on the inside of that sub-polygon I want to extract.

I tried to google for solutions but this is a field of mathematics I am too unfamiliar with to know what to search for. An idea/algorithm of how to tackle this problem would be much appreciated! (I don't need any actual code, just an explanation of how to do this or pseudocode would totally suffice)

Now, unfortunately I don't have some code to show as I need a concept to try out first. I don't know how to tackle this problem and thus can't code anything that could accomplish what I need it to do.

EDIT:

This is just one step in what I am trying to do ultimately, which is polygon triangulation. I have read numerous solutions to that problem and wanted to implement it via trapezoidalization to get monotone polygons and ultimately triangulate those. This is basically the last step of (or I guess the next step after) the trapezoidalization, which is never explained in any resources on the topic that i could find.

The end result of the trapezoidalization are the inner edges which define the split into (in my case vertically monotone) polygons. I just need to split the polygons along those edges "datastructurally" so I can work on the subpolygons individually. I hope that helps to clarify things.


Solution

  • The key of the algorithm that you need, is to know how edges can be ordered:

    Finding out which is the next edge in (counter)clockwise order

    You can calculate the absolute angle, of an edge from node i to node j, with the following formula:

    atan2(jy-iy, jx-ix)
    

    See atan2

    The relative angle between an edge (i, j) and (j, k) is then given by:

    atan2(ky-jy, kx-jx) - atan2(jy-iy, jx-ix)
    

    This expression may yield angles outside of the [-𝜋, 𝜋] range, so you should map the result back into that range by adding or subtracting 2𝜋.

    So when you have traversed an edge (i, j) and need to choose the next edge (j, k), you can then pick the edge that has the smallest relative angle.

    The Algorithm

    The partitioning algorithm does not really need to know upfront which edges are internal edges, so I'll assume you just have an edge list. The algorithm could look like this:

    1. Create an adjacency list, so you have a list of neighboring vertices for every vertex.
      • Add every edge to this adjacency list in both directions, so actually adding two directed edges for each original edge
    2. Pick a directed edge (i, j) from the adjacency list, and remove it from there.
    3. Define a new polygon and add vertex i as its first vertex.
    4. Repeat until you get back to the first vertex in the polygon that's being constructed:
      • add vertex j to the polygon
      • find vertex k among the neighbors of vertex j that is not vertex i, and which minimizes the relative angle with the above given formula.
      • add this angle to a sum
      • Delete this directed edge from the neighbors of vertex j, so it will never be visited again
      • let i = j, j = k
    5. If the sum of angles is positive (it will be 2𝜋) then add the polygon to the list of "positive" polygons, otherwise (it will be -2𝜋) add it to an alternative list.
    6. Keep repeating from step 2 until there are no more directed edges in the adjacency list.
    7. Finally you'll have two lists of polygons. One list will only have one polygon. This will be the original, outer polygon, and can be ignored. The other list will have the partitioning.

    As a demo, here is some code in a runnable JavaScript snippet. It uses one of the examples you pictured in your question (but will consecutive vertex numbering), finds the partitioning according to this algorithm, and shows the result by coloring the polygons that were identified:

    function partition(nodes, edges) {
        // Create an adjacency list
        let adj = [];
        for (let i = 0; i < nodes.length; i++) {
            adj[i] = []; // initialise the list for each node as an empty one
        }
        for (let i = 0; i < edges.length; i++) {
            let a = edges[i][0]; // Get the two nodes (a, b) that this edge connects
            let b = edges[i][1]; 
            adj[a].push(b); // Add as directed edge in both directions
            adj[b].push(a);
        }
        // Traverse the graph to identify polygons, until none are to be found
        let polygons = [[], []]; // two lists of polygons, one per "winding" (clockwise or ccw)
        let more = true;
        while (more) {
            more = false;
            for (let i = 0; i < nodes.length; i++) {
                if (adj[i].length) { // we have unvisited directed edge(s) here
                    let start = i;
                    let polygon = [i]; // collect the vertices on a new polygon
                    let sumAngle = 0;
                    // Take one neighbor out of this node's neighbor list and follow a path
                    for (let j = adj[i].pop(), next; j !== start; i = j, j = next) {
                        polygon.push(j);
                        // Get coordinates of the current edge's end-points
                        let ix = nodes[i][0];
                        let iy = nodes[i][1];
                        let jx = nodes[j][0];
                        let jy = nodes[j][1];
                        let startAngle = Math.atan2(jy-iy, jx-ix);
                        // In the adjacency list of node j, find the next neighboring vertex in counterclockwise order
                        //   relative to node i where we came from.
                        let minAngle = 10; // Larger than any normalised angle
                        for (let neighborIndex = 0; neighborIndex < adj[j].length; neighborIndex++) {
                            let k = adj[j][neighborIndex];
                            if (k === i) continue; // ignore the reverse of the edge we came from
                            let kx = nodes[k][0];
                            let ky = nodes[k][1];
                            let relAngle = Math.atan2(ky-jy, kx-jx) - startAngle; // The "magic"
                            // Normalise the relative angle to the range [-PI, +PI)
                            if (relAngle < -Math.PI) relAngle += 2*Math.PI;
                            if (relAngle >=  Math.PI) relAngle -= 2*Math.PI;
                            if (relAngle < minAngle) { // this one comes earlier in counterclockwise order
                                minAngle = relAngle;
                                nextNeighborIndex = neighborIndex;
                            }
                        }
                        sumAngle += minAngle; // track the sum of all the angles in the polygon
                        next = adj[j][nextNeighborIndex];
                        // delete the chosen directed edge (so it cannot be visited again)
                        adj[j].splice(nextNeighborIndex, 1);
                    }
                    let winding = sumAngle > 0 ? 1 : 0; // sumAngle will be 2*PI or -2*PI. Clockwise or ccw.
                    polygons[winding].push(polygon);
                    more = true;
                }
            }
        }
        // return the largest list of polygons, so to exclude the whole polygon,
        //   which will be the only one with a winding that's different from all the others.
        return polygons[0].length > polygons[1].length ? polygons[0] : polygons[1];
    }
    
    // Sample input:
    let nodes = [[59,25],[26,27],[9,59],[3,99],[30,114],[77,116],[89,102],[102,136],[105,154],[146,157],[181,151],[201,125],[194,83],[155,72],[174,47],[182,24],[153,6],[117,2],[89,9],[97,45]];
    let internalEdges = [[6, 13], [13, 19], [19, 6]];
    // Join outer edges with inner edges to an overall list of edges:
    let edges = nodes.map((a, i) => [i, (i+1)%nodes.length]).concat(internalEdges);
    // Apply algorithm
    let polygons = partition(nodes, edges);
    // Report on results
    document.querySelector("div").innerHTML =
        "input polygon has these points, numbered 0..n<br>" + 
        JSON.stringify(nodes) + "<br>" +
        "resulting polygons, by vertex numbers<br>" +
        JSON.stringify(polygons)
    
    // Graphics handling
    let io = {
        ctx: document.querySelector("canvas").getContext("2d"),
        drawEdges(edges) {
            for (let [a, b] of edges) {
                this.ctx.moveTo(...a);
                this.ctx.lineTo(...b);
                this.ctx.stroke();
            }
        },
        colorPolygon(polygon, color) {
            this.ctx.beginPath();
            this.ctx.moveTo(...polygon[0]);
            for (let p of polygon.slice(1)) {
                this.ctx.lineTo(...p);
            }
            this.ctx.closePath();
            this.ctx.fillStyle = color;
            this.ctx.fill();
        }
    };
    
    // Display original graph
    io.drawEdges(edges.map(([a,b]) => [nodes[a], nodes[b]]));
    // Color the polygons that the algorithm identified
    let colors = ["red", "blue", "silver", "purple", "green", "brown", "orange", "cyan"];
    for (let polygon of polygons) {
        io.colorPolygon(polygon.map(i => nodes[i]), colors.pop());
    }
    <canvas width="400" height="180"></canvas>
    <div></div>