Open schlafly opened 10 months ago
I think that's another way of showing the same thing using Eddie's approach in his gist: for odd oversampling, the Scipy and Astropy convolutions give (visually) the same things, whereas even oversampling yields to sharper differences.
In my head I was guessing that astropy is replacing the even kernel with some half-pixel shifted odd kernel, and so replacing the top-hat with something where the wings ramp linearly down or something, which sounds sane but for a discontinuous kernel like this one would be problematic. But I didn't push hard enough to actually find the code that was doing that and may be off-base. That would make the effective kernel a bit wider than it should be for even oversampling and have an effect like the one we see here.
Right - as @obi-wan76 pointed out, Astropy Boxkernel2D returns an odd-size kernel regardless of the width specified, so that's one difference with the Scipy approach used above.
With that said, I don't measure any difference in the centroids of the images I posted above...
Ah, sorry, I'm not sure I can see obi-wan76's message, thanks. Yes, I didn't test centroid differences for my test image, only for the real Roman images, but you can imagine that I was surprised when I learned that I could reproduce the input centroids well for Roman only for odd-sized convolutions! Note that the systematics were small (~0.01 pix), but we obviously shouldn't be limited by this particular systematic. The real Roman PSFs are of course a bit more complicated than these test images.
I've been doing some Roman simulations and comparing the input simulated positions with the measured positions using photutils and the webbpsf gridded_psf routines. I eventually noticed that I got good performance if I picked odd oversampling factors and bad performance (~0.01 systematic pixel phase errors) if I used even oversampling factors.
I think I understand the root cause to be the convolution over the oversampled pixel around here https://github.com/spacetelescope/webbpsf/blob/develop/webbpsf/gridded_library.py#L260 For the odd case, this is easy and just sums over the N oversampled pixels. For the even case, it's a bit more complicated since if you just sum over the N oversampled pixels you introduce a centroid shift of half a pixel, since you can't center the even array right on a pixel. The astropy convolution machinery is dealing with this somehow, but I haven't looked to see exactly what they're doing. But I don't think they're doing a great job for kernels like this one that are discontinuous.
If I manually replace the PSFs generated by gridded_library with my own versions where I do a different convolution (a la the convolve_scipy here https://gist.github.com/schlafly/a50a27411bdc706c69e093d083a4b6db), then I get good performance for both odd and even oversampling rates.
That gist has a little demo that produces output like
showing that the difference between the astropy convolution and "truth" has this oscillating pattern where oversample 4 is worse than oversample 3, etc.. Meanwhile the scipy version I cooked up is better but starts oscillating around oversample 10. Note that the two versions agree closely for odd oversampling smaller than 7; they're doing the same thing. But the even oversampling is a bit different.
I worry I might not have fully tracked down everything there; I feel like for my intended band-limited 'truth' image I should have seen much better convergence than I did---but I think this is still a good illustration of the issue.
I think I might recommend doing something different for the convolution in the even sampling case, along the lines of my convolve_scipy routine, but you may have better ideas.