How to check if a coordinate lies within a polygon of coordinate points?

How to check if a coordinate lies within a polygon of coordinate points?

Replies

There are some maths involved here.


For convex polygons (as a triangle), I would do this:


First some drawing, starting with a triangle:


A

/ \

/ \

/ P \

B __________ C


P will be inside if it is on the same side of BC than A, on the same side of AB than C, on the same side of AC than B.


How to knwo that a point in on a side ?


BC equation (points on the line) is of the form

(1) a x + b y = c

Hence, if (xA, yA) (xB, yB) (xC, yC) are coordinates of A, B, C

we can compute a, b, c with

(2) a xB + b yB = c

(3) a xC + b yC = c


hence, by subtracting the 2

(4) a (xB - xC) + b (yB - yC) = 0


Let's first consider xB != xC (non vertical), we get:

(5) a = b (yC - yB) / (xB - xC)


Porting into (2)

(2b) b ((yC - yB) / (xB - xC)) xB + b yB = c

which leads by mutilpying everything by (xB - xC) to new equation

(2t) b (yC - yB) xB + b yB (xB - xC) = c (xB - xC)

After simplification


(2q) b (yC xB - yB xC) = c (xB - xC)


If we now port (2q) and (5) into (1)

(1b) b ((yC - yB) / (xB - xC)) x + b y = b (yC xB - yB xC) / (xB - xC)

After multiplying all by (xB - xC)

(1t) b (yC - yB)) x + b y (xB - xC) = b (yC xB - yB xC)

Hence

(1q) (yC - yB) x + y (xB - xC) - (yC xB - yB xC) = 0

And in fact this equation works if xB = xC


With this equation (1q) in hand, we can evaluate for A if the result of equation is positive or negative.


A point inside triangle must be on the same side than A, hence have the same sign for the equation result.


Do the same with AB and BC and you get the 3 equations similar to (1q) to evaluate to know if point is inside.



For a convex polygon, it is exactly the same, but you have as many equations to consider as the number N of sides (or summits) ; and for each side, you have N-2 other summits to test, so for each side you will have N-2 equations similar to (1q) to evaluate.



For non convex, I need to think a bit more about it.


You may read this as well with some interest:

h ttp://www.euclideanspace.com/maths/geometry/elements/intersection/twod/index.htm



EDITED

