maynull
maynull

Reputation: 2046

How can I determine if a point is inside a certain parallelogram in Python?

I'd like to know how to determine if two points d1(920, 52.1), d2(920, 52) are inside of this orange parallelogram, which is comprised of p1~p4.

enter image description here

As you can see, the base of the parallelogram is not parellel to x-axis.

Are there any modules that deal with this kind problems? Or I'd like to get some math help in this problem.

Upvotes: 1

Views: 3250

Answers (3)

Mads Boyd-Madsen
Mads Boyd-Madsen

Reputation: 167

The following solution use the quadrant walk algorithm optimised for parallelograms. It is a derivative of the winding-number algorithm.

That algorithm place the test-point P at the centre of a coordinate system and walks the vertices of the polygon in order, checking for "positive and negative axis-crossings". If the number of axis-crossing is a multiplum of 4 (but not 0), then P is inside the polygon. Intuitively P is inside the polygon, if walking the vertices "cross all four half-axis". Hence the multiplum of 4.

The algorithm can be highly optimised (especially well suited for assembler implementation) and run fast.

Here it is in Delphi, using the same parameter order as the function PtInPlgm given in the answer by MBo

function PtInPlgm2( bx, by, ax, ay, cx, cy, px, py: double ): boolean;
var
  Sum    : integer;
  dx, dy : double;

begin
  Sum := 0;
  dx := cx - bx + ax;
  dy := cy - by + ay;
  if ( px < ax ) <> ( px < bx ) then inc( Sum );
  if ( py < ay ) <> ( py < by ) then inc( Sum );
  if ( px < bx ) <> ( px < cx ) then inc( Sum );
  if ( py < by ) <> ( py < cy ) then inc( Sum );
  if ( px < cx ) <> ( px < dx ) then inc( Sum );
  if ( py < cy ) <> ( py < dy ) then inc( Sum );
  if ( px < dx ) <> ( px < ax ) then inc( Sum );
  if ( py < dy ) <> ( py < ay ) then inc( Sum );
  Result := (Sum mod 4 = 0) AND (Sum <> 0);
end;

I trust it's Ok, that I haven't provided a Python implementation

Upvotes: 1

MBo
MBo

Reputation: 80325

For parallelogram there is approach simpler than for general polygons. Get coordinates of three neighbor vertices in order c, a, b and represent AP vector in basis of AB and AC vectors - coefficients should be in range 0..1.

Note that Delphi and Python function use different argument order. In Delphi base point (a) goes first, then its neighbors (b, c), while Python list contains c,a,b (or b,a,c) points in order.

Delphi code:

function PtInPlgm(ax, ay, bx, by, cx, cy, px, py: Integer): Boolean;
var
  xb, yb, xc, yc, xp, yp, d: Integer;
  bb, cc, oned: Double;
begin
  Result := False;
  xb := bx - ax;
  yb := by - ay;
  xc := cx - ax;
  yc := cy - ay;
  xp := px - ax;
  yp := py - ay;
  d := xb * yc - yb * xc;
  if d <> 0 then begin
    oned := 1 / d;
    bb := (xp * yc - xc * yp) * oned;
    cc := (xb * yp - xp * yb) * oned;
    Result := (bb >= 0) and (cc >= 0) and (bb <= 1) and (cc <= 1);
  end;
end;

Literal Python translation:

def point_inside_prlgm(x,y,poly):
    inside = False
    xb = poly[0][0] - poly[1][0]
    yb = poly[0][1] - poly[1][1]
    xc = poly[2][0] - poly[1][0]
    yc = poly[2][1] - poly[1][1]
    xp = x - poly[1][0]
    yp = y - poly[1][1]
    d = xb * yc - yb * xc;
    if (d <> 0):
        oned = 1.0 / d;
        bb = (xp * yc - xc * yp) * oned
        cc = (xb * yp - xp * yb) * oned
        inside = (bb >= 0) & (cc >= 0) & (bb <= 1) & (cc <= 1)
    return inside

    print(point_inside_prlgm(1, 1, [[1, 2], [0, 0], [2, 0]]))
    print(point_inside_prlgm(-1, 1, [[2, 1], [0, 0], [2, 0]])) 

Upvotes: 2

maynull
maynull

Reputation: 2046

I tried using Shapely, but I read that it's realatively slow. My program needs to use the logic at least 100 times a minute, so I look the keywords up on the internet and found this page: http://www.ariel.com.au/a/python-point-int-poly.html

# determine if a point is inside a given polygon or not
# Polygon is a list of (x,y) pairs.

def point_inside_polygon(x,y,poly):

    n = len(poly)
    inside =False

    p1x,p1y = poly[0]
    for i in range(n+1):
        p2x,p2y = poly[i % n]
        if y > min(p1y,p2y):
            if y <= max(p1y,p2y):
                if x <= max(p1x,p2x):
                    if p1y != p2y:
                        xinters = (y-p1y)*(p2x-p1x)/(p2y-p1y)+p1x
                    if p1x == p2x or x <= xinters:
                        inside = not inside
        p1x,p1y = p2x,p2y

    return inside

It works fast enough and doesn't require additional modules. I think coding this code in C and using it in Python could be a good idea.

Upvotes: 0

Related Questions