# Work the Shell - Calculating the Distance between Two Latitude/Longitude Points

in
Finding your way at the command line.
Grafting It All Together

To get everything to work well, I actually hacked and slashed at the original script to make it a bit more succinct and, of course, invoke the C “distance” program as shown in Listing 1. [Listing 1 also is available on our FTP site at ftp.linuxjournal.com/pub/lj/listings/issue188/10606.tgz.] Ready? It's surprisingly short:

```#!/bin/sh
converter="http://api.maps.yahoo.com/ajax/
↪geocode?appid=onestep&qt=1&id=m&qs="

tmpfile="/tmp/bc.script.\$\$"

# Get lat/long for point 1
addr="\$(echo \$1 | sed 's/ /+/g')"
values="\$(curl -s \$converter\$addr | \
cut -d\" -f13,15 | \
sed 's/[^0-9\.\,\-]//g;s/,\$//')"

lat1=\$(echo \$values | cut -d, -f1)
long1=\$(echo \$values | cut -d, -f2)

# Get lat/long for point 2
addr="\$(echo \$2 | sed 's/ /+/g')"
values="\$(curl -s \$converter\$addr | \
cut -d\" -f13,15 | \
sed 's/[^0-9\.\,\-]//g;s/,\$//')"

lat2=\$(echo \$values | cut -d, -f1)
long2=\$(echo \$values | cut -d, -f2)

# Now we have the lat/long for both points, let's
# figure out the distance between them...
dist=\$(./distance \$lat1 \$long1 \$lat2 \$long2)
echo "\$1 to \$2 is \$dist miles"
exit 0
```

The script would be even shorter if we tweaked the C program to accept x,y location pairs, but I'll leave that one to you. Instead, let's do a few tests:

```\$ farapart.sh \
"union station, denver, co" \
"union station, chicago, il"
union station, denver, co to
union station, chicago, il is 917.984 miles
```

Now, how about something a bit more ambiguous:

```\$ farapart.sh "long beach, ca" "boston, ma"
long beach, ca to boston, ma is 2597.53 miles
```

Well, darn it, that seems way too short. Let's see what Yahoo Maps reports as the distance between those two cities. Sure enough, it reports that the trip should be 3,015 miles, not 2,597 miles.

Debugging the Math Formula

Somewhere there's an error that's giving us poor results. My guess is there's some sort of significant rounding error going on in the C program (because we can verify experimentally that the lat/lon information we're getting is valid, simply by plugging it in to a mapping app and seeing where it places us).

I'm all tapped out on this example, however. It turned out to be far more tricky than I anticipated, and I leave it as an exercise to you, dear reader, to see if you can figure out what's broken in the C program and report your fix to us. We'll publish the best of them next month! Meanwhile, next column, I'll get back to something that's more about the shell and less about mathematics. I mean, heck, I didn't like math when I was working on my computer science degree, so why am I playing with it now?

Dave Taylor has been involved with UNIX since he first logged in to the on-line network in 1980. That means that, yes, he's coming up to the 30-year mark now. You can find him just about everywhere on-line, but start here: www.DaveTaylorOnline.com. In addition to all his other projects, Dave is now a film critic. You can read his reviews at www.DaveOnFilm.com.

______________________

Dave Taylor has been hacking shell scripts for over thirty years. Really. He's the author of the popular "Wicked Cool Shell Scripts" and can be found on Twitter as @DaveTaylor and more generally at www.DaveTaylorOnline.com.

## Comment viewing options

Select your preferred way to display the comments and click "Save settings" to activate your changes.

### Problem of distance variation

I have wrote a bash script (using bc + atan2 definition I found online) to calculate both distances between tree points and angle between two distant points using third as a reference point.

Problem is not your C code, problem is number of decimal places used to calculate.

Using 3 point only 0,5-2Km apart and using 4 to 15 decimal places, I got 15-25% difference in distance calculations. So what your code needs are variables with more decimal points.

Hopes this helps.

### in all the code samples

in all the code samples above, "c = 2 * atan2(sqrt(a), sqrt(1-a));" can be written simpler as "c = 2 * asin(sqrt(a));"

### listing 1 link fix

