API question: Instantiating of time/coord and similar

[Thomas Robitaille]

This is a continuation of the discussion on time/coordinates initially started on the astropy mailing list by mistake (see http://mail.scipy.org/pipermail/astropy/2012-May/001987.html). Please post any future replies here!

Tom

[Thomas Robitaille]

I feel like there are two distinct points in the conversation here - what the *user* API should look like, and what the best framework is to handle coordinate transformations (with this framework being potentially used by higher level tools directly). We don't want users to have to construct/instantiate frames, mappings, frame sets, etc. We want the users to just be able to convert e.g. equatorial to galactic coordinates in the easiest and most intuitive way possible. I think the original question of this post was, what is the easiest and most generic user interface?

Astropysics does it nicely:

galactic = FK5Coords(ra, dec).convert(GalacticCoords)

and this kind of API would allow for new custom coordinate classes for example.

Regarding the similarities between WCS and coordinates, I don't think we want to start merging the two, as this will be confusing from a user point of view (and again, the user is king - Astropy is not just a framework for future tools). For me, WCS means transformations between arbitrary numbers of data and world coordinates. The current astropy.wcs handles arbitrary-dimensional transforms. Marshall, you suggested:

>>> mycoord = myifu.pixel2world(500, 500, 500)
<SpectralCelestialCoord object>

>>>> mycoord.ra
> Value with unit: 200 deg
>>>> mycoord.wavelength
> Value with unit: 6562 Angstrom

But I think this is a bad idea. Do you foresee creating classes for every type of WCS transformation? For example SpectralCelestial, SpectralTime, SpectralCelestialTime, TimeSpectral? And how would you handle arbitrary quantities that WCS allows? Why does mycoord have an ra attribute, and not l/b attributes? And if you do plan to have attributes for different coordinate systems, do you plan to have an attribute for each? (that would be a *lot* of attributes!). Would 'ra' be FK5? FK4? Again, what attributes would be present for arbitrary WCS quantities that aren't time, spectral, or celestial? If you were more specific in the class name, you would need to create SpectralFK5Coord, SpectralGalacticCoord, etc., and as I mentioned before, there is an infinite number of classes you could create like this (IFU is just a very specific case).

To me, it makes more sense to keep them separate, so I would do (for pixel to world transformations):

wx, wy, wz = wcs.wcs_pix2sky(x, y, z)

and I can choose to give more specific names if I want:

ra, dec, nu = wcs.wcs_pix2sky(x, y, z)

(and the order and type is perfectly determined from the order and types in the WCS), and then if I wanted to convert to galactic coordinates, I would take the two coordinates that are celestial, and do, for example (in Astropysics style)

galactic = FK5Coords(ra, dec).convert(GalacticCoords)

Then, galactic being a GalacticCoords object, I would know to expect it has attributes 'l' and 'b'.

Behind the scenes, I have nothing against implementing the mapping between coordinate systems in specific ways, via intermediate classes, which could be exposed for use in tools relying on Astropy.

Finally, regarding units - yes, units are an integral part of transformations, but they are also useful in many other places, so we were planning on developing that separately - it doesn't mean it couldn't be used in the wcs or coords components, but we might want to use units where coordinates and WCS aren't relevant.

Erik T. and Mike D. - since you are officially the coordinators (no pun intended) for astropy.coords and astropy.wcs, do you have thoughts on this?

Cheers,

Tom

[Peter Erwin]

On May 5, 2012, at 10:42 AM, Thomas Robitaille wrote:

> I feel like there are two distinct points in the conversation here - what the *user* API should look like, and what the best framework is to handle coordinate transformations (with this framework being potentially used by higher level tools directly). We don't want users to have to construct/instantiate frames, mappings, frame sets, etc. We want the users to just be able to convert e.g. equatorial to galactic coordinates in the easiest and most intuitive way possible. I think the original question of this post was, what is the easiest and most generic user interface?

I agree with this very much.


To expand a bit on the question of what a user interface might look like, I'd like to suggest that,
from the *input* point of view, it would be very useful if the code could possess a modicum
of flexibility and intelligence in handling input coordinate formats, so that the following :

mycoords = eq_coords(49.3274450, -41.1080650)
mycoords = eq_coords("49.3274450 -41.1080650")
mycoords = eq_coords("49.3274450", "-41.1080650")
mycoords = eq_coords("03h17m18.587s -41d06m29.03s")
mycoords = eq_coords("03:17:18.587 -41:06:29.03")
mycoords = eq_coords("03:17:18.587", "-41:06:29.03")
mycoords = eq_coords("03", "17", "18.587", "-41", "06", "29.03")

etc.

would all work sensibly.

This is pretty close to what the SIMBAD web interface (http://simbad.u-strasbg.fr/simbad/sim-fcoo)
accepts, for example.

Note that I'm assuming there would, of course, be optional keywords with sensible defaults
specifying things like epoch, system, etc. I'd imagine that plausible defaults in the case of equatorial
coordinates would be "J2000" and "ICRS", for example.


cheers,

Peter

=============================================================
Peter Erwin Max-Planck-Insitute for Extraterrestrial
er...@mpe.mpg.de Physics, Giessenbachstrasse
tel. +49 (0)176 2481 7713 85748 Garching, Germany
fax +49 (0)89 30000 3495 http://www.mpe.mpg.de/~erwin

[Wolfgang Kerzendorf]

I'm very glad that this conversation went back to the implementation details ;-)

So to give you a bit of an overview what we did for time:
There's only one time object and it has an internal format resembling JD and using an int and a float.
Once instantiated one can use properties to convert it to other time formats (e.g., mytime.to_utc).Anyways I think that would be great for coordinates as well.

On 2012-05-06, at 7:14 AM, Peter Erwin wrote:

>
> On May 5, 2012, at 10:42 AM, Thomas Robitaille wrote:
>
>> I feel like there are two distinct points in the conversation here - what the *user* API should look like, and what the best framework is to handle coordinate transformations (with this framework being potentially used by higher level tools directly). We don't want users to have to construct/instantiate frames, mappings, frame sets, etc. We want the users to just be able to convert e.g. equatorial to galactic coordinates in the easiest and most intuitive way possible. I think the original question of this post was, what is the easiest and most generic user interface?
>
> I agree with this very much.
>
>
> To expand a bit on the question of what a user interface might look like, I'd like to suggest that,
> from the *input* point of view, it would be very useful if the code could possess a modicum
> of flexibility and intelligence in handling input coordinate formats, so that the following :
>
> mycoords = eq_coords(49.3274450, -41.1080650)
> mycoords = eq_coords("49.3274450 -41.1080650")
> mycoords = eq_coords("49.3274450", "-41.1080650")
> mycoords = eq_coords("03h17m18.587s -41d06m29.03s")
> mycoords = eq_coords("03:17:18.587 -41:06:29.03")
> mycoords = eq_coords("03:17:18.587", "-41:06:29.03")
> mycoords = eq_coords("03", "17", "18.587", "-41", "06", "29.03")
>

Well, I'm not so sure about the greater flexibility, as I always strife for the zen ;-) (explicit is better than implicit; faced with ambiguity resist the temptation to guess).
especially the 03:17:18.587 ... format is ambiguous one doesn't know if these are hours or degrees (and I have seen both these formats in fits headers).
I think we need at least a keyword called format (maybe even two; ra_format, dec_format).

I think the two implementations that I find interesting are (1) keyword based instantiation coord(equatorial=(49, -41)) (2) class method based instantiation coord.from_equatorial(49, -41)
(with me having a slight preference for the second one).

Cheers
W

[Thomas Robitaille]

> I think the two implementations that I find interesting are (1) keyword based instantiation coord(equatorial=(49, -41)) (2) class method based instantiation coord.from_equatorial(49, -41)
> (with me having a slight preference for the second one).

Just a thought - the issue with the examples you are showing is that
they require the coordinate types to be hard-coded. The
astropysics-style example I showed, galactic = FK5Coords(49,
-41).convert(GalacticCoords) has the advantage of being flexible
enough that it is easy to add new classes for new coordinate systems.

Regarding the smart/magic parsing, I think some compromise between all
explicit vs all implicit would be good. Maybe it's ok to have a bit of
magic in cases where it's unambiguous (e.g. 12h23m44s +31d10'44"), but
if there's a risk of ambiguity, we raise an exception?

Cheers,
Tom

[Tim Jenness]

On Saturday, May 5, 2012 1:42:58 AM UTC-7, Thomas Robitaille wrote:

I feel like there are two distinct points in the conversation here - what the *user* API should look like, and what the best framework is to handle coordinate transformations (with this framework being potentially used by higher level tools directly). We don't want users to have to construct/instantiate frames, mappings, frame sets, etc. We want the users to just be able to convert e.g. equatorial to galactic coordinates in the easiest and most intuitive way possible. I think the original question of this post was, what is the easiest and most generic user interface?

Astropysics does it nicely:

galactic = FK5Coords(ra, dec).convert(GalacticCoords)

and this kind of API would allow for new custom coordinate classes for example.

Regarding the similarities between WCS and coordinates, I don't think we want to start merging the two, as this will be confusing from a user point of view (and again, the user is king - Astropy is not just a framework for future tools). For me, WCS means transformations between arbitrary numbers of data and world coordinates. The current astropy.wcs handles arbitrary-dimensional transforms. Marshall, you suggested:

>>> mycoord = myifu.pixel2world(500, 500, 500)
<SpectralCelestialCoord object>

>>>> mycoord.ra
> Value with unit: 200 deg
>>>> mycoord.wavelength
> Value with unit: 6562 Angstrom

But I think this is a bad idea. Do you foresee creating classes for every type of WCS transformation? For example SpectralCelestial, SpectralTime, SpectralCelestialTime, TimeSpectral? And how would you handle arbitrary quantities that WCS allows? Why does mycoord have an ra attribute, and not l/b attributes? And if you do plan to have attributes for different coordinate systems, do you plan to have an attribute for each? (that would be a *lot* of attributes!). Would 'ra' be FK5? FK4? Again, what attributes would be present for arbitrary WCS quantities that aren't time, spectral, or celestial? If you were more specific in the class name, you would need to create SpectralFK5Coord, SpectralGalacticCoord, etc., and as I mentioned before, there is an infinite number of classes you could create like this (IFU is just a very specific case).

AST handles all this by having container objects so you have CmpFrame, CmpMap, CmpRegion and they contain, say, the SkyFrame and the SpecFrame for a data cube or an IFU.

The point about the IFU is that it is a common instrument but can have 2d pixel coordinates and 3d WCS. I don't think the FITS standard can handle a IFU native coordinate frame yet but that's no reason not to support it. AST uses a SwitchMap to enable different mappings to be used depending on the location in the focal plane.

 

To me, it makes more sense to keep them separate, so I would do (for pixel to world transformations):

wx, wy, wz = wcs.wcs_pix2sky(x, y, z)

and I can choose to give more specific names if I want:

ra, dec, nu = wcs.wcs_pix2sky(x, y, z)

(and the order and type is perfectly determined from the order and types in the WCS), and then if I wanted to convert to galactic coordinates, I would take the two coordinates that are celestial, and do, for example (in Astropysics style)

galactic = FK5Coords(ra, dec).convert(GalacticCoords)

ra and dec need to be able to declare that they are a sky projection. It generally doesn't make sense to have them as separate entities until the last moment because the ra and the dec are linked.

 

Then, galactic being a GalacticCoords object, I would know to expect it has attributes 'l' and 'b'.

Behind the scenes, I have nothing against implementing the mapping between coordinate systems in specific ways, via intermediate classes, which could be exposed for use in tools relying on Astropy.

Finally, regarding units - yes, units are an integral part of transformations, but they are also useful in many other places, so we were planning on developing that separately - it doesn't mean it couldn't be used in the wcs or coords components, but we might want to use units where coordinates and WCS aren't relevant.

AST did not have units support initially. It was added later when AST moved from Sky coordinates to spectral and time coordinates where units are more varied.

-- 

Tim Jenness

Joint Astronomy Centre 

[Thomas Robitaille]

> The point about the IFU is that it is a common instrument but can have 2d
> pixel coordinates and 3d WCS. I don't think the FITS standard can handle a
> IFU native coordinate frame yet but that's no reason not to support it. AST
> uses a SwitchMap to enable different mappings to be used depending on the
> location in the focal plane.

Just a quick question about this - is an IFU slice technically not a
3d dataset with third dimension NAXIS3=1? In the original IFU cube,
the frequency will change over the depth of a single voxel, so when
the slice is extracted, would it not remain a 3d dataset in that
respect? (i.e. it has a finite frequency width over the depth of the
slice). So in WCS, I would get different results if I looked at pixel
coordinates (100, 200, 0.5) or (100, 200, 0.7), right?

Tom

[Wolfgang Kerzendorf]

I think that is a sensible way to think about it. Especially when the IFU CUBE has spatial and spectral axes not aligned with the array axes.

Cheers
W


> Tom

[Tim Jenness]

Apologies for replying to myself.

On Mon, May 7, 2012 at 10:25 AM, Tim Jenness <tim.j...@gmail.com> wrote:
>
>
> On Saturday, May 5, 2012 1:42:58 AM UTC-7, Thomas Robitaille wrote:

>>
>> and I can choose to give more specific names if I want:
>>
>> ra, dec, nu = wcs.wcs_pix2sky(x, y, z)
>>

I don't think you can do this because you don't know that RA will be
returned first unless wcs_pix2sky declares that this will always be
the case. You also don't really know if you will get an RA back or a
Galactic longitude.

From my experience with the perl Astro::Coords module this has to
return a SkyCoords object and a frequency/spectral object and those
objects have to be able to declare what they are.

The SkyCoords object would know whether it was FK5 or Galactic or
whatever (in Astro::Coords I internally store in J2000 to make life
easy). It would also know how to convert to another place on the sky.
This could return a simple ra,dec or l,b as floats and that's what I
did originally in Astro::Coords until it became obvious over the years
that it should return angle/hour objects that would be able to
stringify themselves and convert from degrees to radians on demand.

The SpecCoords object needs to be able to convert to velocity as well
(so needs a position on the sky and an observer position) to be useful
(at least for us submm and radio astronomers).

>> (and the order and type is perfectly determined from the order and types
>> in the WCS), and then if I wanted to convert to galactic coordinates, I
>> would take the two coordinates that are celestial, and do, for example (in
>> Astropysics style)
>>
>> galactic = FK5Coords(ra, dec).convert(GalacticCoords)
>>
>

so in my object model

galactic = coords.convert(GalacticCoords)

I think with galactic being a tuple of two angle objects that can
stringify to d:m:s or numify to radians.

[Tim Jenness]

I have no idea about the JWST IFU but for UIST (on UKIRT) you get a 2D
array and no WCS. The pipeline knows where each slice ends up on the
chip so can define the areas that correspond to each wavelength and
spatial slice. Once we have defined that via a SwitchMap and applied
it all our software automatically knows the ra/dec and wavelength of
each pixel (so DS9 can give you feedback as you move your mouse
around) and our rebinning tools automatically know how to regrid to a
3d data cube. All this happens without any special coding in the
applications because they all just use the WCS (which in this case can
not be represented by FITS) and the supporting AST infrastructure.

I'm not really sure how fudging the 3rd pixel coordinate so that you
get different values as you shift by 0.1 pixel works for an IFU. When
you cross a boundary from one slice to the next in 2D how do you know
you are meant to increment the 3rd dimension pixel by a set amount?

Here's a submm example. When people are looking for submm outflows
they have their data cube (ra,dec,velocity) and they want to slice it
up and display chunks that have been collapsed into certain sized
velocity bins. In the past these were generated by a single
application that would split up the plot into separate areas and put a
slice into each area. This worked fine but was not very flexible. We
found the best approach (and this helps in scripting) was to have one
app that generated the channel map (as we call it) and then every
display tool we have can display it like any other images. The channel
maps are standard 2d pixel arrays with each channel being an ra/dec
postage stamp at a particular velocity. As you move from one boundary
to the next the velocity jumps and then stays fixed for the whole
channel. I can probably dig out some example images if there is
interest. The point I'm trying to make is that these kind of datasets
are easy to find.

--
Tim Jenness

[Wolfgang Kerzendorf]

On 2012-05-07, at 1:36 PM, Tim Jenness wrote:

> Apologies for replying to myself.
>
> On Mon, May 7, 2012 at 10:25 AM, Tim Jenness <tim.j...@gmail.com> wrote:
>>
>>
>> On Saturday, May 5, 2012 1:42:58 AM UTC-7, Thomas Robitaille wrote:
>
>>>
>>> and I can choose to give more specific names if I want:
>>>
>>> ra, dec, nu = wcs.wcs_pix2sky(x, y, z)
>>>
>
> I don't think you can do this because you don't know that RA will be
> returned first unless wcs_pix2sky declares that this will always be
> the case.

To channel Obama: Yes we can! Because my wcs will give back three coordinates which are designed (by ME for MY current IFU) to be RA (ICRS), DEC(ICRS), wavelength. The units will be described in wcs.units and will be a tuple of (degree, degree, angstrom). IMHO a wcs object just does a transform without understanding what it means. I mean if one really wants to one can write a wcs that gives back actual sky coordinates, but that can be built on top of the base wcs. What we're trying to do (that's at least my view on it) is provide a base framework to do this all (which needs to be very Vanilla as to cater to the needs of many astronomers; IFU is just one of many 3d datasets).

Cheers
W

[Tim Jenness]

On Mon, May 7, 2012 at 10:53 AM, Wolfgang Kerzendorf
<wkerz...@gmail.com> wrote:
>
> On 2012-05-07, at 1:36 PM, Tim Jenness wrote:
>
>> Apologies for replying to myself.
>>
>> On Mon, May 7, 2012 at 10:25 AM, Tim Jenness <tim.j...@gmail.com> wrote:
>>>
>>>
>>> On Saturday, May 5, 2012 1:42:58 AM UTC-7, Thomas Robitaille wrote:
>>
>>>>
>>>> and I can choose to give more specific names if I want:
>>>>
>>>> ra, dec, nu = wcs.wcs_pix2sky(x, y, z)
>>>>
>>
>> I don't think you can do this because you don't know that RA will be
>> returned first unless wcs_pix2sky declares that this will always be
>> the case.
> To channel Obama: Yes we can! Because my wcs will give back three coordinates which are designed (by ME for MY current IFU) to be RA (ICRS), DEC(ICRS), wavelength. The units will be described in wcs.units and will be a tuple of (degree, degree, angstrom). IMHO a wcs object just does a transform without understanding what it means. I mean if one really wants to one can write a wcs that gives back actual sky coordinates, but that can be built on top of the base wcs. What we're trying to do (that's at least my view on it) is provide a base framework to do this all (which needs to be very Vanilla as to cater to the needs of many astronomers; IFU is just one of many 3d datasets).

