# Why does this numpy construct work and is there a more correct way of doing it

I’m using numpy to extract faces from tetrahdra defined by vertex indices. I have an initial array defining the tehrahedra toplogy

```tetrahedra = np.array([[0, 1, 2, 3], [1, 2, 3, 4], [2, 3, 4, 5]])
```

For each tetrahedra I identify the faces using mask array

```face1 = [True, True, True, False]
face2 = [True, True, False, True]
face3 = [False, True, True, True]
face4 = [True, False, True, True]
```

And I find the following numpy expression yields face defnitions for each tetrahedra

```faces = tetrahedra[:,np.reshape(np.r_[tetrahedra[0][face1 ],tetrahedra[0][face2 ],tetrahedra[0][face3 ],tetrahedra[0][face4 ]], (-1,3))]
```

EDIT : Thanks to @npaulj I now see that this only appears to work because tetrahedra[0] in the indexing notation is actually [0, 1, 2, 3]. This is better expressed by changing the boolean masks to direct index masks as follows,

```mask1_ = np.array([0, 1, 2])
mask2_ = np.array([0, 1, 3])
mask3_ = np.array([1, 2, 3])
mask4_ = np.array([0, 2, 3])
```

and then updating the expression to

```faces = tetrahedra[:,np.reshape(np.r_[face1_, face2_, face3_, face4_], (-1,3))]
```

Now my question is, how is this actually working and is there a preferred/faster way of doing this operation? The output is shown below

THanks in advance for any help with this. I’m tempted to just live with it since it appears to work, but I can’t fugre out why it works which worries me…..

```[[[0 1 2]
[0 1 3]
[1 2 3]
[0 2 3]]
[[1 2 3]
[1 2 4]
[2 3 4]
[1 3 4]]
[[2 3 4]
[2 3 5]
[3 4 5]
[2 4 5]]
```

]

edit Cleaned up version is,

```face_masks = np.array([[0, 1, 2], [0, 1, 3], [1, 2, 3],[0, 2, 3]])
faces = tetrahedra[:,face_masks]
```

## Answer

All the intermediate steps:

```In [77]: tetrahedra = np.array([[0, 1, 2, 3], [1, 2, 3, 4], [2, 3, 4, 5]])
...: face1 = [True, True, True, False]
...: face2 = [True, True, False, True]
...: face3 = [False, True, True, True]
...: face4 = [True, False, True, True]

In [78]: tetrahedra[0]
Out[78]: array([0, 1, 2, 3])

In [79]: tetrahedra[0][face1]
Out[79]: array([0, 1, 2])
```

`r_` concatenates these 4 selections into one array:

```In [80]: np.r_[tetrahedra[0][face1 ],tetrahedra[0][face2 ],tetrahedra[0][face3 ],tetrahedra[0][face4 ]]
Out[80]: array([0, 1, 2, 0, 1, 3, 1, 2, 3, 0, 2, 3])

In [81]: np.reshape(np.r_[tetrahedra[0][face1 ],tetrahedra[0][face2 ],tetrahedra[0][face3
...:  ],tetrahedra[0][face4 ]], (-1,3))
Out[81]:
array([[0, 1, 2],
[0, 1, 3],
[1, 2, 3],
[0, 2, 3]])
```

And finally just use that to index the columns dimension of `tetrahedra`. Index a (3,4) with this (4,3) produces a (3,4,3).

`vstack` could also be used to concatenate the selections, but that’s a minor change:

```In [82]: np.vstack([tetrahedra[0][face1 ],tetrahedra[0][face2 ],tetrahedra[0][face3 ],tet
...: rahedra[0][face4 ]])
Out[82]:
array([[0, 1, 2],
[0, 1, 3],
[1, 2, 3],
[0, 2, 3]])
```

## edit

Or if you don’t want to count on `tetrahedra[0]` being [0,1,2,3], you could just seek the indices of the `True` elements:

```In [106]: np.nonzero([face1,face2,face3,face4])
Out[106]:
(array([0, 0, 0, 1, 1, 1, 2, 2, 2, 3, 3, 3]),
array([0, 1, 2, 0, 1, 3, 1, 2, 3, 0, 2, 3]))
```

The 2nd element of that tuple `[1]` is the same indices you with `r_`.

```In [122]: idx = np.nonzero([face1,face2,face3,face4])
In [123]: idx[1]
Out[123]: array([0, 1, 2, 0, 1, 3, 1, 2, 3, 0, 2, 3])
In [124]: tetrahedra[np.arange(3)[:,None],idx[1]]
Out[124]:
array([[0, 1, 2, 0, 1, 3, 1, 2, 3, 0, 2, 3],
[1, 2, 3, 1, 2, 4, 2, 3, 4, 1, 3, 4],
[2, 3, 4, 2, 3, 5, 3, 4, 5, 2, 4, 5]])
In [125]: tetrahedra[np.arange(3)[:,None],idx[1]].reshape(3,4,3)
Out[125]:
array([[[0, 1, 2],
[0, 1, 3],
[1, 2, 3],
[0, 2, 3]],

[[1, 2, 3],
[1, 2, 4],
[2, 3, 4],
[1, 3, 4]],

[[2, 3, 4],
[2, 3, 5],
[3, 4, 5],
[2, 4, 5]]])
```