Fórmula de Haversine

Desenvolvendo um código para estimar a distância entre dois pontos em uma esfera, deparei-me com alguns conceitos, como a Fórmula de Haversine . Se você está procurando uma forma de converter graus de latitude e longitude em unidades métricas e não se importar com um resultado aproximado, talvez esse post seja de interesse.

/*
	Docs:
	http://andrew.hedges.name/experiments/haversine/
	http://www.movable-type.co.uk/scripts/latlong.html
	http://en.wikipedia.org/wiki/Haversine_formula
*/

#include <iostream>
#include <cmath>

using namespace std;

double calcularDistancia(double lat1, double lon1, double lat2, double lon2) {
	double PI = 4.0*atan(1.0);
	double R = 6371000; //circunfência aprox. da terra em metros
	double dLat = (lat2-lat1)*(PI/180);
	double dLon = (lon2-lon1)*(PI/180);
	lat1 = lat1*(PI/180);
	lat2 = lat2*(PI/180);
	double a = sin(dLat/2) * sin(dLat/2) + sin(dLon/2) * sin(dLon/2) * cos(lat1) * cos(lat2);
	double c = 2 * atan2(sqrt(a), sqrt(1-a));
	double d = R * c;
	return d;
}

int main() {
	cout.precision(15);
	cout << "Distancia entre 1 grau em x direcao negativa: " << calcularDistancia(-23,-51,-24,-51) << endl;
	cout << "Distancia entre 1 grau em y direcao negativa: " << calcularDistancia(-23,-51,-23,-52) << endl;
	cout << "Distancia entre 1 grau em x direcao positiva: " << calcularDistancia(-23,-51,-22,-51) << endl;
	cout << "Distancia entre 1 grau em y direcao positiva: " << calcularDistancia(-23,-51,-23,-50) << endl;
	cout << "Distancia de 1 metro na direcao negativa em x: " << calcularDistancia(-23,-51,-23-8.993210126E-6,-51) << endl;
	cout << "Distancia entre 1 grau em x direcao negativa: " << calcularDistancia( -28.3209481381034, -52.732377, -29.3209481381034, -52.732377) << endl;
	return 0;
}