POJ 2412 /// 空間幾何 經緯度轉三維座標 角度轉弧度 法向量

題目大意：

給定半徑6378km的球上的 多個地點 及其  經緯度

多個詢問 給定三個地點 A B C

A與B的等距點在球上造成一個大圓 

即球面上有一個到兩點距離相等的大圓

且大圓所在平面垂直兩點連線 

求C點到該大圓的最近球面距離 

 

1.特殊狀況即A B爲同一點 此時整個球面上的點都是等距點

則C到等距點的最近球面距離爲0

2.設兩點與球心組成一個平面S1 大圓與球面組成的平面S2

則S1與S2垂直 即兩個平面相差90度

求出球心到C的向量cd與S1的夾角爲ang度

則cd與S2的夾角則爲|ang-90|度

#include <bits/stdc++.h>
using namespace std;
#define INF 0x3f3f3f3f
#define LL long long
#define mem(i,j) memset(i,j,sizeof(i))
#define inc(i,j,k) for(int i=j;i<=k;i++)
#define dec(i,j,k) for(int i=j;i>=k;i--)
const int N=1e5+5;
const double eps=1e-8;
const double PI=acos(-1.0);

int dcmp(double x) {
    if(abs(x)<eps) return 0;
    else return x<0 ? -1:1;
}
struct P {
    double x,y,z;
    P(){}
    P(double x,double y,double z):x(x),y(y),z(z){}
    P operator -(const P& p)const { return P(x-p.x,y-p.y,z-p.z); }
    double dot(const P& p) const { return x*p.x+y*p.y+z*p.z; }
    bool operator ==(const P& p)const {
        return dcmp(x-p.x)==0 && dcmp(y-p.y)==0 && dcmp(z-p.z)==0;
    }
}p[1005];
double Radian(double t) {
    return t*PI/180.0;
}
double lenP(P p) {
    return sqrt(p.dot(p));
}
double Angle(P a,P b) {
    return acos(a.dot(b)/lenP(a)/lenP(b));
}
int tot;
map<string,int>id;
double R=6378.0;

void ptf(string c,int res,string a,string b) {
    cout<<c<<" is ";
    if(res==-1) cout<<"?";
    else cout<<res;
    cout<<" km off "<<a<<"/"<<b<<" equidistance.\n";
}

int main()
{
    ios::sync_with_stdio(false);
    id.clear(); tot=0;

    string s;
    double la,lo;
    while(cin>>s) {
        if(s=="#") break;
        id[s]=++tot;
        cin>>la>>lo;
        p[tot].x=R*cos(Radian(la))*sin(Radian(lo));
        p[tot].y=R*cos(Radian(la))*cos(Radian(lo));
        p[tot].z=R*sin(Radian(la));
    }

    string A,B,C;
    while(cin>>A) {
        if(A=="#") break;
        cin>>B>>C;

        int aid,bid,cid;
        bool flag=0;

        if(!id.count(A)) flag=1; else aid=id[A];
        if(!id.count(B)) flag=1; else bid=id[B];
        if(!id.count(C)) flag=1; else cid=id[C];

        double ans, ang;
        if(flag) ans=-1.0;
        else {
            P a=p[aid], b=p[bid], c=p[cid];
            if(a==b) ans=0.0;
            else {
                ang=Angle(a-b,c); 
                ang=abs(ang-PI/2.0); 
                ans=ang*R+0.5; // 弧長公式
            }
        } 
        ptf(C,(int)ans,A,B);
    }

    return 0;
}

View Code