import numpy as np
from scipy.interpolate import LinearNDInterpolator, griddata

np.set_printoptions(precision=2, suppress=True)

points = np.array([
    [0.0, 0.0],
    [1.0, 0.0],
    [0.0, 1.0],
    [1.0, 1.0],
    [0.5, 0.2],
])
values = points[:, 0] + 2 * points[:, 1]

query_points = np.array([
    [0.25, 0.25],
    [0.75, 0.50],
    [1.20, 0.40],
])

linear = griddata(points, values, query_points, method="linear")
nearest = griddata(points, values, query_points, method="nearest")
with_fill = griddata(points, values, query_points, method="linear", fill_value=-1.0)

reusable_linear = LinearNDInterpolator(points, values)
inside_expected = query_points[:2, 0] + 2 * query_points[:2, 1]

print("query points:")
for row in query_points:
    print(f"  ({row[0]:.2f}, {row[1]:.2f})")
print("linear griddata:", np.round(linear, 2))
print("nearest griddata:", np.round(nearest, 2))
print("linear fill_value:", np.round(with_fill, 2))
print("callable linear:", np.round(reusable_linear(query_points[:2]), 2))
print("expected inside:", np.round(inside_expected, 2))
print("inside match:", np.allclose(linear[:2], inside_expected))
print("outside is nan:", np.isnan(linear[2]))