I'm a submm astronomer by training so I understand that IFU is only
one type of 3d dataset. I was trying to make sure that something
wasn't designed that made it hard to handle them though.

Are you saying that wcs_pix2sky has to be written by each data file
provider? This confuses me since part of the FITS standard was to
define the WCS such that the mapping from pixel to sky was a
documented standard. Or is it that you set attributes in the WCS
object to ensure that you get what you expect when your run pix2sky?
The WCS object in that case sounds very much like an AST FrameSet but
constrained to handle purely WCS. Is there also a wcs_sky2pix method
to go the other way?

Are you also saying that the ra returned above is a simple float
rather than something that knows it is an RA? If it is an object then
I strongly advise that the object combine the ra and the dec into one
entity to make spherical sky coordinate transformations easy.

[Thomas Robitaille]

>> I don't think you can do this because you don't know that RA will be
>> returned first unless wcs_pix2sky declares that this will always be
>> the case.

> To channel Obama: Yes we can! Because my wcs will give back three coordinates which are designed (by ME for MY current IFU) to be RA (ICRS), DEC(ICRS), wavelength. The units will be described in wcs.units and will be a tuple of (degree, degree, angstrom). IMHO a wcs object just does a transform without understanding what it means. I mean if one really wants to one can write a wcs that gives back actual sky coordinates, but that can be built on top of the base wcs. What we're trying to do (that's at least my view on it) is provide a base framework to do this all (which needs to be very Vanilla as to cater to the needs of many astronomers; IFU is just one of many 3d datasets).

