[Numpy-discussion] best way of speeding up a filtering-like algorithm

Chris Barker chris.barker at noaa.gov
Thu Mar 29 13:26:02 EDT 2018


one other note:

As a rule, using numpy array operations from Cython doesn't buy you much,
as you discovered. YOu need to use numpy arrays as n-d containers, and
write the loops yourself.

You may want to check out numba as another alternative -- it DOES optimize
numpy operations.

-CHB



On Wed, Mar 28, 2018 at 5:56 PM, Moroney, Catherine M (398E) <
Catherine.M.Moroney at jpl.nasa.gov> wrote:

> Hello,
>
>
>
> I have the following sample code (pretty simple algorithm that uses a
> rolling filter window) and am wondering what the best way is of speeding it
> up.  I tried rewriting it in Cython by pre-declaring the variables but that
> didn’t buy me a lot of time.  Then I rewrote it in Fortran (and compiled it
> with f2py) and now it’s lightning fast.  But I would still like to know if
> I could rewrite it in pure python/numpy/scipy or in Cython and get a
> similar speedup.
>
>
>
> Here is the raw Python code:
>
>
>
> def mixed_coastline_slow(nsidc, radius, count, mask=None):
>
>
>
>     nsidc_copy = numpy.copy(nsidc)
>
>
>
>     if (mask is None):
>
>         idx_coastline = numpy.where(nsidc_copy == NSIDC_COASTLINE_MIXED)
>
>     else:
>
>         idx_coastline = numpy.where(mask & (nsidc_copy ==
> NSIDC_COASTLINE_MIXED))
>
>
>
>     for (irow0, icol0) in zip(idx_coastline[0], idx_coastline[1]):
>
>
>
>         rows = ( max(irow0-radius, 0), min(irow0+radius+1,
> nsidc_copy.shape[0]) )
>
>         cols = ( max(icol0-radius, 0), min(icol0+radius+1,
> nsidc_copy.shape[1]) )
>
>         window = nsidc[rows[0]:rows[1], cols[0]:cols[1]]
>
>
>
>         npoints = numpy.where(window != NSIDC_COASTLINE_MIXED, True,
> False).sum()
>
>         nsnowice = numpy.where( (window >= NSIDC_SEAICE_LOW) & (window <=
> NSIDC_FRESHSNOW), \
>
>                                 True, False).sum()
>
>
>
>         if (100.0*nsnowice/npoints >= count):
>
>              nsidc_copy[irow0, icol0] = MISR_SEAICE_THRESHOLD
>
>
>
>     return nsidc_copy
>
>
>
> and here is my attempt at Cython-izing it:
>
>
>
> import numpy
>
> cimport numpy as cnumpy
>
> cimport cython
>
>
>
> cdef int NSIDC_SIZE  = 721
>
> cdef int NSIDC_NO_SNOW = 0
>
> cdef int NSIDC_ALL_SNOW = 100
>
> cdef int NSIDC_FRESHSNOW = 103
>
> cdef int NSIDC_PERMSNOW  = 101
>
> cdef int NSIDC_SEAICE_LOW  = 1
>
> cdef int NSIDC_SEAICE_HIGH = 100
>
> cdef int NSIDC_COASTLINE_MIXED = 252
>
> cdef int NSIDC_SUSPECT_ICE = 253
>
>
>
> cdef int MISR_SEAICE_THRESHOLD = 6
>
>
>
> def mixed_coastline(cnumpy.ndarray[cnumpy.uint8_t, ndim=2] nsidc, int
> radius, int count):
>
>
>
>      cdef int irow, icol, irow1, irow2, icol1, icol2, npoints, nsnowice
>
>      cdef cnumpy.ndarray[cnumpy.uint8_t, ndim=2] nsidc2 \
>
>         = numpy.empty(shape=(NSIDC_SIZE, NSIDC_SIZE), dtype=numpy.uint8)
>
>      cdef cnumpy.ndarray[cnumpy.uint8_t, ndim=2] window \
>
>         = numpy.empty(shape=(2*radius+1, 2*radius+1), dtype=numpy.uint8)
>
>
>
>      nsidc2 = numpy.copy(nsidc)
>
>
>
>      idx_coastline = numpy.where(nsidc2 == NSIDC_COASTLINE_MIXED)
>
>
>
>      for (irow, icol) in zip(idx_coastline[0], idx_coastline[1]):
>
>
>
>           irow1 = max(irow-radius, 0)
>
>           irow2 = min(irow+radius+1, NSIDC_SIZE)
>
>           icol1 = max(icol-radius, 0)
>
>           icol2 = min(icol+radius+1, NSIDC_SIZE)
>
>           window = nsidc[irow1:irow2, icol1:icol2]
>
>
>
>           npoints = numpy.where(window != NSIDC_COASTLINE_MIXED, True,
> False).sum()
>
>           nsnowice = numpy.where( (window >= NSIDC_SEAICE_LOW) & (window
> <= NSIDC_FRESHSNOW), \
>
>                                   True, False).sum()
>
>
>
>           if (100.0*nsnowice/npoints >= count):
>
>                nsidc2[irow, icol] = MISR_SEAICE_THRESHOLD
>
>
>
>      return nsidc2
>
>
>
> Thanks in advance for any advice!
>
>
>
> Catherine
>
>
>
> _______________________________________________
> NumPy-Discussion mailing list
> NumPy-Discussion at python.org
> https://mail.python.org/mailman/listinfo/numpy-discussion
>
>


-- 

Christopher Barker, Ph.D.
Oceanographer

Emergency Response Division
NOAA/NOS/OR&R            (206) 526-6959   voice
7600 Sand Point Way NE   (206) 526-6329   fax
Seattle, WA  98115       (206) 526-6317   main reception

Chris.Barker at noaa.gov
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20180329/9f6a1f2a/attachment-0001.html>


More information about the NumPy-Discussion mailing list