I have an array Gamma, Gamma.shape=(20,7,90,144), that represents the dimensions (t,z,y,x) and whose respective indices I denote by [l,k,j,i]. For each (t,y,x), I want to find the lowest value of k such that Gamma[l,k,j,i] > g_crit, where g_crit is some constant. Denote this set of lowest k values as k_low; I've found I can do this nicely with
k_low = np.argmax(Gamma > g_crit, axis=1)
There is another array levs, levs.shape=(7,), corresponding to the z dimension as well. I am ultimately trying to create the array levs_low, such that levs_low.shape=(20,90,144) and levs_low[l,j,i]=levs[k_low[l,j,i]]. I'm stuck on this step. Any ideas please? Thanks.
For your case it seems that np.take() is a good choice:
levs_low = np.take( levs, k_low )
OBS: take seems to give an optimized performance, check this question
This should do the trick:
levs_low=levs[k_low]
>>> Gamma=np.random.rand(20,7,90,144)
>>> k_low = np.argmax(Gamma > .3, axis=1)
>>> levs=np.random.rand(7)
>>> levs_low=levs[k_low]
>>> levs_low.shape
(20, 90, 144)
For a small example:
>>> g=np.random.randint(0,5,(4,4))
>>> g
array([[2, 0, 2, 2],
[2, 0, 1, 0],
[3, 3, 0, 3],
[3, 0, 4, 4]])
>>> k=np.arange(5)*-1
>>> k
array([ 0, -1, -2, -3, -4])
>>> k[g] #Uses indices of g to select values from k. Also same as np.take(k,g)
array([[-2, 0, -2, -2],
[-2, 0, -1, 0],
[-3, -3, 0, -3],
[-3, 0, -4, -4]])
@Saullo Castro's answer is interesting. Strange that there would be a ~5x difference in timings between fancy indexing and np.take.
%timeit levs[k_low]
100 loops, best of 3: 2.3 ms per loop
%timeit np.take( levs, k_low )
1000 loops, best of 3: 439 us per loop
In [33]: np.all(levs[k_low]==np.take(levs,k_low))
Out[33]: True
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With