When you construct the delta matrix here, the coordinate is actually in yx not xy order. However, it does not affect the final results because the sampling region is from -r to r in both x and y dimension. Therefore, the sampled correlation matrix is a transposed version of what it should be.
I just want to raise the issue in case anyone wants to actually use the coordinates to do something else that requires the order to be correct.