Python numpy grid transformation using universal functions

I can't really try the code since I don't have your matrices, but something like this should do the job.

First, instead of declaring conv as a function, get the whole altitude projection for all your data:

conv = np.round(alt / 500.).astype(int)

Using np.round, the numpys version of round, it rounds all the elements of the matrix by vectorizing operations in C, and thus, you get a new array very quickly (at C speed). The following line aligns the altitudes to start in 0, by shifting all the array by its minimum value (in your case, -20):

conv -= conv.min()

the line above would transform your altitude matrix from [-20, 200] to [0, 220] (better for indexing).

With that, interpolation can be done easily by getting multidimensional indices:

t, z, y, x = np.indices(temp.shape)

the vectors above contain all the indices needed to index your original matrix. You can then create the new matrix by doing:

new_matrix[t, conv[t, z, y, x], y, x] = temp[t, z, y, x]

without any loop at all.

Let me know if it works. It might give you some erros since is hard for me to test it without data, but it should do the job.

The following toy example works fine:

A = np.random.randn(3,4,5) # Random 3x4x5 matrix -- your temp matrixB = np.random.randint(0, 10, 3*4*5).reshape(3,4,5) # your conv matrix with altitudes from 0 to 9C = np.zeros((3,10,5)) # your new matrixz, y, x = np.indices(A.shape)C[z, B[z, y, x], x] = A[z, y, x]

C contains your results by altitude.