dmj5:DMJ_FastClosestPoint
|
|||||||||
| PREV CLASS NEXT CLASS | FRAMES NO FRAMES | ||||||||
| SUMMARY: NESTED | FIELD | CONSTR | METHOD | DETAIL: FIELD | CONSTR | METHOD | ||||||||
Object
class
Fast point finder tool class.
Given a set of target points and a test point, find the target
point that is closest to the test point.
The simplest method is brute-force; test against each point
and record the closest. However, for large data sets (e.g.
one million points) this is very slow. Instead, we build a
sorted list of the points and index it with a binary tree.
This reduces the number of points that need to be tested by
locating for us a point that is closest on one particular
axis. From there we walk in each direction looking for close
points.
This was adapted from dmj3.uxf:FastMosaic.
class DMJ_FastClosestPoint {
; Fast point finder tool class.
;
; Given a set of target points and a test point, find the target
; point that is closest to the test point.
;
; The simplest method is brute-force; test against each point
; and record the closest. However, for large data sets (e.g.
; one million points) this is very slow. Instead, we build a
; sorted list of the points and index it with a binary tree.
; This reduces the number of points that need to be tested by
; locating for us a point that is closest on one particular
; axis. From there we walk in each direction looking for close
; points.
;
; This was adapted from dmj3.uxf:FastMosaic.
public:
import "common.ulb"
func DMJ_FastClosestPoint()
m_Points = 0
m_Extra = 0
m_TreeLeft = 0
m_TreeRight = 0
m_OrderLeft = 0
m_OrderRight = 0
m_Order = 0
endfunc
func SetPointData(ComplexArray ppoints, Array pextra)
m_Points = ppoints
m_Extra = pextra
; start with just the first element in the index
int length = m_Points.GetArrayLength()
m_TreeLeft = new IntegerArray(length)
m_TreeRight = new IntegerArray(length)
m_OrderLeft = new IntegerArray(length)
m_OrderRight = new IntegerArray(length)
m_TreeLeft.m_Elements[0] = -1 ; both directions lead nowhere; this node is a leaf
m_TreeRight.m_Elements[0] = -1
m_OrderLeft.m_Elements[0] = -1 ; this node marks the end of the ordered list
m_OrderRight.m_Elements[0] = -1
; insert all other points into the index
int j = 1
int k = 0
int l = 0
complex p = 0
bool leftbranch = false
while (j < length)
; get current point
p = m_Points.m_Elements[j]
; walk the existing tree to find the placement
k = 0
while (k >= 0) ; go until we're at a leaf node
l = k ; save last valid tree node
if (real(p) < real(m_Points.m_Elements[k]) || \
(real(p) == real(m_Points.m_Elements[k]) && imag(p) < imag(m_Points.m_Elements[k]))) ; p is left of k
k = m_TreeLeft.m_Elements[k]
leftbranch = true
else
k = m_TreeRight.m_Elements[k]
leftbranch = false
endif
endwhile
; insert the new node in the tree
m_TreeLeft.m_Elements[j] = -1
m_TreeRight.m_Elements[j] = -1
if (leftbranch)
m_TreeLeft.m_Elements[l] = j
else
m_TreeRight.m_Elements[l] = j
endif
; insert the new node into the sorted list
; this requires updating two links: one on
; each side of the newly-inserted point
if (leftbranch)
if (m_OrderLeft.m_Elements[l] >= 0) ; not the leftmost node so far
m_OrderRight.m_Elements[m_OrderLeft.m_Elements[l]] = j ; make node to left of new point link to it
endif
m_OrderLeft.m_Elements[j] = m_OrderLeft.m_Elements[l] ; link to point left of new point
m_OrderRight.m_Elements[j] = l ; link to point right of new point
m_OrderLeft.m_Elements[l] = j ; make node to right of new point link to it
else
if (m_OrderRight.m_Elements[l] >= 0) ; not the rightmost node so far
m_OrderLeft.m_Elements[m_OrderRight.m_Elements[l]] = j ; make node to right of new point link to it
endif
m_OrderRight.m_Elements[j] = m_OrderRight.m_Elements[l] ; link to point right of new point
m_OrderLeft.m_Elements[j] = l ; link to point left of new point
m_OrderRight.m_Elements[l] = j ; make node to left of new point link to it
endif
j = j + 1
endwhile
; now we have the full order, but the tree is probably
; not balanced; create a single simple index in sorted
; order
; allocate space for the index
m_Order = new IntegerArray(length)
; find the leftmost node
k = 0
while (k >= 0)
l = k
k = m_TreeLeft.m_Elements[k]
endwhile
; copy index numbers
j = 0
while (l >= 0)
m_Order.m_Elements[j] = l
l = m_OrderRight.m_Elements[l]
j = j + 1
endwhile
; release unneeded indices
m_TreeLeft = 0
m_TreeRight = 0
m_OrderLeft = 0
m_OrderRight = 0
endfunc
; given precomputed index, find closest point in our set
int func GetClosestPoint(complex pz)
; start by finding the closest index based solely on the real value
int k = GetClosestRealIndex(pz)
; now scan left and right to find nearby values that are candidates
; for shorter distances
int length = m_Order.GetArrayLength()
int left = k-1
int right = k+1
; complex c = m_Points.m_Elements[m_Order.m_Elements[k]] ; current closest point
float d = GetDistanceIndex(pz, k) ; current closest distance
complex q = 0
float d2 = 0
while (left >= 0 || right < length)
if (left >= 0)
q = m_Points.m_Elements[m_Order.m_Elements[left]] ; coordinate to test
d2 = GetDistanceLimit(real(pz), real(q))
if (d2 > d) ; too far on just the real axis, so nothing in this direction is closer; stop
left = -1
else
d2 = GetDistanceIndex(pz, left)
if (d2 < d) ; this point is closer than our best
d = d2
; c = q
k = left
endif
endif
left = left - 1
endif
if (right < length)
q = m_Points.m_Elements[m_Order.m_Elements[right]] ; coordinate to test
d2 = GetDistanceLimit(real(pz), real(q))
if (d2 > d) ; too far on just the real axis, so nothing in this direction is closer; stop
right = length
else
d2 = GetDistanceIndex(pz, right)
if (d2 < d) ; this point is closer than our best
d = d2
; c = q
k = right
endif
endif
right = right + 1
endif
endwhile
return k
endfunc
; like GetClosestPoint(), but returns a sorted list of
; the closest set of points (as many as requested)
; NOTE: less efficient than GetClosestPoint()
int func GetClosestPointSet(complex pz, int pcount, IntegerArray pindices, FloatArray pdistances, ComplexArray ppoints)
; start by finding the closest point
int k = GetClosestPoint(pz)
; now scan left and right to find nearby values that are candidates
; for shorter distances
int length = m_Order.GetArrayLength()
int left = k-1
int right = k+1
int found = 1
ppoints.m_Elements[0] = m_Points.m_Elements[m_Order.m_Elements[k]] ; current closest point
pdistances.m_Elements[0] = GetDistanceIndex(pz, k) ; current closest distance
complex q = 0
float d2 = 0
int l = 1
int m = 0
while (l < pcount)
pdistances.m_Elements[l] = 1e20
l = l + 1
endwhile
while (left >= 0 || right < length)
if (left >= 0)
q = m_Points.m_Elements[m_Order.m_Elements[left]] ; coordinate to test
d2 = GetDistanceLimit(real(pz), real(q))
if (found >= pcount && d2 > pdistances.m_Elements[found-1]) ; too far on just the real axis, so nothing in this direction is closer; stop
left = -1
else
d2 = GetDistanceIndex(pz, left)
l = 1 ; check every distance to see if this is closer
while (l < pcount)
if (d2 < pdistances.m_Elements[l])
m = found-1 ; distance is closer than previous ones; shift them down
while (m > l)
ppoints.m_Elements[m] = ppoints.m_Elements[m-1]
pdistances.m_Elements[m] = pdistances.m_Elements[m-1]
pindices.m_Elements[m] = pindices.m_Elements[m-1]
m = m - 1
endwhile
ppoints.m_Elements[l] = q
pdistances.m_Elements[l] = d2
pindices.m_Elements[l] = left
if (found < pcount)
found = found + 1
endif
l = pcount + 1
endif
l = l + 1
endwhile
endif
left = left - 1
endif
if (right < length)
q = m_Points.m_Elements[m_Order.m_Elements[right]] ; coordinate to test
d2 = GetDistanceLimit(real(pz), real(q))
if (found >= pcount && d2 > pdistances.m_Elements[found-1]) ; too far on just the real axis, so nothing in this direction is closer; stop
right = length
else
d2 = GetDistanceIndex(pz, right)
l = 1 ; check every distance to see if this is closer
while (l < pcount)
if (d2 < pdistances.m_Elements[l])
m = found-1 ; distance is closer than previous ones; shift them down
while (m > l)
ppoints.m_Elements[m] = ppoints.m_Elements[m-1]
pdistances.m_Elements[m] = pdistances.m_Elements[m-1]
pindices.m_Elements[m] = pindices.m_Elements[m-1]
m = m - 1
endwhile
ppoints.m_Elements[l] = q
pdistances.m_Elements[l] = d2
pindices.m_Elements[l] = right
if (found < pcount)
found = found + 1
endif
l = pcount + 1
endif
l = l + 1
endwhile
endif
right = right + 1
endif
endwhile
return k
endfunc
int func GetClosestRealIndex(complex pz)
; we already have an index that gives us the points in
; sorted real order; just binary search that list
int length = m_Order.GetArrayLength()
int j = 0
int k = 0
int l = 0 ; search range lower bound
int m = length-1 ; search range upper bound
while (l < m)
k = floor((l+m) / 2) ; midpoint between bounds
j = m_Order.m_Elements[k] ; actual point index
if (pz == m_Points.m_Elements[j])
; exact match
l = k
m = k
elseif (real(pz) < real(m_Points.m_Elements[j]) || \
(real(pz) == real(m_Points.m_Elements[j]) && imag(pz) < imag(m_Points.m_Elements[j])))
; lower half
m = k
else
; upper half
if (l == k) ; span of just one; point lies somewhere between
if ((real(pz) < (real(m_Points.m_Elements[m_Order.m_Elements[l]]) + \
real(m_Points.m_Elements[m_Order.m_Elements[m]]))/2) || \
(real(pz) == real(m_Points.m_Elements[m_Order.m_Elements[l]]) && \
real(pz) == real(m_Points.m_Elements[m_Order.m_Elements[m]]) && \
imag(pz) < (imag(m_Points.m_Elements[m_Order.m_Elements[l]]) + \
imag(m_Points.m_Elements[m_Order.m_Elements[m]]))/2))
m = l
else
l = m
k = m
endif
else
l = k
endif
endif
endwhile
return k
endfunc
; distance metric; override to use a different metric
float func GetDistance(complex ppoint1, complex ppoint2)
return |ppoint1-ppoint2|
endfunc
; distance metric, indexed version
float func GetDistanceIndex(complex ppoint, int pindex)
return |ppoint-m_Points.m_Elements[m_Order.m_Elements[pindex]]|
endfunc
; if this distance exceeds the closest distance so far, scanning will stop
float func GetDistanceLimit(float pr1, float pr2)
return sqr(pr1-pr2)
endfunc
; these are public rather than protected
ComplexArray m_Points
Array m_Extra
IntegerArray m_Order
protected:
IntegerArray m_TreeLeft
IntegerArray m_TreeRight
IntegerArray m_OrderLeft
IntegerArray m_OrderRight
default:
}
| Constructor Summary | |
|---|---|
DMJ_FastClosestPoint()
|
|
| Method Summary | |
|---|---|
int |
GetClosestPoint(complex pz)
given precomputed index, find closest point in our set |
int |
GetClosestPointSet(complex pz,
int pcount,
IntegerArray pindices,
FloatArray pdistances,
ComplexArray ppoints)
like GetClosestPoint(), but returns a sorted list of the closest set of points (as many as requested) NOTE: less efficient than GetClosestPoint() |
int |
GetClosestRealIndex(complex pz)
|
float |
GetDistance(complex ppoint1,
complex ppoint2)
distance metric; override to use a different metric |
float |
GetDistanceIndex(complex ppoint,
int pindex)
distance metric, indexed version |
float |
GetDistanceLimit(float pr1,
float pr2)
if this distance exceeds the closest distance so far, scanning will stop |
void |
SetPointData(ComplexArray ppoints,
Array pextra)
|
| Methods inherited from class Object |
|---|
|
| Constructor Detail |
|---|
public DMJ_FastClosestPoint()
| Method Detail |
|---|
public void SetPointData(ComplexArray ppoints,
Array pextra)
public int GetClosestPoint(complex pz)
public int GetClosestPointSet(complex pz,
int pcount,
IntegerArray pindices,
FloatArray pdistances,
ComplexArray ppoints)
public int GetClosestRealIndex(complex pz)
public float GetDistance(complex ppoint1,
complex ppoint2)
public float GetDistanceIndex(complex ppoint,
int pindex)
public float GetDistanceLimit(float pr1,
float pr2)
|
|||||||||
| PREV CLASS NEXT CLASS | FRAMES NO FRAMES | ||||||||
| SUMMARY: NESTED | FIELD | CONSTR | METHOD | DETAIL: FIELD | CONSTR | METHOD | ||||||||