在三维坐标中，给定一个点光源，一个凸多面体，以及一个平面作为地面。
求该凸多面体在地面上阴影的面积。
这三个点共同确定了一个平面，这个平面就是地面。保证这三个点坐标互异且不共线。前三行每行三个实数，每行表示一个点。这三个点共同确定了一个平面，这个平面就是地面。保证这三个点坐标互异且不共线。
接下来一行三个实数，表示一个点。这个点就是点光源。
之后一个整数n，表示凸多面体顶点的数量。
之后n行，每行三个实数，表示凸多面体的一个顶点。
1 #include <iostream> 2 #include <cstdio> 3 #include <cmath> 4 #include <cstring> 5 #include <vector> 6 #include <algorithm> 7 #include <queue> 8 #include <stack> 9 #include <map> 10 #include <set> 11 #include <utility> 12 13 #ifdef DEBUG 14 const int MAXN = 110; 15 #else 16 const int MAXN = 110; 17 #endif 18 19 const double eps = 1e-6; 20 21 using namespace std; 22 23 struct Vector { 24 double x, y, z; 25 Vector(double _x = 0, double _y = 0, double _z = 0): x(_x), y(_y), z(_z) {} 26 Vector operator+(const Vector &rhs) const { 27 return Vector(x + rhs.x, y + rhs.y, z + rhs.z); 28 } 29 Vector operator-(const Vector &rhs) const { 30 return Vector(x - rhs.x, y - rhs.y, z - rhs.z); 31 } 32 Vector operator*(double k) { 33 return Vector(x * k, y * k, z * k); 34 } 35 double len() { 36 return sqrt(x * x + y * y + z * z); 37 } 38 }; 39 typedef Vector Point; 40 41 inline double dot(const Vector &a, const Vector &b) { 42 return a.x * b.x + a.y * b.y + a.z * b.z; 43 } 44 45 inline Vector cross(const Vector &a, const Vector &b) { 46 return Vector(a.y * b.z - b.y * a.z, a.z * b.x - b.z * a.x, a.x * b.y - a.y * b.x); 47 } 48 49 bool flag; 50 inline Vector readv() { 51 double x, y, z; 52 scanf("%lf%lf%lf", &x, &y, &z); 53 if (flag) swap(y, z); 54 return Vector(x, y, z); 55 } 56 57 inline void printv(const Vector &a) { 58 printf("(%lf, %lf, %lf) ", a.x, a.y, a.z); 59 } 60 61 struct Line { 62 Point s; 63 Vector d; 64 Line() {} 65 Line(Point _s, Point _t): s(_s), d(_t - _s) {} 66 }; 67 68 struct Plain { 69 Point s; 70 Vector d; 71 Plain() {} 72 Plain(Point a, Point b, Point c): s(a), d(cross(b - a, c - a)) {} 73 } p; 74 Point s; 75 int n; 76 Point pol[MAXN]; 77 Point sd[MAXN]; 78 79 inline int dcmp(double x) { 80 return x < -eps ? -1 : x > eps ? 1 : 0; 81 } 82 83 inline Point itsct(Line l, Plain p) { 84 double s1 = dot(p.d, l.s - p.s); 85 double s2 = dot(p.d, l.s + l.d - p.s); 86 if (!dcmp(s1 - s2)) { 87 return itsct(Line(l.s, l.s + p.d), p); 88 } 89 double k = s1 / (s1 - s2); 90 return l.s + l.d * k; 91 } 92 93 bool cmp(const Point &a, const Point &b) { 94 return a.x < b.x; 95 } 96 97 Point conv[MAXN]; 98 inline int make_convex() { 99 int cnt = 0; 100 sort(sd, sd + n, cmp); 101 conv[cnt++] = sd[0]; 102 for (int i = 1; i < n; i++) { 103 while (cnt > 1 && dcmp(cross(sd[i] - conv[cnt - 1], conv[cnt - 1] - conv[cnt - 2]).z) < 0) cnt--; 104 if (dcmp((conv[cnt - 1] - sd[i]).len())) conv[cnt++] = sd[i]; 105 } 106 for (int i = n - 2; i >= 0; i--) { 107 while (cnt > 1 && dcmp(cross(sd[i] - conv[cnt - 1], conv[cnt - 1] - conv[cnt - 2]).z) < 0) cnt--; 108 if (dcmp((conv[cnt - 1] - sd[i]).len())) conv[cnt++] = sd[i]; 109 } 110 return cnt; 111 } 112 113 inline void open_file() { 114 freopen("shadow.in", "r", stdin); 115 freopen("shadow.out", "w", stdout); 116 } 117 118 int main() { 119 //open_file(); 120 Point a = readv(), b = readv(), c = readv(); 121 if (!dcmp(cross(b - a, c - a).z)) { 122 swap(a.y, a.z), swap(b.y, b.z), swap(c.y, c.z); 123 flag = true; 124 } 125 p = Plain(a, b, c); 126 s = readv(); 127 scanf("%d", &n); 128 for (int i = 0; i < n; i++) pol[i] = readv(); 129 for (int i = 0; i < n; i++) sd[i] = itsct(Line(s, pol[i]), p); 130 int cnt = make_convex(); 131 double ans = 0; 132 for (int i = 1; i < cnt - 2; i++) { 133 ans += cross(conv[i] - conv[0], conv[i + 1] - conv[0]).len(); 134 } 135 printf("%.2lf ", ans / 2); 136 137 return 0; 138 }
直接给你平面方程
1 #include <algorithm> 2 #include <iostream> 3 #include <cstdlib> 4 #include <cstdio> 5 #include <cmath> 6 7 using namespace std; 8 9 typedef struct p1node 10 { 11 double a,b,c,d; 12 }plane; 13 plane Pl; 14 15 typedef struct p2node 16 { 17 double x,y,z; 18 }point; 19 point temp; 20 point P[ 105 ]; 21 point S[ 105 ]; 22 23 //计算点在平面上的投影 24 int shadow( plane p, point s, int n ) 25 { 26 //求出过s平行于plane的平面 ax+by+c=D 27 double D = p.a*s.x+p.b*s.y+p.c*s.z; 28 if ( D-p.d < 0 ) {//调整方向 29 p.a *= -1;p.b *= -1;p.c *= -1; 30 p.d *= -1; 31 D *= -1; 32 } 33 34 //判断点与面的关系 35 int count = 0; 36 for ( int i = 0 ; i < n ; ++ i ) { 37 double det = p.a*P[i].x+p.b*P[i].y+p.c*P[i].z-D; 38 if ( det < 0 ) count ++; 39 } 40 if ( count == 0 ) return 0; 41 if ( count != n ) return n+1; 42 43 for ( int i = 0 ; i < n ; ++ i ) { 44 //直线方程: (Sx,Sy,Sz) + t(dx,dy,dz) 45 double dx = P[i].x - s.x; 46 double dy = P[i].y - s.y; 47 double dz = P[i].z - s.z; 48 double t = (p.d-p.a*s.x-p.b*s.y-p.c*s.z)/(p.a*dx+p.b*dy+p.c*dz); 49 50 P[i].x = s.x + t*dx; 51 P[i].y = s.y + t*dy; 52 P[i].z = s.z + t*dz; 53 } 54 return n; 55 } 56 57 //坐标系旋转，到x-y平面 58 void change( plane p, int n ) 59 { 60 //平行于x-y平面的不用计算，也不能计算（分母为0） 61 if ( p.a*p.a + p.b*p.b == 0 ) return; 62 for ( int i = 0 ; i < n ; ++ i ) { 63 //绕z轴旋转 64 double cosC = p.b/sqrt(p.a*p.a+p.b*p.b); 65 double sinC = p.a/sqrt(p.a*p.a+p.b*p.b); 66 temp.x = P[i].x*cosC-P[i].y*sinC; 67 temp.y = P[i].x*sinC+P[i].y*cosC; 68 temp.z = P[i].z; 69 P[i] = temp; 70 //绕x轴旋转 71 double cosA = p.c/sqrt(p.a*p.a+p.b*p.b+p.c*p.c); 72 double sinA = sqrt(p.a*p.a+p.b*p.b)/sqrt(p.a*p.a+p.b*p.b+p.c*p.c); 73 temp.x = P[i].x; 74 temp.y = P[i].y*cosA-P[i].z*sinA; 75 temp.z = P[i].y*sinA+P[i].z*cosA; 76 P[i] = temp; 77 } 78 } 79 80 //计算二维凸包面积 81 double crossproduct( point a, point b, point c ) 82 { 83 return (b.x-a.x)*(c.y-a.y)-(c.x-a.x)*(b.y-a.y); 84 } 85 86 bool cmp1( point a, point b ) 87 { 88 if ( a.x == b.x ) 89 return a.y < b.y; 90 return a.x < b.x; 91 } 92 93 bool cmp2( point a, point b ) 94 { 95 return crossproduct( P[0], a, b )>0; 96 } 97 98 void Graham( int n ) 99 { 100 sort( P+0, P+n, cmp1 ); 101 sort( P+1, P+n, cmp2 ); 102 103 int top = -1; 104 if ( n > 0 ) S[++ top] = P[0]; 105 if ( n > 1 ) S[++ top] = P[1]; 106 if ( n > 2 ) { 107 for ( int i = 2 ; i < n ; ++ i ) { 108 while ( crossproduct( S[top-1], S[top], P[i] ) < 0 ) -- top; 109 S[++ top] = P[i]; 110 } 111 } 112 113 double area = 0.0; 114 for ( int i = 2 ; i <= top ; ++ i ) 115 area += crossproduct( S[0], S[i-1], S[i] ); 116 117 printf("%.2lf ",area*0.5); 118 } 119 120 int main() 121 { 122 int n; 123 while ( cin >> Pl.a >> Pl.b >> Pl.c >> Pl.d ) { 124 if ( Pl.a == 0 && Pl.b == 0 && Pl.c == 0 ) break; 125 cin >> n; 126 for ( int i = 0 ; i <= n ; ++ i ) 127 cin >> P[i].x >> P[i].y >> P[i].z; 128 129 int s_count = shadow( Pl, P[n], n ); 130 if ( !s_count ) 131 cout << "0.00" << endl; 132 else if ( s_count > n ) 133 cout << "Infi" << endl; 134 else { 135 change( Pl, s_count ); 136 Graham( s_count ); 137 } 138 } 139 }