p2107.cpp

大概POJ上50道比较难的题的代码

// A B C    D E F    G H I
#include <iostream>
#include <algorithm>
#include <cmath>
using namespace std;
const int MAXN = 50;
struct point{ double x,y,z; };
point a[MAXN*MAXN];
inline double mul(double x1,double y1,double x2,double y2){
    return x1*y2 - x2*y1;
}
int main(){
    freopen("t.txt","r",stdin);
    freopen("r.txt","w",stdout);
    point b,d,e,f,h,ans;
    point vf,hf,vt,ht,N;
    double hu,vu,fac;
    int n,m;
  while(cin >> n >> m){
    cout << n*m << endl;
    for(int i = 0;i < m*n;++i) scanf("%lf%lf%lf",&a[i].x, &a[i].y, &a[i].z);
    for(int i = 0;i < m-2;++i)
        for(int j = 0;j < n-2;++j){
            d = a[i*n+j+1];
            b = a[(i+1)*n+j];
            e = a[(i+1)*n+j+1];
            h = a[(i+1)*n+j+2];
            f = a[(i+2)*n+j+1];
/*            b = a[i*n+j+1];
            d = a[(i+1)*n+j];
            e = a[(i+1)*n+j+1];
            f = a[(i+1)*n+j+2];
            h = a[(i+2)*n+j+1];*/
            vf.x = mul(b.y-h.y, b.z-h.z, e.y-b.y, e.z-b.z);
            vf.y = mul(b.z-h.z, b.x-h.x, e.z-b.z, e.x-b.x);
            vf.z = mul(b.x-h.x, b.y-h.y, e.x-b.x, e.y-b.y);
            
            hf.x = mul(d.y-f.y, d.z-f.z, e.y-d.y, e.z-d.z);
            hf.y = mul(d.z-f.z, d.x-f.x, e.z-d.z, e.x-d.x);
            hf.z = mul(d.x-f.x, d.y-f.y, e.x-d.x, e.y-d.y);
            hu = (e.x-d.x)*(d.x-f.x) + (e.y-d.y)*(d.y-f.y) + (e.z-d.z)*(d.z-f.z);
            vu = (e.x-b.x)*(b.x-h.x) + (e.y-b.y)*(b.y-h.y) + (e.z-b.z)*(b.z-h.z);
            vt.x = ((h.y-e.y)*vf.z - (h.z-e.z)*vf.y) - (h.x-e.x)*vu;
            vt.y = ((h.z-e.z)*vf.x - (h.x-e.x)*vf.z) - (h.y-e.y)*vu;
            vt.z = ((h.x-e.x)*vf.y - (h.y-e.y)*vf.x) - (h.z-e.z)*vu;
            
            ht.x = ((f.y-e.y)*hf.z - (f.z-e.z)*hf.y) - (f.x-e.x)*hu;
            ht.y = ((f.z-e.z)*hf.x - (f.x-e.x)*hf.z) - (f.y-e.y)*hu;
            ht.z = ((f.x-e.x)*hf.y - (f.y-e.y)*hf.x) - (f.z-e.z)*hu;
            N.x = mul(ht.y, ht.z, vt.y, vt.z);
            N.y = mul(ht.z, ht.x, vt.z, vt.x);
            N.z = mul(ht.x, ht.y, vt.x, vt.y);
            fac = sqrt(N.x*N.x + N.y*N.y + N.z*N.z);
            ans.x = N.x / fac;
            ans.y = N.y / fac;
            ans.z = N.z / fac;
            printf("%.4lf %.4lf %.4lf\n",ans.x, ans.y, ans.z);
        }
  }
}