There maybe another method, less analytical, based on ray tracing (I don't have the mathematical confirmation right now, but that seem correct).


If you draw a line starting from a point P inside a polygon (convex or not, but with sides that do not cross each other), it will intersect an odd number of sides of the polygon (1 if convex, 1 or 3 or 5 … if not convex).

If the point P is outside, it will intersect an even number of sides of the polygon.

Then, for a given point P, you could trace rays starting from P every 1 degre, compute the equation of this half line, compute how many segments it intersects. If all numbers are odd, you're inside.


But, at the end, that may be as difficult to compute as the first solution.


Might be some interesting discussion here:


https://stackoverflow.com/questions/217578/how-can-i-determine-whether-a-2d-point-is-within-a-polygon


...as well as you helping to qualify the usual 'depends', beyond what may risk being an over simplified question prompting the usual over complicated replies 😉, until those qualifications are posed and satisfied.

This is not a wheel you want to re-invent. There are a number of geospatial libraries that will perform this type of calculation for you. Unfortunately, they are usually written in languages that you will find extremely difficult to integrate into Swift. Plus, they often carry restrictive open-source licenses.


I suggest clearly identifying the operations you need to perform and then writting a small Objective-C++ wrapper around the Boost Geometry module. It doesn't do all operations, but most of the basic ones. It won't be easy, but I know of no easier method.

If iamsuman has just to deal with simple shapes, it is not a big deal writing the code.


If shapes are triangles, rectagles, lozanges and such simple shapes, computing the regions with equations with

(yC - yB) x + y (xB - xC) - (yC xB - yB xC)

is pretty easy and fast.


Here is a simple code to do it for a triangle ABC, with pt1 inside, pt2 outside


let ptA = CGPoint(x: 5, y: 10)
let ptB = CGPoint(x: 1, y: 2)
let ptC = CGPoint(x: 8, y: 1)
let pt1 = CGPoint(x: 6, y: 4)
let pt2 = CGPoint(x: 9, y: 4)

func positionRelativeToLine(point p: CGPoint, lineStart b: CGPoint, lineEnd c: CGPoint) -> Float {
    let xB = Float(b.x)     // needed to ease compiler work
    let yB = Float(b.y)
    let xC = Float(c.x)
    let yC = Float(c.y)
    let x = Float(p.x)
    let y = Float(p.y)
    return (yC - yB) * x + y * (xB - xC) - (yC * xB - yB * xC)
}

func pointsOnSameSide(first: CGPoint, second: CGPoint, lineStart b: CGPoint, lineEnd c: CGPoint) -> Bool {
   
    let test1 = positionRelativeToLine(point: first, lineStart: b, lineEnd: c)
    let test2 = positionRelativeToLine(point: second, lineStart: b, lineEnd: c)
    return test1 * test2 >= 0
}

func testPointInsideTriangle(point p: CGPoint, summitA: CGPoint, summitB: CGPoint, summitC: CGPoint) -> Bool {
    let testPBC = pointsOnSameSide(first: p, second: summitA, lineStart: summitB, lineEnd: summitC)
    let testPAC = pointsOnSameSide(first: p, second: summitB, lineStart: summitA, lineEnd: summitC)
    let testPAB = pointsOnSameSide(first: p, second: summitC, lineStart: summitA, lineEnd: summitB)
    return testPBC && testPAC && testPAB
}

if testPointInsideTriangle(point: pt1, summitA: ptA, summitB: ptB, summitC: ptC) {
    print("Point 1 is inside")
} else {
    print("Point 1 is outside")
}

if testPointInsideTriangle(point: pt2, summitA: ptA, summitB: ptB, summitC: ptC) {
    print("Point 2 is inside")
} else {
    print("Point 2 is outside")
}


result:

Point 1 is inside

Point 2 is outside


For convexQuad, add the following

let ptD = CGPoint(x: 9, y: 8)

func isQuadConvex(summitA: CGPoint, summitB: CGPoint, summitC: CGPoint, summitD: CGPoint) -> Bool {
   
    let testAB = pointsOnSameSide(first: summitC, second: summitD, lineStart: summitA, lineEnd: summitB)
    let testBC = pointsOnSameSide(first: summitD, second: summitA, lineStart: summitB, lineEnd: summitC)
    let testCD = pointsOnSameSide(first: summitA, second: summitB, lineStart: summitC, lineEnd: summitD)
    let testDA = pointsOnSameSide(first: summitB, second: summitC, lineStart: summitD, lineEnd: summitA)
    return testAB && testBC && testCD && testDA
}

func testPointInsideQuad(point p: CGPoint, summitA: CGPoint, summitB: CGPoint, summitC: CGPoint, summitD: CGPoint) -> Bool {
    let testPAB = pointsOnSameSide(first: p, second: summitC, lineStart: summitA, lineEnd: summitB) && pointsOnSameSide(first: p, second: summitD, lineStart: summitA, lineEnd: summitB)
    let testPBC = pointsOnSameSide(first: p, second: summitA, lineStart: summitB, lineEnd: summitC) && pointsOnSameSide(first: p, second: summitD, lineStart: summitB, lineEnd: summitC)
    let testPCD = pointsOnSameSide(first: p, second: summitB, lineStart: summitC, lineEnd: summitD) && pointsOnSameSide(first: p, second: summitA, lineStart: summitC, lineEnd: summitD)
    let testPDA = pointsOnSameSide(first: p, second: summitB, lineStart: summitA, lineEnd: summitD) && pointsOnSameSide(first: p, second: summitC, lineStart: summitA, lineEnd: summitD)
    return testPAB && testPBC && testPCD && testPDA
}

if isQuadConvex(summitA: ptA, summitB: ptB, summitC: ptC, summitD: ptD) {
    print("Quad is convex")
} else {
    print("Quad is not convex ; use other method")
}

if testPointInsideQuad(point: pt1, summitA: ptA, summitB: ptB, summitC: ptC, summitD: ptD) {
    print("Point 1 is inside Quad")
} else {
    print("Point 1 is outside Quad")
}

if testPointInsideQuad(point: pt2, summitA: ptA, summitB: ptB, summitC: ptC, summitD: ptD) {
    print("Point 2 is inside Quad")
} else {
    print("Point 2 is outside Quad")
}

The second method by tracing rays is more generaland will work also for circular shapes; but need to take care to precision to count how many times a segment crosses shape boundary.

Of course the easy stuff is easy. But how can you even determine that a given polygon is an easy one?

I cannot determine, of course. It all depends what type of use of such detection iamsuman is planning to do.


If he/she plans to do it on triangles or quads, that may be just enough.

If it is for random Mandelbrot curves, that may be another trip !

I'm afraid there are only two realistic cases here. '


1) The polygon is a perfect square, from a mouse selection or viewport;

2) The polygon comes from more arbitrary user input or external input.


If #1 is true and that is all it will ever do, then you don't even have to worry about triangles. Just do the square in a few lines of code.

But if #2 is true, it is then not reasonable to expect the input to stay within the range of simple polygons. It is going to be either self-intersecting or a high-resolution outline of the borders of the Russian Federation.

I don't understand your point.


Any triangle, any convex quad will be handled properly with the test.

And for more than 4 sides, if convex (which is tested), it will work as well by just completing the testPointInPolygon.

