[Numpy-discussion] A Cython apply_along_axis function
Dag Sverre Seljebotn
dagss at student.matnat.uio.no
Thu Dec 2 04:08:12 EST 2010
On 12/02/2010 08:17 AM, Robert Bradshaw wrote:
> On Wed, Dec 1, 2010 at 6:09 PM, John Salvatier
> <jsalvati at u.washington.edu> wrote:
>> On Wed, Dec 1, 2010 at 6:07 PM, Keith Goodman<kwgoodman at gmail.com> wrote:
>>> On Wed, Dec 1, 2010 at 5:53 PM, David<david at silveregg.co.jp> wrote:
>>>> On 12/02/2010 04:47 AM, Keith Goodman wrote:
>>>>> It's hard to write Cython code that can handle all dtypes and
>>>>> arbitrary number of dimensions. The former is typically dealt with
>>>>> using templates, but what do people do about the latter?
>>>> The only way that I know to do that systematically is iterator. There is
>>>> a relatively simple example in scipy/signal (lfilter.c.src).
>>>> I wonder if it would be possible to add better support for numpy
>>>> iterators in cython...
>>> Thanks for the tip. I'm starting to think that for now I should just
>>> template both dtype and ndim.
>>> NumPy-Discussion mailing list
>>> NumPy-Discussion at scipy.org
>> I enthusiastically support better iterator support for cython
> I enthusiastically welcome contributions along this line.
Me too :-)
I guess we're moving into more Cython-list territory, so let's move any
follow-ups there (posting this one both places).
Just in case anybody is wondering what something like this could look
like, here's a rough scetch complete with bugs. The idea would be to a)
add some rudimentary support for using the yield keyword in Cython to
make a generator function, b) inline the generator function if the
generator is used directly in a for-loop. This should result in very
efficient code, and would also be much easier to implement than a
general purpose generator.
cdef array_iter_double(np.ndarray a, int axis=-1):
cdef np.flatiter it
ita = np.PyArray_IterAllButAxis(a, &axis)
cdef Py_ssize_t stride = a.strides[axis], length = a.shape[axis], i
for i in range(length):
yield <double*>(np.PyArray_ITER_DATA(it) + )[i * stride])
# TODO: Probably yield indices as well
# TODO: add faster special-cases for stride == sizeof(double)
# Use NumPy iterator API to sum all values of array with
# arbitrary number of dimensions:
cdef double s = 0, value
for value in array_iter_double(myarray):
s += value
# at this point, the contents of the array_iter_double function is
# and "s += value" simply inserted everywhere "yield" occurs in the
More information about the NumPy-Discussion