I agree with Wolfgang here - I think the base WCS library should just
return values and it is up to the user to figure out what they are.

The approach we've started to take with Astropy is instead to have an
n-dimensional data container (in astropy.nddata) and then plan to have
class-based representations of specific kinds of data (e.g. spectrum,
sky image, IFU, etc.). These will take care of the higher-level type
of interface that you are proposing (for example returning coordinate
and spectral objects in transformations). Of course, this is not set
in stone, so the current discussion is interesting in that respect.

I personally favor the model of having WCS be agnostic of coordinate
type when the user convers pixel and sky coordinates, and leave it up
to higher level classes to deal with specific well-defined cases such
as IFU.

Cheers,
Tom

[Thomas Robitaille]

> Are you saying that wcs_pix2sky has to be written by each data file
> provider? This confuses me since part of the FITS standard was to
> define the WCS such that the mapping from pixel to sky was a
> documented standard. Or is it that you set attributes in the WCS
> object to ensure that you get what you expect when your run pix2sky?
> The WCS object in that case sounds very much like an AST FrameSet but
> constrained to handle purely WCS. Is there also a wcs_sky2pix method
> to go the other way?

pywcs and astropy.wcs implement the WCS transformation as a completely
general (FITS) WCS transformation, and given arrays of pixel values,
will return world coordinates, and vice versa.

