UVALive 4589 Asteroids (三维凸包&&重心)

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <map>
#include <vector>
#include <queue>
#include <set>
#define FRER() freopen("/Users/seven7777777/Desktop/Code/TESTXCODE/TESTXCODE/in.txt","r",stdin)
#define FREW() freopen("/Users/seven7777777/Desktop/Code/TESTXCODE/TESTXCODE/out.txt","w",stdout)
#define go int T;cin>>T;for(int kase=0;kase<T;kase++)
#define debug cout<<"****************"<<endl
#define lowbit(x) x&(-x)
#define eps 1e-6
#define mod 1000000007
using namespace std;
typedef long long ll;
typedef long double ld;
typedef unsigned long long ull;
typedef pair<int,int> pii;
const int maxn = 100 + 7;
struct Point3
{
    double x,y,z;
    Point3(double x = 0 , double y = 0,double z = 0):x(x),y(y),z(z){}
    void print(){
        printf("x:%lf y:%lf z:%lf",x,y,z);
    }
};
typedef Point3 Vector3;
Vector3 operator + (Vector3 A,Vector3 B){
    return Vector3(A.x+B.x,A.y+B.y,A.z+B.z);
}
Vector3 operator - (Vector3 A,Vector3 B){
    return Vector3(A.x-B.x,A.y-B.y,A.z-B.z);
}
Vector3 operator * (Vector3 A,double p){
    return Vector3(A.x*p,A.y*p,A.z*p);
}
Vector3 operator / (Vector3 A,double p){
    return Vector3(A.x/p,A.y/p,A.z/p);
}
double Dot(Vector3 A,Vector3 B){
    return A.x*B.x + A.y*B.y + A.z*B.z;
}
Vector3 Cross(Vector3 A,Vector3 B){
    return Vector3(A.y*B.z - A.z*B.y , A.z*B.x - A.x*B.z ,A.x*B.y - A.y*B.x);
}
double Volume6(Point3 A,Point3 B,Point3 C,Point3 D){
    return Dot(D-A,Cross(B-A, C-A));
}
struct Face{
    int v[3];
    Face(int a=0,int b=0,int c=0){
        v[0] = a;
        v[1] = b;
        v[2] = c;
    }
    Vector3 Normal(Point3* p)const{
        return Cross(p[v[1]]-p[v[0]],p[v[2]]-p[v[0]]);
    }
    int cansee(Point3 *p,int i)const{
        return Dot(p[i]-p[v[0]],Normal(p)) > 0 ? 1 : 0;
    }
    void print(Point3 *p){
        cout<<"a ";p[v[0]].print();cout<<";";
        cout<<"b ";p[v[1]].print();cout<<";";
        cout<<"c ";p[v[2]].print();cout<<";";
        cout<<endl;
    }
};
double rand01(){return rand()/(double)RAND_MAX;}
double randeps(){return (rand01()-0.5)*eps;}
Point3 add_noise(Point3 p){
    return Point3(p.x+randeps(),p.y+randeps(),p.z+randeps());
}
vector<Face> CH3D(Point3* p,int n){
    int vis[maxn][maxn];
    vector<Face>cur;
    cur.push_back(Face(0,1,2));
    cur.push_back(Face(2,1,0));
    for(int i=3;i<n;i++){
        vector<Face>nxt;
        for(int j=0;j<cur.size();j++){
            Face& f = cur[j];
            int res = f.cansee(p,i);
            if(!res) nxt.push_back(f);
            for(int k=0;k<3;k++) vis[f.v[k]][f.v[(k+1)%3]] = res;
        }
        for(int j=0;j<cur.size();j++){
            for(int k=0;k<3;k++){
                int a = cur[j].v[k] , b = cur[j].v[(k+1)%3];
                if(vis[a][b]!=vis[b][a]&&vis[a][b]) nxt.push_back(Face(a,b,i));
            }
        }
        cur = nxt;
    }
    return cur;
}
double length(Vector3 A){
    return sqrt(Dot(A, A));
}
double DistanceToPlane(Point3 p,Point3 p0,Point3 n){
    return fabs(Dot(p-p0,n)/length(n));
}
double solve(Point3 *p,Point3* _p,int n,Point3 avg){
    //avg.print();
    vector<Face>f = CH3D(_p, n);
    Point3 base;
    int len = (int)f.size();
    double totv = 0;
    for(int i=0;i<len;i++){
        //f[i].print(p);
        double v = Volume6(p[f[i].v[0]], p[f[i].v[1]], p[f[i].v[2]], avg);
        base = base + (p[f[i].v[0]]+p[f[i].v[1]]+p[f[i].v[2]]+avg)/4*v;
        totv += v;
    }
    base = base / totv;
    //base.print();
    //debug;
    double ans = 9999999;
    for(int i=0;i<len;i++){
        Vector3 n = f[i].Normal(p);
        ans = min(ans,DistanceToPlane(base, p[f[i].v[0]], n));
        //cout<<ans<<endl;
    }
    //debug;
    return ans;
}
void read(Point3 *p,Point3 *_p,int n,Point3& avg){
    for(int i=0;i<n;i++){
        scanf("%lf%lf%lf",&p[i].x,&p[i].y,&p[i].z);
        _p[i] = add_noise(p[i]);
        avg.x += p[i].x;avg.y += p[i].y;avg.z += p[i].z;
    }
    avg = avg/n;
}
Point3 plant1[maxn],plant2[maxn],_plant1[maxn],_plant2[maxn];
int n,m;
int main(){
    //FRER();
    //FREW();
     while(~scanf("%d",&n)){
         Point3 avg1,avg2;
         read(plant1,_plant1, n, avg1);
         scanf("%d",&m);
         read(plant2,_plant2, m, avg2);
         double ans1 = solve(plant1,_plant1,n,avg1);
         double ans2 = solve(plant2,_plant2,m,avg2);
         printf("%.5f\n",ans1+ans2);
     }
}