It is by no way limited to perfect square, unless I missed some flaw you saw in code.


So, it will work in a lot of general cases (just have to write the testPointIn with arbitrary number of sides, which is reasonnably easy).


For the fun, here is the generic case. As a final point to our interesting discussion, because you are right, it is a question of usage and what type of more complex polynoms may be expected in the app.

extension CGPoint: Hashable {
    public func hash(into hasher: inout Hasher) {
        // Ignore distanceFromOrigin for hashing
        hasher.combine(x)
        hasher.combine(y)
    }
}

typealias Polygon = [CGPoint]

func isConvex(_ p: Polygon) -> Bool {
   
    let n = p.count
    if n <= 2 { return true }
    let pSet = Set(p)       // Possible because Hashable extension
    if pSet.count < n { print("2 summits are the same")}

    var q = p
    q.append(p[0])  // Add first point to get A, B, C, D, A and be able to test AB, BC; CD, DA

    for s in 0..        for other1 in q where other1 != q[s] && other1 != q[s+1] {
            for other2 in q where other2 != other1 && other2 != q[s] && other2 != q[s+1] {
                let test = pointsOnSameSide(first: other1, second: other2, lineStart: q[s], lineEnd: q[s+1])
                if !test { return false }
            }
        }
    }
   
    return true
}

func testPointInside(point pt: CGPoint, polygon p: Polygon) -> Bool? {

    let n = p.count
    if n == 1 { return pt == p[0] }
    if n == 2 {
        let ptA = p[0]
        let ptB = p[1]
        if ptA == ptB { return pt == ptA }
        let x = pt.x ; let y = pt.y
        let xA = ptA.x ; let yA = ptA.y
        let xB = ptB.x ; let yB = ptB.y
        var alpha, beta: CGFloat
        if xB != xA {
            alpha = (x - xA) / (xB - xA)
            beta = (y - yA) / (yB - yA)
            return alpha == beta && alpha >= 0 && alpha <= 1.0
        } else {
            beta = (y - yA) / (yB - yA)
            return x == xA && beta >= 0 && beta <= 1.0
        }
    }
    // Now n >= 3
   
    if !isConvex(p) { return nil }
    var q = p
    q.append(p[0])  // Add first point to get A, B, C, D, A and be able to test AB, BC; CD, DA

    for s in 0..        // Let us test if point pt is on the same side of each side of p as the other summits of p
        for other in q where other != q[s] && other != q[s+1] {
            let test = pointsOnSameSide(first: pt, second: other, lineStart: q[s], lineEnd: q[s+1])
            if !test { return false }
        }
    }
   
    return true
}

let poly = [ptA, ptB, ptC, ptD]
print(poly, "convex ?", isConvex(poly))
let triangle = [ptA, ptB, ptC]
print(triangle, "convex ?", isConvex(triangle))

if let test = testPointInside(point: pt1, polygon: poly) {
    if test {
        print(pt1, "is inside", poly)
    } else {
        print(pt1, "is outside", poly)
    }
} else {
    print("Cannot tell for a non convex polygon")
}

if let test = testPointInside(point: pt2, polygon: poly) {
    if test {
        print(pt2, "is inside", poly)
    } else {
        print(pt2, "is outside", poly)
    }
} else {
    print("Cannot tell for a non convex polygon")
}


[(5.0, 10.0), (1.0, 2.0), (8.0, 1.0), (9.0, 8.0)] convex ? true

[(5.0, 10.0), (1.0, 2.0), (8.0, 1.0)] convex ? true


(6.0, 4.0) is inside [(5.0, 10.0), (1.0, 2.0), (8.0, 1.0), (9.0, 8.0)]

(9.0, 4.0) is outside [(5.0, 10.0), (1.0, 2.0), (8.0, 1.0), (9.0, 8.0)]

It really isn't worth arguing about. But it isn't a coding question. It is a usage question. Simple polygons simply don't exist in the real world. If you have a really simple UI, then that could generate squares. But it won't take long after releasing such a simple UI for someone to ask for a more complicated UI where they can draw their own points. Then you are dealing with self-intersecting polygons, polygons with holes, multi-polygons, and Kenya.


It is better to know ahead of time what the challenges might be down the road. That way, you can limit your input to the data that you can reasonably support. Or, you can integrate a geometry library that can support arbitrary polygons. Integrating a geometry library in Swift isn't going to be easy, but it will be a lot easier than trying to incrementally code your own implementation of the Dimensionally Extended nine-Intersection Model.

make a bezier path (NSBezierPath)

ask it:

https://developer.apple.com/documentation/appkit/nsbezierpath/1520716-contains


another way is to actually make an image, fill the polygon with a color, and then sample the point... if it's that color... then your point is in the polygon.


CGPath also has the ability to check if a point is inside it:

https://developer.apple.com/documentation/coregraphics/cgpath/2427117-contains