> Are you also saying that the ra returned above is a simple float
> rather than something that knows it is an RA? If it is an object then
> I strongly advise that the object combine the ra and the dec into one
> entity to make spherical sky coordinate transformations easy.

Again, I think this should be done at a higher level than the WCS
framework, but I can see that different people will have different
feelings about this.

Tom

[Peter Erwin]

On May 6, 2012, at 4:04 PM, Wolfgang Kerzendorf wrote:
>>
>> To expand a bit on the question of what a user interface might look like, I'd like to suggest that,
>> from the *input* point of view, it would be very useful if the code could possess a modicum
>> of flexibility and intelligence in handling input coordinate formats, so that the following :
>>
>> mycoords = eq_coords(49.3274450, -41.1080650)
>> mycoords = eq_coords("49.3274450 -41.1080650")
>> mycoords = eq_coords("49.3274450", "-41.1080650")
>> mycoords = eq_coords("03h17m18.587s -41d06m29.03s")
>> mycoords = eq_coords("03:17:18.587 -41:06:29.03")
>> mycoords = eq_coords("03:17:18.587", "-41:06:29.03")
>> mycoords = eq_coords("03", "17", "18.587", "-41", "06", "29.03")
>>
> Well, I'm not so sure about the greater flexibility, as I always strife for the zen ;-) (explicit is better than implicit; faced with ambiguity resist the temptation to guess).
> especially the 03:17:18.587 ... format is ambiguous one doesn't know if these are hours or degrees (and I have seen both these formats in fits headers).
> I think we need at least a keyword called format (maybe even two; ra_format, dec_format).

Hmm... I think we have some different conceptions of what "zen" implies ;-)
(I would interpret "zen" as suggesting, among other things, *simplicity* and not having
to constantly worry about fussy details.
Of course, I'm emphasizing zen for the *user*, not for the API designer ;-)

I'll admit I've *never* seen "03:17:18.587" in the RA place to mean anything other
than h:m:s, and I think the latter a reasonable default interpretation. Of course, if
you think that the degrees interpretation is both possible and legitimate, then that's
a job for ra_format and dec_format keywords.

I certainly have no problem with keywords like ra_format and dec_format, as long as
they're *optional*. Forcing users to specify something like
ra_format='sex', ra_units='hour'
when the input is, e.g., "23h12m15.6s" strikes me as willfully perverse.

I think one example of a reasonably good user interface is the matplotlib.pyplot.plot
function. If I want to quickly plot some y values against some x values, it lets me
do just that:
pyplot.plot(x, y)
Then, if I want to start changing things from the defaults, I can do this with
keywords. But I don't *have* to specify color, line thickness, etc., etc., if I don't
need to.

>
> I think the two implementations that I find interesting are (1) keyword based instantiation coord(equatorial=(49, -41)) (2) class method based instantiation coord.from_equatorial(49, -41)
> (with me having a slight preference for the second one).

Of those two, the second is much more preferable, both because it involves less
fussy typing (no unnecessary parentheses and equals sign) and because the concept
of a function which has no standard arguments but *requires* a keyword argument
is, I think, somewhat nonstandard and confusing.

-- Peter

[Erik Bray]

On Mon, May 7, 2012 at 4:01 PM, Peter Erwin <peter...@gmail.com> wrote:
> I certainly have no problem with keywords like ra_format and dec_format, as long as
> they're *optional*. Forcing users to specify something like
>   ra_format='sex', ra_units='hour'
> when the input is, e.g., "23h12m15.6s" strikes me as willfully perverse.

I agree that requiring it is a bit much, but it's still good that it
can be specified. If we're dealing with a lot of coordinates, and we
know ahead of time what their format will be, it's quite excessive to
try to parse out what the format is for each coordinate object
initialized. Perhaps there should be a classmethod that returns a
dict of format keywords for a given coordinate input that can then be
used to instantiate more coordinate objects.

[Erik Tollerud]

Can I propose we split this discussion? If we constantly say "wait,
those two things are interrelated" we'll never get anywhere. As I see
it, there are a few discussions going in one intermingled thread - if
we're ok with splitting the topic, I would recommend doing that if you
reply to this.

* How to go from WCS to coordinate systems. By this I don't
necessarily mean the WCS *standard*, but rather, as Tom says, a scheme
of mapping unitless/data values onto *some* coordinate system. Let me
suggest a different term, so we can discuss it without commingling it
with implementation details of WCS: "DTCT" (for Data-To-Coordinate
Transform). How this maps to coordinate systems both in a physical
and API sense is something we need to think about, and should probably
be a separate thread.

* How to deal with IFU data or similar where spatial and wavelength
(or other physical types) of data are intermixed. This clearly ties
into the above and below point, but I think it needs to be carefully
considered *on its own* (as Vicki and others rightfully pointed out)
with an eye towards practicality for astronomy users.

* Spatial coords API: This is the one thing that can actually be
implemented now, and has use completely independent of the other
issues (sometimes I have an RA/Dec table and just want to know it in
galactic coordinates!). I've made some specific comments on this
topic below, but I think it might be more constructive to work from a
skeleton API. Demitri Muna and I are working on this now, so
hopefully there'll be something a little bit more specific to work
from soon. We then do a relatively simple "it just works"
implementation and we can expand the implementation later to fit with
the above two concerns.

* Should we use AST? This part seems to have lapsed a bit, but as (I
think) I said earlier, I can see no way we can use the implementation
of AST directly as a piece of the astropy core for a variety of
reasons (which I can go into if someone wants me to). What we
*should* do is use implementation lessons their valuable experience
can bring us, and hopefully also make a pyast <-> astropy.coords
mapping possible.

* Time API: this is nominally part of this thread, but the discussion
seems to have lapsed in favor of PR #212, which is fine for me, but it
maybe could use some astropy-dev exposure... but either way it should
be a separate discussion, not explicitly referencing coordinate
systems.

> mycoords = eq_coords(49.3274450, -41.1080650)
> mycoords = eq_coords("49.3274450 -41.1080650")
> mycoords = eq_coords("49.3274450", "-41.1080650")
> mycoords = eq_coords("03h17m18.587s -41d06m29.03s")
> mycoords = eq_coords("03:17:18.587 -41:06:29.03")
> mycoords = eq_coords("03:17:18.587", "-41:06:29.03")

I agree that these should all work. I understand Wolfgang's point
that explicit is better than implicit, but you're removing one of the
most relevant lines for astronomers: "Although practicality beats
purity." astropysics.coords supports this functionality, and use it
literally almost every day. Although note that, I think instead of
`eq_coords`, the class should be more specific - e.g. `ICRSCoords` or
`FK5Coords`, as those are well-defined systems, rather than
"equatorial," which means a lot of different things.

> mycoords = eq_coords("03", "17", "18.587", "-41", "06", "29.03")

This one, on the other hand, is more ambiguous... Its true that in
most cases its likely sexigesimal, but this was the only one that I
looked at and said "I've never seen something like that in a table of
a paper before".

>> when the input is, e.g., "23h12m15.6s" strikes me as willfully perverse.
>
> I agree that requiring it is a bit much, but it's still good that it
> can be specified.  If we're dealing with a lot of coordinates, and we
> know ahead of time what their format will be, it's quite excessive to
> try to parse out what the format is for each coordinate object
> initialized.  Perhaps there should be a classmethod that returns a
> dict of format keywords for a given coordinate input that can then be
> used to instantiate more coordinate objects.

My experience has been that it wouldn't speed things up that much to
pre-determine the parser... Object creation/function call overhead
will probably dominate. *BUT*, if we end up with an array-based
coordinate class (which I think we have to, eventually), this becomes
more important.



These goals can be met with two different interface approaches though,
and I good question is which one should we use. The two choices are:

1. The class initializers do the parsing. e.g.,
ICRSCoords("03h17m18.587s -41d06m29.03s")
2. The class initializers only accept an internal representation, and
some other function is used for the parsing.
ICRSCoords(49.32744583, -41.1080639)
or
ICRSCoords.parse_string("03h17m18.587s -41d06m29.03s")
or
parse_string("03h17m18.587s -41d06m29.03s",coordsys=ICRSCoords)

I favor #1 because it is the easiest to teach to new users, and
there's no obviously compelling reason to separate out the
functionality (it can always be *implemented* internally through a
private parsing function to avoid reinventing the wheel for multiple
coordinate systems).






--
Erik Tollerud

[Wolfgang Kerzendorf]

I very much agree: WCS/IFU to be done later.

So I guess there was some misunderstanding when talking about time and coord. They are also linked physically, but what I meant was that the APIs should be the same.
Along these lines, I think we should talk about: one coordinate class with many initializers vs many coordinate classes with many converters. That should be the same for time and coordinate.

Well I think there's plenty of time for interfaces that do lots of magic once we have a good explicit base API.

>
>> mycoords = eq_coords("03", "17", "18.587", "-41", "06", "29.03")
>
> This one, on the other hand, is more ambiguous... Its true that in
> most cases its likely sexigesimal, but this was the only one that I
> looked at and said "I've never seen something like that in a table of
> a paper before".
>

Well I've seen a few files where the coordinates are space separated (I mean instead of a colon). So this is technically the thing that is described above.

I think we've reached one of the most crucial, but also most difficult parts of astropy.

Some things that (IMHO) need to be addressed:

1. Many classes for an object vs one class with many instantiators (JDTime(242394) vs time(jd=2943240))
2. Arrays of these objects (dtype) vs objects with built-in arrays.
3. How do these things behave in recarrays/tables.
4. How to read with genfromtxt, loadtxt, sql converter functions, ...
5. Separating out parsers vs built-in parsers.
6. All of these questions are adressing both coord and time for a uniform API (In that sense it addresses any construct that has a number and a unit).

Cheers
W
>
>
>
>
> --
> Erik Tollerud

[David Berry]

On 7 May 2012 19:03, Thomas Robitaille <thomas.r...@gmail.com> wrote:
>>> I don't think you can do this because you don't know that RA will be
>>> returned first unless wcs_pix2sky declares that this will always be
>>> the case.
>
>> To channel Obama: Yes we can! Because my wcs will give back three coordinates which are designed (by ME for MY current IFU) to be RA (ICRS), DEC(ICRS), wavelength. The units will be described in wcs.units and will be a tuple of (degree, degree, angstrom). IMHO a wcs object just does a transform without understanding what it means. I mean if one really wants to one can write a wcs that gives back actual sky coordinates, but that can be built on top of the base wcs. What we're trying to do (that's at least my view on it) is provide a base framework to do this all (which needs to be very Vanilla as to cater to the needs of many astronomers; IFU is just one of many 3d datasets).
>
> I agree with Wolfgang here - I think the base WCS library should just
> return values and it is up to the user to figure out what they are.

