-
Notifications
You must be signed in to change notification settings - Fork 0
/
lofi.tri.lisp
219 lines (190 loc) · 8.43 KB
/
lofi.tri.lisp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
;;;; lofi-tri.lisp
(in-package #:lofi.tri)
;;; "lofi-tri" goes here. Hacks and glory await!
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; structs
(defstruct circle
(center nil :type vec)
(radius 0.0 :type float)
(radius-squared 0.0 :type float)
(diameter 0.0 :type float))
(defstruct triangle
(verts #() :type vector)
(circumcircle nil :type circle))
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; Utilities
(defun random-point (&optional (state *random-state*))
"Return an instance of sb-cga:vec initialized with random values."
(apply #'vec (loop repeat 3
collect (random 1.0 state))))
(defun random-point-array (count &optional (state *random-state*))
"Returns an array of <count> random-points."
(let ((result (make-array count :fill-pointer 0)))
(dotimes (i count)
(vector-push (random-point state) result))
result))
(defmacro sort-by (point-set index)
"Sort the input point set by the value in the specified index."
`(sort ,point-set #'< :key (lambda (p) (aref p ,index))))
(defmacro sort-by-x (point-set)
"Sort the input point set by the value at 0"
`(sort-by ,point-set 0))
(defmacro sort-by-y (point-set)
"Sort the input point set by the value at 1"
`(sort-by ,point-set 1))
(defmacro sort-by-z (point-set)
"Sort the input point set by the value at 2"
`(sort-by ,point-set 2))
(defun get-min-max (point-set)
"Return the min and max vectors for the given point set. Effectively the bounding box."
(let* ((first-point (aref point-set 0))
(rest-points (subseq point-set 1))
(minx (aref first-point 0)) (maxx (aref first-point 0))
(miny (aref first-point 1)) (maxy (aref first-point 1))
(minz (aref first-point 2)) (maxz (aref first-point 2)))
(loop :for p :across rest-points :do
(setf minx (min minx (aref p 0)) maxx (max maxx (aref p 0))
miny (min miny (aref p 1)) maxy (max maxy (aref p 1))
minz (min minz (aref p 2)) maxz (max maxz (aref p 2))))
(values (vec minx miny minz) (vec maxx maxy maxz))))
(defun get-bounding-triangle-points (point-set &optional (fudge-factor 10))
(multiple-value-bind (min max) (get-min-max point-set)
(let ((dx (* fudge-factor (- (aref max 0) (aref min 0))))
(dy (* fudge-factor (- (aref max 1) (aref min 1)))))
(make-array 3 :initial-contents
(list (sb-cga:vec (- (aref min 0) dx) (- (aref min 1) (* dy 3)) 0.0)
(sb-cga:vec (- (aref min 0) dx) (+ (aref max 1) dy) 0.0)
(sb-cga:vec (+ (aref max 0) (* dx 3)) (+ (aref max 1) dy) 0.0))))))
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; Meat and potatos
(defun distance (v0 v1 &key 3d squared)
"Calculate the distance between two vectors in 2D or 3D.
Will return the square root of the result unless :squared t"
(declare (type vec v0 v1))
(let* ((diff (vec- v0 v1))
(result (apply #'+ (loop for i from 0 upto (if 3d 2 1)
collect (expt (aref diff i) 2)))))
(if squared
result
(sqrt result))))
(defun midpoint (v0 v1)
"Return a vector representing the midpoint between the two provided vectors."
(declare (type vec v0 v1))
(vec/ (vec+ v0 v1) 2.0))
(defun circumcircle (v0 v1 v2)
"Returns a circle struct representing the circumcircle of the given 3 vertices"
(let* ((v1-v0 (vec- v1 v0))
(v2-v0 (vec- v2 v0))
(v2-v1 (vec- v2 v1))
(v1+v0 (vec+ v1 v0))
(v2+v0 (vec+ v2 v0))
(a (aref v1-v0 0))
(b (aref v1-v0 1))
(c (aref v2-v0 0))
(d (aref v2-v0 1))
(e (+ (* a (aref v1+v0 0))
(* b (aref v1+v0 1))))
(f (+ (* c (aref v2+v0 0))
(* d (aref v2+v0 1))))
(g (* 2.0 (- (* a (aref v2-v1 1))
(* b (aref v2-v1 0)))))
(colinear? (< (abs g) +default-epsilon+))
(cx 0.0) (cy 0.0) (dx 0.0) (dy 0.0))
(if colinear?
(let ((minx (min (aref v0 0) (aref v1 0) (aref v2 0)))
(miny (min (aref v0 1) (aref v1 1) (aref v2 1)))
(maxx (max (aref v0 0) (aref v1 0) (aref v2 0)))
(maxy (max (aref v0 1) (aref v1 1) (aref v2 1))))
(setf cx (/ (+ minx maxx) 2)
cy (/ (+ miny maxy) 2)
dx (- cx minx)
dy (- cy miny)))
;; else
(setf cx (/ (- (* d e) (* b f)) g)
cy (/ (- (* a f) (* c e)) g)
dx (- cx (aref v0 0))
dy (- cy (aref v0 1))))
(let* ((radius-squared (+ (* dx dx)
(* dy dy)))
(radius (sqrt radius-squared)))
(make-circle :center (vec cx cy 0.0)
:radius radius
:radius-squared radius-squared
:diameter (* radius 2)))))
(defun new-triangle (vi0 vi1 vi2 points)
"Returns a new triangle."
(let ((v0 (aref points vi0))
(v1 (aref points vi1))
(v2 (aref points vi2)))
(make-triangle :verts (make-array 3 :initial-contents (list vi0 vi1 vi2))
:circumcircle (circumcircle v0 v1 v2))))
(defun in-circumcircle? (tri p)
"Does point 'p' sit within the circumcircle of 'tri'?"
(declare (type triangle tri) (type vec p))
(let* ((circumcircle (slot-value tri 'circumcircle))
(center (slot-value circumcircle 'center))
(dist-squared (distance center p :squared t)))
(<= dist-squared (slot-value circumcircle 'radius-squared))))
(defmacro edge= (a b)
`(or (and (= (first ,a) (first ,b))
(= (second ,a) (second ,b)))
(and (= (first ,a) (second ,b))
(= (second ,a) (first ,b)))))
(defun unique-edge? (edges a)
(let ((instance-count (length (remove-if-not (lambda (b) (edge= a b)) edges))))
(<= instance-count 1)))
(defun has-shared-verts? (a b)
(declare (type triangle a b))
(let* ((averts (slot-value a 'verts))
(bverts (slot-value b 'verts))
(av0 (aref averts 0))
(av1 (aref averts 1))
(av2 (aref averts 2))
(bv0 (aref bverts 0))
(bv1 (aref bverts 1))
(bv2 (aref bverts 2)))
(or (= bv0 av0) (= bv0 av1) (= bv0 av2)
(= bv1 av0) (= bv1 av1) (= bv1 av2)
(= bv2 av0) (= bv2 av1) (= bv2 av2))))
(defun add-vertex (vi triangles points)
(let* ((edges ())
(unaffected-tris ()))
;; For each triangle in the list we take the edges
;; of any triangle where vert is inside it's circumcircle
;; and append it to the edges list. Otherwise the triangle
;; is collected and stored in unaffected-tris
(setf unaffected-tris
(loop for tri in triangles
if (in-circumcircle? tri (aref points vi))
do (let* ((verts (slot-value tri 'verts))
(e0 (list (aref verts 0) (aref verts 1)))
(e1 (list (aref verts 1) (aref verts 2)))
(e2 (list (aref verts 2) (aref verts 0))))
(setf edges (append edges (list e0 e1 e2))))
else collect tri))
;; Remove any edges that are duplicate so that the edge
;; list only contains the boundary edges.
(setf edges (remove-if-not (curry #'unique-edge? edges) edges))
;; Using the remaining edges and our input vert create
;; new triangles and return them appended to our unaffected-tris list
(append unaffected-tris (loop for edge in edges
collect (let ((vi0 (first edge))
(vi1 (second edge)))
(new-triangle vi0 vi1 vi points))))))
(defun triangulate (points)
(let* (;; Add the coords for a large bounding triangle to the point set
(st-coords (get-bounding-triangle-points points))
(sti0 (length points))
(sti1 (1+ sti0))
(sti2 (1+ sti1))
(points (concatenate 'vector points st-coords))
;; Create the bounding triangle instance
(supertri (new-triangle sti0 sti1 sti2 points))
;; Initialize our triangle list
(triangles (list supertri)))
;; For each point in the list we get an updated set
;; of triangles by retesselating using the new point
(loop for i below (length points)
do (setf triangles (add-vertex i triangles points)))
;; Remove any triangles that share points with the super triangle
(remove-if (curry #'has-shared-verts? supertri) triangles)))