Retrieve image coordinates (X, Y) from geographic coordinates in NetCDF file using Numpy

I’m trying to retrieve image coordinates (X, Y) from geographic coordinates in NetCDF file (Sentinel 3 data) using Numpy. The problem is that the final pixel is too close to the desired pixel but not the closest to it (about one or two pixel far away from the closest in x, y or both x,y). I’m using Numpy’s norm(), min() and where() function. The code is presented below:

from netCDF4 import Dataset
import numpy as np

# Read the NetCDF file with the geographic coordinates information
nc_coord = Dataset('coords.nc', 'r')

lat = nc_coord.variables['latitude'][:]
# Python buffer object pointing to the start of the array’s data
lat = lat.data
lon = nc_coord.variables['longitude'][:]
lon = lon.data

# define your point:
target_lon, target_lat = 23.779675, 35.33786
# lat and lon are x,y matrices, we will add another dimension z in order to stack them
lat = lat[:, :, np.newaxis]
lon = lon[:, :, np.newaxis]
grid = np.concatenate([lat, lon], axis=2)
vector = np.array([target_lat, target_lon]).reshape(1, 1, -1)
subtraction = vector - grid

# vector subtraction
dist = np.linalg.norm(subtraction, axis=2)

# query where in the 2D matrix is the closest point to your points of interest (?)
result =  np.where(dist == dist.min())

# this is your lat/lon geo point translated to a x, y position in the image matrix
target_x_y = result[0][0], result[1][0]

Can you figure out what is happening?

geolocation = Dataset(fname)
points = np.array([geolocation.variables['lon'][:],
                   geolocation.variables['lat'][:]
                  ])
# Exact point
point = np.array([3.60829699, 44.41841389])
dist = np.linalg.norm(point[:, None, None] - points, axis=0)
i, j = np.where(dist == dist.min())
# Should evaluate to True
print(np.allclose(point, points[:, i, j].squeeze()))
print(dist.min())
# Not Exact point, actually closest to neighbouring pixel
point = np.array([3.605, 44.421])
dist = np.linalg.norm(point[:, None, None] - points, axis=0)
i, j = np.where(dist == dist.min())
# Should evaluate to False, as test point not in arrays
print(np.allclose(point, points[:, i, j].squeeze()))
print(dist.min())