# How does ArcGIS compute the distance between two points with a non-equidistant projection?

This is a follow-up question to my previous one, Can you suggest some well-written introductory texts about coordinate system projections?

Let's assume I'm working with the CH1903 map projection, which for all I know is conformal, but not equidistant. Meaning, angles (shape) have been preserved, but not areas, distances, or scale. (At least these have not been preserved exactly). So far so good.

I'm wondering what kind of calculation ArcGIS performs when I now want to calculate the distance between two points. In ArcObjects, I could use the IProximityOperator interface as follows:

IPoint a = ...,
b = ...;

double distance = ((IProximityOperator)a).ReturnDistance(b);

Question: When I'm working with a reference system that does not accurately preserve distances, what would ArcGIS do when I query it for the distance between two points (as shown above)?

• Does it simply do some Pythagorean maths (a2 + b2 = c2) to get the distance, meaning the returned distance will only be as accurate as the projection allows?

• Or will it do something more complicated, like some form of re-projection, to get a more accurate distance?

(The same question, but more generally: Once that geometries have been projected, does ArcGIS perform all calculations simply in Euclidean space, or does the used map projection still influence calculations of distances, angles, areas, etc.?)

• Please create a new question instead of modifying your original. Otherwise, you subvert all the mechanisms on this site: what do ratings mean when two or more questions are in play in one thread? What would it mean to mark one answer as correct? Etc. Commented Aug 23, 2010 at 15:15
• @whuber: Although all that's been written in this thread is still on-topic WRT the original question asked, I agree that there are now really two questions asked here. It's too late to change that now, but will keep your advice in mind for a next time. Commented Aug 23, 2010 at 21:06

If you want a stable method of computing geodesic distances, I recommend Richie Carmichael's wrapper for ESRI's Projection Engine.

Update: I just tried Richie's code with ArcGIS 10.0 on Vista64 and get an exception after calling LoadLibrary. I'll look into that more later.

For now though, here is some code in response to questions in the comments of another answer.

The code compares IProximityOperator for points with and without spatial references. Then it shows how to use an azimuthal equidistant projection (with first point being the point of tangency) to find the great circle distance.

private void Test()
{
IPoint p1 = new PointClass();
p1.PutCoords(-98.0, 28.0);

IPoint p2 = new PointClass();
p2.PutCoords(-78.0, 28.0);

Debug.Print("Euclidian Distance {0}", EuclidianDistance(p1, p2));
Debug.Print("Distance with no spatialref {0}", GetDistance(p1, p2));

ISpatialReferenceFactory srf = new SpatialReferenceEnvironmentClass();
IGeographicCoordinateSystem gcs =
srf.CreateGeographicCoordinateSystem((int)esriSRGeoCSType.esriSRGeoCS_WGS1984);

p1.SpatialReference = gcs;
p2.SpatialReference = gcs;

Debug.Print("Distance with spatialref {0}", GetDistance(p1, p2));
Debug.Print("Great Circle Distance {0}", GreatCircleDist(p1, p2));

}
private double GetDistance(IPoint p1, IPoint p2)
{
return ((IProximityOperator)p1).ReturnDistance(p2);
}

private double EuclidianDistance(IPoint p1, IPoint p2)
{
return Math.Sqrt(Math.Pow((p2.X - p1.X),2.0) + Math.Pow((p2.Y - p1.Y), 2.0));
}

private double GreatCircleDist(IPoint p1, IPoint p2)
{
ISpatialReferenceFactory srf = new SpatialReferenceEnvironmentClass();
IProjectedCoordinateSystem pcs =
srf.CreateProjectedCoordinateSystem((int)esriSRProjCSType.esriSRProjCS_WGS1984N_PoleAziEqui);
pcs.set_CentralMeridian(true, p1.X);
((IProjectedCoordinateSystem2)pcs).LatitudeOfOrigin = p1.Y;
p1.SpatialReference = pcs.GeographicCoordinateSystem;
p1.Project(pcs);
p2.SpatialReference = pcs.GeographicCoordinateSystem;
p2.Project(pcs);
return EuclidianDistance(p1, p2);
}

Here's the output:

Euclidian Distance 20
Distance with no spatialref 20
Distance with spatialref 20
Great Circle Distance 1965015.61318737