It all depends on what precisely you mean by "wcs".

What you're saying sounds very similar to AST which make a clear
distinction between Mappings and Frames. A Mapping is a numerical
recipe for transforming axis values, and knows nothing about what
those axis values represent. This is presumably equivalent to what you
are suggesting.

But in addition you also need another base class to describe a
coordinate system. This is what AST calls a Frame - collection of
attributes that tell you what the axis values represent. Whether they
are spectral, spatial, temporal or whatever, what their units are,
plus all sorts of other useful things. So when you say "it is up to
the user to figure out what they are", I would understand this to mean
that the user (or application code more likely) keeps a Frame object
around that describes what the axis values represent. Without such a
Frame class, it is not possible to write generic application code that
can handle positions and transformations between any arbitrary
coordinate systems.

In this sense , a "WCS" such as you read from a FITS file (what AST
calls a FrameSet") is a collection of Frames telling you what
coordinate systems the FITS header knows about, together with a set of
Mappings telling you how to transform positions between any two of
these Frames.

I really think it's important to start out by choosing some generic
model for the representation of arbitrary coordinate systems and their
inter-relationships, and then instantiate the components of that model
for each required use-case. Otherwise, you end up in morass of
differing requirements that can easily lead to confusion.

David

[David Berry]

If you break the tie between (for instance) RA and Dec, then there is
very little you can do with the resulting values. If I gave you simply
an RA value, just about the only thing you could do with it is to
format it into H:M:S or whatever. Specifically you cannot convert a
single RA value into any other celestial coordinate system. Neither
can you convert it into a corresponding pixel value.

AST and pywcs are in agreement over this - neither of them split up
pairs of celestial (longitude,latitude) values.

David

[David Berry]

On 8 May 2012 06:37, Wolfgang Kerzendorf <wkerz...@gmail.com> wrote:
> I very much agree: WCS/IFU to be done later.

The beauty of the AST-like model is that it allows exactly this sort
of expansion. Do the common stuff first, and the specialist stuff
later. But it relies on having a sufficiently generic model for
describing arbitrary coordinate systems and their inter-connections,
such as the Frame/Mapping/FrameSet approach used in AST.

> So I guess there was some misunderstanding when talking about time and coord. They are also linked physically, but what I meant was that the APIs should be the same.
> Along these lines, I think we should talk about: one coordinate class with many initializers vs many coordinate classes with many converters. That should be the same for time and coordinate.

I would recommend creating a base generic coordinate system class, and
then a set of sub-classes to describe separate physical spaces. So
you would a sub-class to describe "positions on the sky"
, another to describe "position in an electromagnetic spectrum",
another to describe "positions in time", and maybe a fourth to
describe an aggregation of two of the other subclasses. Each physical
space can in general be described by various different coordinate
systems, so sky positions can be given as (RA,Dec), (l,b) etc.,
spectral position as freq, wavelength, velocity, etc, time can be MJD
, epoch etc. The idea is that each sub-class encapsulates all the
information needed to be able to transform between any two such
alternative coordinate systems. So a SkyFrame has an equinox, epoch,
etc, a SpecFRame has a rest frequency, unit, standard of rest, etc.,
and a TimeFrame has a time scale, unit, etc.

Doing it this way means that each class is able to generate Mappings
that transform positions between any two instances of the class, even
if they use different specific coordinate systems.

David

[Tim Jenness]

On Mon, May 7, 2012 at 10:37 PM, Wolfgang Kerzendorf
<wkerz...@gmail.com> wrote:

> Some things that (IMHO) need to be addressed:
>
> 1. Many classes for an object vs one class with many instantiators (JDTime(242394) vs time(jd=2943240))

Is your time model supporting TAI, TDB, TT etc? It has to to support
solar system objects. In which case I assume you want to include
sufficient API to be able to convert from UTC to TAI to TDB and vice
versa. A Julian date on its own doesn't work without context.

Is the JDTime constructor above an alternative to a string
constructor? Or is it meant to be implying UTC JD with other options
being TAITime, TDBTime?

time(jd=2043240, timescale=TAI)
jd = time.tdb

would work and timescale could default to UTC. I have no feel for how
much data is required to store all the known DUT1 values in the class.

--
Tim Jenness

[Perry Greenfield]

In line with the issue of partitioning topics (and perhaps creating
different threads on this), I see the following potential issues
(mostly separable).

1) Terminology. What do we call various things? I've seen the
mathematical function that transforms pixels to coordinates called 3
different things: mappings, transformations, wcs's. It would be good
to adopt a standard name for this since it would be a core element.
Likewise, frame vs coordinate system, etc.

