dstndstn / astrometry.net

Astrometry.net -- automatic recognition of astronomical images
http://astrometry.net
Other
662 stars 186 forks source link

Inconsistent header (both `PC` and `CD` transforms) written in some circumstances #176

Closed mwcraig closed 4 years ago

mwcraig commented 4 years ago

astrometry.net does not seem to correctly recognize a WCS transform in a header if it is expressed using PCi_j instead of CDi_j. It also does not delete PCi_j keywords when replacing the WCS in an image which can result in both PC and CD keywords in the header. That isn't supposed to happen [1] and can result in an inconsistent WCS or a WCS that produces world coordinates that are off by a few arcsec (full disclosure: I've written code that ends up including both myself [2]).

In the PC notation the WCS provides a linear transform PCi_j and a scale CDELTi. The CD notation has a single transform CDi_j that is the product of PCi_j and CDELTi.

Here is one way this can lead to a problem:

  1. Start with an image with an accurate WCS, represented as a PC transform, with CRPIX not at the center; the point is to end up with a different CRPIX for the PC transform than for the CD transform. For the sake of argument assume that CRPIX for the PC is [0, 0].
  2. Run astrometry.net on the image with CRPIX set to center (what matters is that the CRPIX is different than the one for the PC transform) and the --no-verify option. astrometry.net will generate a CD transform.
  3. The result is two transforms in the header:
    • A PC transform that was generated with the reference pixel at [0, 0] even though that is not the value of CRPIX in the header.
    • A CD transform that was generated with the reference pixel at the center.
  4. The FITS file is then read with software (like astropy) that, when presented with both PC and CD in the header, uses the PC transform. The PC transform was generated with the reference pixel at [0, 0] but CRPIX in the file is the center. The result is a WCS that is off by a handful of arcseconds.

To reproduce (using astrometry.net 0.78 or latest master) and astropy 3.2 or newer:

  1. Start with this FITS file TIC_83210867.01-S001-R089-C001-rp-as_pc-crpix-origin.fit, generated by:

    • running on astrometry.net on a file which had an existing WCS and the --crpix-center flag. As reported in #174 this results in an accurate WCS but with CRPIX = [0, 0] instead of at the center.
    • Change the CD transform to a PC transform using astropy and write it out to TIC_83210867.01-S001-R089-C001-rp-as_pc-crpix-origin.fit
  2. Run astrometry.net on that FITS file, ignoring the existing header when doing the plate solve (and other flags to speed things up a bit):

solve-field --obj 100 --scale-low 0.5 --scale-high 0.6 --scale-units arcsecperpix --no-remove-lines --uniformize 0 --use-sextractor --no-plot --corr none --rdls none --match none --wcs none --crpix-center --no-verify TIC_83210867.01-S001-R089-C001-rp-as_pc-crpix-origin.fit
  1. The result, TIC_83210867.01-S001-R089-C001-rp-as_pc-crpix-origin.new, has both PC and CD transforms in it. The PC transform is for a reference pixel of [0, 0] even though though CRPIX is at the center. The CD transform is correct for CRPIX at the center.
  2. Read this new image in astropy, which chooses the PC transform over the CD transform but uses the CRPIX from the new image, which is the image center. The difference in world coordinates for the zero-based pixel point [1000, 1000] is about 2.3 arcsec.

[1] See latest FITS standard paper https://fits.gsfc.nasa.gov/standard40/fits_standard40aa-le.pdf section 8.2, top of p. 31. [2] https://github.com/astropy/astropy/issues/8597

dstndstn commented 4 years ago

Here's the list of image header cards we drop when patching in a new WCS: https://github.com/dstndstn/astrometry.net/blob/master/blind/new-wcs.c#L19

On Wed, Dec 11, 2019 at 2:18 PM Matt Craig notifications@github.com wrote:

astrometry.net does not seem to correctly recognize a WCS transform in a header if it is expressed using PCi_j instead of CDi_j. It also does not delete PCi_j keywords when replacing the WCS in an image which can result in both PC and CD keywords in the header. That isn't supposed to happen [1] and can result in an inconsistent WCS or a WCS that produces world coordinates that are off by a few arcsec (full disclosure: I've written code that ends up including both myself [2]).

In the PC notation the WCS provides a linear transform PCi_j and a scale CDELTi. The CD notation has a single transform CDi_j that is the product of PCi_j and CDELTi.

Here is one way this can lead to a problem:

  1. Start with an image with an accurate WCS, represented as a PC transform, with CRPIX not at the center; the point is to end up with a different CRPIX for the PC transform than for the CD transform. For the sake of argument assume that CRPIX for the PC is [0, 0].
  2. Run astrometry.net on the image with CRPIX set to center (what matters is that the CRPIX is different than the one for the PC transform) and the --no-verify option. astrometry.net will generate a CD transform.
  3. The result is two transforms in the header:
    • A PC transform that was generated with the reference pixel at [0, 0] even though that is not the value of CRPIX in the header.
    • A CD transform that was generated with the reference pixel at the center.
  4. The FITS file is then read with software (like astropy) that, when presented with both PC and CD in the header, uses the PC transform. The PC transform was generated with the reference pixel at [0, 0] but CRPIX in the file is the center. The result is a WCS that is off by a handful of arcseconds.

