getting pointer to fused type function based on np.ndarray dtype

[Thouis (Ray) Jones]

Hello,

I've been working on a Cython implementation of region labeling for
scipy.ndimage. I have a working version, but I've been having trouble
with fused types.

My initial implementation specialized the label() function based on
input and output. Since there are 10 different types for each, the
size of the .c file was quite large. That version is here:
https://github.com/thouis/scipy/blob/048bf06a9098e3c807917d043d3a3a5250d9140f/scipy/ndimage/src/ni_label.pyx

I rewrote it so that the specialized portions were in subfunctions.
However, this required finding a function pointer to the specialized
functions. Unfortunately, I couldn't find a clean way of doing that.
The result is here:
https://github.com/thouis/scipy/blob/048bf06a9098e3c807917d043d3a3a5250d9140f/scipy/ndimage/src/ni_label.pyx

Basically, I ended up writing something like:
if input.dtype == np.int8:
NONZERO = <bint (*)(void *, int) nogil> _NONZERO[np.int8_t]
elif input.dtype == np.uint8:
NONZERO = <bint (*)(void *, int) nogil> _NONZERO[np.uint8_t]
elif input.dtype == np.int16:
NONZERO = <bint (*)(void *, int) nogil> _NONZERO[np.int16_t]
elif input.dtype == np.uint16:
....

Which I find less than optimal. Is there a better way to do this?

Thanks
Thouis Jones

[mark florisson]

Indeed, if both argument types have 10 specializations, you'll get 10
* 10 specializations. There's usually several options, one is to
separate out the main algorithm, but it's still a lot of
specializations. If the specializations for the types aren't dependent
on each other, you could separate out the code dealing with either
argument in two separate functions. Here you don't want to hold the
GIL, so automatic dispatching or python callables aren't going to
work. So I think in this particular case, you're stuck with function
pointers. You could retrieve them in a slightly different way, but it
doesn't win any prize either, but you can create a wrapper def
function that dispatches on the argument and returns a pointer to a
specialized wrapped function:

cimport numpy as np
import numpy as np
cimport cython

cdef extern from *:
ctypedef int Py_intptr_t

ctypedef fused data_t:
np.int8_t
np.uint8_t
np.int16_t
np.uint16_t
np.int32_t
np.uint32_t
np.int64_t
np.uint64_t
np.float32_t
np.float64_t

ctypedef void (*func_t)(np.int32_t)

cdef void func(data_t arg):
print cython.typeof(arg), arg

def get_func(np.ndarray[data_t] a):
return <Py_intptr_t> func[data_t]

cdef func_t f = <func_t> <void *> <Py_intptr_t> get_func(np.arange(10,
dtype=np.int32))
f(20)

[mark florisson]

Hm, I see that this kind of indexing is broken, I'll push a fix.

[mark florisson]

(Fixed now in the latest Cython master on github. Would test on
sage.math, but it's been down for the last couple of weeks.)