CLIM mail archive


[Jeff Morrill: Re: point-in-polygon ]

OK, here's my contribution to the point-in-polygon algorithm.  The issue
comes up every decade or so, as you can see from the references.  I had
lost track of my copy, but JMorrill was able to find it.  Thanks, Jeff.

------- Forwarded Message

To: Ken Anderson <KAnderson@BBN.COM>
Subject: Re: point-in-polygon 
In-reply-to: Ken Anderson's message of Thu, 13 Jan 94 14:58:24 -0500.
Reply-To: jmorrill@BBN.COM
Date: Thu, 13 Jan 94 15:45:41 -0500
Message-ID: <>
From: Jeff Morrill <jmorrill@BBN.COM>

  To: JMorrill@BBN.COM
  Cc: Kanderson@BBN.COM,,, 
  From: Ken Anderson <KAnderson@BBN.COM>
  Subject: point-in-polygon
  Date: Thu, 13 Jan 94 14:58:24 -0500
  Sender: kanderso@BBN.COM
  I know i wrote fortran code a long time ago to decide if a point is in a
  polygon, and i thought i had this in lisp, do you guys recall seening
  anything like this from me?

This is probably what you were looking for.  It was in scigraph.


(defun point-in-polygon (x y map-polygon polygon)
  "Returns T if point X,Y is inside a closed boundary.  Points on the 
  boundary edge are ambiguously either inside or outside.  The
  algorithm counts the number of line segements that intersect the
  semi-infinite line from (x,y) to (+infinity,y).  If there is an odd 
  number of such segements, the point is inside the polygon.  
  The function MAP-POLYGON maps a function over consecutive line 
  segments (X1 Y1 X2 Y2) of the polygon, such that the
  first point is the same as the last.  (The polygon need not be convex.)


  Shimrat, M., 1962, Algorithm 112: position of point relative to
  polygon, CACM,5,434.

  Hall, J.H., 1975, PTLOC - A FORTRAN subroutine for determining the
  position of a point relative to a closed boundary, Math. Geol., 7,

  Anderson, K.R., 1975, Letter to the editor of Math. Geol.

  The final version of this algorithm followed from discussions with
  David Fitterman of the US Geological Survey, Jan 1975."

  (let ((point-in NIL))
    (funcall map-polygon
	     #'(lambda (x1 y1 x2 y2)
		 (if (eq (<= y y1) (>  y y2))	; Segement crosses ray.
		     (when (and (not (= y1 y2))	; Ignore horizontal segement.
				(< (- x x1 (/ (* (- y y1) (- x2 x1))	
					      (- y2 y1)))	; Point is to left.
		       ;; (print (list x1 y1 x2 y2))
		       (setq point-in (not point-in)))))

------- End of Forwarded Message

Main Index | Thread Index