mainGNU Astronomy Utilities - Support: sr #109223, Change default WCS formalism from...


sr #109223: Change default WCS formalism from PCi_j to CDi_j

Submitter:  Lee Kelvin <leeskelvin>
Submitted:  Tue 17 Jan 2017 05:04:40 PM UTC
Category:  All Gnuastro Priority:  5 - Normal
Severity:  1 - Wish Item Group:  * Enhancement
Status:  Need Info Privacy:  Public
Assigned to:  None Open/Closed:  Open
Operating System:  GNU/Linux
* Mandatory Fields

Add a New Comment Rich Markup

Wed 18 Jan 2017 02:58:33 PM UTC, comment #3: 

Unfortunately SWARP doesn't provide the version number of WCSLIB that it is using. (just as a side-note, in Gnuastro, version information is taken very seriously and written to all output headers as described in the book.)

So I guess the only way would be to either dig into its code or ask Emmanuel. But from the formatting of the WCS related keywords, you can guess that SWARP doesn't use WCSLIB's FITS keyword output. In Gnuastro we directly write WCSLIB's produced headers into the FITS file. You can see that both the ordering and the comments on the keywords differs between SWARP and WCSLIB (Gnuastro's outputs). I think SWARP writes the keywords itsself directly from the structure.

Maybe we can make a feature request to WSCLIB. Currently `wcshdo' does take some flags for telling it how to print the header information. So we might be able to convince Mark Calabretta to add another flag for writing the `CDi_j' elements.

I just changed the category of this issue to "All Gnuastro" because all Gnuastro's programs ultimately use a single function to read and write FITS files, thus this discussion applies to all of them that deal with FITS WCS information (the bug #50073 and bug #50073 also applied to all Gnuastro's programs).

Mohammad Akhlaghi <makhlaghi>
Group administrator
Wed 18 Jan 2017 02:36:42 PM UTC, comment #2: 

Thanks for your reply. I agree, manually converting from PC/CDELT to CD is possible, but adds another layer of complexity that may break in the future. This is therefore not an ideal solution long term.

I am curious how other software packages operate in this regard. For example, the image warping and mosaicking package SWarp ( written by Emmanuel Bertin produces output files natively in the CD formalism. As far as I am aware, this package also requires WCSLIB, and I don't think they use a converter to convert from PC to CD after-the-fact.

Ideally, one would want to find a way to force WCSLIB to output WCS information natively in the CD formalism, rather than add on a converter to enforce this after-the-fact. If this is not possible, then this should remain as a feature request for now.

Lee Kelvin <leeskelvin>
Wed 18 Jan 2017 02:39:57 AM UTC, comment #1: 

Thanks for this very intriguing suggestion. Indeed the `CDi_j' formalism is the most common in most imaging archives.

This is what the WCSLIB manual says under `wcsprm::altin':

"These alternate specifications of the linear transformation matrix are translated immediately to PCi_ja by wcsset() and are invisible to the lower-level WCSLIB routines. In particular, wcsset() resets wcsprm::cdelt to unity if CDi_ja is present (and no PCi_ja )."

I couldn't find any better explanation. I would guess that since WCSLIB is intended for any kind of data (not just 2D images), for example spectra, 3D data cubes, wide full-sky surveys, and etc, it has chosen this formalism since the scale is separate and more clearly defined. I should just add that with the fix to bug #50073, the `CDELTi's in the `PCi_j' formalism are corrected.

Ultimately, if its not possible to tell WCSLIB to print the CD matrix, it is not hard to manually change the `PCi_j' and `CDELTi' key words to `CDi_j' immediately before saving then in the FITS. I am just slightly concerned because WCSLIB has already been written to be the most general library for the job, taking the most general applications in mind. But here, we are only focused on 2D imaging surveys. This is also the type of data I deal with most, but we don't want to limit Gnuastro to this type of data.

Just for reference the conversion from `wcsprm' to text is done by the call to `wcshdo' in `lib/fits.c' and `gal_fits_add_wcs_to_header' is in charge of writing that text into the output FITS header. Some corrections can be done in either of these two steps (before writing into the FITS or after it).

So, let's keep this discussion going on here (not in the task tracker, so everyone can easily post comments). Please share your thoughts on what the best way forward can be.

Mohammad Akhlaghi <makhlaghi>
Group administrator
Tue 17 Jan 2017 05:04:40 PM UTC, original submission:  

The FITS standard has two means by which to describe the translation from pixel position to the adopted world coordinate system (WCS): PCi_j and CDi_j, where i and j are the WCS and pixel coordinate axis numbers respectively. (N.B.: there exists a third standard involving rotational information using the header keyword CROTA, but usage of this is now deprecated).

Both the PCi_j and CDi_j formalisms provide linear transformation matrices between pixel axes j and coordinate axes i.

Pixel scale in the PCi_j formalism is encoded in another header keyword, CDELTi. Pixel scale is already natively encoded in the CDi_j formalism.

ImageWarp applies the PCi_j formalism to warped (rotated) output files. This occurs even if the input file adheres to the CDi_j formalism. I believe that it would make more sense to adopt the CDi_j formalism for describing WCS information rather than PCi_j, if at all possible.

The reasons for considering this may be:
1) The CDi_j formalism entirely describes the linear transformation matrix and pixel scale, and is therefore simpler.
2) The majority of modern surveys adopt the CDi_j formalism as its native WCS encoding. By way of example, I attach a tarball containing small 3x3 FITS files with data produced at these facilities: VST, SDSS, Chandra, Wise, Hubble, Spitzer, VLA, DSS. Each of these datasets adopts the CDi_j formalism.
3) As the majority of modern data uses the CDi_j formalism, adopting it within ImageWarp necessitates less translation of the input data.

I'm aware that WCS manipulation within ImageWarp is handled by WCSLIB. Perhaps there is a flag which allows ImageWarp to select the preferred output formalism when invoking WCSLIB?

Alternatives to enforcing use of the CDi_j formalism always could be:
1) ImageWarp adopts the same formalism as the input dataset,
2) ImageWarp allows the user to specify on the command line their preferred formalism.

Lee Kelvin

Lee Kelvin <leeskelvin>


(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:

Attached Files
file #39496:  wcsexamples.tar.gz added by leeskelvin (4KiB - application/gzip)


Depends on the following items: None found

Items that depend on this one: None found


Carbon-Copy List
  • -email is unavailable- added by makhlaghi (Posted a comment)
  • -email is unavailable- added by leeskelvin (Submitted the item)

    There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

    Only logged-in users can vote.


    Follow 3 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2017-01-18 makhlaghi StatusNone Need Info
    2017-01-18 makhlaghi CategoryWarp All Gnuastro
    2017-01-17 leeskelvin Attached File- Added wcsexamples.tar.gz, #39496

    Back to the top

    Powered by Savane 3.13-3cd0.
    Corresponding source code