CODEWITHSUNDEEP
pythonmatlabnumpyinterpolation
Wkjung
Wkjung

Reputation: 23

Translating float index interpolation from MATLAB to Python

For example, I have a index array

ax = [0, 0.2, 2] #start from index 0: python

and matrix I

I= 
10    20    30    40    50
10    20    30    40    50
10    20    30    40    50
10    20    30    40    50
10    20    30    40    50

In MATLAB, by running this code

[gx, gy] = meshgrid([1,1.2,3], [1,1.2,3]);
I = [10:10:50];
I = vertcat(I,I,I,I,I)
SI = interp2(I,gx,gy,'bilinear');

The resulting SI is

SI =
10    12    30
10    12    30
10    12    30

I tried to do the same interpolation in Python, using NumPy. I first interpolate row-wise, then column-wise

import numpy as np
ax = np.array([0.0, 0.2, 2.0])
ay = np.array([0.0, 0.2, 2.0])
I = np.array([[10,20,30,40,50]])
I = np.concatenate((I,I,I,I,I), axis=0)
r_idx = np.arange(1, I.shape[0]+1)
c_idx = np.arange(1, I.shape[1]+1)
  
I_row = np.transpose(np.array([np.interp(ax, r_idx, I[:,x]) for x in range(0,I.shape[0])]))
I_col = np.array([np.interp(ay, c_idx, I_row[y,:]) for y in range(0, I_row.shape[0])])

SI = I_col

However, the resulting SI is

SI =
10    10    20
10    10    20
10    10    20

Why are my results using Python different from those using MATLAB?

Upvotes: 1

Views: 112

Answers (1)

lpeak
lpeak

Reputation: 181

It seems that you over-corrected yourself by passing from MATLAB to Python, as shown by your first code excerpt.

ax = [0, 0.2, 2] #start from index 0: python

In numpy logic this sequence does not represents the indexes but the coordinate for the function to interpolate. Since you already take care of incrementing the coordinate to be compatible with matlab here:

r_idx = np.arange(1, I.shape[0]+1)
c_idx = np.arange(1, I.shape[1]+1)

You can reuse the same interpolation coordinate that you used in Matlab:

ax = [1,1.2,3]

Full code:

import numpy as np
ax = np.array([1.0, 1.2, 3.0])
ay = np.array([1.0, 1.2, 3.0])
I = np.array([[10,20,30,40,50]])
I = np.concatenate((I,I,I,I,I), axis=0)
r_idx = np.arange(1, I.shape[0]+1)
c_idx = np.arange(1, I.shape[1]+1)

I_row = np.transpose(np.array([np.interp(ax, r_idx, I[:,x]) for x in range(0,I.shape[
0])]))
I_col = np.array([np.interp(ay, c_idx, I_row[y,:]) for y in range(0, I_row.shape[0])]
)

SI = I_col

and result:

array([[10., 12., 30.],
       [10., 12., 30.],
       [10., 12., 30.]])

Explanation about the bug

Since ax represented coordinates the first two values 0.0 and 0.2 were before the first coordinate of r_idx. According to the documentation, the interpolation will default to I[:,x][0].

Upvotes: 1

Related Questions