cython

1
# tag: numpy
2
# You can ignore the previous line.
3
# It's for internal testing of the cython documentation.
4

5
import numpy as np
6

7
# "cimport" is used to import special compile-time information
8
# about the numpy module (this is stored in a file numpy.pxd which is
9
# distributed with Numpy).
10
# Here we've used the name "cnp" to make it easier to understand what
11
# comes from the cimported module and what comes from the imported module,
12
# however you can use the same name for both if you wish.
13
cimport numpy as cnp
14

15
# It's necessary to call "import_array" if you use any part of the
16
# numpy PyArray_* API. From Cython 3, accessing attributes like
17
# ".shape" on a typed Numpy array use this API. Therefore we recommend
18
# always calling "import_array" whenever you "cimport numpy"
19
cnp.import_array()
20

21
# We now need to fix a datatype for our arrays. I've used the variable
22
# DTYPE for this, which is assigned to the usual NumPy runtime
23
# type info object.
24
DTYPE = np.int64
25

26
# "ctypedef" assigns a corresponding compile-time type to DTYPE_t. For
27
# every type in the numpy module there's a corresponding compile-time
28
# type with a _t-suffix.
29
ctypedef cnp.int64_t DTYPE_t
30

31
# "def" can type its arguments but not have a return type. The type of the
32
# arguments for a "def" function is checked at run-time when entering the
33
# function.
34
#
35
# The arrays f, g and h is typed as "np.ndarray" instances. The only effect
36
# this has is to a) insert checks that the function arguments really are
37
# NumPy arrays, and b) make some attribute access like f.shape[0] much
38
# more efficient. (In this example this doesn't matter though.)
39
def naive_convolve(cnp.ndarray f, cnp.ndarray g):
40
    if g.shape[0] % 2 != 1 or g.shape[1] % 2 != 1:
41
        raise ValueError("Only odd dimensions on filter supported")
42
    assert f.dtype == DTYPE and g.dtype == DTYPE
43

44
    # The "cdef" keyword is also used within functions to type variables. It
45
    # can only be used at the top indentation level (there are non-trivial
46
    # problems with allowing them in other places, though we'd love to see
47
    # good and thought out proposals for it).
48
    #
49
    # For the indices, the "int" type is used. This corresponds to a C int,
50
    # other C types (like "unsigned int") could have been used instead.
51
    # Purists could use "Py_ssize_t" which is the proper Python type for
52
    # array indices.
53
    cdef int vmax = f.shape[0]
54
    cdef int wmax = f.shape[1]
55
    cdef int smax = g.shape[0]
56
    cdef int tmax = g.shape[1]
57
    cdef int smid = smax // 2
58
    cdef int tmid = tmax // 2
59
    cdef int xmax = vmax + 2 * smid
60
    cdef int ymax = wmax + 2 * tmid
61
    cdef cnp.ndarray h = np.zeros([xmax, ymax], dtype=DTYPE)
62
    cdef int x, y, s, t, v, w
63

64
    # It is very important to type ALL your variables. You do not get any
65
    # warnings if not, only much slower code (they are implicitly typed as
66
    # Python objects).
67
    cdef int s_from, s_to, t_from, t_to
68

69
    # For the value variable, we want to use the same data type as is
70
    # stored in the array, so we use "DTYPE_t" as defined above.
71
    # NB! An important side-effect of this is that if "value" overflows its
72
    # datatype size, it will simply wrap around like in C, rather than raise
73
    # an error like in Python.
74
    cdef DTYPE_t value
75
    for x in range(xmax):
76
        for y in range(ymax):
77
            s_from = max(smid - x, -smid)
78
            s_to = min((xmax - x) - smid, smid + 1)
79
            t_from = max(tmid - y, -tmid)
80
            t_to = min((ymax - y) - tmid, tmid + 1)
81
            value = 0
82
            for s in range(s_from, s_to):
83
                for t in range(t_from, t_to):
84
                    v = x - smid + s
85
                    w = y - tmid + t
86
                    value += g[smid - s, tmid - t] * f[v, w]
87
            h[x, y] = value
88
    return h
89