;;; This software is being provided to you, the LICENSEE, by the
;;; Massachusetts Institute of Technology (M.I.T.), and the NEC Research
;;; Institute (NECI) under the following license.  By obtaining, using
;;; and/or copying this software, you agree that you have read,
;;; understood, and will comply with these terms and conditions:
;;; 
;;; Permission to use, copy, modify and distribute, including the right to
;;; grant others rights to distribute at any tier, this software and its
;;; documentation for research purposes only and without fee or royalty is
;;; hereby granted, provided that you agree to comply with the following
;;; copyright notice and statements, including the disclaimer, and that
;;; the same appear on ALL copies of the software and documentation,
;;; including modifications that you make for internal use or for
;;; distribution:
;;; 
;;; Copyright 1995 by the Massachusetts Institute of Technology.  All
;;; rights reserved.  Copyright 1995 by NEC Research Institute Inc. All
;;; rights reserved
;;; 
;;; THIS SOFTWARE IS PROVIDED "AS IS", AND M.I.T. and NECI MAKE NO
;;; REPRESENTATIONS OR WARRANTIES, EXPRESS OR IMPLIED.  By way of example,
;;; but not limitation, M.I.T. and NECI MAKE NO REPRESENTATIONS OR
;;; WARRANTIES OF MERCHANTABILITY OR FITNESS FOR ANY PARTICULAR PURPOSE OR
;;; THAT THE USE OF THE LICENSED SOFTWARE OR DOCUMENTATION WILL NOT
;;; INFRINGE ANY THIRD PARTY PATENTS, COPYRIGHTS, TRADEMARKS OR OTHER
;;; RIGHTS.
;;; 
;;; The name of the Massachusetts Institute of Technology or M.I.T., or
;;; NEC Research Institute, or NECI may NOT be used in advertising or
;;; publicity pertaining to distribution of the software.  Title to
;;; copyright in this software and any associated documentation shall at
;;; all times remain with M.I.T. and NECI, and USER agrees to preserve
;;; same.
;;; 

;;; -*- Mode: LISP; Base: 10; Package: 2DG; Syntax: Common-Lisp -*-
(in-package '2DG)

(export '(inside? point-inside-rectangle? point-inside-polygon? point-inside-convex-polygon? polygon-diameter
	  intersect-rectangle-and-convex-polygon? line-segment-crosses-rectangle? intersect-line-segment-and-rectangle?))

(defun INSIDE? (object container)
  (case (geometric-type object)
    (point (case (geometric-type container)
	     (polygon (point-inside-polygon? object container))
	     (rectangle (point-inside-rectangle? object container))
	     (t (error "Don't know how to find out if a point is inside ~A" container))))
    (line-segment (case (geometric-type container)
		    (rectangle (and (inside? (line-segment-point-1 object) container)
				    (inside? (line-segment-point-1 object) container)))
		    (t (error "Don't know how to find out if a line segment is inside ~A" container))))
    (line nil)
    (polygon (case (geometric-type container)
	       (rectangle (loop for point in (polygon-points object)
				always
				(inside? point container)))
	       (t (error "Don't know how to find out if a polygon is inside ~A" container))))
    (rectangle (case (geometric-type container)
		 (rectangle (and (point-inside-rectangle? (first-point object) container)
				 (point-inside-rectangle? (second-point object) container)))
		 (t (error "Don't know how to find out if a rectangle is inside ~A" container))))))

(defun POINT-INSIDE-RECTANGLE? (pt rect)
  (let ((x (point-x pt))
	(y (point-y pt)))
    (and (~<= (rectangle-min-x rect) x)
	 (~<= x (rectangle-max-x rect))
	 (~<= (rectangle-min-y rect) y)
	 (~<= y (rectangle-max-y rect)))))

