Calculating distance between 2 points on Earth using the Haversine formula

2016 Update
The code and the Xcode project for this article is now also available on GitHub.

A while ago, I worked on a project where I had to implement code in C/C++ that could calculate the distance between 2 locations within the United States. It would be nice if the Earth was flat, since this would make my job a little easier, but the Earth is round and I needed to find a way that I could calculate distance. From the research that I have done, I found that plenty of people have worked on this problem before and one popular article is sited over and over again is this: R. W. Sinnott, “Virtues of the Haversine”, Sky and Telescope 68 (2), 159 (1984).

The Math used in all of the articles I read refer to the Haversine formula. Every example I found of this formula were either written in plain Math or shown in JavaScript. Although it was nice finding examples of this in JavaScript, there were a couple of problems with this:

1. C datatypes aren’t the same as JavaScript datatypes and
2. All examples used Latitude and Longitude, but I had to use Zip codes

So I had some more hurdles to overcome. They didn’t seem like big ones, but ones that I had to overcome nonetheless.

Now the funny thing about this project was that this was one of those endeavors where a spec wasn’t available before or during development. As a consequence, I would have designed the code differently than what I had have delivered. Not that there is anything wrong with what I ended up with, it’s just that I like neat code and what I wrote for this article is what I would have done if I had known all of the variables.

Example Code

What I am going to show you is how to calculate distance between two points within the United States, using the Haversine formula, implemented in C. The code will show you how to do that using Latitude and Longitude coordinates. The example application that is posted will show you how to use Zip codes in place of Latitude and Longitude.

The example code is demonstrated in an application written for Mac OS X. Win32 developers shouldn’t have to worry, since the Math code is done in a separate C file. Objective-C is used only for the UI file so you won’t have to read my Objective-C code to understand how to implement this formula in your own application. For the zip code implementation, I use the standard template library to contain the information I need in map so that I can extract the Latitude and Longitude for a particular postal code.

Creating the example code required these steps:

1. Creating a Cocoa Application
2. Creating the Haversine C code and
3. Creating code that can handle Zip codes and finally
4. Creating a routine to check the validity of the output

Creating a Cocoa Application

This is pretty simple step that most Apple developers are familiar with. I’m using Xcode 3.0, however, I set the project settings to be compatible with Xcode 2.4. Once the Cocoa application was created, I needed to create an interface that would take my Longitude and Latitude parameters and display the resulting output.

I also created a section that would take in Zip codes and get the result that way too! The logic is pretty simple and contained within the Objective-C code.

Creating the Haversine C code

This part was a little tough to work with. Every example I saw on the web used JavaScript. Here is the Math and JavaScript code I found at http://www.movable-type.co.uk/scripts/latlong.html :

Haversine Formula:
R = earth’s radius (mean radius = 6,371km)
?lat = lat₂? lat₁
?long = long₂? long₁

a = sin²(?lat/2) + cos(lat₁) * cos(lat₂) * sin²(?long/2)
c = 2 * atan2(?a, ?(1?a))

d = R * c

(Note that angles need to be in radians to pass to trig functions).

JavaScript Implementation:

var R = 6371; // km
var dLat = (lat2-lat1).toRad();
var dLon = (lon2-lon1).toRad();
var a = Math.sin(dLat/2) * Math.sin(dLat/2) +

Math.cos(lat1.toRad()) * Math.cos(lat2.toRad()) * Math.sin(dLon/2) * Math.sin(dLon/2);

var c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1-a));
var d = R * c;

Now for the C code!

To convert this code, there were three things I had to do to implement this as C code: A. Include the <math.h> file into my C file so that I get pi. B. I had to convert our Latitude and Longitude differences into radians, by multiplying the result by ( pi / 180 ) and C. I used the pow function to get the squared value instead of doing multiplying my value twice. Here is what I got:

double CalculateDistance( double nLat1, double nLon1, double nLat2, double nLon2 )
{
    int nRadius = 6371; // Earth's radius in Kilometers
    // Get the difference between our two points
    // then convert the difference into radians
    double nDLat = (nLat2 - nLat1) * (M_PI/180);
    double nDLon = (nLon2 - nLon1) * (M_PI/180);
    double nA = pow ( sin(nDLat/2), 2 ) + cos(nLat1) * cos(nLat2) * pow ( sin(nDLon/2), 2 );

    double nC = 2 * atan2( sqrt(nA), sqrt( 1 - nA ));
    double nD = nRadius * nC;

    return nD; // Return our calculated distance
}

Check out the application and you will see that it does NOT do the conversion properly. Here is where I had to do some digging around to see what was happening. The code on the website said one thing but the JavaScript code implemented on the website said another thing! I found that there were two lines omitted from the webpage: lat1 = lat1.toRad(), lat2 = lat2.toRad();. Of course, I had to add this to my code:

double CalculateDistance( double nLat1, double nLon1, double nLat2, double nLon2 )
{
    double nRadius = 6371; // Earth's radius in Kilometers
    // Get the difference between our two points
    // then convert the difference into radians

    double nDLat = ToRad(nLat2 - nLat1);
    double nDLon = ToRad(nLon2 - nLon1);

    // Here is the new line
    nLat1 =  ToRad(nLat1);
    nLat2 =  ToRad(nLat2);

    double nA = pow ( sin(nDLat/2), 2 ) + cos(nLat1) * cos(nLat2) * pow ( sin(nDLon/2), 2 );

    double nC = 2 * atan2( sqrt(nA), sqrt( 1 - nA ));
    double nD = nRadius * nC;

    return nD; // Return our calculated distance
}

Anyone who is looking the JavaScript code will wonder where to get the toRad function. This was actually declared within  the webpage, and I implemented my own C version called ToRad.

Now let’s move onto incorporating Zip codes into our example!

Creating code that can handle Zip codes

Now, to be able to use Zip codes in our example means that we have to be able to map our zip code to a Latitude/Longitude pair. Fortunately for us, there is a complete listing of postal codes at http://sourceforge.net/projects/zips/ , complete with town, state and postal code. All we have to do, is to read the values into our application and then create a routine for getting our Latitude/Longitude with a given Zip code.

For those who don’t know, Zip is actually an acronym for Zoning Improvement Plan. If you do a little research, you will find that there is an interesting history behind this 5 digit number!

In order for us to use our text file, there are a couple of things we have to do:

1. Read the text file contents into memory.
2. Parse the file contents into a container variable.
3. Create a function that allows us to easy obtain the Latitude and Longitude for a given Zip code.

Since the above listed steps can be accomplished using cross platform friendly POSIX code, I will implement the C code in a C++ file to keep things nice and neat. The reason for putting my code in a C++ file is because I wanted to create a C++ struct that would hold my Latitude Longitude pair. Using a C++ struct allows me to create private and public struct members and even a constructor.

Create a routine to check the validity of the output

This part was pretty easy. Since the JavaScript implementation works, I created a simpler version of the JavaScript code and posted that page at the end of this article. Using that code, and Google, I can see that the application is outputting the right distance.

I hoped you enjoyed reading this article as much as I did writing it 😉 Check out the sample code and happy coding!

Haversine application for Mac OS X

JavaScript Implementation of Haversine Formula

C Source Code and Xcode Project