Go: // Great-circular distance
// Reference: https://ja.wikipedia.org/wiki/大円距離

package main

import (
	"bufio"
	"os"
	"strings"
	"strconv"
	"math"
	"fmt"
)

var scanner = bufio.NewScanner(os.Stdin)
const (
	R = 6378.1
	M = math.Pi/180
)

func main() {
	for scanner.Scan() {
		buf := strings.Split(scanner.Text(), " ")
		input := make([]float64, len(buf))
		for i := range input {
			input[i], _ = strconv.ParseFloat(buf[i], 64)
		}
		if input[0] == -1 && input[1]== -1&&
			input[2] == -1 && input[3] == -1 {
			break
		}
		y1, x1, y2, x2 := input[0], input[1], input[2], input[3]
		a := math.Cos(y1*M)*math.Cos(y2*M)*math.Cos((x1-x2)*M) +
			math.Sin(y1*M)*math.Sin(y2*M)
		fmt.Println(int(R*math.Acos(a)+0.5))
	}
}

C++: #include <cstdio>
#include <cmath>
#define r 6378.1

int main(void){
	double a,b,c,d;
	while(1){
		scanf("%lf %lf %lf %lf",&a,&b,&c,&d);
		if(a+b+c+d<-2.0)break;
		a=90.0-a;
		c=90.0-c;
		double x1=cos(b*M_PI/180.0)*sin(a*M_PI/180.0);
		double x2=cos(d*M_PI/180.0)*sin(c*M_PI/180.0);
		double y1=sin(b*M_PI/180.0)*sin(a*M_PI/180.0);
		double y2=sin(d*M_PI/180.0)*sin(c*M_PI/180.0);
		double z1=cos(a*M_PI/180.0);
		double z2=cos(c*M_PI/180.0);
		double e=sqrt(pow(x1-x2,2.0)+pow(y1-y2,2.0)+pow(z1-z2,2.0));
		double ang=acos(1.0-e*e/2.0);
		printf("%.0f\n",r*ang+0.5);
	}
	return 0;
}

C++: #include <iostream>
#include <algorithm>
#include <cmath>
using namespace std;

const double pi = acos(-1);

int main() {
  double a, b, c, d;
  while(cin>>a>>b>>c>>d && !(a==-1&&b==-1&&c==-1&&d==-1)) {
    double x1,y1,z1,x2,y2,z2,rad;
    a = pi*a/180.0;
    b = pi*b/180.0;
    c = pi*c/180.0;
    d = pi*d/180.0;
    x1 = cos(b)*cos(a);
    z1 = sin(b)*cos(a);
    y1 = sqrt(1.0-x1*x1-z1*z1)*(a>0.0?1:-1);
    x2 = cos(d)*cos(c);
    z2 = sin(d)*cos(c);
    y2 = sqrt(1.0-x2*x2-z2*z2)*(c>0.0?1:-1);
    rad = acos(x1*x2+y1*y2+z1*z2);
    cout << (int)(6378.1*rad+0.5) << endl;
  }
  return 0;
}

C++: #include <cstdio>
#include <cmath>
#define r 6378.1

int main(void){
	double a,b,c,d;
	while(1){
		scanf("%lf %lf %lf %lf",&a,&b,&c,&d);
		if(b+d<-1.0)break;
		a=90.0-a;
		c=90.0-c;
		double x1=cos(b*M_PI/180.0)*sin(a*M_PI/180.0);
		double x2=cos(d*M_PI/180.0)*sin(c*M_PI/180.0);
		double y1=sin(b*M_PI/180.0)*sin(a*M_PI/180.0);
		double y2=sin(d*M_PI/180.0)*sin(c*M_PI/180.0);
		double z1=cos(a*M_PI/180.0);
		double z2=cos(c*M_PI/180.0);
		double e=sqrt(pow(x1-x2,2.0)+pow(y1-y2,2.0)+pow(z1-z2,2.0));
		double ang=acos(1.0-e*e/2.0);
		printf("%.0f\n",r*ang+0.5);
	}
	return 0;
}

C++: #include <iostream>
#include <cmath>

using namespace std;

#define R 6378.1
#define PI 3.1415926535

int main() {
	double a, b, c, d;
	while (cin >> a >> b >> c >> d) {
		if (a == -1 && b == -1 && c == -1 && d == -1) {
			break;
		}
		double x1, y1, z1, x2, y2, z2;
		x1 = R * cos((a/360.0)*2.0*PI) * cos((b/360.0)*2.0*PI);
		y1 = R * cos((a/360.0)*2.0*PI) * sin((b/360.0)*2.0*PI);
		z1 = R * sin((a/360.0)*2.0*PI);
		x2 = R * cos((c/360.0)*2.0*PI) * cos((d/360.0)*2.0*PI);
		y2 = R * cos((c/360.0)*2.0*PI) * sin((d/360.0)*2.0*PI);
		z2 = R * sin((c/360.0)*2.0*PI);
		double r = sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1)+(z2-z1)*(z2-z1));
		double t = R * acos((R*R+R*R-r*r)/(2.0*R*R));
		double tt = t - (int)t;
		if (tt >= 0.5) {
			t += 1.0;
		}
		cout << (int)t << endl;
	}
	return 0;
}

C++: #include <bits/stdc++.h>
using namespace std ;

#define pb(n) push_back(n)
#define fi first
#define se second
#define all(r) (r).begin(),(r).end()
#define gsort(st,en) sort((st),(en),greater<int>())
#define vmax(ary) *max_element(all(ary))
#define vmin(ary) *min_element(all(ary))
#define debug(x) cout<<#x<<": "<<x<<endl
#define fcout(n) cout<<fixed<<setprecision((n))
#define scout(n) cout<<setw(n)
#define vary(type,name,size,init) vector< type> name(size,init)