To reproduce (using astrometry.net 0.78 or latest master) and astropy 3.2 or newer:

1.

Start with this FITS file TIC_83210867.01-S001-R089-C001-rp-as_pc-crpix-origin.fit https://www.dropbox.com/s/f0rx4dtoryxjuzb/TIC_83210867.01-S001-R089-C001-rp-as_pc-crpix-origin.fit?dl=0, generated by:

  • running on astrometry.net on a file which had an existing WCS and the --crpix-center flag. As reported in #174 https://github.com/dstndstn/astrometry.net/issues/174 this results in an accurate WCS but with CRPIX = [0, 0] instead of at the center.

    • Change the CD transform to a PC transform using astropy and write it out to TIC_83210867.01-S001-R089-C001-rp-as_pc-crpix-origin.fit 2.

    Run astrometry.net on that FITS file, ignoring the existing header when doing the plate solve (and other flags to speed things up a bit):

solve-field --obj 100 --scale-low 0.5 --scale-high 0.6 --scale-units arcsecperpix --no-remove-lines --uniformize 0 --use-sextractor --no-plot --corr none --rdls none --match none --wcs none --crpix-center --no-verify TIC_83210867.01-S001-R089-C001-rp-as_pc-crpix-origin.fit

  1. The result, TIC_83210867.01-S001-R089-C001-rp-as_pc-crpix-origin.new https://www.dropbox.com/s/uha9qure3n2xrbh/TIC_83210867.01-S001-R089-C001-rp-as_pc-crpix-origin.new?dl=0, has both PC and CD transforms in it. The PC transform is for a reference pixel of [0, 0] even though though CRPIX is at the center. The CD transform is correct for CRPIX at the center.
  2. Read this new image in astropy, which chooses the PC transform over the CD transform but uses the CRPIX from the new image, which is the image center. The difference in world coordinates for the zero-based pixel point [1000, 1000] is about 2.3 arcsec.

[1] See latest FITS standard paper https://fits.gsfc.nasa.gov/standard40/fits_standard40aa-le.pdf section 8.2, top of p. 31. [2] astropy/astropy#8597 https://github.com/astropy/astropy/issues/8597

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/dstndstn/astrometry.net/issues/176?email_source=notifications&email_token=AAIEH7KYYXBHPWUCDEYZ6ILQYE4H5A5CNFSM4JZT4OQKYY3PNVWWK3TUL52HS4DFUVEXG43VMWVGG33NNVSW45C7NFSM4H73KO5Q, or unsubscribe https://github.com/notifications/unsubscribe-auth/AAIEH7OAYFBF4PXFVK3QCBTQYE4H5ANCNFSM4JZT4OQA .

mwcraig commented 4 years ago

Excellent, I'll try to get a PR done today or tomorrow.

dstndstn commented 4 years ago

That's awesome, thank you!

On Wed, Dec 11, 2019 at 2:45 PM Matt Craig notifications@github.com wrote:

Excellent, I'll try to get a PR done today or tomorrow.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/dstndstn/astrometry.net/issues/176?email_source=notifications&email_token=AAIEH7I7N2BXQBBRUWKD7B3QYE7OVA5CNFSM4JZT4OQKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOEGUK35Q#issuecomment-564702710, or unsubscribe https://github.com/notifications/unsubscribe-auth/AAIEH7NKO7SJD5QBZQT45VTQYE7OVANCNFSM4JZT4OQA .

mwcraig commented 4 years ago

Can you point me towards the part of the code that reads in existing FITS headers? There is a separate problem if there is a pre-existing WCS in PC form then currently astrometry.net recognizes there is a WCS header, but isn't processing it correctly.

It should just be a couple additional lines to calculate CD from PC and CDELT and then continue on.

(I can open a separate issue for this if you want)

mwcraig commented 4 years ago

Example of the failure to verify an existing header in PC form; clearly the PC matrix isn't being used to choose the quad size:

Solving...
Reading file "/Users/mcraig/Dropbox/MSUM/Research/for-astrometry.net-bug-report/TIC_83210867.01-S001-R089-C001-rp-as_pc-crpix-origin.axy"...
Verifying WCS using indices with quads of size [24576, 246540] arcmin
Got 0 solutions.
dstndstn commented 4 years ago

Here, I think: https://github.com/dstndstn/astrometry.net/blob/master/blind/augment-xylist.c#L799

On Thu, Dec 12, 2019 at 1:23 PM Matt Craig notifications@github.com wrote:

Example of the failure to verify an existing header in PC form; clearly the PC matrix isn't being used to choose the quad size:

Solving... Reading file "/Users/mcraig/Dropbox/MSUM/Research/for-astrometry.net-bug-report/TIC_83210867.01-S001-R089-C001-rp-as_pc-crpix-origin.axy"... Verifying WCS using indices with quads of size [24576, 246540] arcmin Got 0 solutions.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/dstndstn/astrometry.net/issues/176?email_source=notifications&email_token=AAIEH7JNID4CGJWRJTLGS43QYJ6QRA5CNFSM4JZT4OQKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOEGXSCJI#issuecomment-565125413, or unsubscribe https://github.com/notifications/unsubscribe-auth/AAIEH7OHJHAGRABM7PELEHDQYJ6QRANCNFSM4JZT4OQA .