Added fix for #16245

[manodeep]

Description

This pull request fixes #16245

• By checking this box, the PR author has requested that maintainers do NOT use the "Squash and Merge" button. Maintainers should respect this when possible; however, the final decision is at the discretion of the maintainer that merges the PR.

[github-actions]

Thank you for your contribution to Astropy! 🌌 This checklist is meant to remind the package maintainers who will review this pull request of some common things to look for.

• Do the proposed changes actually accomplish desired goals?
• Do the proposed changes follow the Astropy coding guidelines?
• Are tests added/updated as required? If so, do they follow the Astropy testing guidelines?
• Are docs added/updated as required? If so, do they follow the Astropy documentation guidelines?
• Is rebase and/or squash necessary? If so, please provide the author with appropriate instructions. Also see instructions for rebase and squash.
• Did the CI pass? If no, are the failures related? If you need to run daily and weekly cron jobs as part of the PR, please apply the "Extra CI" label. Codestyle issues can be fixed by the bot.
• Is a change log needed? If yes, did the change log check pass? If no, add the "no-changelog-entry-needed" label. If this is a manual backport, use the "skip-changelog-checks" label unless special changelog handling is necessary.
• Is this a big PR that makes a "What's new?" entry worthwhile and if so, is (1) a "what's new" entry included in this PR and (2) the "whatsnew-needed" label applied?
• At the time of adding the milestone, if the milestone set requires a backport to release branch(es), apply the appropriate "backport-X.Y.x" label(s) before merge.

[manodeep]

Ideally the test python code (in #16245) should also be included with this fix.

[astrofrog]

This fix looks good but I can investigate the test failure tomorrow. Looks puzzling at first glance!

[astrofrog]

The failure is related, it seems to be related to the fact if e.g. crval is modified in place then flag doesn't get reset. I need to look more closely to understand why set() would change results in that case though.

[manodeep]

Huh - that's weird. set_crval() calls note_change(), and the getter-setters are set to set_crval()

[astrofrog]

I think this specific example changes one element of the array in-place so I'm not sure if the setter gets called?

[manodeep]

Based on this SO post (and a couple of others I read) seems that the setter does not get called if an element of the array is modified in-place. The only way that functionality could happen is through adding __setitems__ for each "property" - which seems unnecessarily complex and overkill in this case.

One option is to move the change back to wcs_types() - both locally within astropy and upstream (if Mark agrees)

[mcara]

In comments:

Therefore, if x.flag == WCSSET, there is no need to call wcsset() and we can return immediately.

I really would like to know Mark's opinion on this. Why wasn't this check implemented in wcsset() itself as this function is quite expensive. I wonder whether there is a reason not to do it. Could someone reach out to Mark (on this and the other - wcs_types issue)?

[manodeep]

In comments:

Therefore, if x.flag == WCSSET, there is no need to call wcsset() and we can return immediately.

I really would like to know Mark's opinion on this. Why wasn't this check implemented in wcsset() itself as this function is quite expensive. I wonder whether there is a reason not to do it. Could someone reach out to Mark (on this and the other - wcs_types issue)?

Yeah agreed. In an email with Tom and me, Mark said that such if conditions exist within other functions that call wcsset. However, he also said that wcsset can be called by other functions without setting wcs->flag = 0. Here's the relevant part of Mark's email from last night (responding to a broader query about visibility of the WCSSET symbol and the proposed fix):

On the face of it, if wcsset() looked at wcsprm::flag and returned early
if set (to WCSSET), I believe that would fix the problem.  However, I
can't do that because there may be applications that use wcsset() to do
force majeur resetting of the wcsprm struct (i.e. without setting
wcsprm::flag = 0) - indeed, such code exists within WCSLIB.

There are also many instances within WCSLIB code of the following:

  if (wcs->flag != WCSSET) {
    if ((status = wcsset(wcs))) return status;
  }

which is pretty much what you want.  However, rather than expose WCSSET
(and the others), I think it may be better to provide an enquiry
function, wcsisset(), which just returns wcs->flag == WCSSET.  Would
that do the trick for you?  Alternatively, just testing wcsprm::flag
for a non-zero value should be quite reliable.

Moving the fix to wcs_types() would solve the immediate issue but would still (potentially) leave other memory race issues for the various lin/prj/ etc structs within struct wcsprm.

[mcara]

BTW, the unit test that's failing, test_get_coord_range_nan_regression, is an example of a bad way of setting CRVAL.

[manodeep]

BTW, the unit test that's failing, test_get_coord_range_nan_regression, is an example of a bad way of setting CRVAL.

That's good to know! I do not have the context to know that :)

[mcara]

context is exactly the same as in #16279 (comment): "setter does not get called if an element of the array is modified in-place". This was brought up numerous times in the past although I cannot find past discussions.

Related to the failing test, I am not sure that changing

wcs.wcs.crval[0] = 0

to

wcs.wcs.crval[0] = 0
wcs.wcs.set()

OR

wcs.wcs.crval = [0, *wcs.wcs.crval[1:]]

would solve the failure.

[mcara]

Probably if this PR is merged, something like this will not work:

w = WCS()
w.wcs.ctype = ['RA---TAN', 'DEC--TAN']
w.wcs.ctype[1] = 'GLON--TAN'
w.wcs.ctype[0] = 'GLAT--TAN'
assert w.wcs.lattyp == 'GLAT'

Now one will have to call w.wcs.set() before assert for this to work.

[manodeep]

context is exactly the same as in #16279 (comment): "setter does not get called if an element of the array is modified in-place". This was brought up numerous times in the past although I cannot find past discussions.

I meant I lacked the context for what would be considered reasonable user behaviour (and not the code context).

Probably if this PR is merged, something like this will not work:

w = WCS()
w.wcs.ctype = ['RA---TAN', 'DEC--TAN']
w.wcs.ctype[1] = 'GLON--TAN'
w.wcs.ctype[0] = 'GLAT--TAN'
assert w.wcs.lattyp == 'GLAT'

Now one will have to call w.wcs.set() before assert for this to work.

This PR should not be merged at this stage. We need to develop a better fix that allows the user to update individual items.

[manodeep]

After some digging around, I see that these lines within wcslib state that wcs->flag needs to be set to 0 whenever any of the named wcsprm members are modified.

*   int flag
*     (Given and returned) This flag must be set to zero whenever any of the
*     following wcsprm struct members are set or changed:
*
*       - wcsprm::naxis (q.v., not normally set by the user),
*       - wcsprm::crpix,
*       - wcsprm::pc,
*       - wcsprm::cdelt,
*       - wcsprm::crval,
*       - wcsprm::cunit,
*       - wcsprm::ctype,
*       - wcsprm::lonpole,
*       - wcsprm::latpole,
*       - wcsprm::restfrq,
*       - wcsprm::restwav,
*       - wcsprm::npv,
*       - wcsprm::pv,
*       - wcsprm::nps,
*       - wcsprm::ps,
*       - wcsprm::cd,
*       - wcsprm::crota,
*       - wcsprm::altlin,
*       - wcsprm::ntab,
*       - wcsprm::nwtb,
*       - wcsprm::tab,
*       - wcsprm::wtb.

Since the setters are not invoked when changing a single element, wcs->flag can't be reset to 0 with the existing call to note_change() from the setters. It seems to me that setting individual elements within these specified struct elements violates the software contract of wcslib. At the same time, there is no obvious/easy way of fixing the single-element change and resetting wcs->flag (outside of keeping an internal (deep?)copy of wcs and then comparing all the relevant fields).

[manodeep]

One way to detect whether any field has changed could be to store a hash value of the contents of the relevant fields (and ignoring memory addresses). If a change is detected, then wcs->flag could be reset to 0, and the fields would be reset/updated correctly.

Another option could be to allow user access to wcs->flag, where the user could manually update the flag to 0 after changing any relevant single element (but that puts the onus on the user to get the correct behaviour, something I tend to avoid doing as much as possible)

[astrofrog]

Another option on the astropy side might be to use a numpy array subclass to store the wcsprm elements so that we can correctly set the flag if individual items are set?

[manodeep]

Noting the reproducer here as well: https://gist.github.com/manodeep/23e5a2cc037752ce861bdb2998314137

[manodeep]

Another option on the astropy side might be to use a numpy array subclass to store the wcsprm elements so that we can correctly set the flag if individual items are set?

That is beyond my python / astropy expertise.