3

Consider example matrix array:

[[0 1 2 1 0]
 [1 1 2 1 0]
 [0 1 0 0 0]
 [1 2 1 0 0]
 [1 2 2 3 2]]

What I need to do:

  1. find maxima in every row
  2. select smaller surrounding of the maxima from every row (3 values in this case)
  3. paste the surrounding of the maxima into new array (narrower)

For the example above, the result is:

[[ 1.  2.  1.]
 [ 1.  2.  1.]
 [ 0.  1.  0.]
 [ 1.  2.  1.]
 [ 2.  3.  2.]]

My current working code:

import numpy as np

A = np.array([
    [0, 1, 2, 1, 0],
    [1, 1, 2, 1, 0],
    [0, 1, 0, 0, 0],
    [1, 2, 1, 0, 0],
    [1, 2, 2, 3, 2],
])

b = A.argmax(axis=1)

C = np.zeros((len(A), 3))
for idx, loc, row in zip(range(len(A)), b, A):
    print(idx, loc, row)
    C[idx] = row[loc-1:loc+2]

print(C)

My question:

How to get rid of the for loop and replace it with some cheaper numpy operation?

Note:

This algorithm is for straightening broken "lines" in video stream frames with thousands of rows.

matousc
  • 3,698
  • 10
  • 40
  • 65

2 Answers2

1

Approach #1

We can have a vectorized solution based on setting up sliding windows and then indexing into those with b-offsetted indices to get desired output. We can leverage np.lib.stride_tricks.as_strided based scikit-image's view_as_windows to get sliding windows. More info on use of as_strided based view_as_windows.

The implementation would be -

from skimage.util.shape import view_as_windows

L = 3 # window length
w = view_as_windows(A,(1,L))[...,0,:]
Cout = w[np.arange(len(b)),b-L//2]

Being a view-based method, this has the advantage of being memory-efficient and hence good on performance too.

Approach #2

Alternatively, a one-liner by creating all those indices with outer-addition would be -

A[np.arange(len(b))[:,None],b[:,None] + np.arange(-(L//2),L//2+1)]
Divakar
  • 218,885
  • 19
  • 262
  • 358
  • Will this also work if the maximum of a row is at index 0 or -1? – Joe Jul 04 '19 at 09:15
  • 1
    @Joe OP need to clarify what's the intended output as it would go out of bounds and we might end up with variable number of elements per row. – Divakar Jul 04 '19 at 09:21
1

This works by making and array with all the desired indices, but somehow using that directly on A results in a 3D array, hence the subsequent indexing... Probably not optimal, but definitely another way of doing it!

import numpy as np

A = np.array([
    [0, 1, 2, 1, 0],
    [1, 1, 2, 1, 0],
    [0, 1, 0, 0, 0],
    [1, 2, 1, 0, 0],
    [1, 2, 2, 3, 2],
])

b = A.argmax(axis = 1).reshape(-1, 1)
index = b + np.arange(-1,2,1).reshape(1, -1)
A[:,index][np.arange(b.size),np.arange(b.size)]
holmrenser
  • 425
  • 3
  • 11
  • This will not work if the maximum of a row is at index 0 or -1. Not sure how the OP wants to handle this case. – Joe Jul 04 '19 at 09:23