Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Pixel oversampling when nxy is even produces asymmetric PSF #7

Open
ptbrown1729 opened this issue Dec 12, 2022 · 2 comments
Open

Pixel oversampling when nxy is even produces asymmetric PSF #7

ptbrown1729 opened this issue Dec 12, 2022 · 2 comments

Comments

@ptbrown1729
Copy link

  • psfmodels version: '0.3.2'
  • Python version: 3.9
  • Operating System: Windows 10 Enterpise

Description

I have a use case where it is helpful to work with even-sized PSF arrays. The documentation in PSFmodels says not to do this, but I thought it might be useful to look into what the problem is. I am expecting the output PSF's to always be symmetric about the center, which is at pixel (nxy-1)//2. I find that this is the case except when nxy is an even number.

It looks like part of the issue is in this code

  R = new double[npx_];
  int idx = 0;
  double xi, yi;
  for (int y = -xymax_; y <= xymax_; ++y)
  {
    for (int x = -xymax_; x <= xymax_; ++x)
    {
      xi = (double)x - xp_;
      yi = (double)y - yp_;
      R[idx] = sqrt(xi * xi + yi * yi);
      ++idx;
    }
  }
}

where xymax_ = ((nx_)*p.sf - 1) / 2;

If, for example nx_ = 4 and p.sf = 3 then the coordinates would get grouped (-5, -4, -3); (-2, -1, 0); (1, 2, 3); (4, 5, ?), and the asymmetry between the 0th and 2nd grouping is clear

It seems like this should work for the even case also if the loop went from

xymin_  = -((nx_ - 1) // 2) * p.sf - (p.sf//2) ;
xymax_ = ((nx_ + 2) // 2) * p.sf - (p.sf + 1)//2;

for the even case there should be p.sf more coordinate points which are greater than zero versus less than zero.

nx_=4 and p.sf = 3 would then get grouped (-4, -3, -2); (-1, 0, 1); (2, 3, 4); (5, 6, 7)

What I Did

import psfmodels
import matplotlib.pyplot as plt

dxy = 0.065
wvl = 0.488
na = 1.3
neven = 4
nodd = 3

# without oversampling
sf = 1
psf_even1 = psfmodels.vectorial_psf(0, neven, dxy, wvl=wvl, params={"NA": na, "sf": sf}, normalize=False)
psf_odd1 = psfmodels.vectorial_psf(0, nodd, dxy, wvl=wvl, params={"NA": na, "sf": sf}, normalize=False)

# with oversampling
sf = 3
psf_even3 = psfmodels.vectorial_psf(0, neven, dxy, wvl=wvl, params={"NA": na, "sf": sf}, normalize=False)
psf_odd3 = psfmodels.vectorial_psf(0, nodd, dxy, wvl=wvl, params={"NA": na, "sf": sf}, normalize=False)


figh = plt.figure()
ax = figh.add_subplot(2, 2, 1)
ax.imshow(psf_even1[0])
ax.set_title("PSF, even sf = 1")

ax = figh.add_subplot(2, 2, 2)
ax.imshow(psf_odd1[0])
ax.set_title("PSF, odd sf = 1")

ax = figh.add_subplot(2, 2, 3)
ax.imshow(psf_even3[0])
ax.set_title("PSF, even sf = 3")

ax = figh.add_subplot(2, 2, 4)
ax.imshow(psf_odd3[0])
ax.set_title("PSF, odd sf = 3")

Figure_1

@tlambert03
Copy link
Owner

hey @ptbrown1729, I'm so sorry this one slipped through the cracks.

I suspect you've long moved on to something that works for you. You're absolutely correct that using even sized arrays wouldn't be a difficult change. that "warning" was largely inherited from some legacy code and I hadn't looked into recently.

was this issue essentially just a suggestion for modification? If so, thanks very much, I appreciate you looking into it and I can incorporate these changes

@ptbrown1729
Copy link
Author

@tlambert03 - no worries, and yes just a suggestion for modification/improvement. This is an issue that took me a long time to track down, and still seems like unexpected behavior.

Where I use these functions, I wrap them with code to throw an error if oversampling is requested with an even sized grid.

I guess there is a bigger issue with supporting even-sized PSF's at all, which is there is not a natural choice of center location. So the PSF starts to depend on what routine you use for convolving an image with it.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants