Repository navigation
Function or method to compute the RA and Dec offsets between two coordinates? #4209
Description
This might be a naive question, but I could not find the answer in the documentation. Is there a function or method to compute the RA and Dec offsets between two coordinates? I'm looking for a more direct way of doing this:
import numpy as np
from astropy.coordinates import SkyCoord
mm1 = SkyCoord('03h28m55.58s', '+31d14m37.06s')
mm2 = SkyCoord('03h28m55.71s', '+31d14m35.33s')
ra_offset = (mm2.ra - mm1.ra) * np.cos(mm1.dec.to('radian'))
dec_offset = (mm2.dec - mm1.dec)
print (ra_offset.to('arcsec'), dec_offset.to('arcsec'))EDIT by @eteq:
The discussion below has lead to two follow on issues, which when closed/merged will mean we can close this:
- Spherical_offsets_to method #4338 (
spherical_offsetmethod) - Create AstrometricFrame #4339 (
AstrometricFrame)
Activity
PBarmby commentedon Oct 1, 2015
I think you could get this from the
separationandposition_anglemethods, seeSkyCoord convenience functions -- but that wouldn't be any simpler than what you suggest above.
pllim commentedon Oct 1, 2015
@PBarmby , I saw that as well, but the documentation only gives answer in angle, and @smaret wants RA and DEC offsets separately.
PBarmby commentedon Oct 1, 2015
With both separation and position angle, you could compute the separate RA and dec offsets. But it wouldn't be any simpler than what @smaret is proposing above.
smaret commentedon Oct 1, 2015
Thanks for your answers. Wouldn't it make sense to have a convenience function for this?
embray commentedon Oct 1, 2015
I could swear I've seen this exact same question asked recently, like a few weeks ago, but I can't recall where....
eteq commentedon Oct 2, 2015
This has come up a few times in different contexts (although I'm not entirely sure if it has come up on-list or just in conversations I've had)... This certainly could be implemented easily (basically just place the snippet you gave here, @smaret, but as a method on
SkyCoord). The question is whether it should be.My main reasoning for questioning that is because it's not really clear how to implement it in a frame-independent way. In the case of your example (where I presume both coordinates are in ICRS), it's pretty straightforward, because both frames have the same set of coordinates. But what if, instead, I were to do:
should the
drabe ICRS delta-RA, or FK4 delta-RA? In fact, those are intrinsically different (ICRS is non-rotating, while FK4 is, so the deltas depend on equinox) or what if mm2 is Galactic coordinates? Then it's a radical difference between the two frames. This is in contrast to e.g.separation(which @PBarmby rightly pointed out), which is in principle the same regardless of what the two frames were that were giving the separation. Even worse, what if one of the coordinates is in a frame that prefers to be expressed in a cartesian representation? Then isoffset_tosupposed to return delta-RA cos(dec)/delta-Dec/delta-distance, or delta-x/delta-y/delta-z?All that said, there are possible solutions. I see a few straight away (in order of increasing difficulty:
mm1is ICRS. Then just assume "offset" always means (lon2 - lon1)*cos(lat), lat2-lat1 and don't worry about cartesian or cylindrical definitions of "offset"Thoughts (incl. @smaret @PBarmby @pllim @adrn @astrofrog)? I certainly consider this open for debate, but the ambiguities (plus the ease of sticking with status quo) have thus far tilted me in the direction of # 1.
smaret commentedon Oct 2, 2015
Angular offsets with respect to a given reference may be seen as a user-defined coordinate frame. So perhaps another option would be implement such user-defined frames? Converting any coordinates to offsets would then look like:
The various options (cartesian/spherical/cylindrical, ICRS/FK4/Galactic) could be defined when creating the
offset_mm1frame. I don't know how difficult it would be to implement this though.astrofrog commentedon Oct 2, 2015
I'd like to suggest option 4 (though to be honest I'm also ok with 1): raise an error unless the frames are the same, i.e. let the user deal with transformations but provide a convenience method that finds the differences and returns a tuple of quantities for example, if the frames match.
mhvk commentedon Oct 2, 2015
I think it may be good to think a bit more clearly about the cases where this is useful. When I have needed offsets in RA, Dec, often it was simply for getting telescopes to point right (acquire locally bright star, offset to target). For this purpose, I think a convenience method is useful, as long as it is obvious that one gets the offset relative to a given source in that source's system. I don't see any reason not to implement this; the same presumably already holds for
position_angle(which is also system dependent).The other use case I can see is doing astrometry. For that purpose, I've always found @smaret's view the easiest to work with: when you do relative positions, you are effectively defining a new coordinate system, and it is useful for that system to be at latitude = longitude = 0, with latitude and longitude aligned at the centre with lines of constant RA, Dec (or other system). In this system, offsets in arcseconds are then automatically the right offsets in RA, Dec. But with it comes of course that the system is chosen around a reference point, and this reference point also defines the system in terms of rotation, etc.
Note that there is some precedent in
Time: when one subtracts twoTimeinstances, the difference gets thescaleof the first if that's possible (ortaiif not; e.g., time differences cannot have UTC scale).pllim commentedon Oct 2, 2015
I'm okay with Option 1.
embray commentedon Oct 2, 2015
Then maybe this could go in the observation planning affiliated package, if it's not appropriate to go directly on the SkyCoord class (which doesn't necessarily have to do everything).
eteq commentedon Nov 22, 2015
@mhvk - I really like the way you've split these concepts up. I absolutely see your point that the use case that I had in mind ("offsets" like you need to give a telescope operator) are subtly different from "astrometric offsets" like what you and @smaret are discussion.
So here's a new proposal (partly motivated by #4323 and an in-person conversation I had where someone asked how to do offsets):
spherical_offsetsfunction and/or method (it could have other names if there are suggestions). This would take two coordinate objects, and would only run if the two objects are in the same frame (otherwise it would raise an error with a helpful message). It would return(ΔRA * cos(Dec), ΔDec), both asAngleobjects. This is similar to @astrofrog's option 4, but by naming itspherical_something we make it clear that it's specifically for computing offsets in the spherical system, which avoids one of my concerns above.AstrometricFrame? This would be a frame with one frame attribute: the "origin" coordinate. This would along the lines of what @mhvk and @smaret are proposing, although it has to be a full rotation (rather than justcos(dec)flattening) to account for positions that might be farther away than the small angle approximation allows. Note that the code to do this is already in theGalactocentricFrametransformation - that's actually even more complicated than this because it allows an additional "roll" angle (whereas here we're saying the new y-axis has to point in the same direction as the old y-axis at the astrometric point). We could also add aSkyCoordconvenience method to get the astrometric frame of a coordinate, allowing use likecoord1.transform_to(coord2.astrometric_frame()), which would then trivially give us the angles we'd need forspherical_offsets.If this sounds good, I can create issues for each of those two concrete sub-pieces.
taldcroft commentedon Nov 22, 2015
I'm not sure about difficulty of implementation, but I definitely like the idea of expressing offsets via a new frame as described by @mhvk. This has another advantage that it doesn't rely on small angle approximations and is essentially exact and transitive (at least in the fiducial frame, presumably that of the first coordinate(s)). I do this kind of thing all the time with spacecraft attitudes, where we are using quaternions but still need the concept of an offset.
With this idea you also get all the general representation power for free.
Downside is that it would definitely be slow compared to the simple
spherical_offsetsmethod at least in the scalar case.mhvk commentedon Nov 22, 2015
@eteq - sounds very good. As you write, all the pieces would seem to in place; the main work may be to write the documentation and test cases...
The one question I have is whether the
spherical_offsetsshould return the approximate answers, or just do the full transformation and return the longitude and latitude in the new system (which are equal to those offsets). My tendency is towards the latter, as I think our aim should be to give the right result, even if it is substantially slower than the quick calculation for small offsets (this is also what we decided about angular separation, after long discussion...).eteq commentedon Nov 23, 2015
OK, I've just created two new issues (#4338 and #4339 ) - when those are both done we can close this issue. I'd suggest moving specifics of implementation to those two issues. (Note I've already responded to one of @mhvk's comments here in #4338)
r-owen commentedon Mar 31, 2016
I agree with @mhvk that you should go with the full transformation. It reduces to the approximate answer for small offsets far from the pole, and gives the right answer everywhere -- at least with an explanation of what the values mean, e.g. length of arc along a great circle = sqrt(dlong^2 + dlat^2)