;; Much of this code was written by Francisco Javier Blanco-Silva <fbs@math.purdue.edu>

;; A Polygon is a (very!) abstract data type, defined as follows: if you're on
;; a vertex of a polygon, you can follow next-polygon-vertex until you get back
;; to  the original vertex; along the way, all the vertices you touch are Polygon-vertices

;; the invariants are  (and (eq (Edge-vertex-2 (Polygon-vertex-backward-edge vertex)) vertex)
;;                          (eq (Edge-vertex-1 (Polygon-vertex-forward-edge vertex))  vertex))
;; for all vertices on the polygon.

(define-class Polygon-vertex Boundary-vertex
  ((= forward-edge maybe-uninitialized:)
   (= backward-edge maybe-uninitialized:)))

(define (next-polygon-vertex v)
  (Edge-vertex-2 (Polygon-vertex-forward-edge v)))

(define (previous-polygon-vertex v)
  (Edge-vertex-1 (Polygon-vertex-backward-edge v)))

(define (Triangulation<-points l)   ;; l is list of Points
  
  (define (convex-corner? v)      ;Decides if a given vertex in a given polygon is in a convex corner.
    
    (define determ                       ; 3x3 determinant (I'm sorry) 
      (lambda (a b c d e f g h i)
	(declare (flonum))
	(- (+ (* a e i) 
	      (* b f g) 
	      (* c d h))
	   (+ (* c e g) 
	      (* b d i) 
	      (* a f h)))))
    
    (define (convex-angle? ff ss tt)
      (declare (flonum))
      (let ((point-1 (Vertex-point ff))
	    (point-2 (Vertex-point ss))
	    (point-3 (Vertex-point tt)))
	(let ((a (Point-x point-1)) 
	      (b (Point-y point-1)) 
	      (c (Point-x point-2)) 
	      (d (Point-y point-2)) 
	      (e (Point-x point-3))
	      (f (Point-y point-3)))
	  (> (determ 1.0 1.0 1.0 a c e b d f) 0.0))))
    
    (convex-angle? (previous-polygon-vertex v) v (next-polygon-vertex v)))
  
  (define (other-vertices-outside-triangle? v)
    (let ((previous-vertex (previous-polygon-vertex v))
	  (next-vertex     (next-polygon-vertex v)))
      (let ((point-1 (Vertex-point previous-vertex))
	    (point-2 (Vertex-point v))
	    (point-3 (Vertex-point next-vertex)))
	(let ((x1 (Point-x point-1))   ; really dumb naming convention
	      (x2 (Point-y point-1))
	      (y1 (Point-x point-2))
	      (y2 (Point-y point-2))
	      (z1 (Point-x point-3))
	      (z2 (Point-y point-3)))
	  (let ((vertex-in-triangle?
		 (lambda (v)
		   (let ((p (Vertex-point v)))
		     (let ((aa (Point-x p))
			   (bb (Point-y p)))
		       (declare (flonum))
		       (let ((det (- (* (- x1 z1) (- y2 z2)) (* (- y1 z1) (- x2 z2))))
			     (det2 (- (* (- aa z1) (- y2 z2)) (* (- y1 z1) (- bb z2))))
			     (det3 (- (* (- x1 z1) (- bb z2)) (* (- aa z1) (- x2 z2)))))
			 (and (<= 0.0 det2) (<= det2 det)
			      (<= 0.0 det3) (<= det3 det))))))))
	    (let loop ((v (next-polygon-vertex next-vertex)))
	      (if (eq? v previous-vertex)
		  #t
		  (if (vertex-in-triangle? v)
		      #f
		      (loop (next-polygon-vertex v))))))))))
  
  ;; for each point on the polygon, make a vertex
  
  (let ((vertices (map1 (lambda (p)
			 (instantiate Polygon-vertex
				      point: p))
		       l)))
    
    ;; make the edges connecting the vertices; a forward edge always has
    ;; vertex-1 the edge it's leaving and vertex-2 the connecting edge
    
    (do ((l vertices (cdr l)))
	((null? (cdr l)) ; at last vertex
	 (Polygon-vertex-forward-edge-set! (car l) 
					   (instantiate Boundary-edge
							vertex-1: (car l)
							vertex-2: (car vertices))))
      (Polygon-vertex-forward-edge-set! (car l) 
					(instantiate Boundary-edge
						     vertex-1: (car l)
						     vertex-2: (cadr l))))
    
    ;; put in the backward edges
    
    (do ((l vertices (cdr l)))
	((null? (cdr l)) ; at last vertex
	 (Polygon-vertex-backward-edge-set! (car vertices)
					    (Polygon-vertex-forward-edge (car l))))
      (Polygon-vertex-backward-edge-set! (cadr l)
					 (Polygon-vertex-forward-edge (car l))))
    
    ;; add the existing edges to Vertex-edges.
    
    (for-each (lambda (v)
		(Vertex-edges-set! v (list (Polygon-vertex-forward-edge v)
					   (Polygon-vertex-backward-edge v)))
		(Boundary-vertex-boundary-edges-set! 
		 v 
		 (Boundary-vertex-boundary-edges<-list (list (Polygon-vertex-forward-edge v)
							     (Polygon-vertex-backward-edge v)))))
	      vertices)

    ;; the main loop.  We will be constructing new edges and triangles, so we have lists to
    ;; add them to.  We also need a vertex on the remaining polygon whenever we start this loop

    (let loop ((edges (map1 (lambda (v) (Polygon-vertex-forward-edge v))
			   vertices))
	       (triangles '())
	       (vertex-on-polygon (car vertices)))
      (if (eq? (next-polygon-vertex (next-polygon-vertex vertex-on-polygon)) ; one triangle left
	       (previous-polygon-vertex vertex-on-polygon))

	  ;; Here was have only three vertices, so we make the final triangle, convert all the
	  ;; lists in Vertex-edges and return the triangulation.

	  (let ((final-triangle
		 (instantiate Triangle
			      vertex-1: vertex-on-polygon
			      vertex-2: (next-polygon-vertex vertex-on-polygon)
			      vertex-3: (previous-polygon-vertex vertex-on-polygon)
			      edge-1: (Polygon-vertex-forward-edge vertex-on-polygon)
			      edge-2: (Polygon-vertex-forward-edge (next-polygon-vertex vertex-on-polygon))
			      edge-3: (Polygon-vertex-backward-edge vertex-on-polygon))))
	    (for-each (lambda (v)
			(Vertex-edges-set! v (Vertex-edges<-list (Vertex-edges v))))
		      vertices)
	    (instantiate Triangulation
			 vertices:  (Triangulation-vertices<-list vertices)
			 edges:     (Triangulation-edges<-list edges)
			 triangles: (Triangulation-triangles<-list (cons final-triangle triangles))
			 has-polygonal-boundary?: #t))
	  
	  ;; Here we examine each vertex on the remaining polygon in turn.  When we find a convex
	  ;; corner for which the triangle made up of three adjacent vertices does not contain any
	  ;; other vertices, we construct that triangle, remove it from the polygon, and go through
	  ;; it again.  Unfortunately, it seems to be an O(N^3) algorithm.  Is it?

	  (let inner ((v vertex-on-polygon))
	    (if (and (convex-corner? v)
		     (other-vertices-outside-triangle? v))
		(let ((v-1 (previous-polygon-vertex v))
		      (v+1 (next-polygon-vertex v))
		      (e-1 (Polygon-vertex-backward-edge v))
		      (e+1 (Polygon-vertex-forward-edge v)))
		  (with-co-instantiation
		   ((new-edge Edge 
			      vertex-1: v-1 
			      vertex-2: v+1 
			      )
		    (new-triangle Triangle 
				  vertex-1: v-1 
				  vertex-2: v 
				  vertex-3: v+1
				  edge-1: e-1 
				  edge-2: e+1 
				  edge-3: new-edge))
		   (Edge-for-each-vertex (lambda (v)
					   (Vertex-edges-set! v (cons new-edge (Vertex-edges v))))
					 new-edge)
		   (Polygon-vertex-backward-edge-set! v+1 new-edge)
		   (Polygon-vertex-forward-edge-set!  v-1 new-edge)
		   (loop (cons new-edge edges)
			 (cons new-triangle triangles)
			 v-1)))
		(inner (next-polygon-vertex v))))))))