• [报告]ACMICPC 2012 Regionals Asia Changchun I Polaris of Pandora / ZJU 3663


    Abstract

    ACM-ICPC 2012 Regionals Asia - Changchun I Polaris of Pandora / ZJU 3663

    三维计算几何

    Body

    Source

    http://acm.zju.edu.cn/onlinejudge/showProblem.do?problemCode=3663

    Description

    一个半径为R的球,给定起点终点的经纬度和飞行高度H。平行光从正北方向往正南方向照射,问飞行路程中被光找到的比例。

    Solution

    本场悲剧的我们大概也就这水平了……

    本质上就是求飞行路径所在的大圆和半径为R的圆柱的交点。建立直角坐标系硬算就好了……

    本题精度坑爹(h和R可能相差一万倍),不可以把高度和半径归一化(注意我中间注释的那两行),但是不归一化我算出来的点甚至和起点终点不共面……当然也可能是我姿势不优美……总之给1y的hl大神跪了。

    Code

    #include <iomanip>
    #include <cassert>
    #include <cstdio>
    #include <cstring>
    #include <cmath>
    #include <iostream>
    #include <algorithm>
    using namespace std;
    
    const double EPS = 1e-8;
    const double PI = acos(-1.0);
    
    inline int sgn(double d) {
        if (fabs(d)<EPS) return 0;
        return d>0?1:-1;
    }
    
    struct sp {
        double x, y, z;
        sp() {}
        sp(double a, double b, double c): x(a), y(b), z(c) {}
        void write() {
            printf("%f %f %f\n", x, y, z);
        }
    };
    
    sp operator-(const sp &u, const sp &v) {
        return sp(u.x-v.x, u.y-v.y, u.z-v.z);
    }
    
    double dot(const sp &u, const sp &v) {
        return u.x*v.x+u.y*v.y+u.z*v.z;
    }
    
    sp det(const sp &u, const sp &v) {
        return sp(u.y*v.z-v.y*u.z,
                u.z*v.x-u.x*v.z,
                u.x*v.y-u.y*v.x);
    }
    
    double nrm(const sp &u) {
        return sqrt(u.x*u.x+u.y*u.y+u.z*u.z);
    }
    
    double dis(const sp &u, const sp &v) {
        double dx = u.x-v.x;
        double dy = u.y-v.y;
        double dz = u.z-v.z;
        return sqrt(dx*dx+dy*dy+dz*dz);
    }
    
    double ang(const sp &u, const sp &v) {
        double d = dot(u, v);
        return acos(d/(nrm(u)*nrm(v)));
    }
    
    double r, h, as, ns, at, nt, R;
    sp s, t, hmr, mdk;
    double a, b, c, tmp;
    
    bool quad(double a, double b, double c, double &x1, double &x2) {
        double delta = b*b-4*a*c;
        if (sgn(delta)<=0) return 0;
        x1 = (-b+sqrt(delta))/(2*a);
        x2 = (-b-sqrt(delta))/(2*a);
        return 1;
    }
    
    bool coplane(const sp &a, const sp &b, const sp &c, const sp &d) {
        return sgn(dot(det(c-a, b-a), d-a))==0;
    }
    
    bool inside(const sp &s, const sp &p, const sp &t) {
        return sgn(dot(det(s,p),det(s,t)))>=0 && sgn(dot(det(p,t),det(s,t)))>=0;
    }
    
    int main() {
        while (cin>>r>>h>>as>>ns>>at>>nt) {
            if (sgn(as)==0 && sgn(at)==0) {
                puts("100.000");
                continue;
            }
            as += PI/2; at += PI/2;
            ns += PI; nt += PI;
            //h /= r; r = 1;
            R = r+h;
            s = sp(R*sin(as)*cos(ns), R*sin(as)*sin(ns), R*cos(as));
            t = sp(R*sin(at)*cos(nt), R*sin(at)*sin(nt), R*cos(at));
            a = s.y*t.z-s.z*t.y; b = s.z*t.x-s.x*t.z; c = s.x*t.y-s.y*t.x;
            tmp = -sqrt(R*R-r*r);
            hmr.z = mdk.z = tmp;
            assert(sgn(a)!=0 || sgn(b)!=0);
            if (sgn(b)) {
                if (quad(a*a+b*b, 2*a*c*tmp, c*c*(R*R-r*r)-b*b*r*r, hmr.x, mdk.x)) {
                    hmr.y = (-c*tmp-a*hmr.x)/b;
                    mdk.y = (-c*tmp-a*mdk.x)/b;
                }
                else {
                    puts("100.000");
                    continue;
                }
            }
            else {
                if (quad(a*a+b*b, 2*b*c*tmp, c*c*(R*R-r*r)-a*a*r*r, hmr.y, mdk.y)) {
                    hmr.x = (-c*tmp-b*hmr.y)/a;
                    mdk.x = (-c*tmp-b*mdk.y)/a;
                }
                else {
                    puts("100.000");
                    continue;
                }
            }
            //assert(coplane(s, t, hmr, mdk));
            int cnt = 0;
            if (inside(s, hmr, t)) ++cnt;
            if (inside(s, mdk, t)) ++cnt;
            double tot=ang(s, t), dark;
            if (cnt==0) {
                if (inside(hmr, s, mdk) && inside(hmr, t, mdk)) dark = tot;
                else dark = 0;
            }
            else if (cnt==1) {
                if (inside(s, mdk, t)) swap(hmr, mdk);
                assert(inside(s, hmr, t)&&!inside(s, mdk, t));
                if (inside(hmr, s, mdk)) dark = ang(hmr, s);
                else dark = ang(hmr, t);
            }
            else dark = ang(hmr, mdk);
            double ans = (1-dark/tot)*100;
            printf("%.3f\n", ans);
            //cout<<fixed<<setprecision(3)<<ans<<endl;
        }
        return 0;
    }
  • 相关阅读:
    生活网站
    input 输入值的监听 禁止输入特殊字符
    jq 插件分享
    css3 特效分享
    sharepoint---RBS回收站清空设置
    &&和||
    DataTable得到某行某列的值
    后台刷新当前页面和弹出对话框跳转页面
    c#.netGr idView1在div不局中
    DIV UL LI
  • 原文地址:https://www.cnblogs.com/jffifa/p/2725035.html
Copyright © 2020-2023  润新知