I think it would be interesting to test this against the projection engine dll (pe.dll). Will post results if I ever get Richie's code to work.

Update: Once I changed Richies code to compile for x86, I got it to run. Interesting ... the great circle distance it give me is 1960273.80162999 - a significant difference from that returned from the azimuthal equidistant method above.

• The reason for the discrepancy is likely because the line segment (in the PCS) connecting the points is not a projected geodesic, which would be curvilinear when projected. Accordingly you get a smaller value than you should. A test of this theory is simple to do: take a simple geodesic (such as the equator) and compare two calculations of the distance between two widely separated points on the geodesic. One is the direct calculation, as in your code; the other breaks the geodesic into segments, directly calculates the segment lengths, and adds them up. The latter should be more accurate. Commented Aug 23, 2010 at 15:13

In ArcGIS 10, check out IGeometryServer2 which now has GetDistanceGeodesic (geodesic distance between two geometries), GetLengthsGeodesic (return the geodesic lengths of each polyline), and DensifyGeodesic (densify a polyline by plotting points along the geodesic lines connecting the vertices, uses IPolycurve4::GeodesicDensify) methods.

As mentioned in the other answers, ArcGIS still uses mostly planar calculations.

Melita Kennedy

Some commentary on the other answers (not enough rep yet to comment directly!).

Esri's azimuthal equidistant projection supports ellipsoids. The GreatCircleDist code is creating a PCS that uses an ellipsoid/spheroid-based GCS, thus distances from the center/origin point will be geodesic distances, not great circle distances. It could also be simplified. We know the projected coords of the first point because it's the center of the projection: 0,0. So only the 2nd point needs to be projected. A simplified EuclidianDistance function could then be used.

I checked the results against the pe.dll's geodesic functions and it matched. It looks like Richie's app is using a sphere, so it is returning great circle distances/coordinates in its test application. That's why the results don't match. I didn't recognize the radius values; I think I need to talk to him about it!

• Melita - Great to see you here! Commented Oct 7, 2010 at 3:00
• I agree, welcome aboard! Commented Oct 13, 2010 at 21:48

The accuracy of any answer about ArcGIS is subject to change at any time--for all we know, new procedures will be introduced in the very next service pack without warning or documentation. That being said, ESRI software has for a long time used Euclidean calculations (e.g., the Pythagorean formula for distances) whenever projected coordinates are used. Often, in calculations like the one you illustrate, the software doesn't even have access to the projection information, so what else can be done?

Your question itself seems to suggest that Euclidean distance calculations for an equidistant projection are correct. Nothing could be further from the truth. For a one-point equidistant projection, the Euclidean distance to the base point is guaranteed to equal the geodesic distance; for a two-point equidistant projection, the Euclidean distance to either base point is guaranteed to equal the geodesic distance. In return for those guarantees, the metric distortion between all other pairs of points is typically greatly increased compared to other projections one might choose.

• @whuber: Thanks for answering. Concerning the first paragraph: I figured ArcGIS might see that the CH1903 map projection (which uses the Bessel 1841 ellipsoid) is used, and then project the points back onto that ellipsoid via the datum, and then do distance calculations on the ellipsoid. From your answer I take that ArcGIS won't do all that and will remain in Euclidean X-Y space to do calculations. (How about other GIS software?) -- 2nd paragraph: You're right of course, thanks for clarifying this point. Commented Aug 20, 2010 at 18:32
• A hidden reprojection mechanism is possible only if the point objects maintain references to a projection. I don't believe they do. Commented Aug 20, 2010 at 18:37
• @whuber: Would it be sufficient (for more accurate calculations) to know the ellipsoid used for the projection? AFAIK, ArcGIS stores a reference to the used projection with each feature class (data layer). Commented Aug 21, 2010 at 16:00
• Actually IPoint, which derives from IGeometry, has SpatialReference as a property. help.arcgis.com/en/sdk/10.0/arcobjects_net/componenthelp/… However, I don't think ReturnDistance uses it. It might be worth testing to see if that's changed though. Commented Aug 21, 2010 at 20:49
• @stakx I've updated my answer to include code that shows that setting spatialref has no impact on ReturnDistance. Commented Aug 23, 2010 at 3:36