[AstroPy] WCS query
Perry Greenfield
perry at stsci.edu
Wed Jun 15 16:02:54 EDT 2005
On Jun 15, 2005, at 1:02 PM, Stephen Walton wrote:
> But I still find it somewhat awkward, and users here are going to find
> it very confusing, that what we old Fortran hands think of as FITS
> pixel (i,j) winds up displayed in matplotlib at Cartesian coordinates
> (i-1,j-1) even though that pixel would have to be addressed as
> data[j-1,i-1] in software. This is because John Hunter, not
> unreasonably, adopted the MATLAB convention for imshow(), contour(),
> and their relatives that the first array coordinate is vertical and
> the second horizontal. Effectively, then, imshow() 'undoes' the
> transpose done by PyFITS so that the displayed image actually comes
> out in the same orientation on the screen as it does in IRAF (if
> origin='lower' is used in imshow, as I do) even though the actual
> array has to be addressed with the indices in reverse order from
> Fortran.
>
This whole area hinges on the definitions of what the conventions are
(which can be confusing). I'm not sure I'd classify it as John Hunter's
convention but more properly a Numeric/numarray convention. That's
because the data as read from the FITS file retain the same order they
did in the FITS file (PyFITS doesn't reorder the data). Since the last
index for Numeric/numarray represents the adjacent data values, that's
where the issue really arises, not in matplotlib other than it is
adopting the convention that adjacent data values are associated with x
as far as images go (but having two different conventions for which way
y should be displayed, i.e., up or down).
As I've mentioned, this is probably the biggest sore point astronomers
are going to have with numarray (or Numeric). Different approaches
could have been used, but all have their drawbacks:
1) numarray should have adopted the other convention. It could have,
but then it would have been incompatible with Numeric in this respect
and that would have been bad. [Numeric chose (I believe) to do it that
way since it means that arr[i,j] == arr[i][j] instead of
arr[i,j]==arr[j][i] and that was thought important enough to be
consistent for]
2) PyFITS could actually reorder the data on reading. I suppose, but
this can be an expensive operation for large datasets and it renders
memory mapping useless. Even after doing so, the presumption that x
corresponds to adjacent data values is no longer true (so if algorithms
are written that presume that, large inefficiencies will result)
3) PyFITS could read the data in the normal order and fiddle with the
array stride attributes to have the same effect. This would allow
continued use of memory mapping. Although now the first index
represents adjacent data, all other Numeric/numarray libraries that
presume the last index usually represents adjacent data will be wrong.
Any of these choices results in some sort of mismatch. We concluded
that, although initially painful, the best thing was to remain as
consistent with the existing language and array package semantics
rather than introduce tricks to make it look like another language's
semantics. The former is much more painful up front, the latter
introduces long-term problems that never go away. Perhaps that's a bad
marketing strategy. It's a little like the 0-based vs 1-based indexing
issue. One can try to set things up so that you can use 1-based
indexing in a 0-based language, but in the end it's just more trouble
than learning to use it the way the language works (in my opinion
anyway). But I'll agree that this is probably the single biggest
irritation that astronomers will face.
> Incidentally, MATLAB's provided fitsread routine transposes FITS data
> on input, so in MATLAB data(i,j) refers to the same pixel as in
> Fortran (and image displays *look* transposed).
>
> I suppose long-time CFITSIO and PyFITS users are just used to all of
> this.
>
That depends. There was never consensus even in our group about which
approach was best (though there did seem to be an age correlation, but
maybe I imagined that).
Perry
More information about the AstroPy
mailing list