The shortest distance between two locations on the surface of Earth (or any planet) is known as the Great Circle Distance. Except for relatively short distances these cannot be measured on a map due to the distortion and flattening necessary to represent a sphere on a flat plane. However, the calculation of the Great Circle Distance from two coordinates is simple, although I suspect generations of midshipmen might not have agreed. In this post I will write a short program in C to calculate the distance from London to various other cities round the world.
As always, Wikipedia provides a comprehensive reference on the subject of calculating Great Circle distances. The principle is very simple, if you understand radians.
Everyone is familiar with using degrees to measure angles - we all gain an understanding of them when we are given our first school protractor. Slightly less well understood is the radian. (If you already know all about radians skip to the next paragraph!) Imagine a circle (for the purposes of this project, imagine the circle is the flat part of a hemisphere, a sphere cut in half) and now imagine drawing a line from the centre to the edge. The length of that line is of course the radius. Now draw round the circumference of the circle by the same distance as the radius, finally drawing another line back to the centre of the circle. The angle between the two lines meeting at the centre is 1 radian. The size of 1 radian in degrees is 180/π which is an irrational number but 57.29577951° to 8dp.
If you skipped my scholarly discourse on radians, welcome back. If you think about it for a moment, it should become clear that if the angle between two points on the surface of our little blue planet is 1 radian, then the great circle distance between them is the same as the radius of the planet. If the angle is 2 radians, the distance is twice the radius; if the angle is 0.5 radians the distance is half the radius and so in. In short, the Great Circle Distance between two points in the angle between them multiplied by the radius. The units of distance are irrelevant as long as they are the same (don't mix miles and kilometres for example) but we must use radians rather than degrees. This is the reason we are using radians - they reduce the problem to one simple multiplication.
We're getting ahead of ourselves a bit here, as we first need to calculate the angle between two points from their latitudes and longitudes. I'll "borrow" the formula for this from Wikipedia, but first let's define some terms.
- Φ1 and λ1 - latitude and longitude in radians of first location
- Φ2 and λ2 - latitude and longitude in radians of second location
- ΔΦ and Δλ - absolute differences between latitudes and longitudes
- Δσ - angle between locations
Don't worry, I won't be using Greek in the code. Now for the formula - note use of dot notation for multiplication.
So basically we need to apply this formula to a pair of coordinates, then multiply the result by the radius of the Earth.
The Code
Create a new folder and within it create these empty files. You can download the code as a zip or pull it from Github if you prefer.
- greatcircle.h
- greatcircle.c
- main.c
Source Code Links
Open greatcircle.h and enter the following.
greatcircle.h
#include<stdbool.h> //-------------------------------------------------------- // STRUCT greatcircle //-------------------------------------------------------- typedef struct greatcircle { char* name1; double latitude1_degrees; double longitude1_degrees; double latitude1_radians; double longitude1_radians; char* name2; double latitude2_degrees; double longitude2_degrees; double latitude2_radians; double longitude2_radians; double central_angle_radians; double central_angle_degrees; double distance_kilometres; double distance_miles; bool valid; } greatcircle; //-------------------------------------------------------- // FUNCTION PROTOTYPES //-------------------------------------------------------- greatcircle* greatcircle_create(char* name1, double latitude1_degrees, double longitude1_degrees, char* name2, double latitude2_degrees, double longitude2_degrees); void greatcircle_print(greatcircle* gc); void greatcircle_calculate(greatcircle* gc); void greatcircle_free(greatcircle* gc);
In greatcircle.h we create a struct to hold the various variables needed. Most pieces of information exist twice: the latitudes, longitudes and angle in both radians and degrees, and the distance in miles and kilometres. There is also a boolean valid flag. The header files also contains four function prototypes which I'll describe as we implement them. Let's move on to doing that now by opening greatcircle.c and entering the following.
greatcircle.c
#include<stdlib.h> #include<stdio.h> #include<string.h> #include<math.h> #include"greatcircle.h" #define DEGREES_IN_RADIAN 57.29577951 #define MEAN_EARTH_RADIUS_KM 6371 #define KILOMETRES_IN_MILE 1.60934 //-------------------------------------------------------- // FUNCTION PROTOTYPES //-------------------------------------------------------- static void validate_degrees(greatcircle* gc); static void calculate_radians(greatcircle* gc); static void calculate_central_angle(greatcircle* gc); static void calculate_distance(greatcircle* gc); //-------------------------------------------------------- // FUNCTION greatcircle_create //-------------------------------------------------------- greatcircle* greatcircle_create(char* name1, double latitude1_degrees, double longitude1_degrees, char* name2, double latitude2_degrees, double longitude2_degrees) { greatcircle* gc = malloc(sizeof(greatcircle)); gc->name1 = malloc(strlen(name1) + 1); strcpy(gc->name1, name1); gc->name2 = malloc(strlen(name2) + 1); strcpy(gc->name2, name2); gc->latitude1_degrees = latitude1_degrees; gc->longitude1_degrees = longitude1_degrees; gc->latitude2_degrees = latitude2_degrees; gc->longitude2_degrees = longitude2_degrees; return gc; } //-------------------------------------------------------- // FUNCTION PROTOTYPES //-------------------------------------------------------- void greatcircle_calculate(greatcircle* gc) { validate_degrees(gc); if(gc->valid) { calculate_radians(gc); calculate_central_angle(gc); calculate_distance(gc); } } //-------------------------------------------------------- // FUNCTION greatcircle_print //-------------------------------------------------------- void greatcircle_print(greatcircle* gc) { if(gc->valid) { printf("%s\n", gc->name1); printf("Latitude: %lf degrees / %lf radians \n", gc->latitude1_degrees, gc->latitude1_radians); printf("Longitude: %lf degrees / %lf radians \n", gc->longitude1_degrees, gc->longitude1_radians); printf("%s\n", gc->name2); printf("Latitude: %lf degrees / %lf radians\n", gc->latitude2_degrees, gc->latitude2_radians); printf("Longitude: %lf degrees / %lf radians\n", gc->longitude2_degrees, gc->longitude2_radians); printf("Valid: %s\n", gc->valid ? "Yes" : "No"); printf("Central Angle %lf radians / %lf degrees\n", gc->central_angle_radians, gc->central_angle_degrees); printf("Distance %lf kilometres / %lf miles\n", gc->distance_kilometres, gc->distance_miles); } else { printf("Latitude and longitude for %s to %s are invalid\n", gc->name1, gc->name2); } puts("---------------------------------------------------------------"); } //-------------------------------------------------------- // FUNCTION greatcircle_free //-------------------------------------------------------- void greatcircle_free(greatcircle* gc) { free(gc->name1); free(gc->name2); free(gc); } //-------------------------------------------------------- // STATIC FUNCTION validate_degrees //-------------------------------------------------------- static void validate_degrees(greatcircle* gc) { gc->valid = true; if(gc->latitude1_degrees < -90.0 || gc->latitude1_degrees > 90.0) { gc->valid = false; } if(gc->longitude1_degrees < -180.0 || gc->longitude1_degrees > 180.0) { gc->valid = false; } if(gc->latitude2_degrees < -90.0 || gc->latitude2_degrees > 90.0) { gc->valid = false; } if(gc->longitude2_degrees < -180.0 || gc->longitude2_degrees > 180.0) { gc->valid = false; } } //-------------------------------------------------------- // STATIC FUNCTION calculate_radians //-------------------------------------------------------- static void calculate_radians(greatcircle* gc) { gc->latitude1_radians = gc->latitude1_degrees / DEGREES_IN_RADIAN; gc->longitude1_radians = gc->longitude1_degrees / DEGREES_IN_RADIAN; gc->latitude2_radians = gc->latitude2_degrees / DEGREES_IN_RADIAN; gc->longitude2_radians = gc->longitude2_degrees / DEGREES_IN_RADIAN; } //-------------------------------------------------------- // STATIC FUNCTION calculate_central_angle //-------------------------------------------------------- static void calculate_central_angle(greatcircle* gc) { double longitudes_abs_diff; if(gc->longitude1_radians > gc->longitude2_radians) { longitudes_abs_diff = gc->longitude1_radians - gc->longitude2_radians; } else { longitudes_abs_diff = gc->longitude2_radians - gc->longitude1_radians; } gc->central_angle_radians = acos(sin(gc->latitude1_radians) * sin(gc->latitude2_radians) + cos(gc->latitude1_radians) * cos(gc->latitude2_radians) * cos(longitudes_abs_diff)); gc->central_angle_degrees = gc->central_angle_radians * DEGREES_IN_RADIAN; } //-------------------------------------------------------- // STATIC FUNCTION calculate_distance //-------------------------------------------------------- static void calculate_distance(greatcircle* gc) { gc->distance_kilometres = MEAN_EARTH_RADIUS_KM * gc->central_angle_radians; gc->distance_miles = gc->distance_kilometres / KILOMETRES_IN_MILE; }
In greatcircle.c the first thing to note is that we #define a few values, DEGREES_IN_RADIAN, MEAN_EARTH_RADIUS_KM and KILOMETRES_IN_MILE, and then prototype four static functions which I'll discuss as we come to them.
The first function is greatcircle_create which takes as arguments the names, latitudes and longitudes of our two locations. It then uses malloc to get us memory for a new greatcircle struct and sets some of the properties from the arguments before returning a pointer to the new struct.
The next function is greatcircle_calculate which is central to the whole project as it calculates the missing values from the latitudes and longitudes in degrees. Before it does this it checks the degrees are valid with the validate_degrees function, and if so calls calculate_radians, calculate_central_angle and calculate_distance.
Next we have greatcircle_print, which simply prints out the various members of the greatcircle struct.
We then have the greatcircle_free function which must be called after we have finished with the greatcircle struct, as it and its two name fields are created using malloc.
Now let's move on to the four static functions called by greatcircle_calculate. The first is validate_degrees which simply checks the latitudes are between -90° and +90°, and that the longitudes are between -180° and +180°.
The calculate_radians function does its stuff by dividing the degrees by DEGREES_IN_RADIAN which was #defined at the top of the file.
The calculate_central_angle function calculates the absolute difference between the longitudes, and then implements this formula
which we saw above to actually calculate the angle in radians. Finally it calculates the degree equivalent.
The last function in this file is calculate_distance, which as I have mentioned is simply a matter of multiplying the angle in radians by the radius of the earth. We also do the kilometres to miles conversion here.
We can now move on to the main function which uses the code we have written for seven different cities. For each city pair we need to do the following:
- Create a new greatcircle struct with greatcircle_create
- Calculate the radians and distances with greatcircle_calculate
- Print out the results with greatcircle_print
- Free up the memory with greatcircle_free
The last set of values includes the invalid 91° north just to see what happens...!
main.c
#include<stdio.h> #include<math.h> #include<time.h> #include<locale.h> #include<stdlib.h> #include "greatcircle.h" //-------------------------------------------------------- // FUNCTION main //-------------------------------------------------------- int main(int argc, char* argv[]) { puts("--------------------------"); puts("| codedrome.com |"); puts("| Great Circle Distances |"); puts("--------------------------\n"); greatcircle* gc_london_tokyo = greatcircle_create("London", 51.507222, -0.1275, "Tokyo", 35.683333, 139.683333); greatcircle_calculate(gc_london_tokyo); greatcircle_print(gc_london_tokyo); greatcircle_free(gc_london_tokyo); greatcircle* gc_london_new_york = greatcircle_create("London", 51.507222, -0.1275, "New York", 40.7127, -74.0059); greatcircle_calculate(gc_london_new_york); greatcircle_print(gc_london_new_york); greatcircle_free(gc_london_new_york); greatcircle* gc_london_new_delhi = greatcircle_create("London", 51.507222, -0.1275, "New Delhi", 28.613889, 77.208889); greatcircle_calculate(gc_london_new_delhi); greatcircle_print(gc_london_new_delhi); greatcircle_free(gc_london_new_delhi); greatcircle* gc_london_sydney = greatcircle_create("London", 51.507222, -0.1275, "Sydney", -33.865, 151.209444); greatcircle_calculate(gc_london_sydney); greatcircle_print(gc_london_sydney); greatcircle_free(gc_london_sydney); greatcircle* gc_london_cape_town = greatcircle_create("London", 51.507222, -0.1275, "Cape Town", -33.925278, 18.423889); greatcircle_calculate(gc_london_cape_town); greatcircle_print(gc_london_cape_town); greatcircle_free(gc_london_cape_town); greatcircle* gc_london_rio_de_janeiro = greatcircle_create("London", 51.507222, -0.1275, "Rio de Janeiro", -22.908333, -43.196389); greatcircle_calculate(gc_london_rio_de_janeiro); greatcircle_print(gc_london_rio_de_janeiro); greatcircle_free(gc_london_rio_de_janeiro); greatcircle* gc_london_oblivion = greatcircle_create("London", 51.507222, -0.1275, "Oblivion", 91, 360); greatcircle_calculate(gc_london_oblivion); greatcircle_print(gc_london_oblivion); greatcircle_free(gc_london_oblivion); return EXIT_SUCCESS; }
The code is now finished so we can compile and run it with the following commands.
Compile and Run
gcc main.c greatcircle.c -std=c11 -lm -o main ./main
Program Output (partial)
-------------------------- | codedrome.com | | Great Circle Distances | -------------------------- London Latitude: 51.507222 degrees / 0.898971 radians Longitude: -0.127500 degrees / -0.002225 radians Tokyo Latitude: 35.683333 degrees / 0.622792 radians Longitude: 139.683333 degrees / 2.437934 radians Valid: Yes Central Angle 1.500399 radians / 85.966540 degrees Distance 9559.043090 kilometres / 5939.728765 miles
As well as the names, latitudes and longitudes we supplied we now have the latitudes and longitudes in radians, the central angle in both degrees and radians, and the distance in kilometres and miles.