UVa 11817 - Tunnelling the Earth

将平面上的亮点坐标转化成空间坐标，为球面上的两点。求所对应的弧长与弦长的差值，转化为空间坐标的函数为

double torad(double deg) {
  return deg / 180 * M_PI;
}
void get_coord(double lat, double lng, double& x, double& y, double& z) {
  lat = torad(lat);
  lng = torad(lng);
  x = R * cos(lat) * cos(lng);
  y = R * cos(lat) * sin(lng);
  z = R * sin(lat);
}

/*************************************************************************
 > File Name: 11817.cpp
 > Author: AcToy
 > Mail: ycsgldy@163.com 
 > Created Time: 2013年07月19日 星期五 19时45分23秒
 ************************************************************************/

#include <algorithm>
#include <iostream>
#include <cstring>
#include <cstdlib>
#include <climits>
#include <cstdio>
#include <string>
#include <vector>
#include <queue>
#include <cmath>
#include <stack>
#include <map>
#include <set>

using namespace std;
typedef unsigned int u32;
typedef long long i64;
typedef unsigned long long u64;
typedef vector<int> IV;
typedef vector<bool> BV;
typedef pair<int,int> II;
typedef vector<II> IIV;
#define For(t,v,c) for(t::const_iterator v=c.begin(); v!=c.end(); ++v)
const int INF = 0x7FFFFFFF;
const double eps = 1E-10;
const double PI = acos(-1);
const double R = 6371009;
double torad(double deg) {
	return deg / 180 * M_PI;
}
void get_coord(double lat, double lng, double& x, double& y, double& z) {
	lat = torad(lat);
	lng = torad(lng);
	x = R * cos(lat) * cos(lng);
	y = R * cos(lat) * sin(lng);
	z = R * sin(lat);
}
double sqr(double x) { return x * x; }
double dist(double x1, double y1, double z1, double x2, double y2, double z2) {
	return sqrt(sqr(x1 - x2) + sqr(y1 - y2) + sqr(z1 - z2));
}
int main() {
	int Case;
	double lat1, lng1, lat2, lng2, d1, d2;
	double x1, y1, z1, x2, y2, z2;
	scanf("%d", &Case);
	while(Case--) {
		scanf("%lf%lf%lf%lf", &lat1, &lng1, &lat2, &lng2);
    get_coord(lat1, lng1, x1, y1, z1);
		get_coord(lat2, lng2, x2, y2, z2);
		d1 = dist(x1, y1, z1, x2, y2, z2);
		d2 = 2 * asin(d1 / (2 * R)) * R;
		printf("%.0lf\n", d2 - d1);
	}
  return 0;
}