Splitting polygons with split lines?

[Thomas]

Hi,

Just started using this library, and so far it looks awesome! One thing that I
haven't figured out yet, is how to split polygons with a line, or polyline. Is
that even possible? It crossed my mind to make a closed loop of the split
line, big enough to guarantee inclusion of half the polygon that needs
splitting, but it seemed overly complicated. Is there an easier way?

Thanks

[Angus Johnson]

It's not possible to use lines or polylines to clip polygons. However, as a
workaround, you could convert a line/polyline into a polygon before clipping.

[Lord Helmchen]

This would indeed be a nice feature useful for GIS applications. It should be possible to do this after the intersection calculation by walking the closed polygons clockwise starting and ending at an intersection with a line. If you come a cross another intersection with a line, follow it right/left (depending on polygon orientation) and mark it as visited. Then back at the start, continue with the next unvisited intersection between a polygon and a line.

If this shall not be a part of clipper, then it should be possible to get the nescessary intersection information from the ZFillFuntion. Exposing the Edges and Path information directly would be a bit easier though saving reimplementation.

[Lord Helmchen]

Here is some code that works for my simple test cases. It can split a polygon with a line. It should work with more than one line, but I did not test it. It is not efficient and could be optimized easily with more access to Clippers internals, but I just needed a quick hack to split Google Earth KML-Polygons. (The program will be open source in https://sourceforge.net/p/dereglobus/ as soon as it is somewhat usable.)
Updated it with better working code.

using ClipperLib;
using System;
using System.Collections.Generic;
using System.Linq;
using System.Windows;
using System.Windows.Media;

namespace KmlEdit
{
    public class Clipping
    {
        public IList<KmlShape> Union(IEnumerable<KmlShape> list)
        {
            Clipper clipper = new Clipper();
            IList<KmlShape> ret = new List<KmlShape>();
            foreach (KmlShape ks in list)
            {
                var path = ks.Points.Select(p => new IntPoint(p.X * 10000, p.Y * 10000)).ToList();
                clipper.AddPath(path, (ks.ShapeType == ShapeType.Polyline) ? PolyType.ptSubject : PolyType.ptClip, (ks.ShapeType == ShapeType.Polyline) ? false : true);
            }
            PolyTree result = new PolyTree();
            clipper.Execute(ClipType.ctUnion, result, PolyFillType.pftNonZero, PolyFillType.pftNonZero);
            PolyNode n = result;
            while (n != null)
            {
                var path = n.Contour.Select(p => new Point(p.X / 10000.0, p.Y / 10000.0));
                var ks = new KmlShape(new PointCollection(path), n.IsOpen ? ShapeType.Polyline : ShapeType.Polygon, "Unionresult");
                ret.Add(ks);
                n = n.GetNext();
            }
            return ret;
        }

        List<bool> isclosed;
        List<List<IntPoint>> paths;

        class Intersection
        {
            public IntPoint bot1;
            public IntPoint top1;
            public IntPoint bot2;
            public IntPoint top2;
            public IntPoint pt;
            //public int visited = 0;

            public IntPoint[] GetEdgeWithZ(long pathZ)
            {
                if (top1.Z == pathZ)
                    return new IntPoint[] { bot1, top1 };
                if (top2.Z == pathZ)
                    return new IntPoint[] { bot2, top2 };
                return null;
            }
            public IntPoint[] GetEdgeWithOtherZ(long pathZ)
            {
                if (top1.Z != pathZ)
                    return new IntPoint[] { bot1, top1 };
                if (top2.Z != pathZ)
                    return new IntPoint[] { bot2, top2 };
                return null;
            }
        }

        private List<Intersection> intersections = new List<Intersection>();

        public void ZFillCallback(IntPoint bot1, IntPoint top1,
        IntPoint bot2, IntPoint top2, ref IntPoint pt)
        {
            intersections.Add(new Intersection() { bot1 = bot1, top1 = top1, bot2 = bot2, top2 = top2, pt = pt });
        }

        public IList<KmlShape> SplitPolygonsWithLines(IEnumerable<KmlShape> list)
        {
            Clipper clipper = new Clipper();
            clipper.ZFillFunction = ZFillCallback;
            paths = new List<List<IntPoint>>();
            isclosed = new List<bool>();
            List<bool> isclockwise = new List<bool>();
            int i = 0;
            foreach (KmlShape ks in list)
            {
                i++;
                var path = ks.Points.Select(p => new IntPoint(p.X * 10000, p.Y * 10000, i)).ToList();
                paths.Add(path);
                isclosed.Add(ks.ShapeType == ShapeType.Polygon);
                if (!isclosed[i - 1])
                    isclockwise.Add(false);
                else
                    isclockwise.Add(IsClockwise(path));
                clipper.AddPath(path, (ks.ShapeType == ShapeType.Polyline) ? PolyType.ptSubject : PolyType.ptClip, (ks.ShapeType == ShapeType.Polyline) ? false : true);
            }
            PolyTree result = new PolyTree();
            intersections.Clear();
            clipper.Execute(ClipType.ctIntersection, result, PolyFillType.pftNonZero, PolyFillType.pftNonZero);
            List<List<IntPoint>> results = new List<List<IntPoint>>();
            for (int direction = 1; direction <= 2; direction++)
                foreach (var p in intersections)
                {
                    // für jeden Kreuzungspunkt das zugehörige geschlossene polygon umlaufen, an jeder intersection abbiegen, bis wir wieder hier angekommen sind.
                    var pathIndex1 = (int)p.top1.Z - 1;
                    var pathIndex2 = (int)p.top2.Z - 1;
                    long pathZ = 0;
                    if (isclosed[pathIndex1])
                    {
                        pathZ = p.top1.Z;
                    }
                    else if (isclosed[pathIndex2])
                    {
                        pathZ = p.top2.Z;
                    }
                    else
                        //skip when neither is a poylgon
                        continue;
                    //skip when neither is a line
                    if (isclosed[pathIndex2] && isclosed[pathIndex1])
                        continue;
                    //starting and end point
                    IntPoint start = p.pt;

                    var closedEdge = p.GetEdgeWithZ(pathZ);
                    if (closedEdge == null)
                        continue;
                    IntPoint current = closedEdge[0];
                    var path = paths[(int)pathZ - 1];
                    int ci = path.IndexOf(current), ni = path.IndexOf(closedEdge[1]);
                    if ((ci > ni || (ci == path.Count-1 && ni == 0)) == (direction == 1))
                        current = closedEdge[1];

                    var nextPathZ = p.GetEdgeWithOtherZ(pathZ)[0].Z;

                    var r = WalkPath(p, current, nextPathZ, direction, start, isclockwise[(int)pathZ - 1] & (direction == 1));
                    if (r != null)
                    {
                        r.Insert(0, start);
                        results.Add(r);
                        continue;
                    }
                }

            IList<KmlShape> ret = new List<KmlShape>();
            i = 0;
            foreach (var c in results)
            {
                i++;
                KmlShape ks = new KmlShape(
                    new PointCollection(c.Select(p => new Point(p.X / 10000.0, p.Y / 10000.0))),
                    ShapeType.Polygon,
                    String.Format("Result {0}", i));
                ret.Add(ks);
            }
            return ret;
        }

        List<IntPoint> WalkPath(Intersection intersection, IntPoint lastPoint, long pathZ, int direction, IntPoint startEndPoint, bool walkclockwise)
        {
            var currentEdge = intersection.GetEdgeWithZ(pathZ);
            var prevEdge = intersection.GetEdgeWithOtherZ(pathZ);
            //use last point to get the edges orientation
            if (prevEdge[0] != lastPoint)
            {
                prevEdge[1] = prevEdge[0];
                prevEdge[0] = lastPoint;
            }
            var path = paths[(int)pathZ - 1];

            //product between edges + isclockwise + direction to determine which way is correct
            IntPoint current, next;
            int incr;

            var rad = AngleBetween(prevEdge, currentEdge);
            if ((rad < 0 == walkclockwise)) // rad < 0 : path goes left, rad > 0 path goes right (follow it right, when clockwise, else go left)
            {
                current = currentEdge[0]; //dont' switch
                next = currentEdge[1];
            }
            else
            {
                current = currentEdge[1]; //switch
                next = currentEdge[0];
            }
            var ci = path.IndexOf(current);
            {
                var ni = path.IndexOf(next);
                if (ci == path.Count - 1 && ni == 0)
                    incr = 1;
                else if (ni == path.Count - 1 && ci == 0)
                    incr = -1;
                else
                    incr = ni - ci;
                //discard ni from here on and just use ci as counter
                ci = ni;
            }

            //assume the intersection as lastPoint, since we want to know the distance from there
            lastPoint = intersection.pt;

            List<IntPoint> result = new List<IntPoint>();
            //while not at path end (or start depending on direction) ...
            while (true)
            {
                var intersects = GetIntersections(lastPoint, current, next);
                foreach (var ei in intersects)
                {
                    //finished, we are back at the start
                    if (ei.pt == startEndPoint)
                        return result;
                    if (ei == intersection)
                        continue;
                    //just walk those paths which are either lines or from the original polygon
                    if ((startEndPoint.Z != ei.top1.Z) && (startEndPoint.Z != ei.top2.Z) && isclosed[(int)ei.top1.Z - 1] && isclosed[(int)ei.top2.Z - 1])
                        continue;
                    //walk intersection until we know wether it leads back to the polygon
                    var nextEdge = ei.GetEdgeWithOtherZ(pathZ);
                    var nextZ = nextEdge[0].Z;
                    var r = WalkPath(ei, current, nextZ, direction, startEndPoint, walkclockwise);
                    if (r == null)
                        continue;
                    else
                    {
                        result.Add(ei.pt);
                        result.AddRange(r);
                        return result;
                    }
                }
                //finished, we are back at the start
                if (next == startEndPoint)
                    return result;

                //continue with next vertex
                result.Add(next);
                current = next;
                lastPoint = current;

                //increase/decrease ci
                ci += incr;
                if (!isclosed[(int)pathZ - 1])
                {
                    //dead end
                    if (ci < 0 || ci >= path.Count)
                        return null;
                }
                else
                {
                    if (ci < 0)
                        ci = path.Count - 1;
                    if (ci >= path.Count)
                        ci = 0;
                }

                next = path[ci];
                //loop ...
            }
        }

        public double AngleBetween(IntPoint[] e1, IntPoint[] e2)
        {
            double rad = Math.Atan2(e1[1].Y - e1[0].Y, e1[1].X - e1[0].X) - Math.Atan2(e2[1].Y - e2[0].Y, e2[1].X - e2[0].X);
            if (rad > Math.PI)
                rad -= 2 * Math.PI;
            return rad;
        }

        public bool IsClockwise(List<IntPoint> points)
        {
            return SignedPolygonArea(points) < 0;
        }

        public float SignedPolygonArea(List<IntPoint> points)
        {
            // Add the first point to the end.
            int num_points = points.Count;
            IntPoint[] pts = new IntPoint[num_points + 1];
            points.CopyTo(pts, 0);
            pts[num_points] = points[0];

            // Get the areas.
            float area = 0;
            for (int i = 0; i < num_points; i++)
            {
                area +=
                    (pts[i].X * pts[i + 1].Y) -
                    (pts[i + 1].X * pts[i].Y) / 2;
            }

            // Return the result.
            return area;
        }

        IOrderedEnumerable<Intersection> GetIntersections(IntPoint intersection, IntPoint current, IntPoint next)
        {
            var l = intersections.Where(i => (i.top1 == current && i.bot1 == next || i.top1 == next && i.bot1 == current || i.top2 == current && i.bot2 == next || i.top2 == next && i.bot2 == current));
            var m = l.Where(i => GetDirectedDistance(intersection, next, i.pt) > 0).OrderBy(i => GetDirectedDistance(intersection, next, i.pt));
            return m;
        }

        public double GetDirectedDistance(IntPoint source, IntPoint direction, IntPoint point)
        {
            double dirx = direction.X - source.X;
            double diry = direction.Y - source.Y;
            var dy = (point.Y - source.Y) / diry;
            var dx = (point.X - source.X) / dirx;
            return (dx == 0.0) ? dy : dx;

        }

        public double GetSqrDistance(IntPoint p1, IntPoint p2)
        {
            return ((p1.X - p2.X) * (p1.X - p2.X) + (p1.Y - p2.Y) * (p1.Y - p2.Y));
        }
    }
}