(defun POINT-INSIDE-POLYGON? (point polygon)
  "
  Purpose:  Determine if a point is inside of or on one of the edges of a simple polygon.
	Works for degenerate polygons as long as no point is repeated.
  Returns:  T if POINT is inside or on an edge of POLYGON.
  "
  (declare (values boolean))
  (let ((vertices (polygon-points polygon)))
    (cond ((= 2 (length vertices))
	   (intersect-point-and-line-segment? point (make-line-segment (car vertices) (cadr vertices))))
	  ((member point vertices :test #'point-equal)
	   t)					;POINT was a POLYGON vertex
	  (t
	   (multiple-value-bind (extreme-vertex extreme-position) (extreme-point vertices)
	     ;; Reorder the polygon from the vertex with min y-coordinate (and min x-coordinate if there is a tie)
	     (let ((reordered-vertices (append (nthcdr extreme-position vertices)
					       (butlast vertices (- (length vertices) extreme-position 1))))
		   (point-x (point-x point))
		   (point-y (point-y point)))
	       (when (~<= (point-y extreme-vertex) (point-y point))	;POLYGON doesn't sit totally above POINT
		   (loop for (x y) in (mapcar #'point-coordinates-list reordered-vertices)
			 with count = 0		
			 with last-x = (point-x extreme-vertex)
			 with last-y = (point-y extreme-vertex)
			 with horizontal-intersect
			 do
		     ;; Count the number of intersections of the polygon with a line drawn horizontally to the left from POINT.
		     ;; An odd number of intersections means we are inside.
		     (cond ((~= last-y point-y y)
			    ;; The current polygon side is horizontal and in line with POINT.  Does POINT lie on it?
			    (and (monotonic? last-x point-x x) t))
			   ((~= last-x point-x x)
			    ;; The current polygon side is vertical and in line with POINT.  Does POINT lie on it?
			    (and (monotonic? last-y point-y y) t))
			   ((~= y point-y)
			    ;; POINT is horizontal with one endpoint of the current segment.  Count a crossing if it is the upper endpoint and it is to the left POINT.
			    (and (> y last-y) (< x point-x) (incf count)))
			   ((~= last-y point-y)
			    ;; POINT is horizontal with one endpoint of the current segment.  Count a crossing if it is the upper endpoint and it is to the left POINT.
			    (and (> last-y y) (< last-x point-x) (incf count)))
			   ((monotonic? last-y point-y y)
			    ;; We crossed the horizontal of POINT
			    (setf horizontal-intersect (horizontal-intersect point-y last-x last-y x y))
			    (cond ((~= horizontal-intersect point-x)
				   ;; POINT lies on this edge
				   (return t))
				  ((< horizontal-intersect point-x)
				   ;; The polygon edge is to the left of POINT
				   (incf count)))))
		     ;; Keep track of last vertex
		     (setf last-x x
			   last-y y)
			 finally
			   (return (oddp count))))))))))

(defun POINT-INSIDE-CONVEX-POLYGON? (point convex-polygon)
  "
  Purpose:  Determine if a point is inside of or on one of the edges of a convex polygon.
	Works for degenerate polygons as long as no point is repeated.
  Doesn't actually check that polygon is convex.
  Returns:  T if POINT is inside or on an edge of CONVEX-POLYGON.
  "
  (declare (values boolean))
  (let ((vertices (polygon-points convex-polygon)))
    (cond ((= 2 (length convex-polygon))
	   (intersect-point-and-line-segment? point (make-line-segment (car vertices) (cadr vertices))))
	  ((member point convex-polygon :test #'point-equal)
	   t)					;POINT was a CONVEX-POLYGON vertex
	  (t
	   (multiple-value-bind (extreme-vertex extreme-position) (extreme-point vertices)
	     ;; Reorder the convex-polygon from the vertex with min y-coordinate (and min x-coordinate if there is a tie)
	     (let ((reordered-vertices (append (nthcdr extreme-position convex-polygon)
					       (butlast vertices (- (length vertices) extreme-position 1))))
		   (point-x (point-x point))
		   (point-y (point-y point)))
	       (when (~<= (point-y extreme-vertex) (point-y point))
		   nil				;CONVEX-POLYGON sits totally above POINT
		   (loop for (x y) in (mapcar #'point-coordinates-list reordered-vertices)
			 with count = 0
			 with last-x = (point-x extreme-vertex)
			 with last-y = (point-y extreme-vertex)
			 with horizontal-intersect
			 do
		     ;; Count the number of intersections of the polygon with a line drawn horizontally to the left from POINT.  A single intersection means we are inside.
		     (cond ((~= last-y point-y y)
			    ;; The current polygon side is horizontal and in line with POINT.  Does POINT lie on it?
			    (and (monotonic? last-x point-x x) (return t)))
			   ((~= last-x point-x x)
			    ;; The current polygon side is vertical and in line with POINT.  Does POINT lie on it?
			    (and (monotonic? last-y point-y y) (return t)))
			   ((~= y point-y)
			    ;; POINT is horizontal with one endpoint of the current segment.  Count a crossing if it is the upper endpoint and it is to the left POINT.
			    (and (> y last-y) (< x point-x) (incf count)))
			   ((~= last-y point-y)
			    ;; POINT is horizontal with one endpoint of the current segment.  Count a crossing if it is the upper endpoint and it is to the left POINT.
			    (and (> last-y y) (< last-x point-x) (incf count)))
			   ((monotonic? last-y point-y y)
			    ;; We crossed the horizontal of POINT
			    (setf horizontal-intersect (horizontal-intersect point-y last-x last-y x y))
			    (cond ((~= horizontal-intersect point-x)
				   ;; POINT lies on this edge
				   (return t))
				  ((< horizontal-intersect point-x)
				   ;; The polygon edge is to the left of POINT
				   (incf count)))))
		     ;; If we intersect the CONVEX-POLYGON twice, we can't be inside it
		     (when (= count 2) (return nil))
		     ;; Keep track of last vertex
		     (setf last-x x
			   last-y y)
			 finally
			   (return (= count 1))))))))))

(defun MONOTONIC? (a b c)
  "
  Purpose:  Determine if three numbers are monotonically increasing or decreasing.
  Args:	A	-- a number
	B	-- a number
	C	-- a number
  Returns:  T if the sequence A B C is monotonic
  "
  (declare (values boolean))
  (or (<= a b c)
      (>= a b c)))

(defun HORIZONTAL-INTERSECT (y x1 y1 x2 y2)
  "
  Purpose:  Assuming that y1 /= y2, find the intersection of the line through (X1 Y1) and (X2 Y2) with the horizontal line through Y.
  Args:	Y	-- y-coordinate of horizontal line of interest
	X1	-- x-coordinate of one point on line
	Y1	-- y-coordinate of one point on line
	X2	-- x-coordinate of other point on line
	Y2	-- y-coordinate of other point on line
  Returns:  x-coordinate of intersection point
  "
  (declare (values x-coordinate))
  (if (= x1 x2)
      x1
      ;; Assume y /= y1 /= y2
      (+ x2 (* (/ (- y y2) (- y1 y2)) (- x1 x2)))))

(defun POLYGON-DIAMETER (poly)
  (let ((points (polygon-points poly)))
    (loop for pt1 in points
	  for rest-points on (cdr points)
	  maximize
	  (loop for pt2 in rest-points
		maximize
		(distance pt1 pt2)))))

(defun INTERSECT-RECTANGLE-AND-CONVEX-POLYGON? (rect poly)
  (or (loop for pt in (polygon-points poly)
	    thereis
	    (inside? pt rect))
      (loop for line in (polygon-line-segments poly)
	    thereis
	    (line-segment-crosses-rectangle? line rect))))

(defun INTERSECT-LINE-SEGMENT-AND-RECTANGLE? (ls rect)
  (or (inside? (line-segment-point-1 ls) rect)
      (inside? (line-segment-point-2 ls) rect)
      (line-segment-crosses-rectangle? ls rect)))

(defun LINE-SEGMENT-CROSSES-RECTANGLE? (ls rect)
  (let* ((pt1 (line-segment-point-1 ls))
	 (pt2 (line-segment-point-2 ls))
	 (x1 (point-x pt1))
	 (y1 (point-y pt1))
	 (x2 (point-x pt2))
	 (y2 (point-y pt2))
	 (m (/ (- y2 y1) (- x2 x1)))
	 (b (- y1 (* m x1)))
	 (bottom-cross (/ (- (rectangle-min-y rect) b) m))
	 (top-cross  (/ (- (rectangle-max-y rect) b) m))
	 (left-cross (+ b (* (rectangle-min-x rect) m)))
	 (right-cross (+ b (* (rectangle-max-x rect) m))))
    (or (and (~<= (min x1 x2) (rectangle-min-x rect))
	     (~<= (rectangle-min-x rect) (max x1 x2))
	     (~<= (rectangle-min-y rect) left-cross)
	     (~<= left-cross (rectangle-max-y rect)))
	(and (~<= (min x1 x2) (rectangle-max-x rect))
	     (~<= (rectangle-max-x rect) (max x1 x2))
	     (~<= (rectangle-min-y rect) right-cross)
	     (~<= right-cross (rectangle-max-y rect)))
	(and (~<= (min y1 y2) (rectangle-min-y rect))
	     (~<= (rectangle-min-y rect) (max y1 y2))
	     (~<= (rectangle-min-x rect) bottom-cross)
	     (~<= bottom-cross (rectangle-max-x rect)))
	(and (~<= (min y1 y2) (rectangle-max-y rect))
	     (~<= (rectangle-max-y rect) (max y1 y2))
	     (~<= (rectangle-min-x rect) top-cross)
	     (~<= top-cross (rectangle-max-x rect))))))
	