2) Is there any objection to clearly distinguishing frame/coordinate
system issues from mapping/transform issues? It sure seems to me that
the mapping/transform functionality can be used quite separately from
frames/coordinate systems. Frames may need mappings, but the converse
isn't true.

3) Are units simply a input/output only issue? Will the internal
representations use a consistent unit for coordinates? (that's what I
tend to think, but maybe I'm off here).

4) How many classes are needed to represent frames/coordinate systems?
I think there is a middle ground here perhaps between the two extremes
I've seen suggested.

5) Lower level details (keywords vs methods, etc)

I may have left a few out, but maybe this would help focus discussion
points.

[Perry Greenfield]

What if we had a few basic types of frames/coordinate systems like

Sky
Spectral
Time

and some way of indicating combinations of these where the order
indicates what order the coordinates come in.

(Sky, Spectral) rather than SpectralCelestral. And that each Sky
coordinate variant is not a new class, but something the class know
how to represent itself in, perhaps using labeled arrays. Say:

coords = skyinstance.as("FK5")

Where coords would have each axis labeled "ra", "dec" but

coords = skyinstance.as("Galactic") would have labels "l", "b", or even

coords = skyinstance.as("altaz", earthcoords=<xxx>, time=<yyy>)

Perry

[David Berry]

On 8 May 2012 15:59, Perry Greenfield <pe...@stsci.edu> wrote:
> In line with the issue of partitioning topics (and perhaps creating
> different threads on this), I see the following potential issues (mostly
> separable).

A good list.

> 1) Terminology. What do we call various things? I've seen the mathematical
> function that transforms pixels to coordinates called 3 different things:
> mappings, transformations, wcs's. It would be good to adopt a standard name
> for this since it would be a core element. Likewise, frame vs coordinate
> system, etc.

I would extend this item to include defining as well as naming the
"various things".

Taking the risk of talking too much, in AST a "Mapping" defines two
"transformations" - a "forward" transformation and an "inverse"
transformation. That is, if a Mapping goes from Frame A to Frame B,
then its "forward transformation" transform Frame A positions into
Frame B positions, and it's "inverse transformation" transforms Frame
B positions into Frame A positions.

Also, in AST a "Frame" in general describes a whole set of related
"coordinate systems" - all the coordinate systems that can be used to
describe positions in physical domain represented by the Frame. Of
course, a Frame only represents one of these coordinate systems at any
one time. Examples of such domains include "the celestial sphere", "an
EM spectrum", "time", "telescope focal plane", "pixel array", etc. So
the SkyFrame class can be used to represent positions on the sky in a
variety of supported coordinate system, and also knows how to create
Mappings that transform between any pair of these coordinate systems.
And the same goes for TimeFrame and SpecFrame.

A "WCS" in the FITS sense is a collection of coordinate systems,
including pixel coordinates and potentially several different forms of
"world" coordinates. Don''t forget FITS-WCS allows you to have several
different "alternate" world coordinate systems defined by one header.
This equates directly with the idea of a FrameSet in AST.

> 2) Is there any objection to clearly distinguishing frame/coordinate system
> issues from mapping/transform issues? It sure seems to me that the
> mapping/transform functionality can be used quite separately from
> frames/coordinate systems. Frames may need mappings, but the converse isn't
> true.

Agreed.

> 3) Are units simply a input/output only issue? Will the internal
> representations use a consistent unit for coordinates? (that's what I tend
> to think, but maybe I'm off here).

Even if it's just an input/output issue, you still need to know how to
convert between different units in order to be able to perform that
input/output. And if you have this unit conversion infrastructure
available, why not use it directly within the Frame class, and so get
the extra flexibility that provides?

> 4) How many classes are needed to represent frames/coordinate systems? I
> think there is a middle ground here perhaps between the two extremes I've
> seen suggested.

A single class for everything is one extreme, and the other extreme is
separate classes for every possible coordinate system that could be
used to define positions in any of the supported physical domains. AST
is in the middle in that it has a separate Frame sub-class for each
physical domain - sky,spectrum,time. It seems to work well in
practice, in that each such class encapsulates in a single place all
your knowledge about the corresponding physical domain.

David

[Perry Greenfield]

On May 8, 2012, at 11:30 AM, David Berry wrote:
>
>
>> 1) Terminology. What do we call various things? I've seen the
>> mathematical
>> function that transforms pixels to coordinates called 3 different
>> things:
>> mappings, transformations, wcs's. It would be good to adopt a
>> standard name
>> for this since it would be a core element. Likewise, frame vs
>> coordinate
>> system, etc.
>
> I would extend this item to include defining as well as naming the
> "various things".
>
> Taking the risk of talking too much, in AST a "Mapping" defines two
> "transformations" - a "forward" transformation and an "inverse"
> transformation. That is, if a Mapping goes from Frame A to Frame B,
> then its "forward transformation" transform Frame A positions into
> Frame B positions, and it's "inverse transformation" transforms Frame
> B positions into Frame A positions.
>

I basically take this as implicit (that transforms are invertable).
Yes, we can define ones that are not one-to-one, but it seems that for
the normal case, that invertability can be presumed (and that if an
analytic model isn't provided, an iterative approach on the forward
transform is used to invert by default).

> Also, in AST a "Frame" in general describes a whole set of related
> "coordinate systems" - all the coordinate systems that can be used to
> describe positions in physical domain represented by the Frame. Of
> course, a Frame only represents one of these coordinate systems at any
> one time. Examples of such domains include "the celestial sphere", "an
> EM spectrum", "time", "telescope focal plane", "pixel array", etc. So
> the SkyFrame class can be used to represent positions on the sky in a
> variety of supported coordinate system, and also knows how to create
> Mappings that transform between any pair of these coordinate systems.
> And the same goes for TimeFrame and SpecFrame.
>
> A "WCS" in the FITS sense is a collection of coordinate systems,
> including pixel coordinates and potentially several different forms of
> "world" coordinates. Don''t forget FITS-WCS allows you to have several
> different "alternate" world coordinate systems defined by one header.
> This equates directly with the idea of a FrameSet in AST.

