Skip to content
This repository has been archived by the owner on Feb 23, 2018. It is now read-only.

Commit

Permalink
Imtermediate step. KD-traversal is going to be replaced.
Browse files Browse the repository at this point in the history
  • Loading branch information
pavlo-alkhimov committed Feb 24, 2011
1 parent 759bda0 commit a5d59ed
Show file tree
Hide file tree
Showing 3 changed files with 184 additions and 131 deletions.
5 changes: 4 additions & 1 deletion README
Original file line number Diff line number Diff line change
@@ -1,4 +1,7 @@
;; -*- mode: org; indent-tabs-mode: nil -*-

* Tasks
** TODO Implement kd-tree traversal
** TODO KD-Tree traversal
*** DONE AABB-Ray intersection. :traversal:
*** DONE Ray-Split intersection :traversal:
*** TODO TA-B-rec algorithm
282 changes: 152 additions & 130 deletions traverse/traverse-kd-tree.lisp
Original file line number Diff line number Diff line change
@@ -1,130 +1,150 @@
(in-package #:kd)

(defun method-parser (x0 y0 z0 x1 y1 z1 x y z dx dy dz method)
"[xyz]0-[xyz]1 is a AABB, where [xyz]0 <= [xyz]1;
[xyz] is origin of RAY, d[xyz] is direction of the RAY;
METHOD is one of: :[xyz]-line | :(xy|yz|xz)-plane | :arbitrary | :zero-direction;
Return value: entry point of RAY to AABB."
(case method
(:x-line (when (and (<= y0 y y1) (<= z0 z z1) (<= x x0))
(list x0
y
z
(/ (- x0 x) dx)
(/ (- x1 x) dx))))
(:y-line (when (and (<= x0 x x1) (<= z0 z z1) (<= y y0))
(list x
y0
z
(/ (- y0 y) dy)
(/ (- y1 y) dy))))
(:z-line (when (and (<= x0 x x1) (<= y0 y y1) (<= z z0))
(list x
y
z0
(/ (- z0 z) dz)
(/ (- z1 z) dz))))
(:xy-plane (let* ((t-from-x (/ (- x0 x) dx))
(t-from-y (/ (- y0 y) dy)))
(when (and (<= 0.0 t-from-x) (<= 0.0 t-from-y) (<= z0 z z1))
(let* ((x-hit (+ x (* dx t-from-y)))
(y-hit (+ y (* dy t-from-x))))
(if (<= x0 x-hit x1)
(list x-hit y0 z)
(when (<= y0 y-hit y1)
(list x0 y-hit z)))))))
(:xz-plane (let* ((t-from-x (/ (- x0 x) dx))
(t-from-z (/ (- z0 z) dz))
(x-hit (+ x (* dx t-from-z)))
(z-hit (+ z (* dz t-from-x))))
(when (and (<= 0.0 t-from-x) (<= 0.0 t-from-z) (<= y0 y y1))
(if (<= x0 x-hit x1)
(list x-hit y z0)
(when (<= z0 z-hit z1)
(list x0 y z-hit))))))
(:yz-plane (let* ((t-from-y (/ (- y0 y) dy))
(t-from-z (/ (- z0 z) dz))
(y-hit (+ y (* dy t-from-z)))
(z-hit (+ z (* dz t-from-y))))
(when (and (<= 0.0 t-from-y) (<= 0.0 t-from-z) (<= x0 x x1))
(if (<= y0 y-hit y1)
(list x y-hit z0)
(when (<= z0 z-hit z1)
(list x y0 z-hit))))))
(:arbitrary (let* ((te (/ (- x0 x))) ;; hit YZ plane at x = x0 dx
(y-res (+ y (* dy te)))
(z-res (+ z (* dz te))))
(if (and (<= 0.0 te) (<= y0 y-res y1) (<= z0 z-res z1))
(list x0 y-res z-res)
(let* ((te (/ (- y0 y) dy)) ;; hit XZ plane at y = y0
(x-res (+ x (* dx te)))
(z-res (+ z (* dz te))))
(if (and (<= 0.0 te) (<= x0 x-res x1) (<= z0 z-res z1))
(list x-res y0 z-res)
(let* ((te (/ (- z0 z) dz)) ;; hit XY plane at z = z0
(x-res (+ x (* dx te)))
(y-res (+ y (* dy te))))
(when (and (<= 0.0 te) (<= x0 x-res x1) (<= y0 y-res y1))
(list x-res y-res z0))))))))))

(defun intersection-of-aabb-with-ray (aabb ray)
(declare (optimize (debug 3)))
(let* ((x (aref ray 0 0))
(y (aref ray 0 1))
(z (aref ray 0 2))
(dx (aref ray 1 0))
(dy (aref ray 1 1))
(dz (aref ray 1 2))
(aabb-data (corners aabb))
(x0 (car aabb-data)) (aabb-data (cdr aabb-data))
(y0 (car aabb-data)) (aabb-data (cdr aabb-data))
(z0 (car aabb-data)) (aabb-data (cdr aabb-data))
(x1 (car aabb-data)) (aabb-data (cdr aabb-data))
(y1 (car aabb-data)) (aabb-data (cdr aabb-data))
(z1 (car aabb-data))
(method
(cond ((not (= 0.0 (* dx dy dz))) :arbitrary) ;; the most frequent case is tested first
((= 0.0 dx dy dz) :zero-direction)
((= 0.0 dy dz) :x-line)
((= 0.0 dx dz) :y-line)
((= 0.0 dx dy) :z-line)
((= 0.0 dz) :xy-plane)
((= 0.0 dx) :yz-plane)
((= 0.0 dy) :xz-plane)
(t (error "Failed to test the (~a,~a,~a) vector." dx dy dz))))
;; To make dx, dy and dz >=0, the space can be flipped:
(flipped-x (when (< dx 0.0)
(setf dx (- dx))
(setf x (- x))
(setf x0 (- x1)) ;; !!! 0 -> 1
(setf x1 (- x0)))) ;; !!! 1 -> 0
(flipped-y (when (< dy 0.0)
(setf dy (- dy))
(setf y (- y))
(setf y0 (- y1)) ;; !!! 0 -> 1
(setf y1 (- y0)))) ;; !!! 1 -> 0
(flipped-z (when (< dz 0.0)
(setf dz (- dz))
(setf z (- z))
(setf z0 (- z1)) ;; !!! 0 -> 1
(setf z1 (- z0)))) ;; !!! 1 -> 0
(result (case method
(:x-line (when (and (<= y0 y y1) (<= z0 z z1) (<= x x0))
(list x0 y z)))
(:y-line (when (and (<= x0 x x1) (<= z0 z z1) (<= y y0))
(list x y0 z)))
(:z-line (when (and (<= x0 x x1) (<= y0 y y1) (<= z z0))
(list x y z0)))
(:xy-plane (let* ((t-from-x (/ (- x0 x) dx))
(t-from-y (/ (- y0 y) dy)))
(when (and (<= 0.0 t-from-x) (<= 0.0 t-from-y) (<= z0 z z1))
(let* ((x-hit (+ x (* dx t-from-y)))
(y-hit (+ y (* dy t-from-x))))
(if (<= x0 x-hit x1)
(list x-hit y0 z)
(when (<= y0 y-hit y1)
(list x0 y-hit z)))))))
(:xz-plane (let* ((t-from-x (/ (- x0 x) dx))
(t-from-z (/ (- z0 z) dz))
(x-hit (+ x (* dx t-from-z)))
(z-hit (+ z (* dz t-from-x))))
(when (and (<= 0.0 t-from-x) (<= 0.0 t-from-z) (<= y0 y y1))
(if (<= x0 x-hit x1)
(list x-hit y z0)
(when (<= z0 z-hit z1)
(list x0 y z-hit))))))
(:yz-plane (let* ((t-from-y (/ (- y0 y) dy))
(t-from-z (/ (- z0 z) dz))
(y-hit (+ y (* dy t-from-z)))
(z-hit (+ z (* dz t-from-y))))
(when (and (<= 0.0 t-from-y) (<= 0.0 t-from-z) (<= x0 x x1))
(if (<= y0 y-hit y1)
(list x y-hit z0)
(when (<= z0 z-hit z1)
(list x y0 z-hit))))))
(:arbitrary (let* ((te (/ (- x0 x))) ;; hit YZ plane at x = x0 dx
(y-res (+ y (* dy te)))
(z-res (+ z (* dz te))))
(if (and (<= 0.0 te) (<= y0 y-res y1) (<= z0 z-res z1))
(list x0 y-res z-res)
(let* ((te (/ (- y0 y) dy)) ;; hit XZ plane at y = y0
(x-res (+ x (* dx te)))
(z-res (+ z (* dz te))))
(if (and (<= 0.0 te) (<= x0 x-res x1) (<= z0 z-res z1))
(list x-res y0 z-res)
(let* ((te (/ (- z0 z) dz)) ;; hit XY plane at z = z0
(x-res (+ x (* dx te)))
(y-res (+ y (* dy te))))
(when (and (<= 0.0 te) (<= x0 x-res x1) (<= y0 y-res y1))
(list x-res y-res z0))))))))))))
(when result
(when flipped-x
(setf (nth 0 result) (- (nth 0 result))))
(when flipped-y
(setf (nth 1 result) (- (nth 1 result))))
(when flipped-z
(setf (nth 2 result) (- (nth 2 result)))))
(values result
method
(when flipped-x :flipped-x)
(when flipped-y :flipped-y)
(when flipped-z :flipped-z)))
(macrolet ((swap-inverse (a b) `(let ((c ,a))
(setf ,a (- ,b))
(setf ,b (- c))))
(invert (a b) `(progn (setf ,a (- ,a))
(setf ,b (- ,b)))))
(let* ((x (aref ray 0 0)) (y (aref ray 0 1)) (z (aref ray 0 2))
(dx (aref ray 1 0)) (dy (aref ray 1 1)) (dz (aref ray 1 2))
(aabb-data (corners aabb))
(x0 (car aabb-data)) (aabb-data (cdr aabb-data))
(y0 (car aabb-data)) (aabb-data (cdr aabb-data))
(z0 (car aabb-data)) (aabb-data (cdr aabb-data))
(x1 (car aabb-data)) (aabb-data (cdr aabb-data))
(y1 (car aabb-data)) (aabb-data (cdr aabb-data))
(z1 (car aabb-data))
(method
(cond ((not (= 0.0 (* dx dy dz))) :arbitrary) ;; the most frequent case is tested first
((= 0.0 dx dy dz) :zero-direction)
((= 0.0 dy dz) :x-line)
((= 0.0 dx dz) :y-line)
((= 0.0 dx dy) :z-line)
((= 0.0 dz) :xy-plane)
((= 0.0 dx) :yz-plane)
((= 0.0 dy) :xz-plane)
(t (error "Failed to test the (~a,~a,~a) vector." dx dy dz))))
;; To make d[xyz]>=0, flip the space
(flipped-x (when (< dx 0.0)
(invert x dx)
(swap-inverse x0 x1)
-1.0))
(flipped-y (when (< dy 0.0)
(invert y dy)
(swap-inverse y0 y1)
-1.0))
(flipped-z (when (< dz 0.0)
(invert z dz)
(swap-inverse z0 z1)
-1.0))
(result (method-parser x0 y0 z0 x1 y1 z1 x y z dx dy dz method)))
(when result
(when flipped-x
(setf (nth 0 result)
(- (nth 0 result))))
(when flipped-y
(setf (nth 1 result)
(- (nth 1 result))))
(when flipped-z
(setf (nth 2 result)
(- (nth 2 result)))))
(values result
method
(list (or flipped-x 1.0)
(or flipped-y 1.0)
(or flipped-z 1.0))))))

