#include <stdio.h>
#include <math.h>

const double earth_R = 6378.137;    // 地球の平均半径
const double PI = 3.141592653589793238462643383279; // 円周率

typedef struct {                // 地点データ
    double ido;
    double keido;
    double v[3];                // 単位ベクトル x,y,z 成分
} POS;

int main()
{
    int i;
    POS pos[2];                 // 地点の情報
    double angle;               // 単位ベクトルのなす角度
    double dv[3];               // 2点を結ぶベクトル (x,y,z)

    // 入力
    printf("地球上の2点間の距離を求める(地球の平均半径=%fkmと仮定)\n\n", earth_R);
    printf("  2地点の緯度経度をDegree(-180.0～180.0)で入力\n");
    printf("    (※ 北緯 35度30分・西経 100度の場合、'35.5 -100' と入力\n");
    printf("        南緯 35度30分・東経 100度の場合、'-35.5 100' と入力 )\n");
    for (i = 0; i < 2; i++) {
        printf("・地点 %d の緯度経度：", i + 1);
        scanf("%lf %lf", &pos[i].ido, &pos[i].keido);
    }
    // それぞれの単位ベクトルをを求める
    for (i = 0; i < 2; i++) {
        pos[i].v[0] = cos(pos[i].keido * PI / 180.0) * cos(pos[i].ido * PI / 180.0);
        pos[i].v[1] = sin(pos[i].keido * PI / 180.0) * cos(pos[i].ido * PI / 180.0);
        pos[i].v[2] = sin(pos[i].ido * PI / 180.0);
    }
    // 2地点を結ぶベクトル dv[] (x,y,z) を求める
    for (i = 0; i < 3; i++) {
        dv[i] = pos[1].v[i] - pos[0].v[i];
    }
    // 単位ベクトルのなす角度 angle[rad] を求める
    angle = 0.0;                // 内積計算temp
    for (i = 0; i < 3; i++) {
        angle += pos[0].v[i] * pos[1].v[i];
    }
    angle = acos(angle);
    // 表示
    printf("\n地点 1→2の\n");
    printf("・球面上の最短距離       = %fkm\n", earth_R * angle);
    printf("・最短距離(地下トンネル) = %fkm\n", earth_R * sqrt(dv[0] * dv[0] + dv[1] * dv[1] + dv[2] * dv[2]));

    // 終了
    return 0;
}