function g = resample(f, sx, sy)
%% Scales an image using nearest neighbor interpolation.
% f - the input image
% sx - the scaling factor in x-direction
% sy - the scaling factor in y-direction
g = zeros(size(f, 1) * sx, size(f, 2) * sy, class(f));
for i = 1:size(g, 1)
for j = 1:size(g, 2)
fi = 1 + round((i - 0.5) / sx - 0.5);
fj = 1 + round((j - 0.5) / sy - 0.5);
g(i, j) = f(fi, fj);
end
end
end
How to find this arithmetic expression? Step by step.
Convert a 1-based index i to a 0-based index i’
$$i' = Z(i) = i - 1$$$$i = Z^{-1}(i') = 1 + i'$$Convert a 0-based index i’ to a Cartesian coordinate
The first pixel has coordinate 0.5, the next pixel 1.5 and so on. This provides consistent spacing between pixels and follows the standard image representation on screens and other media.
$$x = C(i') = 0.5 + i'$$$$i' = C^{-1}(x) = \mathrm{round}(x - 0.5)$$Scale
$$x' = S(x) = s\cdot x$$$$x = S(x') = x' / s$$Chaining the transformations together
Now, given an 1-based index i_f from the input image, how to calculate the corresponding 1-based index i_g in the output image? By chaining the transformations together.
$$i_g = (Z^{-1} \circ C^{-1} \circ S \circ C \circ Z)(i_f)$$By inverting the composite function, we obtain the sought-after function
$$i_f = (Z^{-1} \circ C^{-1} \circ S^{-1} \circ C \circ Z)(i_g)$$We can now plug in the function definitions
$$ \begin{aligned} i_f &= (Z^{-1} \circ C^{-1} \circ S^{-1} \circ C \circ Z)(i_g) \\ &= (Z^{-1} \circ C^{-1} \circ S^{-1} \circ C)(i_g - 1) \\ &= (Z^{-1} \circ C^{-1} \circ S^{-1})(i_g - 0.5) \\ &= (Z^{-1} \circ C^{-1})(i_g / s - 0.5) \\ &= Z^{-1}(\mathrm{round}(i_g / s - 0.5) - 0.5) \\ &= 1 + \mathrm{round}(i_g / s - 0.5) - 0.5 \end{aligned} $$Compare this expression with the algorithm.