(defun intersection-of-split-with-ray (&key aabb ray axis-index split point)
(defun intersection-of-split-with-ray (&key aabb ray axis-index split point flips)
(declare (optimize (debug 3))
(type fixnum axis-index)
(type coordinate split)
(type list point)
(type (simple-array coordinate (2 3)) ray))
(let* ((d (aref ray 1 axis-index));; i.e. DX
(p (nth axis-index point)));; i.e. X
;; RAY touches the POINT?
(if (= p split)
;; RAY lies in plane?
(if (= 0.0 d)
(values point
:in-plane)
(values point
(if (> 0.0 d) :R :L)))
;; is RAY outgoing?
(if (< 0.0 (* (- p split) d))
(let* ((d (aref ray 1 axis-index)) ;; i.e. DX
(p (nth axis-index point))
(flip (nth axis-index flips)))
(if (= 0.0 d)
(if (= p split)
(values point :in-plane)
(if (< (* flip p) (* flip split))
(values point :left-only)
(values point :right-only)))
(if (< 0.0 (* d (- p split)))
(values point
(if (> 0.0 d) :R :L))
(let* ((hit (let* ((te (/ (- split
(aref ray 0 axis-index))
(if (< (* flip p) (* flip split))
:left-only
:right-only))
;; Here we potentially cross the SPLIT
(let* ((hit (let* ((te (/ (- split (aref ray 0 axis-index))
d)))
(loop for i from 0 below 3
collect (+ (aref ray 0 i)
Expand All @@ -151,29 +171,32 @@
10.0 10.0 10.0)))
(ray1 (make-array '(2 3)
:element-type 'coordinate
:initial-contents '(( -1.0 -1.0 11.0)
( 1.0 1.0 -3.0)))))
:initial-contents '((11.0 11.0 5.0)
(-1.0 -1.0 0.0)))))
(multiple-value-bind
(point method fx fy fz)
(point method flips)
(intersection-of-aabb-with-ray aabb1 ray1)
(multiple-value-bind
(a b)
(intersection-of-split-with-ray :aabb aabb1 :ray ray1
:axis-index 2 :split 5.0
:point point)
(values a b point method fx fy fz)))))
(intersection-of-split-with-ray :aabb aabb1
:ray ray1
:axis-index 2
:split 6.0
:point point
:flips flips)
(values a b point method)))))

(defun test2 ()
#|(defun test2 ()
(let* ((p (load-patch "d:/Paul.revised//git.repos//github//srt//data//dodecahedron.obj"))
(k (build-tree p :aabb (corners (aabb p))))
(r (make-array '(2 3)
:element-type 'coordinate
:initial-contents '((-10.0 -10.0 -10.0)
( 1.0 1.0 1.0)))))
(setf (tree p) k)
(intersection-of-patch-with-ray p r)))
(intersection-of-patch-with-ray p r)))|#

(defun intersection-of-patch-with-ray (patch ray)
#|(defun intersection-of-patch-with-ray (patch ray)
(declare (type tri-patch patch))
(multiple-value-bind
(point method fx fy fz)
Expand All @@ -192,9 +215,9 @@
(corners (aabb patch))
ray
:axis-index 0
:point point)))))
:point point)))))|#

