Welcome, guest | Sign In | My Account | Store | Cart
"""
N-D Bresenham line algo
"""
import numpy as np
def _bresenhamline_nslope(slope):
    """
    Normalize slope for Bresenham's line algorithm.

    >>> s = np.array([[-2, -2, -2, 0]])
    >>> _bresenhamline_nslope(s)
    array([[-1., -1., -1.,  0.]])

    >>> s = np.array([[0, 0, 0, 0]])
    >>> _bresenhamline_nslope(s)
    array([[ 0.,  0.,  0.,  0.]])

    >>> s = np.array([[0, 0, 9, 0]])
    >>> _bresenhamline_nslope(s)
    array([[ 0.,  0.,  1.,  0.]])
    """
    scale = np.amax(np.abs(slope), axis=1).reshape(-1, 1)
    zeroslope = (scale == 0).all(1)
    scale[zeroslope] = np.ones(1)
    normalizedslope = np.array(slope, dtype=np.double) / scale
    normalizedslope[zeroslope] = np.zeros(slope[0].shape)
    return normalizedslope

def _bresenhamlines(start, end, max_iter):
    """
    Returns npts lines of length max_iter each. (npts x max_iter x dimension) 

    >>> s = np.array([[3, 1, 9, 0],[0, 0, 3, 0]])
    >>> _bresenhamlines(s, np.zeros(s.shape[1]), max_iter=-1)
    array([[[ 3,  1,  8,  0],
            [ 2,  1,  7,  0],
            [ 2,  1,  6,  0],
            [ 2,  1,  5,  0],
            [ 1,  0,  4,  0],
            [ 1,  0,  3,  0],
            [ 1,  0,  2,  0],
            [ 0,  0,  1,  0],
            [ 0,  0,  0,  0]],
    <BLANKLINE>
           [[ 0,  0,  2,  0],
            [ 0,  0,  1,  0],
            [ 0,  0,  0,  0],
            [ 0,  0, -1,  0],
            [ 0,  0, -2,  0],
            [ 0,  0, -3,  0],
            [ 0,  0, -4,  0],
            [ 0,  0, -5,  0],
            [ 0,  0, -6,  0]]])
    """
    if max_iter == -1:
        max_iter = np.amax(np.amax(np.abs(end - start), axis=1))
    npts, dim = start.shape
    nslope = _bresenhamline_nslope(end - start)

    # steps to iterate on
    stepseq = np.arange(1, max_iter + 1)
    stepmat = np.tile(stepseq, (dim, 1)).T

    # some hacks for broadcasting properly
    bline = start[:, np.newaxis, :] + nslope[:, np.newaxis, :] * stepmat

    # Approximate to nearest int
    return np.array(np.rint(bline), dtype=start.dtype)

def bresenhamline(start, end, max_iter=5):
    """
    Returns a list of points from (start, end] by ray tracing a line b/w the
    points.
    Parameters:
        start: An array of start points (number of points x dimension)
        end:   An end points (1 x dimension)
            or An array of end point corresponding to each start point
                (number of points x dimension)
        max_iter: Max points to traverse. if -1, maximum number of required
                  points are traversed

    Returns:
        linevox (n x dimension) A cumulative array of all points traversed by
        all the lines so far.

    >>> s = np.array([[3, 1, 9, 0],[0, 0, 3, 0]])
    >>> bresenhamline(s, np.zeros(s.shape[1]), max_iter=-1)
    array([[ 3,  1,  8,  0],
           [ 2,  1,  7,  0],
           [ 2,  1,  6,  0],
           [ 2,  1,  5,  0],
           [ 1,  0,  4,  0],
           [ 1,  0,  3,  0],
           [ 1,  0,  2,  0],
           [ 0,  0,  1,  0],
           [ 0,  0,  0,  0],
           [ 0,  0,  2,  0],
           [ 0,  0,  1,  0],
           [ 0,  0,  0,  0],
           [ 0,  0, -1,  0],
           [ 0,  0, -2,  0],
           [ 0,  0, -3,  0],
           [ 0,  0, -4,  0],
           [ 0,  0, -5,  0],
           [ 0,  0, -6,  0]])
    """
    # Return the points as a single array
    return _bresenhamlines(start, end, max_iter).reshape(-1, start.shape[-1])


