i-am-spartacus i-am-spartacus - 5 months ago 9
Python Question

numpy ndarray indexing - retrieving indexes from tuple

I have asked a similar question before, but I'm still not completely sure how numpy organises its indices.

I am working with many 3D arrays, all of which are the same size. due to later operations (view as window with scipy and others) I need to slice the arrays which I am doing with a series of operations looking like this:

imFrag.append(Padded[:100, :100, :100)

which splits the arrays into 8 pieces. I am trying to obtain the original indices for each of the slices. I can do this for the whole 3d array using:

np.where(Mat == Mat)

which gives me a tuple containing x, y, and z components. is there a way I can do this?



Lets try something small and 2d:

In [102]: Mat = np.arange(16).reshape(4,4)    
In [103]: sub = Mat[2:,:2]

In [104]: Mat
array([[ 0,  1,  2,  3],
       [ 4,  5,  6,  7],
       [ 8,  9, 10, 11],
       [12, 13, 14, 15]])

In [105]: sub
array([[ 8,  9],
       [12, 13]])

Are you wondering how to tell were sub occurred in Mat?

In general practice the best thing is to hang on to the indexing tuple

In [106]: ind=(slice(2,None),slice(None,2))

In [107]: Mat[ind]
array([[ 8,  9],
       [12, 13]])

In [108]: Mat[ind] += sub   # duplicate the sub values

In [109]: Mat
array([[ 0,  1,  2,  3],
       [ 4,  5,  6,  7],
       [16, 18, 10, 11],
       [24, 26, 14, 15]])

(np.s_ is a nice way of constructing such a tuple if you prefer the colon syntax).

I think I could deduce this slicing from information in Mat and sub, e.g.

In [120]: Mat.__array_interface__['data'][0],Mat.shape,Mat.strides,Mat.itemsize
Out[120]: (169022104, (4, 4), (16, 4), 4)

In [121]: sub.__array_interface__['data'][0],sub.shape,sub.strides
Out[121]: (169022136, (2, 2), (16, 4))

but it requires some understanding of how data is stored and accessed. For simple slices like this is shouldn't be too hard. For more general ones, such as ones with steps and transposes, it would be harder.

But I never had need to do this. Hanging on to the original slicing tuple is easier. And if you do advanced indexing (which creates a copy rather than a view) hanging on to the indexing or masking is the only way.


following on your comment:

In [140]: I,J=np.where(sub==sub) 

In [141]: ind
Out[141]: (slice(2, None, None), slice(None, 2, None))

In [142]: Mat[2+I,0+J]    # 0 inplace of None for J
Out[142]: array([16, 18, 24, 26])

So yes you can use indices from sub to find corresponding elements in Mat. The use of where(sub==sub) to get all the indices bugs me a bit. meshgrid and mgrid will work but they require generating ranges. I can't off hand think of a function that takes a shape and gives corresponding I,J.


will do, but isn't pretty.