Repository navigation
Right Ascension Angular Separation cos(dec) concern #4323
Description
Hi, I am checking some RA angular separations with astropy.
has some examples which I started from:
c1 = SkyCoord('12h00m01.0s', '-00d00m00s', frame='icrs')
c2 = SkyCoord('12h00m00.0s', '-00d00m00s', frame='icrs')
print()
print('dRA: ', (c1.ra-c2.ra).arcsec)
print('dDec: ', (c1.dec-c2.dec).arcsec)
print('AngSep: ', c1.separation(c2).arcsec)
dRA: 15.0
dDec: 0.0
AngSep: 15.0
as expected:
Now moving to Dec=-60 where Cosine(Dec)=0.5
c1 = SkyCoord('12h00m01.0s', '-60d00m00s', frame='icrs')
c2 = SkyCoord('12h00m00.0s', '-60d00m00s', frame='icrs')
print()
print('dRA: ', (c1.ra-c2.ra).arcsec)
print('dDec: ', (c1.dec-c2.dec).arcsec)
print('AngSep: ', c1.separation(c2).arcsec)
dRA: 15.0
dDec: 0.0
AngSep: 7.49999999875
Angsep is OK, but I was expecting for dRA, I was expecting dRA = 7.5arcsec at dec=-60.0
Activity
mhvk commentedon Nov 19, 2015
@richardgmcmahon -- if you just subtract the
RA, you get the difference inRAunits. This is more naturally seen if you keep RA in time units: there is a one second difference in time between the two objects crossing the same location, independent of declination. Or, by analogy with longitude and latitude, the difference between two lines of longitude in angular units does not depend on latitude.mhvk commentedon Nov 19, 2015
p.s. This is why proper motion in
RAis often explicitly labelledmu_alpha cos(delta). But you are just calculating Delta alpha.richardgmcmahon commentedon Nov 19, 2015
@mhvk I am not sure whether you are agreeing with me that there is an issue.
In my opinion arc seconds means a great circle second of arc angle and this is the unit I would expect. Tt is dangerous in my opinion c1.ra-c2.ra).arcsec is not returning seconds of arc as expected.
c1, c2 'know' what the declination is so I would expect that this be used.
My two cents
mhvk commentedon Nov 19, 2015
I actually do not think that there is an issue; rather, the code does what I would hope it would do. Here, my thinking is informed by the fact that it is useful to be able to calculate the difference between longitudes, in a way that is independent of one's location on a sphere, and that this is most logically done via just subtracting the two longitudes.
Indeed, thinking a bit more, I'm not even sure that for a general case I understand the use of getting the answer in any different way. For instance, what would you have expected to get for
If you would like
cos(delta)to be used here, which declination should one pick? And what physically would the answer mean?richardgmcmahon commentedon Nov 20, 2015
I understand your thinking. However the difference in longitude is the same in seconds of time and not the same in seconds of arc at different declinations. IF one was to implement a solution, I would use the dec of the first position. Maybe there is an IAU/Geo recommendation on this for long/lat since we do not want to reinvent a recommendation. Alternatively document the potential for confusion.
Consider the example if I wanted to get from (ra1, dec1) to a point (r2, dec2) by three routes.
(i) create circle route (which is the shortest)
(ii) travel/fly/draw from ra1(long1) to ra2(long2) along dec1(lat1) and then travel from dec1(lat1) to dec2(lat2) at ra2(long2)
(iii) ravel/fly/draw from ra1(long1) to ra2(long2) along dec2(lat2) and then travel from dec2(lat2) to dec1(lat1) at ra2(long2)
A real case is where I want to publish the offsets in both ra, dec from a central object and I would use the dec of this primary object.
taldcroft commentedon Nov 20, 2015
@richardgmcmahon - in astropy coordinates the
rainput value is fundamentally anAngle, which is simply a pure (mathematical) angle in radians. Theraanddecinputs are thealphaandbetaparameters of an Euler angle transform, and again these values are pure angles. This definition is consistent with formal standards in astronomy.If we renamed
raanddecto bealphaandbeta(which is what they really are) then perhaps you would be less surprised at the result when subtractingc1.alpha - c2.alpha. Fundamentally speaking that subtraction operation is not a physically meaningful quantity (i.e. it's unrelated to separation). To attach a physical meaning (separation) you must include the correspondingbetavalues in the computation , butc.alphaas a standalone value does not know whatbetais.To put a more "going-forward" spin on this, what about a separation method that returns a 2-tuple of values which are what you originally wanted, name the local separation of two nearby coordinates as
Angleobjects. I think this would be generally useful for anyone computing proper motions or specifying the 2-d offset of two nearby objects. One needs to be aware that this parametrization of offsets basically fails at the poles, but that is something astronomers have accepted from what I see.richardgmcmahon commentedon Nov 20, 2015
Thanks, the going forward approach looks good. My main concern is the innocent new graduate student making an unforeseen error; some of this is education and documentation.
mhvk commentedon Nov 20, 2015
Agreed that the case of offsets from a location is a common operation -- it is usually what telescope operators ask for! -- but also that in practice the only solution is to have a method. We'd also have to be quite clear about the order in which the offsets would have to be applied for them to be valid; if we use the declination of the first object, then I guess the offset in ra has to be applied first (this would seem to make most sense). Maybe
skycoord1.offset_to(skycoord2)?eteq commentedon Nov 22, 2015
I think the discussion in #4209 is about nearly the same thing as this. @richardgmcmahon, do you agree that the thing the original issue there is asking about (and the subsequent discussion focuses on) would address this? I.e., an "offset" method? (See my thoughts there on the topic, but of course I can certainly be outvoted!)
If yes, then I think this specific issue can be closed, but we should continue the discussion in #4209.
richardgmcmahon commentedon Nov 23, 2015
Yes, an offset method would cover it. However I am also concerned that there is an opportunity for user error if they
follow the catalog matching example here:
http://docs.astropy.org/en/stable/coordinates/matchsep.html#matching-catalogs
dra = (matches.ra - c.ra).arcmin
It is a common practice to analyze the dra and ddec distributions of offsets when catalog matching to look for systematic
offsets so the example should explicitly have cos(dec) specified for dra.
In this case I would suggest that the example is reversed and become:
dra = (c.ra - matches.ra ).arcmin * cos(c.dec)
eteq commentedon Nov 23, 2015
@richardgmcmahon - I see your point there, but that particular docs section is focused on the matching functionality, and I'm a bit concerned that change purely on its own will make it more confusing what's being taught in that example.
Perhaps if we follow the route I suggested in #4209, the solution will just be to replace the
dra =andddec =lines withdra, ddec = spherical_offset(matches, c)? That would seem to be the best of both worlds. Then only if #4209 is instead not addressed in that manner, the example can get updated with a bit more detail including thecosfactor?eteq commentedon Nov 23, 2015
Oh, and in #4338 I made the change I just dicussed above.
richardgmcmahon commentedon Jul 7, 2016
#5163 is an issue that relates to the implemented solution in astropy 1.2; so it looks like this can be closed