It is actually ftp://ftp.linuxjournal.com/pub/lj/listings/issue188/10606l1.txt and not ftp://ftp.linuxjournal.com/pub/lj/listings/issue188/10606.tgz as mentioned above :)
great article! Thanks! - Israel Torres

### No bugs I believe, just works as it should

Hello,
I was surprised but the answer is: there is no bug in C program code, the mileage is CORRECT.

The shortest distance between 2 points on earth surface of course is not a straight line, that would be impossible unless you are able to dig a tunnel going underneath!

Remember, when talking about earth you have to think 3D!

This interesting web page shows how the haversine formula works (http://www.movable-type.co.uk/scripts/latlong.html) and a javascript implements the formula on the web page...
Of course the results are shown in Km and not in Miles because Science is polite with numbers...
Anyway using the coordinates Long Beach and Boston the result is 4180 km (that is ~2597 Mi) that is correct.
The hypothetically shortest way to go from East to West on the Northern Earth surface wouldn't be walking on a Parallel (Meridians and Parallels speaking), but would be walking on an arch pushed towards North (for example it would be a route similar to the following Boston-Michigan-NorthDakota-Utha-California-LongBeach).

Just think about what aiplanes do to save as much fuel as possible. Their routes are very different from what you can imagine...

This is php code I use to calculate (in Km) the shortest distance between 2 points on Earth surface (it is raw but it works):

``` <?php function getdistance(\$lat1,\$lon1,\$lat2,\$lon2) { \$R = 6371; \$dLat = (\$lat2 - \$lat1) * 3.14159265359 / 180 ; \$dLon = (\$lon2 - \$lon1) * 3.14159265359 / 180 ; \$lat1 = \$lat1 * 3.14159265359 / 180 ; \$lat2 = \$lat2 * 3.14159265359 / 180 ; \$a = sin(\$dLat/2) * sin(\$dLat/2) + cos(\$lat1) * cos(\$lat2) * sin(\$dLon/2) * sin(\$dLon/2); \$c = 2 * atan2(sqrt(\$a),sqrt(1-\$a)); \$d = \$R * \$c; return \$d; } \$LAT1 = \$_REQUEST["lat1"]; \$LON1 = \$_REQUEST["lon1"]; \$LAT2 = \$_REQUEST["lat2"]; \$LON2 = \$_REQUEST["lon2"]; print (getdistance(\$LAT1,\$LON1,\$LAT2,\$LON2)); ?> ```

You can either call it via POST or via GET.
In fact _REQUEST is magic. :)

take care,
Claudio

### code

hi, I am gayatri i am writing linux c code for distance calculation using haversine formula i am not getting correct if there is any errors plz correct it and mail me plz.

```#include
#include
main()
{
//convert degrees to radians float dtor(float degrees)
{
return(degrees *pi/180);
}

//convert radians to degrees float rtod(float radians)
{
}
//calculate distance form lat1/lon1 to lat2/log2 using haversine formula

//Note lat1/lon1/lat2/lon2 must be in radians

//returns float distance in feet
float calcDistance(float lat1, float lon1, float lat2, float lon2)
{

float dlon,dlat,a,c;
float dist=0.0;

lon1=lon1*-1.0;//make west=positive

lon2=lon2*-1.0;
dlon=lon2-lon1;
dlat=lat2-lat1;
a=pow(sin(dlat/2),2)+cos(lat1)*cos(lat2)*pow(sin(dlon/2),2);
c=2*atan(sqrt(a),sqrt(1-a));
dist=20925656.2*c; //radius of the earth (6378140 meters) in feet 20925656.2
return(dist);
}

```
Webcast
How to Build an Optimal Hadoop Cluster to Store and Maintain Unlimited Amounts of Data Using Microservers

Realizing the promise of Apache® Hadoop® requires the effective deployment of compute, memory, storage and networking to achieve optimal results. With its flexibility and multitude of options, it is easy to over or under provision the server infrastructure, resulting in poor performance and high TCO. Join us for an in depth, technical discussion with industry experts from leading Hadoop and server companies who will provide insights into the key considerations for designing and deploying an optimal Hadoop cluster.