(defun intersection-of-kd-with-ray (node aabb ray &key (axis-index 0) (point nil))
#|(defun intersection-of-kd-with-ray (node aabb ray &key (axis-index 0) (point nil))
(and node
(if (and (l node)
(not (r node)))
Expand All @@ -207,5 +230,4 @@
(format t "Entering the leaf split by ~a~%" (case (mod axis-index 3) (0 'x) (1 'y) (2 'z)))
(or (intersection-of-kd-with-ray (l node) l-aabb ray :axis-index next-axis :point point)
(intersection-of-kd-with-ray (r node) r-aabb ray :axis-index next-axis :point point))
(format t "Leaving the leaf split by ~a~%" (case (mod axis-index 3) (0 'x) (1 'y) (2 'z))))))))

(format t "Leaving the leaf split by ~a~%" (case (mod axis-index 3) (0 'x) (1 'y) (2 'z))))))))|#
28 changes: 28 additions & 0 deletions types/aabb.lisp
Original file line number Diff line number Diff line change
Expand Up @@ -5,3 +5,31 @@
:type list))
(:documentation "Axis-Aligned Bounding Box.
CORNERS are represented as list (x0 y0 z0 x1 y1 z1)"))

(defmethod initialize-instance :after ((box aabb) &key from-list)
(setf (slot-value box 'corners)
(copy-tree from-list)))

(defun calc-aabb (vertexes &key (start -1) (end -1))
(let* ((start (if (= -1 start) 0 start))
(end (if (= -1 end)
(array-dimension vertexes 0)
end)))
(make-instance 'aabb
:from-list (loop for i from start below end
minimizing (aref vertexes i 0) into x0
minimizing (aref vertexes i 1) into y0
minimizing (aref vertexes i 2) into z0
maximizing (aref vertexes i 0) into x1
maximizing (aref vertexes i 1) into y1
maximizing (aref vertexes i 2) into z1
finally (return (list x0 y0 z0 x1 y1 z1))))))

(defun split-aabb (aabb &key axis position)
(let* ((l (copy-tree aabb))
(r (copy-tree l)))
(setf (nth (+ axis 3) l)
position)
(setf (nth axis r)
position)
(values l r)))

0 comments on commit a5d59ed

Please sign in to comment.