if __name__ == "__main__":
    import doctest
    doctest.testmod()

Diff to Previous Revision

--- revision 1 2012-04-25 21:44:27
+++ revision 2 2012-05-16 16:15:04
@@ -2,20 +2,20 @@
 N-D Bresenham line algo
 """
 import numpy as np
-def bresenhamline_nslope(slope):
+def _bresenhamline_nslope(slope):
     """
     Normalize slope for Bresenham's line algorithm.
 
     >>> s = np.array([[-2, -2, -2, 0]])
-    >>> bresenhamline_nslope(s)
+    >>> _bresenhamline_nslope(s)
     array([[-1., -1., -1.,  0.]])
 
     >>> s = np.array([[0, 0, 0, 0]])
-    >>> bresenhamline_nslope(s)
+    >>> _bresenhamline_nslope(s)
     array([[ 0.,  0.,  0.,  0.]])
 
     >>> s = np.array([[0, 0, 9, 0]])
-    >>> bresenhamline_nslope(s)
+    >>> _bresenhamline_nslope(s)
     array([[ 0.,  0.,  1.,  0.]])
     """
     scale = np.amax(np.abs(slope), axis=1).reshape(-1, 1)
@@ -25,52 +25,46 @@
     normalizedslope[zeroslope] = np.zeros(slope[0].shape)
     return normalizedslope
 
-def bresenhamline_step(start, normalizedslope, err):
+def _bresenhamlines(start, end, max_iter):
     """
-    Performs a single incremental step according to Bresenham's algorithm.
+    Returns npts lines of length max_iter each. (npts x max_iter x dimension) 
 
-    >>> s = np.array([[4, 4, 4, 0]])
-    >>> ns = bresenhamline_nslope(np.array([[-2, -2, -2, 0]]))
-    >>> e = np.zeros((1, 4))
-    >>> bresenhamline_step(s, ns, e)
-    (array([[ 3.,  3.,  3.,  0.]]), array([[ 0.,  0.,  0.,  0.]]))
+    >>> s = np.array([[3, 1, 9, 0],[0, 0, 3, 0]])
+    >>> _bresenhamlines(s, np.zeros(s.shape[1]), max_iter=-1)
+    array([[[ 3,  1,  8,  0],
+            [ 2,  1,  7,  0],
+            [ 2,  1,  6,  0],
+            [ 2,  1,  5,  0],
+            [ 1,  0,  4,  0],
+            [ 1,  0,  3,  0],
+            [ 1,  0,  2,  0],
+            [ 0,  0,  1,  0],
+            [ 0,  0,  0,  0]],
+    <BLANKLINE>
+           [[ 0,  0,  2,  0],
+            [ 0,  0,  1,  0],
+            [ 0,  0,  0,  0],
+            [ 0,  0, -1,  0],
+            [ 0,  0, -2,  0],
+            [ 0,  0, -3,  0],
+            [ 0,  0, -4,  0],
+            [ 0,  0, -5,  0],
+            [ 0,  0, -6,  0]]])
+    """
+    if max_iter == -1:
+        max_iter = np.amax(np.amax(np.abs(end - start), axis=1))
+    npts, dim = start.shape
+    nslope = _bresenhamline_nslope(end - start)
 
-    >>> s = np.array([[0, 0, 9, 0]])
-    >>> ns = bresenhamline_nslope(-1 * s)
-    >>> e = np.zeros((1, 4))
-    """
-    nextpoint = start + normalizedslope + err
-    nextvox = np.rint(nextpoint)
-    err = nextpoint - nextvox
-    return nextvox, err
+    # steps to iterate on
+    stepseq = np.arange(1, max_iter + 1)
+    stepmat = np.tile(stepseq, (dim, 1)).T
 
