diff --git a/scikits/image/graph/__init__.py b/scikits/image/graph/__init__.py index da648dff..1b55961b 100644 --- a/scikits/image/graph/__init__.py +++ b/scikits/image/graph/__init__.py @@ -1,10 +1,10 @@ try: from spath import shortest_path - from trace_path import trace_path + from mcp import MCP, MCP_Geometric, route_through_array except ImportError: - print """*** The shortest path extension has not been compiled. Run + print """*** The cython extensions have not been compiled. Run python setup.py build_ext -i in the source directory to build in-place. Please refer to INSTALL.txt -for further detail.""" +for further detail.""" \ No newline at end of file diff --git a/scikits/image/graph/setup.py b/scikits/image/graph/setup.py index 6b9950d6..214f36c9 100644 --- a/scikits/image/graph/setup.py +++ b/scikits/image/graph/setup.py @@ -1,7 +1,6 @@ #!/usr/bin/env python from scikits.image._build import cython - import os.path base_path = os.path.abspath(os.path.dirname(__file__)) @@ -14,12 +13,15 @@ def configuration(parent_package='', top_path=None): # This function tries to create C files from the given .pyx files. If # it fails, we build the checked-in .c files. - cython(['spath.pyx'], working_path=base_path) - cython(['trace_path.pyx'], working_path=base_path) + cython(['_spath.pyx'], working_path=base_path) + cython(['_mcp.pyx'], working_path=base_path) + cython(['heap.pyx'], working_path=base_path) - config.add_extension('spath', sources=['spath.c'], + config.add_extension('_spath', sources=['_spath.c'], include_dirs=[get_numpy_include_dirs()]) - config.add_extension('trace_path', sources=['trace_path.c'], + config.add_extension('_mcp', sources=['_mcp.c'], + include_dirs=[get_numpy_include_dirs()]) + config.add_extension('heap', sources=['heap.c'], include_dirs=[get_numpy_include_dirs()]) return config @@ -32,4 +34,4 @@ if __name__ == '__main__': url = 'http://stefanv.github.com/scikits.image/', license = 'Modified BSD', **(configuration(top_path='').todict()) - ) + ) \ No newline at end of file diff --git a/scikits/image/graph/spath.pxd b/scikits/image/graph/spath.pxd deleted file mode 100644 index 35fbae27..00000000 --- a/scikits/image/graph/spath.pxd +++ /dev/null @@ -1,3 +0,0 @@ -cimport numpy as np - -cpdef shortest_path(np.ndarray, int reach=?) \ No newline at end of file diff --git a/scikits/image/graph/spath.pyx b/scikits/image/graph/spath.pyx deleted file mode 100644 index 5b8d8f13..00000000 --- a/scikits/image/graph/spath.pyx +++ /dev/null @@ -1,82 +0,0 @@ -# -*- python -*- - -import numpy as np -cimport numpy as np - -cdef extern from "math.h": - double fabs(double f) - -cpdef shortest_path(np.ndarray arr, int reach=1): - """Find the shortest left-to-right path through an array. - - Parameters - ---------- - arr : (M, N) ndarray of float64 - reach : int, optional - By default (``reach = 1``), the shortest path can only move - one row up or down for every column it moves forward (i.e., - the path gradient is limited to 1). `reach` defines the - number of rows that can be skipped at each step. - - Returns - ------- - p : ndarray of int - For each column, give the row-coordinate of the - shortest path. - cost : float - Cost of path. This is the absolute sum of all the - differences along the path. - - """ - if arr.ndim != 2: - raise ValueError("Expected 2-D array as input") - - cdef np.ndarray[np.double_t, ndim=2] data = \ - np.ascontiguousarray(arr, dtype=np.double) - - cdef int M = arr.shape[0] - cdef int N = arr.shape[1] - - cdef np.ndarray[np.int_t, ndim=2] node = \ - np.empty((M, N), dtype=int) - - cdef np.ndarray[np.double_t, ndim=2] cost = \ - np.empty((M, N), dtype=np.double) - - cdef np.ndarray[np.int_t] out = np.empty((N,), dtype=int) - - cdef int c, r, rb, r_min_node - cdef int r_bracket_min = 0, r_bracket_max = 0 - cdef double delta0 = 0, delta1 = 0 - - cost[:, 0] = 0 - - for c in range(1, N): - for r in range(M): - r_bracket_min = r - reach - r_bracket_max = r + reach - - if r_bracket_min < 0: - r_bracket_min = 0 - if r_bracket_max > M - 1: - r_bracket_max = M - 1 - - node[r, c] = r_bracket_min - for rb in range(r_bracket_min, r_bracket_max + 1): - delta0 = fabs(data[r, c] - data[rb, c - 1]) - delta1 = fabs(data[r, c] - data[node[r, c], c - 1]) - if delta0 < delta1: - node[r, c] = rb - - cost[r, c] = cost[node[r, c], c - 1] + \ - fabs(data[r, c] - data[node[r, c], c - 1]) - - # Find minimum cost path - r_min_node = cost[:,-1].argmin() - - # Backtrack - out[N - 1] = r_min_node - for c in range(N - 1, 0, -1): - out[c - 1] = node[out[c], c] - - return out, cost[r_min_node, N - 1] diff --git a/scikits/image/graph/tests/test_trace_path.py b/scikits/image/graph/tests/test_trace_path.py deleted file mode 100644 index a10770cb..00000000 --- a/scikits/image/graph/tests/test_trace_path.py +++ /dev/null @@ -1,71 +0,0 @@ -import numpy as np -from numpy.testing import * - -from scikits.image.graph import trace_path - -a = np.ones((8,8), dtype=np.float32) -a[1:-1, 1] = 0 -a[1, 1:-1] = 0 - -## array([[ 1., 1., 1., 1., 1., 1., 1., 1.], -## [ 1., 0., 0., 0., 0., 0., 0., 1.], -## [ 1., 0., 1., 1., 1., 1., 1., 1.], -## [ 1., 0., 1., 1., 1., 1., 1., 1.], -## [ 1., 0., 1., 1., 1., 1., 1., 1.], -## [ 1., 0., 1., 1., 1., 1., 1., 1.], -## [ 1., 0., 1., 1., 1., 1., 1., 1.], -## [ 1., 1., 1., 1., 1., 1., 1., 1.]], dtype=float32) - -def test_basic(): - costs, return_path = trace_path(a, (1, 6), [(7, 2)]) - assert_array_equal(costs, - [[ 1., 1., 1., 1., 1., 1., 1., 1.], - [ 1., 0., 0., 0., 0., 0., 0., 1.], - [ 1., 0., 1., 1., 1., 1., 1., 1.], - [ 1., 0., 1., 2., 2., 2., 2., 2.], - [ 1., 0., 1., 2., 3., 3., 3., 3.], - [ 1., 0., 1., 2., 3., 4., 4., 4.], - [ 1., 0., 1., 2., 3., 4., 5., 5.], - [ 1., 1., 1., 2., 3., 4., 5., 6.]]) - - assert_array_equal(return_path, - [[(1, 6), - (1, 5), - (1, 4), - (1, 3), - (1, 2), - (2, 1), - (3, 1), - (4, 1), - (5, 1), - (6, 1), - (7, 2)]]) - -def test_no_diagonal(): - costs, path = trace_path(a, (1, 6), [(7, 2)], diagonal_steps=False) - assert_array_equal(costs, - [[ 2., 1., 1., 1., 1., 1., 1., 2.], - [ 1., 0., 0., 0., 0., 0., 0., 1.], - [ 1., 0., 1., 1., 1., 1., 1., 2.], - [ 1., 0., 1., 2., 2., 2., 2., 3.], - [ 1., 0., 1., 2., 3., 3., 3., 4.], - [ 1., 0., 1., 2., 3., 4., 4., 5.], - [ 1., 0., 1., 2., 3., 4., 5., 6.], - [ 2., 1., 2., 3., 4., 5., 6., 7.]]) - assert_array_equal(path, - [[(1, 6), - (1, 5), - (1, 4), - (1, 3), - (1, 2), - (1, 1), - (2, 1), - (3, 1), - (4, 1), - (5, 1), - (6, 1), - (6, 2), - (7, 2)]]) - -if __name__ == "__main__": - run_module_suite() diff --git a/scikits/image/graph/trace_path.pyx b/scikits/image/graph/trace_path.pyx deleted file mode 100644 index 0e362650..00000000 --- a/scikits/image/graph/trace_path.pyx +++ /dev/null @@ -1,151 +0,0 @@ -# -*- python -*- - -import numpy as numpy -cimport numpy as numpy -cimport cython - -@cython.boundscheck(False) -def trace_path(numpy.ndarray[numpy.float32_t, ndim=2] costs not None, - start, ends, diagonal_steps=True): - """Find the lowest-cost path from the start point to each given end point. - - Costs are given by the input array: a move onto any given position in the - costs array adds that cost to the path. Paths may be constrained to - vertical and horizontal moves only by passing False for the diagonal_steps - parameter. Costs must be non-negative! - - The array of cumulative costs from the starting point, and a list of paths - from the start to each end point are returned. - - Parameters - ---------- - costs : ndarray - start : tuple of ints - ``(x, y)`` position (i.e., ``(column, row)``) of starting position. - ends : list of tuple of ints - ``[(x1, y1), (x2, y2), ...]`` List of end points. - diagonal_steps : bool - Whether to allow diagonal steps (True, by default). - - Notes - ----- - Paths are found by (more or less) breadth-first search outward from the - starting point: each time a lower-cost route to a given pixel is found, that - pixel is marked "active"; the neighbors of all active pixels are then - examined to see if their costs can be lowered as well. This continues until - no pixels are marked active. - - """ - if costs.min() < 0: - raise ValueError("All costs must be non-negative.") - - try: - a, b = start - except: - raise ValueError("The start point must be an (x, y) pair") - - if not (0 <= a < costs.shape[0] and 0 <= b < costs.shape[1]): - raise ValueError("The start point must fall within the array") - - for end in ends: - try: - a, b = end - except: - raise ValueError("All end points must be (x, y) pairs") - if not (0 <= a < costs.shape[0] and 0 <= b < costs.shape[1]): - raise ValueError("The end points must fall within the array") - - cdef numpy.ndarray[numpy.float32_t, ndim=2] cumulative_costs = \ - numpy.empty_like(costs) - - cumulative_costs.fill(numpy.inf) - cumulative_costs[start] = 0 - costs_shape = (costs.shape[0], costs.shape[1]) - cdef numpy.ndarray[numpy.uint8_t, ndim=2] active_nodes = \ - numpy.zeros(costs_shape, dtype=numpy.uint8) - - active_nodes[start] = 1 - cdef numpy.ndarray[numpy.uint8_t, ndim=2] parent_nodes = \ - numpy.empty(costs_shape, dtype=numpy.uint8) - - parent_nodes.fill(255) - cdef numpy.ndarray[numpy.int8_t, ndim=2] offsets - if diagonal_steps: - offsets = numpy.array([[-1, -1], - [-1, 0], - [-1, 1], - [ 0, -1], - [ 0, 1], - [ 1, -1], - [ 1, 0], - [ 1, 1]], dtype=numpy.int8) - else: - offsets = numpy.array([[-1, 0], - [0, -1], - [0, 1], - [1, 0]], dtype=numpy.int8) - - cdef Py_ssize_t x, y, ox, oy, xo, yo, i - cdef Py_ssize_t a_xmax, a_xmin, a_ymax, a_ymin, tmp_xmax, \ - tmp_xmin, tmp_ymax, tmp_ymin - cdef unsigned int xmax, ymax, active, num_steps - xmax = costs.shape[0] - 1 - ymax = costs.shape[1] - 1 - num_steps = 0 - tmp_xmax = tmp_xmin = start[0] - tmp_ymax = tmp_ymin = start[1] - cdef float current_cost, current_cumulative_cost, cumulative_cost, new_cost - - while True: - active = 0 - # iterate over array - for x in range(0, xmax + 1): - for y in range(0, ymax + 1): - if active_nodes[x, y]: - active_nodes[x, y] = 0 - active = 1 - current_cumulative_cost = cumulative_costs[x, y] - - # iterate over offsets - for i in range(8): - ox = offsets[i, 0] - oy = offsets[i, 1] - xo = x + ox - yo = y + oy - if xo < 0 or xo > xmax or yo < 0 or yo > ymax: - continue - - current_cost = costs[xo, yo] - new_cost = current_cost + current_cumulative_cost - - # if a cheaper path to a given point is found, - # activate that point - if cumulative_costs[xo, yo] > new_cost: - cumulative_costs[xo, yo] = new_cost - parent_nodes[xo, yo] = i - active_nodes[xo, yo] = 1 - - if not active: - break - - cdef unsigned int startx, starty - startx = start[0] - starty = start[1] - return_paths = [] - # Trace the paths from the endpoints to the start - for end in ends: - path = None - x = end[0] - y = end[1] - if cumulative_costs[x, y] != numpy.inf: - path = [(x, y)] - while not (x == startx and y == starty): - i = parent_nodes[x, y] - ox = offsets[i, 0] - oy = offsets[i, 1] - x -= ox - y -= oy - path.append((x, y)) - path.reverse() - return_paths.append(path) - return cumulative_costs, return_paths