Ray-casting algorithm: Difference between revisions

{{task|Classic CS problems and programs}}

{{Wikipedia|Point_in_polygon}}

Given a point and a polygon, check if the point is inside or outside the polygon using the [[wp:Point in polygon#Ray casting algorithm|ray-casting algorithm]].

 

An intuitive explanation of why it works is that every time we cross

a border, we change "country" (inside-outside, or outside-inside), but

 

So the main part of the algorithm is how we determine if a ray intersects a segment. The following text explain one of the possible ways.

[[Image:posslope.png|128px|thumb|right]]

[[Image:negslope.png|128px|thumb|right]]

Let us take into account a segment AB (the point A having y coordinate always smaller than B's y coordinate, i.e. point A is always below point B) and a point P. Let us use the cumbersome notation PAX to denote the angle between segment AP and AX, where X is always a point on the horizontal line passing by A with x coordinate bigger than the maximum between the x coordinate of A and the x coordinate of B. As explained graphically by the figures on the right, if PAX is greater than the angle BAX, then the ray starting from P intersects the segment AB. (In the images, the ray starting from P_A does not intersect the segment, while the ray starting from P_B in the second picture, intersects the segment).

 

'''if''' Py < Ay '''or''' Py > By '''then'''

'''return''' false

'''return''' false

'''else'''

'''end''' '''if'''

 

 

=={{header|AutoHotkey}}==

}

 

 

 

#include <stdlib.h>

#include <math.h>

 

 

 

 

 

{

 

 

 

 

{

 

 

{

 

 

{

 

 

 

int main()

{

 

 

 

 

 

(do ((in-p nil)) ((endp polygon) in-p)

(when (ray-intersects-segment point (pop polygon))

((null m-blue) t)

((null m-red) nil)

 

Testing code

#((0 . 0) (10 . 0) (10 . 10) (0 . 10)

(2.5 . 2.5) (7.5 . 2.5) (7.5 . 7.5) (2.5 . 7.5)

do (format t "~&~w ~:[outside~;inside ~]."

test-point

 

=={{header|Fortran}}==

{{works with|Fortran|95 and later}}

The following code uses the <tt>Points_Module</tt> defined [[Closest pair problem#Fortran|here]].

 

This module ''defines'' "polygons".

use Points_Module

implicit none

end subroutine free_polygon

 

 

The ray casting algorithm module:

use Polygons

implicit none

end function point_is_inside

 

 

'''Testing'''

use Points_Module

use Ray_Casting_Algo

polys(4) = create_polygon(pts, (/ 11,12, 12,10, 10,13, 13,14, 14,9, 9,11 /) )

 

 

p = (/ point(5,5), point(5, 8), point(-10, 5), point(0,5), point(10,5), &

end do

 

=={{header|Haskell}}==

 

type Point = (Rational, Rational)

(py /= by || ay < py)

((ax, ay), (bx, by)) = side

 

=={{header|J}}==

NB. bool=. points crossPnP polygon

crossPnP=: 4 : 0"2

p2=. (x0-/X) < (x0-x1) * (y0-/Y) % (y0 - y1)

2|+/ p1*.p2

 

Sample data:

SQUAREV=: SQUAREV, 2.5 2.5 , 7.5 0.1 , 7.5 7.5 ,: 2.5 7.5

 

STRANGE=: (0 4,4 3,3 7,7 6,6 2,2 1,1 5,:5 0) , .{ SQUAREV

 

 

Testing:

1 1 1 0 0 1 1 1 0 1 0

1 1 1 0 0 1 1 1 0 1 0

0 1 1 0 0 1 1 1 0 1 0

 

=={{header|OCaml}}==

{{Trans|C}}

type polygon = {

(s.y +. _eps) < (min a.y b.y) ||

s.x > (max a.x b.x) then (false)

else

let ca = if a.x <> b.x then (coeff_ang a b) else (huge_val) in

let my = min_y ~a ~b in

let cp = if (s.x -. my.x) <> 0.0 then (coeff_ang my s) else (huge_val) in

;;

 

print_newline()

) test_points;

 

=={{header|Perl}}==

use List::Util qw(max min);

 

 

return ($m_blue >= $m_red) ? 1 : 0;

 

Testing:

sub point

{

push @poly, [ $pts->[$sides->[$i]-1], $pts->[$sides->[$i+1]-1] ];

}

}

 

point(3,0), point(7,0), point(7,10), point(3,10) );

 

 

my @p = ( point(5,5), point(5, 8), point(-10, 5), point(0,5), point(10,5), &

point(8,5), point(10,10) );

 

no strict 'refs';

print "$pol\n";

foreach my $tp ( @p ) {

print "\t($tp->[0],$tp->[1]) " .

}

 

=={{header|Python}}==

from pprint import pprint as pp

import sys

 

def polypp(poly):

 

if __name__ == '__main__':

Pt(x=10, y=10))

for poly in polys:

polypp(poly)

 

'''Sample output'''

 

'''Helper routine to convert Fortran Polygons and points to Python'''

point = Pt

pts = (point(0,0), point(10,0), point(10,10), point(0,10),

print ' ', ',\n '.join(str(e) for e in p.edges) + '\n )),'

print ' ]'

 

=={{header|R}}==

count <- 0

for(side in polygon) {

}

}

 

 

point <- function(x,y) list(x=x, y=y)

}

pol

 

 

pts <- list(point(0,0), point(10,0), point(10,10), point(0,10),

p$x, p$y, point_in_polygon(polygons[[polysi]], p), names(polygons[polysi])))

}

 

=={{header|Smalltalk}}==

{{works with|GNU Smalltalk}}

The class Segment holds the code to test if a ray starting from a point (and going towards "right") intersects the segment (method <tt>doesIntersectRayFrom</tt>); while the class Polygon hosts the code to test if a point is inside the polygon (method <tt>pointInside</tt>).

|pts|

Segment class >> new: points [ |a|

^ ( cnt \\ 2 = 0 ) not

]

 

'''Testing'''

 

points := {

].

' ' displayNl.

 

=={{header|Tcl}}==

proc point_in_polygon {point polygon} {

} {

puts "$point in $poly = [point_in_polygon $point $poly]"

 

=={{header|Ursala}}==

value if it's outside, using the algorithm described above.

For points on the boundary the result is unspecified.

 

in =

@lrzyCipPX ~|afatPRZaq ~&EZ+fleq~~lrPrbr2G&& ~&B+fleq~~lrPrbl2G!| -&

~&Y+ ~~lrPrbl2G fleq,

This test program tries it on a couple of examples.

 

examples =

in* <

((0.5,0.6),<(0.,0.),(1.,2.),(1.,0.)>),

output:

<pre><true,false></pre>