-def bresenhamline_gen(start, end):
-    """
-    Returns a generator of points from (start, end] by ray tracing a line b/w
-    the points.
-    Parameters:
-        start: An array of start points (number of points x dimension)
-        end:   An end points (1 x dimension)
-            or An array of end point corresponding to each start point
-                (number of points x dimension)
+    # some hacks for broadcasting properly
+    bline = start[:, np.newaxis, :] + nslope[:, np.newaxis, :] * stepmat
 
-    Returns:
-        A generator that returns a point being traversed at each step
-        corresponding to each start point.
-
-    >>> [p for p in bresenhamline_gen(np.array([[3, 1, 3, 0]]), np.zeros(4))]
-    [array([[2, 1, 2, 0]]), array([[1, 0, 1, 0]]), array([[0, 0, 0, 0]])]
-    """
-    dim = start.shape[1]
-    err = np.zeros(start.shape)
-    slope = np.array(end - start, dtype=np.float32)
-    nslope = bresenhamline_nslope(slope)
-    cur_vox = start
-    while (np.sum(nslope != 0)):
-        cur_vox, err = bresenhamline_step(cur_vox, nslope, err)
-        # reached end ?
-        nslope[(np.abs(end - cur_vox) < 1).all(1)] = np.zeros(dim) 
-        yield np.array(cur_vox, dtype=start.dtype)
+    # Approximate to nearest int
+    return np.array(np.rint(bline), dtype=start.dtype)
 
 def bresenhamline(start, end, max_iter=5):
     """
@@ -81,40 +75,37 @@
         end:   An end points (1 x dimension)
             or An array of end point corresponding to each start point
                 (number of points x dimension)
-        max_iter: Max points to traverse
+        max_iter: Max points to traverse. if -1, maximum number of required
+                  points are traversed
 
     Returns:
         linevox (n x dimension) A cumulative array of all points traversed by
         all the lines so far.
 
-    >>> s = np.array([[0, 0, 3, 0]])
-    >>> bresenhamline(s, np.zeros(s.shape[1]), max_iter=9)
-    array([[0, 0, 2, 0],
-           [0, 0, 1, 0],
-           [0, 0, 0, 0]])
+    >>> s = np.array([[3, 1, 9, 0],[0, 0, 3, 0]])
+    >>> bresenhamline(s, np.zeros(s.shape[1]), max_iter=-1)
+    array([[ 3,  1,  8,  0],
+           [ 2,  1,  7,  0],
+           [ 2,  1,  6,  0],
+           [ 2,  1,  5,  0],
+           [ 1,  0,  4,  0],
+           [ 1,  0,  3,  0],
+           [ 1,  0,  2,  0],
+           [ 0,  0,  1,  0],
+           [ 0,  0,  0,  0],
+           [ 0,  0,  2,  0],
+           [ 0,  0,  1,  0],
+           [ 0,  0,  0,  0],
+           [ 0,  0, -1,  0],
+           [ 0,  0, -2,  0],
+           [ 0,  0, -3,  0],
+           [ 0,  0, -4,  0],
+           [ 0,  0, -5,  0],
+           [ 0,  0, -6,  0]])
+    """
+    # Return the points as a single array
+    return _bresenhamlines(start, end, max_iter).reshape(-1, start.shape[-1])
 
-    >>> s = np.array([[3, 1, 9, 0]])
-    >>> bresenhamline(s, np.zeros(s.shape[1]), max_iter=9)
-    array([[3, 1, 8, 0],
-           [2, 1, 7, 0],
-           [2, 1, 6, 0],
-           [2, 1, 5, 0],
-           [1, 0, 4, 0],
-           [1, 0, 3, 0],
-           [1, 0, 2, 0],
-           [0, 0, 1, 0],
-           [0, 0, 0, 0]])
-    """
-    dim = start.shape[1]
-    linevox = np.zeros((0, dim), dtype=start.dtype)
-    linegen = bresenhamline_gen(start, end)
-    try:
-        for i in range(max_iter):
-            cur_vox = linegen.next()
-            linevox = np.vstack((linevox, cur_vox))
-    except StopIteration:
-        pass
-    return np.array(linevox, dtype=start.dtype)
 
 if __name__ == "__main__":
     import doctest

History