Yes, I understand. But I don't want our discussion to be framed by
what FITS does for WCS. In my view when we use the term WCS, it
shouldn't carry FITS limitations.

>
>
>> 3) Are units simply a input/output only issue? Will the internal
>> representations use a consistent unit for coordinates? (that's what
>> I tend
>> to think, but maybe I'm off here).
>
> Even if it's just an input/output issue, you still need to know how to
> convert between different units in order to be able to perform that
> input/output. And if you have this unit conversion infrastructure
> available, why not use it directly within the Frame class, and so get
> the extra flexibility that provides?

Certainly, and I didn't mean to imply that the frame class couldn't
use it. But it's use can still be limited to certain areas.

Perry

[David Berry]

On 8 May 2012 16:44, Perry Greenfield <pe...@stsci.edu> wrote:
>
> On May 8, 2012, at 11:30 AM, David Berry wrote:
>>
>>
>>
>>> 1) Terminology. What do we call various things? I've seen the
>>> mathematical
>>> function that transforms pixels to coordinates called 3 different things:
>>> mappings, transformations, wcs's. It would be good to adopt a standard
>>> name
>>> for this since it would be a core element. Likewise, frame vs coordinate
>>> system, etc.
>>
>>
>> I would extend this item to include defining as well as naming the
>> "various things".
>>
>> Taking the risk of talking too much, in AST a "Mapping" defines two
>> "transformations" - a "forward" transformation and an "inverse"
>> transformation. That is, if a Mapping goes from Frame A to Frame B,
>> then its "forward transformation" transform Frame A positions into
>> Frame B positions, and it's "inverse transformation" transforms Frame
>> B positions into Frame A positions.
>>
> I basically take this as implicit (that transforms are invertable). Yes, we
> can define ones that are not one-to-one, but it seems that for the normal
> case, that invertability can be presumed (and that if an analytic model
> isn't provided, an iterative approach on the forward transform is used to
> invert by default).

I was mainly addressing your point about how to assign meaning to the
various names we are throwing around - just suggesting a distinction
between the words "Mapping" and "transformation".

>
>> Also, in AST a "Frame" in general describes a whole set of related
>> "coordinate systems" - all the coordinate systems that can be used to
>> describe positions in physical domain represented by the Frame. Of
>> course, a Frame only represents one of these coordinate systems at any
>> one time. Examples of such domains include "the celestial sphere", "an
>> EM spectrum", "time", "telescope focal plane",  "pixel array", etc. So
>> the SkyFrame class can be used to represent positions on the sky in a
>> variety of supported coordinate system, and also knows how to create
>> Mappings that transform between any pair of these coordinate systems.
>> And the same goes for TimeFrame and SpecFrame.
>>
>> A "WCS" in the FITS sense is a collection of coordinate systems,
>> including pixel coordinates and potentially several different forms of
>> "world" coordinates. Don''t forget FITS-WCS allows you to have several
>> different "alternate" world coordinate systems defined by one header.
>> This equates directly with the idea of a FrameSet in AST.
>
>
> Yes, I understand. But I don't want our discussion to be framed by what FITS
> does for WCS. In my view when we use the term WCS, it shouldn't carry FITS
> limitations.

I agree absolutely. FITS-WCS is unnecessarily prescriptive about the
nature of the mappings between the different available coordinate
systems. But the basic idea of a set of coordinate system objects
connected by mapping objects is of much wider relevance than just
FITS-WCS.

>>
>>
>>> 3) Are units simply a input/output only issue? Will the internal
>>> representations use a consistent unit for coordinates? (that's what I
>>> tend
>>> to think, but maybe I'm off here).
>>
>>
>> Even if it's just an input/output issue, you still need to know how to
>> convert between different units in order to be able to perform that
>> input/output. And if you have this unit conversion infrastructure
>> available, why not use it directly within the Frame class, and so get
>> the extra flexibility that provides?
>
>
> Certainly, and I didn't mean to imply that the frame class couldn't use it.
> But it's use can still be limited to certain areas.

Yes. The AST SkyFrame is limited to using radians as the numerical
representation of longitude and latitude. With hind-sight, there was
no real reason for such a limitation though. We could certainly have
supported degrees as well as radians.

David

[Wolfgang Kerzendorf]

So inspired by Perry's quest for clarity, I put a poll up at: http://poll.pollcode.com/0gs6

I think everyone agrees that we want some class that can convert between pixel and physical coordinates (not necessary sky and in some way or another). Currently this class name is set as wcs in the nddata structure, but we can change this. I think it doesn't really matter, let's just pick one and move on.

different "alternate" world coordinate systems defined by one header.

This equates directly with the idea of a FrameSet in AST.

Yes, I understand. But I don't want our discussion to be framed by what FITS does for WCS. In my view when we use the term WCS, it shouldn't carry FITS limitations.

Agree very very much!!!

Cheers

   Wolfgang

[Perry Greenfield]

before we go too far with this, I do distinguish between two things
here:

A mathematical function that has no intrinsic units (and I think that
is a useful class to have here), and one that provides a function, and
a tie to physical units and a frame of reference. For me the first is
a transformation and carries no standard frame of reference or units,
and the second is a WCS

Perry

[David Berry]

On 8 May 2012 17:02, Wolfgang Kerzendorf <wkerz...@gmail.com> wrote:
> So inspired by Perry's quest for clarity, I put a poll up
> at: http://poll.pollcode.com/0gs6
>
> I think everyone agrees that we want some class that can convert between
> pixel and physical coordinates (not necessary sky and in some way or
> another). Currently this class name is set as wcs in the nddata structure,
> but we can change this. I think it doesn't really matter, let's just pick
> one and move on.

What *exactly* is it that we are giving a name to? Is it the thing
that describes a numerical recipe for converting one set of numbers
into another? Or is it the thing that knows how to create such a
recipe?

It's surely more important to decide how many base classes we need and
what each represents, than what they are called. Names are secondary.

David

[James Turner]

Yes, I'm not quite sure how general the thing being named is meant
to be. I had envisaged a "transform" that would convert from any
numeric co-ordinate system to any other kind (which may or may not
be a WCS). It could be from pixels to pixels (eg. some array
subsection or a distortion correction) or it could be from frequency
to wavelength. Is that the scope Wolfgang is asking about?

I'd just like to say *+1* to both of Perry's points, viz.:

- I basically take this as implicit (that transforms are invertable).
- But I don't want our discussion to be framed by what FITS does for WCS.

Cheers,

James.

PS. Sorry, can't keep up with all of the list traffic right now :-(.

[Tim Jenness]

On Tue, May 8, 2012 at 10:08 AM, James Turner <jtu...@gemini.edu> wrote:

> I'd just like to say *+1* to both of Perry's points, viz.:
>
> - I basically take this as implicit (that transforms are invertable).
>

I think it's good to be explicit here since we've had a lot of
experience of people defining things without supplying an inverse. I
might be mis-remembering but the Spitzer WCS definition didn't bother
with an inverse and it took a lot of work to provide the inverse
inside AST.

> - But I don't want our discussion to be framed by what FITS does for WCS.

That is a very good idea. Otherwise you become overly constrained in
your design. You do have to worry about how to store your transforms
to a FITS file if the FITS WCS standard doesn't support your
particular form. (AST provides a FitsChan class for doing all the FITS
header parsing and FITS header generation and has a mode to generate a
compatible header that only AST can understand as a fallback).

--
Tim Jenness

[James Turner]

> Just a quick question about this - is an IFU slice technically not a
> 3d dataset with third dimension NAXIS3=1? In the original IFU cube,
> the frequency will change over the depth of a single voxel, so when
> the slice is extracted, would it not remain a 3d dataset in that
> respect? (i.e. it has a finite frequency width over the depth of the
> slice). So in WCS, I would get different results if I looked at pixel
> coordinates (100, 200, 0.5) or (100, 200, 0.7), right?

I'm not sure whether I've understood this question correctly (what is
special about frequency?) but I think the answer is no -- each pixel
or voxel or whatever you call it is equivalent to a infinitesimal
point sample of a 3D "image" that has been convolved with the pixel
size. However, if you were to extract planes at different wavelengths,
0.5 & 0.7, the results would be different because you're interpolating
and resampling at a slightly different position. (Hope I didn't miss
the point.)

Cheers,

James.

[James Turner]

> * How to deal with IFU data or similar where spatial and wavelength
> (or other physical types) of data are intermixed. This clearly ties
> into the above and below point, but I think it needs to be carefully
> considered *on its own* (as Vicki and others rightfully pointed out)
> with an eye towards practicality for astronomy users.

Of course this also applies to common long slit spectra, where each
pixel may also generally have 3 world (ra, dec, lambda) co-ordinates,
not just to IFU data.

Cheers,

JAmes.

[James Turner]

On 08/05/12 13:15, Tim Jenness wrote:
>> - I basically take this as implicit (that transforms are invertable).
>
> I think it's good to be explicit here since we've had a lot of
> experience of people defining things without supplying an inverse. I
> might be mis-remembering but the Spitzer WCS definition didn't bother
> with an inverse and it took a lot of work to provide the inverse
> inside AST.

Considering the general transform case, not only do I forsee
compelling use cases for this, but it seems a bit of a "no brainer",
since the inverse transform is what one normally needs for resampling,
isn't it? (ie. for each output pixel, find the corresponding input
co-ordinates to sample at). For world co-ordinates (and various things
like subtracting sky from the raw data format, a la Kelson 2003), it
has to be able to work the other way around.

