HDU-3692 三维几何，二维凸包

Problem Description
On the planet Pandora, huge rocks are floating in the sky. Despite the beautiful scenery, these rocks bring some critical problem: they block the sunshine, and shade on the ground. Na’vi (“people” living on the planet) cannot grow plants in the shade because of the lack of sunshine. In order to predicate the amount of food they can get, they have to calculate the area of the shade.

Let’s assume the “sun” as a point-source of light, and the ground as a infinite flat plane. Also, to simplify this problem, assume all rocks are convex polyhedrons.

Now here is a mathematical problem. Given the position of the point-source of light, the position of a convex polyhedron and the position of an infinite plane, you are required to calculate the area of the shade. Please note that light travels in a strictly straight line and the light source and the convex polyhedron are on the same side of the plane. The light source is not placed in the convex polyhedron, nor on the polyhedron or on the plane.

Input
The input contains no more than 100 cases.
Each case contains several lines which are formatted as follows.

a b c d
n
x1 y1 z1
x2 y2 z2
…
xn yn zn
x0 y0 z0

a, b, c and d means that the plane is ax+by+cz=d. n means the number of vertex of the convex polyhedron. It is guaranteed that n is no more than 100. Following it are n lines. Each line contains three float numbers, indicating the position of a point which is a vertex of the polyhedron. The last line also contains three float numbers, which means the position of the point-source of light.

The input is ended by a=0, b=0, c=0 and d=0.

Output
For each case, if there is no shade on the plane, print “0.00”.
If the area of the shade is infinite, print “Infi”.
Otherwise, print out the area of the shade. Please round the result to two digits after the decimal point.

Sample Input
0 0 1 0
8
1 1 1
1 -1 1
-1 1 1
-1 -1 1
1 1 0
1 -1 0
-1 1 0
-1 -1 0
2 2 1
1 0 0 -2
8
1 1 1
1 1 -1
1 -1 1
1 -1 -1
-1 1 1
-1 1 -1
-1 -1 1
-1 -1 -1
2 0 0
1 1 1 1
0
0 0 0
0 0 0 0

Sample Output
Infi
64.00
0.00

多面体每个点与点光源连一条直线，算出与平面的交点cross
若交点数为0： 0.00
若交点数<多面体点数: Infi
else:
把这些三维点的坐标转换为二维坐标
方法：
在平面内找到一对单位正交向量 i,j 和一个原点O
cross.x - > (cross.x-O.x) * i
cross.y - > (cross.y-O.y) * j
（* 为点乘）

