https://docs.xarray.dev/en/stable/user-guide/indexing.html#vectorized-indexing / https://stackoverflow.com/questions/66027917/python-xarray-vectorized-indexing and show how this can be used for non-orthogonal indexing (orthogonal = each dimension separately, and you get intersection), like nearest point sampling