Jump to content

K-d tree: Difference between revisions

K-d tree (view source)

Revision as of 11:03, 10 November 2021

842 bytes added ,  2 years ago
→‎{{header|Common Lisp}}: make code organization more idiomatic, remove redefinition spamming
m (Minor edit to C++ code)
(→‎{{header|Common Lisp}}: make code organization more idiomatic, remove redefinition spamming)
Line 495:
The 3D data set is the set of coordinates from (0,0,0) to (9,9,9) (ie. uniformly distributed), while the ordinates of the random target are positive real < 10.
<lang lisp>
(in-package cl-user)
(defun main ()
(let ((dims 0) (target nil) (hits 0))
;;; distance node to target:
;;; returns squared euclidean distance, or squared semi distance if option set
(defun distance (n &optional (semi nil))
(if semi (expt (- (nth (first n) (second n)) (nth (first n) target)) 2)
(reduce #'+ (mapcar (lambda (x y) (* (- x y) (- x y))) (second n) target))))
 
(defvar *random-target*
;;; returns true if target is to its left in axis dim
(list (float (/ (random 1000) 100))
(defun target< (n)
(< (nth (first n) target)float (nth/ (firstrandom n1000) (second n100))))
(float (/ (random 1000) 100))))
 
(defun distance (n target &key (semi nil))
;;; return the next child when nn searching, return opposing child if option oppose set
"distance node to target:
(defun next-node (n &optional (oppose nil))
returns squared euclidean distance, or squared semi distance if option set"
(if (or (and (target< n) (not oppose)) (and (not (target< n)) oppose)) (third n) (fourth n)))
(if semi
(expt (- (nth (first n) (second n))
(nth (first n) target))
2)
(reduce #'+
(mapcar (lambda (x y) (* (- x y) (- x y)))
(second n)
target))))
 
(defun target< (n target)
;;; a kdtree is a binary tree where nodes are:
"returns true if target is to its left in axis dim"
;;; terminal: (axis data-point), or
(< (nth (first n) target)
;;; branch: (axis split-point (left-kdtree) (right-kdtree))
(defunnth make-kdtree(axisfirst datan) (second n))))
(if (null data) nil
(if (eql (length data) 1) ; singleton?
(list axis (first data)) ;; terminal node
;; else branch node:
;; #pts=odd splits list into 2 even parts with sp in middle
;; #pts=even splits list into 2 uneven parts with shorter length first (but never nil)
(let ((sd (sort (copy-list data) #'< :key (lambda (x) (nth axis x)))) ;; sort the axis ordinates
(sp (truncate (/ (length data) 2))) ;; get mid pt
(nxta (mod (1+ axis) dims)))
(list axis (nth sp sd) (make-kdtree nxta (subseq sd 0 sp)) (make-kdtree nxta (subseq sd (1+ sp))))))))
 
(defun next-node (n target &key (opposite nil))
;;; depth first visit all nodes in kdtree and optionally apply a function to each node visited
"return the next child when nn searching, return opposing child if option
(defun visit-kdtree (kdt &key (node-function null))
oppose set"
(when kdt
(if (or (and (target< n target) (not opposite))
(when node-function (funcall node-function kdt))
(visit-kdtreeand (thirdnot kdt)(target< :node-functionn node-functiontarget)) opposite))
(third n)
(visit-kdtree (fourth kdt) :node-function node-function)))
(fourth n)))
;;; count of the terminal nodes
(defun count-nodes (kdt)
(if kdt
(if (eql (length kdt) 2) 1
(+ 1 (count-nodes (third kdt)) (count-nodes (fourth kdt))))
0))
;;; nearest neighbour search
(defun nn-kdtree (kdt node-stack)
(when kdt
;; stage 1 - find the 'closest' terminal node using insertion logic
(let ((best (do ((node kdt (next-node node))) ((not (next-node node)) (incf hits) node) ;; return first best est.
(push node node-stack) (incf hits)))) ; iteration
;; stage 2 - unwind the path, at each node if node is closer then make it best
(do ((node (pop node-stack) (pop node-stack))) ((null node) best) ;; return nearest pt
;; iteration: update best if node is closer
(when (< (distance node) (distance best))
(setf best node))
 
(defun make-kdtree (axis data dims)
;; venture down opposing side if split point is inside HS
"a kdtree is a binary tree where nodes are:
(let ((opposing-best
terminal: (axis data-point), or
(if (< (distance node 'semi) (distance best)) ; use semi dist here
branch: (axis split-point (nnleft-kdtree) (nextright-node node 'opposite) (listkdtree))"
(cond ((null data) nil)
nil))) ;; otherwise ignore this subtree
((eql (length data) 1) ; singleton?
(list axis (first data))) ; terminal node
(when (and opposing-best (< (distance opposing-best) (distance best)))
;; else branch node:
(setf best opposing-best)))))))
;; #pts=odd splits list into 2 even parts with sp in middle
;; #pts=even splits list into 2 uneven parts with shorter length first (but never nil)
(t
(let ((sd (sort (copy-list data)
#'<
:key (lambda (x) (nth axis x)))) ; sort the axis ordinates
(sp (truncate (length data) 2)) ; get mid pt
(nxta (mod (1+ axis) dims)))
(list axis
(nth sp sd)
(make-kdtree nxta
(subseq sd 0 sp)
dims)
(make-kdtree nxta
(subseq sd (1+ sp))
dims))))))
 
(defun visit-kdtree (kdt &key (node-function nil))
"depth first visit all nodes in kdtree and optionally apply a function to
each node visited"
(when kdt
(when node-function (funcall node-function kdt))
(visit-kdtree (third kdt) :node-function node-function)
(visit-kdtree (fourth kdt) :node-function node-function)))
 
(defun count-nodes (kdt)
"count of the terminal nodes"
(cond ((null kdt) 0)
((= (length kdt) 2) 1)
(t (+ 1
(count-nodes (third kdt))
(count-nodes (fourth kdt))))))
 
(defvar *hits*)
 
(defun nn-kdtree (kdt node-stack target)
"nearest neighbour search"
(when kdt
;; stage 1 - find the 'closest' terminal node using insertion logic
(let ((best (do ((node kdt (next-node node target)))
((not (next-node node target))
(incf *hits*)
node) ; return first best est.
(push node node-stack)
(incf *hits*))))
;; stage 2 - unwind the path, at each node if node is closer then make it best
(do ((node (pop node-stack) (pop node-stack)))
((null node) best) ; return nearest pt
;; iteration: update best if node is closer
(when (< (distance node target) (distance best target))
(setf best node))
;; venture down opposing side if split point is inside HS
(let ((opposing-best
(if (< (distance node target :semi t) (distance best target))
(nn-kdtree (next-node node target :opposite t)
(list)
target)
nil))) ; otherwise ignore this subtree
(when (and opposing-best
(< (distance opposing-best target) (distance best target)))
(setf best opposing-best)))))))
 
(defun process (data target
&key (render nil)
&aux (dims (length target)))
"process one set of data & optionally display tree"
(let ((*hits* 0))
(let* ((kdt (make-kdtree 0 data dims))
(nn (nn-kdtree kdt (list) target)))
(when render
(visit-kdtree kdt
:node-function
(lambda (n)
(format t
"~A node: axis:~A point: ~A target:~A semi-distance-sqd:~A euclidean-distance-sqd:~A~%"
(if (not (next-node n target))
"TERMINAL"
"BRANCH")
(first n)
(second n)
target
(distance n target :semi t)
(distance n target)))))
(format t "~%NN to ~A is ~A, distance ~A [tree has ~A nodes, ~A were visited.]~%"
target
(second nn)
(sqrt (distance nn target))
(count-nodes kdt)
*hits*))))
 
(defun main ()
;; TASK 1 - nn search small set of 2D points
(process '((2 3) (5 4) (9 6) (4 7) (8 1) (7 2))
'(9 2)
:render t)
 
;; TASK 2 - nn search 1000 coordinate points in 3D space
;;; process one set of data & optionally display tree
(defun process (datalet tgt &optional((ll (render nillist)))
(setfdotimes target(i tgt10)
(setf dims (dotimes (lengthj target)10)
(setfdotimes hits(k 010)
(let* ((kdt (make-kdtree 0 data)) (nn (nn-kdtree kdt (push (list i j k) ll))))
(when render ll)
(visit-kdtree kdt*random-target*))
:node-function (lambda (n)
(format t "~A node: axis:~A point: ~A target:~A semi-distance-sqd:~A euclidean-distance-sqd:~A~%"
(if (not (next-node n)) "TERMINAL" "BRANCH") (first n) (second n)
target (distance n 'semi) (distance n)))))
(format t "~%NN to ~A is ~A, distance ~A [tree has ~A nodes, ~A were visited.]~%" target (second nn) (sqrt (distance nn)) (count-nodes kdt) hits)))
;; MAIN: TASK 1 - nn search small set of 2D points
(process '((2 3) (5 4) (9 6) (4 7) (8 1) (7 2)) '(9 2) 'render)
;; TASK 2 - nn search 1000 coordinate points in 3D space
(process
(progn (let ((ll (list))) (dotimes (i 10) (dotimes (j 10) (dotimes (k 10) (push (list i j k) ll)))) ll))
(list (float (/ (random 1000) 100)) (float (/ (random 1000) 100)) (float (/ (random 1000) 100))))))
 
</lang>
Anonymous user
Retrieved from "https://rosettacode.org/wiki/Special:MobileDiff/213125"
Cookies help us deliver our services. By using our services, you agree to our use of cookies.
More information