> (AST provides a FitsChan class for doing all the FITS header parsing
> and FITS header generation and has a mode to generate a compatible
> header that only AST can understand as a fallback).

This may be a good model for figuring out how to store the info. in
FITS in as compatible a way as possible. Not everything can be in the
standard FITS header but presumably we want a standard FITS header
to reflect it as closely as practical??

Cheers,

James.

[David Berry]

If an AST FrameSet can be converted to FITS-WCS then it does the
conversion and produces a standard FITS-WCS header. Otherwise, it
gives you the choice of using AST own native format for storing the
WCS information (which is completely unlike FITS-WCS), or producing no
header at all. Producing a standard FITS_WCS header that describes a
rough approximation to the WCS sounds dangerous. How would subsequent
users know that it's only a rough approximation, and not the real
thing?

David

[Perry Greenfield]

On May 8, 2012, at 2:27 PM, James Turner wrote:
>
>
>> (AST provides a FitsChan class for doing all the FITS header parsing
>> and FITS header generation and has a mode to generate a compatible
>> header that only AST can understand as a fallback).
>
> This may be a good model for figuring out how to store the info. in
> FITS in as compatible a way as possible. Not everything can be in the
> standard FITS header but presumably we want a standard FITS header
> to reflect it as closely as practical??

Dealing with FITS is going to introduce warts at some level. In
thinking about this I've come to the following alternatives and
conclusions.

1) There is no good general way of mapping what we want into FITS WCS.
Some sort of compromises must be made.
2) One approach is to approximate the astropy model with a FITS WCS
model. The advantage of doing this is that it will likely use an
approach that tools support (e.g., DS9). The obvious disadvantages are
that it loses some accuracy, and it isn't rount-trippable (information
is lost about the original model)
3) One sticks the exact model into FITS by
a) sticking into some opaque FITS item. Only our tool can
understand it. It's round-trippable. But it's opaque.
b) sticking it into some complex system of FITS keywords. The
advantage is that at least FITS tools can access the information in
pieces, if not as a whole. Why, maybe some alternative FITS libraries
will grow up around it. But I see this a bit like keeping a dead horse
on life support :-)
4) trying to preserve accuracy by using pixel mapping tables in FITS
WCS. The advantage is that it is a FITS WCS standard, and little or no
accuracy is lost. The disadvantage is that it is very inefficient by
putting large arrays into the file, and it loses the original
information (it's not round-trippable). And I can think of cases where
even this doesn't necessarily support all WCS models we may want.

These are not mutually exclusive. Personally, I would do both 2) and
3a) (besides using something other than FITS ;-)

Perry

[David Berry]

In practice, most WCS information is likely going to conform to the
FITS-WCS model simply because it was read from a FITS file that uses
FITS-WCS. The most common operations such as rotation, scaling, etc,
on such WCS will not change this. So why not another option:

5) Represent in FITS-WCS if doing so looses no information, and use
some other way of encoding the WCS into a FITS header if not.

David

[Perry Greenfield]

On May 8, 2012, at 2:57 PM, David Berry wrote:

> In practice, most WCS information is likely going to conform to the
> FITS-WCS model simply because it was read from a FITS file that uses
> FITS-WCS. The most common operations such as rotation, scaling, etc,
> on such WCS will not change this. So why not another option:
>
> 5) Represent in FITS-WCS if doing so looses no information, and use
> some other way of encoding the WCS into a FITS header if not.
>
> David

But if you are a data creator (e.g., an observatory), this isn't
always true for resampled data. It's basically not the case for most
HST files (yes, we put them in FITS WCS now, but wish we had something
better), and I suspect the same will be true for JWST.

Perry

[James Turner]

I think it's also not the case that FITS WCSs generally aren't a
rough approximation; they come with all kinds of telescope pointing
errors and the like; our raw spectra even have an imaging WCS. I
suppose the danger is more when the results are supposed to be
science-quality though. One could add a warning in a comments.
However, I wasn't trying to prescribe a solution so much as raise
the question. Perry's 2 & 3a suggestion sounds sensible but those
who implement it and/or have the most relevant experience can decide.

One day it would be good for what we do to turn into some kind of
standard but first things first.

James.

[David Berry]

I sympathise. These are exactly the issues we faced with AST and
FITS-WCS. We didn't really have any hope of introducing an alternative
WCS format that would have any hope of being widely adopted outside
the UK, so went for the pragmatic approach of using fits-wcs if
possible and our own format otherwise.

David

[Victoria G. Laidler]

The Astrolib Coords package defined separate classes for the
Position, which located a point on the sky, and the Coord,
which was a representation of that point.

From the docstrings:

"The Position class is the basic class in the coords
library. It is designed to permit users to define a
position and then access many representations of it."

"A Coord is distinct from a Position by being intrinsically
expressed in
a particular set of units.

Each Coord subclass knows how to parse its own input,
and convert itself
into the internal representation (decimal degrees) used
by the package."