Unexpected behavior from `zoom` spline interpolation + fix

[bommaritom]

Hi Diyu,

I'm currently looking into downsampling methods for the multires feature, including with arrays of odd size.

While doing that, I noticed an odd behavior with the scipy.ndimage zoom function. Namely, if zoom=0.5, the function does not behave as expected. I ran the zoom function on a series of 1D arrays with zoom=0.5 and recorded the input and output array lengths. On the left side are listed the original array sizes, and on the right are listed the output array sizes:

• 3 --> 2
• 4 --> 2
• 5 --> 2
• 6 --> 3
• 7 --> 4
• 8 --> 4
• 9 --> 4
• 10 --> 5
• 11 --> 6
• 12 --> 6
• 13 --> 6
• 14 --> 7

and so on. This is unexpected; I would expect

• m --> ceil(m/2)

but this is not the case for m = 5, 9, 13, 17, etc as seen above.

I'm not exactly sure how the zoom function works, but this probably has something to do with the way floating-point values are handled. In any case, it seems like we don't want the function to be doing this.

I was able to "hack" around this by setting zoom=0.5000001. In any case if we want to continue using the zoom method to downsample, it looks like this is a viable fix for this issue.

[bommaritom]

Here is the Python code I ran; I couldn't figure out how to attach a .py file:

import numpy as np
from scipy.ndimage import zoom

print('original size --> size after zoom')
for size in np.arange(50)+100.:
  data = np.arange(size)
  zoom_data_grid = zoom(data, 0.5000001)
  zoom_data_no_grid = zoom(data, 0.5000001, grid_mode=True)

  print(str(np.size(data)) + '-->' + str(np.size(zoom_data_no_grid)))

[dyang37]

Hi Marco, I'm not sure about the exact behavior of zoom=0.5 when array size contains odd number, but in the case of multi-resolution, the zoom function always deals with array size with all even numbers, so this should be fine for multi-resolution. In the case where array size contains odd numbers, we simply stop going to the lower resolution space, and the case zoom=0.5 is never used in that case: https://github.com/cabouman/mbircone/blob/1441b039aba62b957e0ce7f1c2c3434ac973ba2d/mbircone/interface_cy_c.pyx#L282-L286

Please also see this PR for more detailed discussion: https://github.com/cabouman/mbircone/pull/102#issue-1466804750

[dyang37]

Closing this issue since we already discussed the zoom function in a previous PR.