#define rep(i,n) for(int i = 0; i < (int)(n);++i)
#define REP(i,a,b) for(int i = (a);i < (int)(b);++i)
#define repi(it,array) for(auto it = array.begin(),end = array.end(); it != end;++it)
#define repa(n,array) for(auto &n :(array))

using ll = long long;
using vi = vector<int>;
using vl = vector<ll>;
using dict = map<string,int>;
using pii = pair<int,int> ;

const int mod = 1000000007;
constexpr int imax = ((1<<30)-1)*2+1 ;
constexpr int inf = 100000000;
constexpr double PI = acos(-1.0) ;
double eps = 1e-10 ;
const int dy[] = {-1,0,1,0};
const int dx[] = {0,-1,0,1};

inline bool value(int x,int y,int w,int h){
  return (x >= 0 && x < w && y >= 0 && y < h);
}

template<typename T>
void Unique(vector<T> &v){
  sort(all(v));
  v.erase(unique(all(v)),v.end());
}

template<typename T>
T ston(string& str, T n){
  istringstream sin(str) ;
  T num ;
  sin >> num ;
  return num ;
}

void Ans(bool f){
  if(f) cout << "YES"<<endl;
  else cout << "NO"<<endl;
}

int main(){
  cin.tie(0);
  ios::sync_with_stdio(false);
  const int r = 6378.1;
  double x1,y1,x2,y2;
  while(cin >> y1 >> x1 >> y2 >> x2){
    if(y1 == -1 && x1 == -1 && y2 == -1 && x2 == -1) break;
    if(x1 > 180){
      x1 = x1 - 360;
    }
    if(x2 > 180){
      x2 = x2 - 360;
    }
    y1 = y1 * M_PI / 180;
    y2 = y2 * M_PI / 180;
    double d = r * acos(sin(y1) * sin(y2) + cos(y1) * cos(y2) * cos((x2 - x1) * M_PI/180));
    cout << round(d) << endl;
  }
  return 0;
}

C++: #include<iostream>
#include<cmath>
using namespace std;
#define REP(i,b,n) for(int i=b;i<n;i++)
#define rep(i,n)   REP(i,0,n)
#define PI acos(-1)
const double R=6378.01;

double dist(double x,double y,double z){
  return sqrt(x*x+y*y+z*z+1e-10);
}

void calc(double  &x,double &y,double &z,
	  double a,double b){
  double phi  = PI/2*a/90.0;
  double theta= PI*b/180.0;
  z=R*sin(phi);
  x=R*cos(phi)*cos(theta);
  y=R*cos(phi)*sin(theta);
}

double solve(double *a){
  double x[2],y[2],z[2];
  rep(i,2){
    calc(x[i],y[i],z[i],
	 a[2*i],a[2*i+1]);
  }
  
  double C=dist(x[0]-x[1],y[0]-y[1],z[0]-z[1]);
  

  double theta = acos((2*R*R-C*C)/(2*R*R));
  if ( !(-2*R*R<=2*R*R-C*C))theta=acos(-1);
  else if ( !(2*R*R-C*C<=2*R*R))theta=acos(1);
  return R*theta;
}

main(){
  double in[4];
  while(cin>>in[0]>>in[1]>>in[2]>>in[3]){
    if ( in[0] == -1 && in[1] == -1 && in[2] == -1 && in[3]==-1)break;
    printf("%.0lf\n",solve(in));
  }
}

C++: #include <iostream>
#include <math.h>

using namespace std;

int main() {
	int i;
	double px[2],py[2],pz[2],x[2],y[2],len,s,pi=3.141592;
	while(cin >> y[0]) {
		if (y[0]==-1) break;
		cin >> x[0] >> y[1] >> x[1];
		for (i=0;i<2;i++) {
			px[i]=cos(pi*y[i]/180); pz[i]=sin(pi*y[i]/180);
			py[i]=sin(pi*x[i]/180)*px[i]; px[i]=cos(pi*x[i]/180)*px[i];
		}
		len=sqrt((px[0]-px[1])*(px[0]-px[1])+(py[0]-py[1])*(py[0]-py[1])+(pz[0]-pz[1])*(pz[0]-pz[1]));
		s=acos(1-len*len/2);
		cout << (int)(6378.1*s+0.5) << endl;
	} 
	return 0;
}

C++: #include <iostream>
#include <vector>
#include <algorithm>
#include <cmath>

#define REP(i,k,n) for(int i=k;i<n;i++)
#define rep(i,n) for(int i=0;i<n;i++)
#define EPS 1e-10
#define PI acos(-1)
#define R 6378.1

using namespace std;

struct P
{
	double x,y,z;

	double distance(const P& o)
	{
		return sqrt( (x-o.x)*(x-o.x) + (y-o.y)*(y-o.y) + (z-o.z)*(z-o.z) );
	}
};

double toRad(double a)
{
	return (a/180)*PI;
}

int main()
{
	double a,b,c,d;
	while(cin >> a >> b >> c >> d)
	{
		if(a == -1 && b == -1 && c == -1 && d == -1) break;

		P p[2];

		p[0].x = R*cos(toRad(a))*cos(toRad(b));
		p[0].y = R*cos(toRad(a))*sin(toRad(b));
		p[0].z = R*sin(toRad(a));

		p[1].x = R*cos(toRad(c))*cos(toRad(d));
		p[1].y = R*cos(toRad(c))*sin(toRad(d));;
		p[1].z = R*sin(toRad(c));

		double angle = acos((R*R + R*R - p[0].distance(p[1])*p[0].distance(p[1])) / (2*R*R));
		cout << (int)(PI*R*(angle/PI) + 0.5) << endl;

	}
	return 0;
}