The task:

I am trying to interpolate a vector field on a regular grid, i.e.:

The issue:

I am using the RegularGridInterpolator from scipy to do this. However, it seems that the resulting vector field is turned with respect to the original:

Anyone knows why?

Python code to reproduce example:

from scipy.interpolate import RegularGridInterpolator
import matplotlib.pyplot as plt
import numpy as np

# ORIGINAL 
# Number of points (NxN)
N = 50
# Boundaries
ymin = -2.; ymax = 2.
xmin = -2.; xmax = 2.
# Create Meshgrid
x = np.linspace(xmin,xmax, N)
y = np.linspace(ymin,ymax, N)
xx, yy = np.meshgrid(x, y)
# Vector Field
Fx  = np.cos(xx + 2*yy)
Fy  = np.sin(xx - 2*yy)
# Plot vector field
fig, ax = plt.subplots()
ax.quiver(x, y, Fx, Fy)
plt.title("Original")
plt.show()

# REDUCED 
# Number of points (NxN)
N = 10
# Boundaries
ymin = -2.; ymax = 2.
xmin = -2.; xmax = 2.
# Create Meshgrid
x = np.linspace(xmin,xmax, N)
y = np.linspace(ymin,ymax, N)
xx, yy = np.meshgrid(x, y)
# Vector Field
Fx  = np.cos(xx + 2*yy)
Fy  = np.sin(xx - 2*yy)
# Plot vector field
fig, ax = plt.subplots()
ax.quiver(x, y, Fx, Fy)
plt.title("Reduced")
plt.show()

# INTERPOLATED VERSION BASED ON REDUCED

# Iterpolate
my_interpolating_function_x = RegularGridInterpolator((x, y), Fx)
my_interpolating_function_y = RegularGridInterpolator((x, y), Fy)
# Create Meshgrid
N = 50
x = np.linspace(xmin,xmax, N)
y = np.linspace(ymin,ymax, N)
grid = np.meshgrid(x, y)
new_points = np.vstack(list(map(np.ravel, grid))).T
# Interpolate
F_x_inter = my_interpolating_function_x(new_points)
F_y_inter = my_interpolating_function_y(new_points)
# reshape
F_x_inter = np.reshape(F_x_inter,(50,50))
F_y_inter = np.reshape(F_y_inter,(50,50))
#plot
fig, ax = plt.subplots()
ax.quiver(x, y, F_x_inter, F_y_inter)
plt.title("Interpolated")
plt.show()