#include<iostream>
#include<string.h>
#include<math.h>
#include<stdio.h>
#include<algorithm>
#include<vector>
using namespace std;
typedef long long ll;
const int N=5e4+7;
const int maxp=1e3+7;
const double inf=1e200;
const double eps=1e-8;
const double pi =acos(-1.0);
int sgn(double x){
	if(fabs(x)<eps)return 0;
	return x>0?1:-1;
}
struct Point3{//三维几何
	double x,y,z;
	Point3(double _x=0,double _y=0,double _z=0){x=_x;y=_y;z=_z;}
	void input(){scanf("%lf%lf%lf",&x,&y,&z);}
	void show(){printf("%.2lf %.2lf %.2lf\n",x,y,z);}
	bool operator ==(const Point3 &b){return !sgn(x-b.x)&&sgn(y-b.y)&&!sgn(z-b.z);}
	bool operator  <(const Point3 &b){return !sgn(x-b.x)?(!sgn(y-b.y)?sgn(z-b.z)<0:y<b.y):x<b.x;}
	double len(){return sqrt(x*x+y*y+z*z);}
	double len2(){return x*x+y*y+z*z;}
	double distance(const Point3 &b){return sqrt((x-b.x)*(x-b.x)+(y-b.y)*(y-b.y)+(z-b.z)*(z-b.z));}
	Point3 operator -(const Point3 &b){return Point3(x-b.x,y-b.y,z-b.z);}
	Point3 operator +(const Point3 &b){return Point3(x+b.x,y+b.y,z+b.z);}
	Point3 operator *(const double &k){return Point3(x*k,y*k,z*k);}
	Point3 operator /(const double &k){return Point3(x/k,y/k,z/k);}
	double operator *(const Point3 &b)const{return x*b.x+y*b.y+z*b.z;}//点乘
	Point3 operator ^(const Point3 &b)const{return Point3(y*b.z-z*b.y,z*b.x-x*b.z,x*b.y-y*b.x);}//叉乘
	double rad(Point3 a,Point3 b){
		Point3 p=(*this);
		return acos(((a-p)*(b-p))/(a.distance(p)*b.distance(p)));
	}
	Point3 trunc(double r){//变换长度
		double l=len();
		if(!sgn(l))return *this;
		r/=l;
		return Point3(x*r,y*r,z*r);
	}
};
struct Line3{
	Point3 s,e;
	Line3(){}
	Line3(Point3 _s,Point3 _e){s=_s;e=_e;}
	bool operator ==(const Line3 v){return (s==v.s)&&(e==v.e);}
	void input(){s.input();e.input();}
	double length(){return s.distance(e);}
	double dispointtoline(Point3 p){return ((e-s)^(p-s)).len()/s.distance(e);}//点到直线距离
	double dispointtoseg(Point3 p){//点到线段距离
		if(sgn((p-s)*(e-s))<0||sgn((p-e)*(s-e))<0)return min(p.distance(s),e.distance(p));
		return dispointtoline(p);
	}
	Point3 lineprog(Point3 p){return s+(((e-s)*((e-s)*(p-s)))/((e-s).len2()));}//返回点p在直线上的投影
	Point3 rotate(Point3 p,double ang){//p 绕此向量逆时针arg 角度
		if(!sgn(((s-p)^(e-p)).len()))return p;
		Point3 f1=(e-s)^(p-s),f2=(e-s)^(f1);
		double len=((s-p)^(e-p)).len()/s.distance(e);
		f1=f1.trunc(len),f2=f2.trunc(len);
		Point3 h=p+f2,pp=h+f1;
		return h+((p-h)*cos(ang))+((pp-h)*sin(ang));
	}
	bool pointonline(Point3 p){return !sgn(((s-p)^(e-p)).len());}//点在直线上
	bool pointonseg(Point3 p){return !sgn(((s-p)^(e-p)).len())&&sgn((s-p)*(e-p))<=0;}//点在线段上
	double linetoline(Line3 v){//两直线距离
		Point3 n=(s-e)^(v.s-v.e);
		if(!sgn(n.len()))return dispointtoline(v.s);//平行
		return fabs((s-v.s)*n)/n.len();
	}
	Point3 crosspoint(Line3 v){//求两直线的交点//要保证两直线不平行或重合
		double a1=((v.e-v.s)^(s-v.s)).len();
		double a2=((v.e-v.s)^(e-v.s)).len();
		return Point3((s.x*a2-e.x*a1)/(a2-a1),(s.y*a2-e.y*a1)/(a2-a1),(s.z*a2-e.z*a1)/(a2-a1));
	}
};
struct Plane{
	Point3 a,b,c,o;//平面上的三个点，以及法向量
	Plane(){}
	Plane(Point3 _a,Point3 _b,Point3 _c){a=_a;b=_b;c=_c;o=pvec();}
	Point3 pvec(){return (b-a)^(c-a);}
	Plane(double _a,double _b,double _c,double _d){//ax+by+cz+d=0
		o=Point3(_a,_b,_c);
		if(sgn(_a)!=0)a=Point3((-_d-_c-_b)/_a,1,1);
		else if(sgn(_b)!=0)a=Point3(1,(-_d-_c-_a)/_b,1);
		else if(sgn(_c)!=0)a=Point3(1,1,(-_d-_a-_b)/_c);
	}
	bool pointonplane(Point3 p){return !sgn((p-a)*o);}//点在平面上的判断
	double angleplane(Plane f){return acos(o*f.o)/(o.len()*f.o.len());}//两平面夹角
	int crossline(Line3 u,Point3 &p){//平面和直线的交点，返回值是交点个数
		double x=o*(u.e-a);
		double y=o*(u.s-a);
		double d=x-y;
		if(!sgn(d))return 0;
		p=((u.s*x)-(u.e*y))/d;
		return 1;
	}
	Point3 pointtoplane(Point3 p){//点到平面最近点(也就是投影)
		Line3 u=Line3(p,p+o);
		crossline(u,p);
		return p;
	}
	int crossplane(Plane f,Line3 &u){//平面和平面的交线
		Point3 oo=o^f.o;
		Point3 v=o^oo;
		double d=fabs(f.o*v);
		if(!sgn(d))return 0;
		Point3 q=a+(v*(f.o*(f.a-a))/d);
		u=Line3(q,q+oo);return 1;
	}
};
struct Point{
	double x,y;
	Point(double x=0,double y=0):x(x),y(y){}
	bool operator ==(const Point& b){return !sgn(x-b.x)&&!sgn(y-b.y);}
	bool operator  <(const Point& b)const{return !sgn(x-b.x)?y<b.y:x<b.x;}
	double operator^(const Point& b){return x*b.y-y*b.x;}//叉积
	double operator*(const Point& b){return x*b.x+y*b.y;}//点积
	Point operator +(const Point& b){return Point(x+b.x,y+b.y);}
	Point operator -(const Point& b){return Point(x-b.x,y-b.y);}
	Point operator *(const double& k){return Point(x*k,y*k);}
	Point operator /(const double& k){return Point(x/k,y/k);}
	double len(){return hypot(x,y);}
	double len2(){return x*x+y*y;}
	double rad(Point a,Point b){Point p=*this;return fabs(atan2(fabs((a-p)^(b-p)),(a-p)*(b-p)));}//pa和pb夹角
	double angle(){return atan2(y,x);}//倾斜角
	double angle(Point B){//夹角
		Point A=(*this);
		return acos(A*B/A.len()/B.len());
	}
	double distance(Point b){return hypot(x-b.x,y-b.y);}
	Point trunc(double r){//化为长度为r的向量
		double l=len();
		if(!sgn(l))return *this;
		r/=l;return Point(x*r,y*r);
	}
	Point rotleft(){return Point(-y,x);}//逆时针旋转90 度
	Point rotright(){return Point(y,-x);}//顺时针旋转90 度
	Point rotate(double rad){//rad为弧度 逆时针旋转
		Point A=(*this);
		double c=cos(rad),s=sin(rad);
		return Point(A.x*c-A.y*s,A.x*s+A.y*c);
	}
	Point rotate(Point p,double angle){//绕p逆时针旋转后点
		Point v=((*this)-p).rotate(angle);
		return Point(p.x+v.x,p.y+v.y);
	}
	void input(){scanf("%lf%lf",&x,&y);}
	void show(){printf("(%.2lf %.2lf)\n",x,y);}
	void zero(){
		if(!sgn(x))x=0.0;
		if(!sgn(y))y=0.0;
	}
};
struct polygon{
	int n;
	Point p[maxp];
	void input(int _n){n=_n;for(int i=0;i<n;++i)p[i].input();}
	void add(Point q){p[n++]=q;}
	struct cmp{
		Point p;
		cmp(const Point &p0){p=p0;}
		bool operator()(const Point &aa,const Point &bb){
			Point a=aa,b=bb;
			int d=sgn((a-p)^(b-p));
			if(!d)return sgn(a.distance(p)-b.distance(p))<0;
			return d>0;
		}
	};//进行极角排序//首先需要找到最左下角的点//需要重载号好Point 的< 操作符(min 函数要用)
	void norm(){
		Point mi=p[0];
		for(int i=1;i<n;++i)mi=min(mi,p[i]);
		sort(p,p+n,cmp(mi));
	}//得到凸包,点编号是0-n-1//两种凸包的方法
	void getconvex(polygon &convex){//注意如果有影响 要特判下所有点共点 或者共线的特殊情况//测试LightOJ1203 LightOJ1239
		sort(p,p+n),convex.n=n;
		for(int i=0;i<min(n,2);++i)convex.p[i]=p[i];
		if(convex.n==2&&(convex.p[0]==convex.p[1]))--convex.n;//特判
		if(n<=2)return;
		int &top=convex.n;top=1;
		for(int i=2;i<n;++i){
			while(top&&sgn((convex.p[top]-p[i])^(convex.p[top-1]-p[i]))<=0)--top;
			convex.p[++top]=p[i];
		}int temp=top;
		convex.p[++top]=p[n-2];
		for(int i=n-3;i>=0;--i){
			while(top!=temp&&sgn((convex.p[top]-p[i])^(convex.p[top-1]-p[i]))<=0)--top;
			convex.p[++top]=p[i];
		}if(convex.n==2&&(convex.p[0]==convex.p[1]))--convex.n;//特判
		convex.norm();//原来得到的是顺时针的点 排序后逆时针
	}
	double getarea(){//得到面积
		double sum=0;
		for(int i=0;i<n;++i)sum+=(p[i]^p[(i+1)%n]);
		return fabs(sum)/2;
	}
};
Point3 p[maxp],sum;
Plane sp;
polygon con,ans;
int n;
int main(){
	double a,b,c,d;
	while(scanf("%lf%lf%lf%lf",&a,&b,&c,&d),sgn(a)||sgn(b)||sgn(c)||sgn(d)){
		sp=Plane(a,b,c,-d);
		con.n=0;
		scanf("%d",&n);
		for(int i=0;i<n;++i)p[i].input();
		sum.input();
		int cnt=0;
		Point3 cross;
		for(int i=0;i<n;++i){
			Line3 L(sum,p[i]);
			if(sp.crossline(L,cross))
				if((p[i]-sum)*(cross-sum)>0)p[cnt++]=cross;
		}
		if(cnt==0)puts("0.00");
		else if(cnt<n)puts("Infi");
		else{
			Point3 o=sp.a;
			Point3 s=p[0];
			Point3 x=(s-o).trunc(1);
			Point3 y=((s-o)^sp.o).trunc(1);
			for(int i=0;i<cnt;++i)con.add(Point(x*(p[i]-o),y*(p[i]-o)));
			con.getconvex(ans);
			printf("%.2lf\n",ans.getarea());
		}
	}
	return 0;
}