# Calculating distance between two points using ArcPy

I have 4 coordinates on map.

``````lat1 = radians(47.434756)
``````

I want calculate distance between this.

``````from math import sin, cos, sqrt, atan2, radians
R = 6373.0

dlon = lon2 - lon1
dlat = lat2 - lat1

a = sin(dlat / 2)**2 + cos(lat1) * cos(lat2) * sin(dlon / 2)**2
c = 2 * atan2(sqrt(a), sqrt(1 - a))

distance = R * c/100

print("Result:", distance)
``````

My code say that distance is `0.1212`. ArcGIS says `0,11062` (shape length).

But when I measured with a ruler - its shows 9.23km.

How do I correctly calculate distance with ArcPy? There's no need to use spherical haversine computation when ArcPy already has the full Inverse (aka Reverse) Problem of Geodesy implemented via the Esri Projection Engine.

The key here is to use the arcpy.PointGeometry type (`arcpy.Point` is a simple helper class to store double precision values), with an explicit `SpatialReference` parameter, set to your preferred geographic coordinate system. This example uses WGS84:

``````import arcpy

sr = arcpy.SpatialReference(4326)               # WGS84
sr.setFalseOriginAndUnits(-400,-400,10000000)   # ten millionth of a degree ~= 1.1112 cm resolution

p1 = arcpy.PointGeometry(arcpy.Point(42.971468,47.434756),sr)  # always use spatial_reference!
p2 = arcpy.PointGeometry(arcpy.Point(42.946061,47.542419),sr)
a1,d1 = p1.angleAndDistanceTo(p2,'GEODESIC')
a2,d2 = p2.angleAndDistanceTo(p1,'GEODESIC')

print("Distance = {:8.3f} meters, Bearing1->2 = {:9.5f} degrees".format(d1,a1))
print("Distance = {:8.3f} meters, Bearing2->1 = {:9.5f} degrees".format(d2,a2))
``````

The results are in degrees and meters (unless the SpatialReference is projected, in which case it's in the units of the projection):

``````Distance = 12122.175 meters, Bearing1->2 =  -9.07843 degrees
Distance = 12122.175 meters, Bearing2->1 = 170.90285 degrees
``````

Note that, even for such close points, the angle of 2->1 is not a full complement of the angle from 1->2, since the path follows the curvature of the Earth.

If the lat and lon are flipped from the values you provided (`Point` takes x,y parameters, which is lon,lat)...

``````p3 = arcpy.PointGeometry(arcpy.Point(47.434756,42.971468),sr)
p4 = arcpy.PointGeometry(arcpy.Point(47.542419,42.946061),sr)
a3,d3 = p3.angleAndDistanceTo(p4,'GEODESIC')
a4,d4 = p4.angleAndDistanceTo(p3,'GEODESIC')

print("Distance = {:8.3f} meters, Bearing3->4 = {:9.5f} degrees".format(d3,a3))
print("Distance = {:8.3f} meters, Bearing4->3 = {:9.5f} degrees".format(d4,a4))
``````

then the distance is 9.23 kilometers.

``````Distance = 9227.107 meters, Bearing3->4 = 107.77524 degrees
Distance = 9227.107 meters, Bearing4->3 = -72.15139 degrees
``````