i-am-spartacus - 10 months ago 34

Python Question

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?

Thanks!

Answer

Lets try something small and 2d:

```
In [102]: Mat = np.arange(16).reshape(4,4)
In [103]: sub = Mat[2:,:2]
In [104]: Mat
Out[104]:
array([[ 0, 1, 2, 3],
[ 4, 5, 6, 7],
[ 8, 9, 10, 11],
[12, 13, 14, 15]])
In [105]: sub
Out[105]:
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]
Out[107]:
array([[ 8, 9],
[12, 13]])
In [108]: Mat[ind] += sub # duplicate the sub values
In [109]: Mat
Out[109]:
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`

.

```
I,J=np.array(list(np.ndindex(sub.shape))).T
```

will do, but isn't pretty.