[SciPy-user] Efficient "windowing"

[Robert Kern]

Jose Luis Gomez Dans wrote:
> Hi,
> I am trying to apply a simple algorithm to a 2D matrix (an image). What I want to do is, for each pixel, choose the highest (... lowest) value in its 8- or 4-connected neighbours. I have done this using weave.inline, and using a couple of loops, but I was curious if there's some way of doing this using numpy slice syntax? My (allegedly, unelegant) attempts have been versions of the following:
> 
> b[1:-1,1:-1] = scipy.array([a[0:-2,1:-1] , a[2:,1:-1] , a[1:-1,0:-2] ,\
>       a[1:-1,2:],a[0:-2,0:-2], a[0:-2,2:], a[2:,0:-2], a[2:,2:]],'f').max()
> 
> They don't work, because the max() call at the end refers to the whole array, so you are given a constant value array, equal to the max. value of a. Using for loops is very slow when dealing with large arrays.

.max() takes an `axis` argument. It also takes an `out` argument that will help
you save you from making some large temporaries.

In [1]: from numpy import *

In [2]: n = 256

In [5]: a = random.random((n,n)).astype(float32)

In [7]: b = zeros([n,n], dtype=float32)

In [8]: neighbors = array([a[0:-2,1:-1], a[2:,1:-1], a[1:-1,0:-2], a[1:-1,2:],
a[0:-2,0:-2], a[0:-2,2:], a[2:,0:-2], a[2:,2:]])

In [9]: neighbors.shape
Out[9]: (8, 254, 254)

In [10]: neighbors.max(axis=0, out=b[1:-1,1:-1])
Out[10]:
array([[ 0.94350582,  0.94350582,  0.94350582, ...,  0.98800218,
         0.97231197,  0.61088812],
       [ 0.94350582,  0.82977545,  0.94350582, ...,  0.61088812,
         0.97231197,  0.97231197],
       [ 0.94350582,  0.94350582,  0.94350582, ...,  0.95455241,
         0.95455241,  0.95455241],
       ...,
       [ 0.92812324,  0.92812324,  0.95367759, ...,  0.92305821,
         0.92305821,  0.92305821],
       [ 0.95591497,  0.95591497,  0.95367759, ...,  0.96583498,
         0.85454881,  0.92305821],
       [ 0.92812324,  0.95591497,  0.92812324, ...,  0.96583498,
         0.92305821,  0.92305821]], dtype=float32)

In [11]: b
Out[11]:
array([[ 0.        ,  0.        ,  0.        , ...,  0.        ,
         0.        ,  0.        ],
       [ 0.        ,  0.94350582,  0.94350582, ...,  0.97231197,
         0.61088812,  0.        ],
       [ 0.        ,  0.94350582,  0.82977545, ...,  0.97231197,
         0.97231197,  0.        ],
       ...,
       [ 0.        ,  0.95591497,  0.95591497, ...,  0.85454881,
         0.92305821,  0.        ],
       [ 0.        ,  0.92812324,  0.95591497, ...,  0.92305821,
         0.92305821,  0.        ],
       [ 0.        ,  0.        ,  0.        , ...,  0.        ,
         0.        ,  0.        ]], dtype=float32)


-- 
Robert Kern

"I have come to believe that the whole world is an enigma, a harmless enigma
 that is made terrible by our own mad attempt to interpret it as though it had
 an underlying truth."
  -- Umberto Eco