Trying to calculate winding numbers - a better submission added by delmarre on Mon Nov 20 21:13:33 2017

;;; Trying to implement the code at
;;; http://www.geomalgorithms.com/a03-_inclusion.html#wn_PnPoly()

(define is-left?
  (lambda (p0 p1 p2)
    (- (* (- (vector-ref p1 0) (vector-ref p0 0))
	  (- (vector-ref p2 1) (vector-ref p0 1)))
       
       (* (- (vector-ref p2 0) (vector-ref p0 0))
	  (- (vector-ref p1 1) (vector-ref p0 1))))))


;; Tests

(define v1 (vector 0 0))
(define v2 (vector 5 5))

(display "Should get >0, 0, <0")
(newline)
(display (is-left? v1 v2 (vector -5 0)))
(newline)
(display (is-left? v1 v2 (vector -50 -50)))
(newline)
(display (is-left? v1 v2 (vector 12 -50)))

(define find-wn
  (lambda (p ply #!optional (wn 0))
    ; Check if we have tested all of the edges yet
    (if (> (length ply) 1)
      ; If we haven't
      ; Do these checks with the edge made 
      ; of points ply[i] and ply[i + 1]
      (cond ((and (<= (vector-ref (car ply) 1)
		      (vector-ref p 1))
		  (> (vector-ref (cadr ply) 1)
		     (vector-ref p 1))
		  (> (is-left? (car ply)
			       (cadr ply)
			       p)
		     0))

	     (find-wn p (cdr ply) (+ wn 1)))

	    ((and (<= (vector-ref (cadr ply) 1)
		      (vector-ref p 1))
		  (< (is-left? (car ply)
			       (cadr ply)
			       p)
		     0))
	     (find-wn p (cdr ply) (- wn 1)))

	    (else (find-wn p (cdr ply) wn)))
      wn)))

;; Tests

; Something like the diagram for this method

(define g 
  (list
    (vector 0 0)
    (vector 0 10)
    (vector 6 10)
    (vector 6 2)
    (vector 10 2)
    (vector 10 6)
    (vector 3 6)
    (vector 3 8)
    (vector 14 8)
    (vector 14 0)
    (vector 0 0)))

; Test points chosen to correspond to the areas
; in the diagram so that we have known expected values

(newline)
(newline)
(display "Should get 0, 0, 1 2")
(newline)
; A point in the empty space that is bordered
; by the shape
(display (find-wn (vector 8 5) g))
(newline)
; A point very far outside
(display (find-wn (vector 25 25) g))
(newline)
; A point in the blue wn=1 section
; near the origin
(display (find-wn (vector 2 1) g))
(newline)
; A point in the green wn=2 section
; in the 'overlapping' portion of the diagram
(display (find-wn (vector 4 7) g))

;;; is-left? now seems to be working however find-wn is wildly wrong