summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--doc/source/reference/c-api.array.rst66
-rw-r--r--doc/source/reference/c-api.types-and-structures.rst12
-rw-r--r--numpy/core/code_generators/generate_numpy_api.py2
-rw-r--r--numpy/core/code_generators/numpy_api_order.txt1
-rw-r--r--numpy/core/include/numpy/_neighborhood_iterator_imp.h133
-rw-r--r--numpy/core/include/numpy/ndarrayobject.h51
-rw-r--r--numpy/core/src/multiarray/iterators.c113
-rw-r--r--numpy/core/src/multiarray/multiarraymodule.c5
8 files changed, 383 insertions, 0 deletions
diff --git a/doc/source/reference/c-api.array.rst b/doc/source/reference/c-api.array.rst
index 87e2f8c96..fe41dac8d 100644
--- a/doc/source/reference/c-api.array.rst
+++ b/doc/source/reference/c-api.array.rst
@@ -1946,6 +1946,72 @@ Broadcasting (multi-iterators)
loop should be performed over the axis that won't require large
stride jumps.
+Neighborhood iterator
+---------------------
+
+Neighborhood iterators are subclasses of the iterator object, and can be used
+to iter over a neighborhood of a point. For example, you may want to iterate
+over every voxel of a 3d image, and for every such voxel, iterate over an
+hypercube. Neighborhood iterator automatically handle boundaries, thus making
+this kind of code much easier to write than manual boundaries handling, at the
+cost of a slight overhead.
+
+.. cfunction:: PyObject* PyArray_NeighborhoodIterNew(PyArrayIterObject* iter, npy_intp bounds)
+
+ This function creates a new neighborhood iterator from an existing
+ iterator. The neighborhood will be computed relatively to the position
+ currently pointed by *iter*. The *bounds* argument is expected to be a (2
+ * iter->ao->nd) arrays, such as the range bound[2*i]->bounds[2*i+1] defines
+ the range where to walk for dimension i (both bounds are included in the
+ walked coordinates). The bounds should be ordered for each dimension
+ (bounds[2*i] <= bounds[2*i+1]).
+
+ - The iterator holds a reference to iter
+ - Return NULL on failure (in which case the reference count of iter is not
+ changed)
+ - iter itself can be a Neighborhood iterator: this can be useful for .e.g
+ automatic boundaries handling
+ - the object returned by this function should be safe to use as a normal
+ iterator
+ - If the position of iter is changed, any subsequent call to
+ PyArrayNeighborhoodIter_Next is undefined behavior, and
+ PyArrayNeighborhoodIter_Reset must be called.
+ - If the coordinates point to a point outside the array, the iterator
+ points to an item which contains 0 of the same type as iter. More
+ elaborate schemes (constants, mirroring, repeat) may be implemented later.
+
+ .. code-block:: c
+
+ PyArrayIterObject *iter;
+ PyArrayNeighborhoodIterObject *neigh_iter;
+ iter = PyArray_IterNew(x);
+
+ //For a 3x3 kernel
+ bounds = {-1, 1, -1, 1};
+ neigh_iter = (PyArrayNeighborhoodIterObject*)PyArrayNeighborhoodIter_New(iter, bounds);
+
+ for(i = 0; i < iter->size; ++i) {
+ for (j = 0; j < neigh_iter->size; ++j) {
+ // Walk around the item currently pointed by iter->dataptr
+ PyArrayNeighborhoodIter_Next(neigh_iter);
+ }
+
+ // Move to the next point of iter
+ PyArrayIter_Next(iter);
+ PyArrayNeighborhoodIter_Reset(neigh_iter);
+ }
+
+.. cfunction:: int PyArrayNeighborhoodIter_Reset(PyArrayNeighborhoodIterObject* iter)
+
+ Reset the iterator position to the first point of the neighborhood. This
+ should be called whenever the iter argument given at
+ PyArray_NeighborhoodIterObject is changed (see example)
+
+.. cfunction:: int PyArrayNeighborhoodIter_Next(PyArrayNeighborhoodIterObject* iter)
+
+ After this call, iter->dataptr points to the next point of the
+ neighborhood. Calling this function after every point of the
+ neighborhood has been visited is undefined.
Array Scalars
-------------
diff --git a/doc/source/reference/c-api.types-and-structures.rst b/doc/source/reference/c-api.types-and-structures.rst
index b99702e11..ece7f341b 100644
--- a/doc/source/reference/c-api.types-and-structures.rst
+++ b/doc/source/reference/c-api.types-and-structures.rst
@@ -917,6 +917,18 @@ PyArrayMultiIter_Type
to be broadcast together. On return, the iterators are adjusted for
broadcasting.
+PyArrayNeighborhoodIter_Type
+----------------------------
+
+.. cvar:: PyArrayNeighborhoodIter_Type
+
+ This is an iterator object that makes it easy to loop over an N-dimensional
+ neighborhood.
+
+.. ctype:: PyArrayNeighborhoodIterObject
+
+ The C-structure corresponding to an object of :cdata:`PyArrayNeighborhoodIter_Type` is
+ the :ctype:`PyArrayNeighborhoodIterObject`.
PyArrayFlags_Type
-----------------
diff --git a/numpy/core/code_generators/generate_numpy_api.py b/numpy/core/code_generators/generate_numpy_api.py
index 4d639c3e5..33e3aeaa4 100644
--- a/numpy/core/code_generators/generate_numpy_api.py
+++ b/numpy/core/code_generators/generate_numpy_api.py
@@ -28,6 +28,7 @@ extern NPY_NO_EXPORT PyTypeObject PyArrayFlags_Type;
extern NPY_NO_EXPORT PyTypeObject PyArrayIter_Type;
extern NPY_NO_EXPORT PyTypeObject PyArrayMapIter_Type;
extern NPY_NO_EXPORT PyTypeObject PyArrayMultiIter_Type;
+extern NPY_NO_EXPORT PyTypeObject PyArrayNeighborhoodIter_Type;
extern NPY_NO_EXPORT PyTypeObject PyBoolArrType_Type;
extern NPY_NO_EXPORT PyBoolScalarObject _PyArrayScalar_BoolValues[2];
#else
@@ -39,6 +40,7 @@ NPY_NO_EXPORT PyTypeObject PyArrayFlags_Type;
NPY_NO_EXPORT PyTypeObject PyArrayIter_Type;
NPY_NO_EXPORT PyTypeObject PyArrayMapIter_Type;
NPY_NO_EXPORT PyTypeObject PyArrayMultiIter_Type;
+NPY_NO_EXPORT PyTypeObject PyArrayNeighborhoodIter_Type;
NPY_NO_EXPORT PyTypeObject PyBoolArrType_Type;
NPY_NO_EXPORT PyBoolScalarObject _PyArrayScalar_BoolValues[2];
#endif
diff --git a/numpy/core/code_generators/numpy_api_order.txt b/numpy/core/code_generators/numpy_api_order.txt
index e111f0fd3..3d1e40add 100644
--- a/numpy/core/code_generators/numpy_api_order.txt
+++ b/numpy/core/code_generators/numpy_api_order.txt
@@ -174,3 +174,4 @@ PyArray_MultiIterFromObjects
PyArray_GetEndianness
PyArray_GetNDArrayCFeatureVersion
PyArray_Acorrelate
+PyArray_NeighborhoodIterNew
diff --git a/numpy/core/include/numpy/_neighborhood_iterator_imp.h b/numpy/core/include/numpy/_neighborhood_iterator_imp.h
new file mode 100644
index 000000000..b110f1c15
--- /dev/null
+++ b/numpy/core/include/numpy/_neighborhood_iterator_imp.h
@@ -0,0 +1,133 @@
+#ifndef _NPY_INCLUDE_NEIGHBORHOOD_IMP
+#error You should not include this header directly
+#endif
+/*
+ * Private API (here for inline)
+ */
+static NPY_INLINE int
+_PyArrayNeighborhoodIter_IncrCoord(PyArrayNeighborhoodIterObject* iter);
+static NPY_INLINE int
+_PyArrayNeighborhoodIter_SetPtr(PyArrayNeighborhoodIterObject* iter);
+
+/*
+ * Inline implementations
+ */
+static NPY_INLINE int PyArrayNeighborhoodIter_Reset(PyArrayNeighborhoodIterObject* iter)
+{
+ int i;
+
+ for (i = 0; i < iter->nd; ++i) {
+ iter->coordinates[i] = iter->bounds[i][0];
+ }
+ _PyArrayNeighborhoodIter_SetPtr(iter);
+
+ return 0;
+}
+
+/*
+ * Update to next item of the iterator
+ *
+ * Note: this simply increment the coordinates vector, last dimension
+ * incremented first , i.e, for dimension 3
+ * ...
+ * -1, -1, -1
+ * -1, -1, 0
+ * -1, -1, 1
+ * ....
+ * -1, 0, -1
+ * -1, 0, 0
+ * ....
+ * 0, -1, -1
+ * 0, -1, 0
+ * ....
+ */
+#define _UPDATE_COORD_ITER(c) \
+ wb = iter->coordinates[c] < iter->bounds[c][1]; \
+ if (wb) { \
+ iter->coordinates[c] += 1; \
+ return 0; \
+ } \
+ else { \
+ iter->coordinates[c] = iter->bounds[c][0]; \
+ }
+
+static NPY_INLINE int _PyArrayNeighborhoodIter_IncrCoord(PyArrayNeighborhoodIterObject* iter)
+{
+ int i, wb;
+
+ for (i = iter->nd - 1; i >= 0; --i) {
+ _UPDATE_COORD_ITER(i)
+ }
+
+ return 0;
+}
+
+/*
+ * Version optimized for 2d arrays, manual loop unrolling
+ */
+static NPY_INLINE int _PyArrayNeighborhoodIter_IncrCoord2D(PyArrayNeighborhoodIterObject* iter)
+{
+ int wb;
+
+ _UPDATE_COORD_ITER(1)
+ _UPDATE_COORD_ITER(0)
+
+ return 0;
+}
+#undef _UPDATE_COORD_ITER
+
+#define _INF_SET_PTR(c) \
+ bd = iter->coordinates[c] + iter->_internal_iter->coordinates[c]; \
+ if (bd < 0 || bd > iter->dimensions[c]) { \
+ iter->dataptr = iter->zero; \
+ return 1; \
+ } \
+ offset = iter->coordinates[c] * iter->strides[c]; \
+ iter->dataptr += offset;
+
+/* set the dataptr from its current coordinates */
+static NPY_INLINE int _PyArrayNeighborhoodIter_SetPtr(PyArrayNeighborhoodIterObject* iter)
+{
+ int i;
+ npy_intp offset, bd;
+
+ iter->dataptr = iter->_internal_iter->dataptr;
+
+ for(i = 0; i < iter->nd; ++i) {
+ _INF_SET_PTR(i)
+ }
+
+ return 0;
+}
+
+static NPY_INLINE int _PyArrayNeighborhoodIter_SetPtr2D(PyArrayNeighborhoodIterObject* iter)
+{
+ npy_intp offset, bd;
+
+ iter->dataptr = iter->_internal_iter->dataptr;
+
+ _INF_SET_PTR(0)
+ _INF_SET_PTR(1)
+
+ return 0;
+}
+#undef _INF_SET_PTR
+
+/*
+ * Advance to the next neighbour
+ */
+static NPY_INLINE int PyArrayNeighborhoodIter_Next2D(PyArrayNeighborhoodIterObject* iter)
+{
+ _PyArrayNeighborhoodIter_IncrCoord2D(iter);
+ _PyArrayNeighborhoodIter_SetPtr2D(iter);
+
+ return 0;
+}
+
+static NPY_INLINE int PyArrayNeighborhoodIter_Next(PyArrayNeighborhoodIterObject* iter)
+{
+ _PyArrayNeighborhoodIter_IncrCoord (iter);
+ _PyArrayNeighborhoodIter_SetPtr(iter);
+
+ return 0;
+}
diff --git a/numpy/core/include/numpy/ndarrayobject.h b/numpy/core/include/numpy/ndarrayobject.h
index 0826f59a6..57c4daab5 100644
--- a/numpy/core/include/numpy/ndarrayobject.h
+++ b/numpy/core/include/numpy/ndarrayobject.h
@@ -909,6 +909,57 @@ typedef struct {
} PyArrayMapIterObject;
+typedef struct {
+ PyObject_HEAD
+
+ /*
+ * PyArrayIterObject part: keep this in this exact order
+ */
+ int nd_m1; /* number of dimensions - 1 */
+ npy_intp index, size;
+ npy_intp coordinates[NPY_MAXDIMS];/* N-dimensional loop */
+ npy_intp dims_m1[NPY_MAXDIMS]; /* ao->dimensions - 1 */
+ npy_intp strides[NPY_MAXDIMS]; /* ao->strides or fake */
+ npy_intp backstrides[NPY_MAXDIMS];/* how far to jump back */
+ npy_intp factors[NPY_MAXDIMS]; /* shape factors */
+ PyArrayObject *ao;
+ char *dataptr; /* pointer to current item*/
+ npy_bool contiguous;
+
+ /*
+ * New members
+ */
+ npy_intp nd;
+
+ /* Dimensions is the dimension of the array */
+ npy_intp dimensions[NPY_MAXDIMS];
+ /* Bounds of the neighborhood to iterate over */
+ npy_intp bounds[NPY_MAXDIMS][2];
+
+ /* Neighborhood points coordinates are computed relatively to the point pointed
+ * by _internal_iter */
+ PyArrayIterObject* _internal_iter;
+ /* To keep a reference to the zero representation correponding to the dtype
+ * of the array we iterate over */
+ char* zero;
+} PyArrayNeighborhoodIterObject;
+
+/*
+ * Neighborhood iterator API
+ */
+static NPY_INLINE int
+PyArrayNeighborhoodIter_Reset(PyArrayNeighborhoodIterObject* iter);
+static NPY_INLINE int
+PyArrayNeighborhoodIter_Next(PyArrayNeighborhoodIterObject* iter);
+static NPY_INLINE int
+PyArrayNeighborhoodIter_Next2D(PyArrayNeighborhoodIterObject* iter);
+
+/* Include inline implementations - functions defined there are not considered
+ * public API */
+#define _NPY_INCLUDE_NEIGHBORHOOD_IMP
+#include "_neighborhood_iterator_imp.h"
+#undef _NPY_INCLUDE_NEIGHBORHOOD_IMP
+
/* The default array type
*/
#define NPY_DEFAULT_TYPE NPY_DOUBLE
diff --git a/numpy/core/src/multiarray/iterators.c b/numpy/core/src/multiarray/iterators.c
index ca5518de5..842522d73 100644
--- a/numpy/core/src/multiarray/iterators.c
+++ b/numpy/core/src/multiarray/iterators.c
@@ -1731,3 +1731,116 @@ NPY_NO_EXPORT PyTypeObject PyArrayMultiIter_Type = {
0, /* *tp_next */
#endif
};
+
+/*========================= Neighborhood iterator ======================*/
+
+/*NUMPY_API*/
+NPY_NO_EXPORT PyObject*
+PyArray_NeighborhoodIterNew(PyArrayIterObject *x, intp *bounds)
+{
+ int i;
+ PyArrayNeighborhoodIterObject *ret;
+
+ ret = _pya_malloc(sizeof(*ret));
+ if (ret == NULL) {
+ return NULL;
+ }
+ PyObject_Init((PyObject *)ret, &PyArrayNeighborhoodIter_Type);
+
+ array_iter_base_init((PyArrayIterObject*)ret, x->ao);
+ Py_INCREF(x);
+ ret->_internal_iter = x;
+
+ ret->nd = x->ao->nd;
+
+ /* Compute the neighborhood size and copy the shape */
+ ret->size = 1;
+ for (i = 0; i < ret->nd; ++i) {
+ ret->bounds[i][0] = bounds[2 * i];
+ ret->bounds[i][1] = bounds[2 * i + 1];
+ ret->size *= (bounds[2*i+1] - bounds[2*i]) + 1;
+ }
+
+ for (i = 0; i < ret->nd; ++i) {
+ ret->dimensions[i] = x->ao->dimensions[i];
+ }
+ ret->zero = PyArray_Zero(x->ao);
+
+ /*
+ * XXX: we force x iterator to be non contiguous because we need
+ * coordinates... Modifying the iterator here is not great
+ */
+ x->contiguous = 0;
+
+ PyArrayNeighborhoodIter_Reset(ret);
+
+ return (PyObject*)ret;
+}
+
+static void neighiter_dealloc(PyArrayNeighborhoodIterObject* iter)
+{
+ PyDataMem_FREE(iter->zero);
+ Py_DECREF(iter->_internal_iter);
+
+ array_iter_base_dealloc((PyArrayIterObject*)iter);
+ _pya_free((PyArrayObject*)iter);
+}
+
+NPY_NO_EXPORT PyTypeObject PyArrayNeighborhoodIter_Type = {
+ PyObject_HEAD_INIT(NULL)
+ 0, /*ob_size*/
+ "numpy.neigh_internal_iter", /*tp_name*/
+ sizeof(PyArrayNeighborhoodIterObject), /*tp_basicsize*/
+ 0, /*tp_itemsize*/
+ (destructor)neighiter_dealloc, /*tp_dealloc*/
+ 0, /*tp_print*/
+ 0, /*tp_getattr*/
+ 0, /*tp_setattr*/
+ 0, /*tp_compare*/
+ 0, /*tp_repr*/
+ 0, /*tp_as_number*/
+ 0, /*tp_as_sequence*/
+ 0, /*tp_as_mapping*/
+ 0, /*tp_hash */
+ 0, /*tp_call*/
+ 0, /*tp_str*/
+ 0, /*tp_getattro*/
+ 0, /*tp_setattro*/
+ 0, /*tp_as_buffer*/
+ Py_TPFLAGS_DEFAULT, /*tp_flags*/
+ 0, /* tp_doc */
+ 0, /* tp_traverse */
+ 0, /* tp_clear */
+ 0, /* tp_richcompare */
+ 0, /* tp_weaklistoffset */
+ 0, /* tp_iter */
+ (iternextfunc)0, /* tp_iternext */
+ 0, /* tp_methods */
+ 0, /* tp_members */
+ 0, /* tp_getset */
+ 0, /* tp_base */
+ 0, /* tp_dict */
+ 0, /* tp_descr_get */
+ 0, /* tp_descr_set */
+ 0, /* tp_dictoffset */
+ (initproc)0, /* tp_init */
+ 0, /* tp_alloc */
+ 0, /* tp_new */
+ 0, /* tp_free */
+ 0, /* tp_is_gc */
+ 0, /* tp_bases */
+ 0, /* tp_mro */
+ 0, /* tp_cache */
+ 0, /* tp_subclasses */
+ 0, /* tp_weaklist */
+ 0, /* tp_del */
+
+#ifdef COUNT_ALLOCS
+ /* these must be last and never explicitly initialized */
+ 0, /* tp_allocs */
+ 0, /* tp_frees */
+ 0, /* tp_maxalloc */
+ 0, /* tp_prev */
+ 0, /* *tp_next */
+#endif
+};
diff --git a/numpy/core/src/multiarray/multiarraymodule.c b/numpy/core/src/multiarray/multiarraymodule.c
index a10636f16..084a3dd02 100644
--- a/numpy/core/src/multiarray/multiarraymodule.c
+++ b/numpy/core/src/multiarray/multiarraymodule.c
@@ -2674,6 +2674,11 @@ PyMODINIT_FUNC initmultiarray(void) {
if (PyType_Ready(&PyArrayMultiIter_Type) < 0) {
return;
}
+ PyArrayNeighborhoodIter_Type.tp_new = PyType_GenericNew;
+ if (PyType_Ready(&PyArrayNeighborhoodIter_Type) < 0) {
+ return;
+ }
+
PyArrayDescr_Type.tp_hash = PyArray_DescrHash;
if (PyType_Ready(&PyArrayDescr_Type) < 0) {
return;