bzoj-1502 月下柠檬树

题意：

一颗树由n个圆台组成，现在有倾斜角为alpha的光；

不计树干阴影，光线沿直线传播，求这个树在水平地面投影的面积；

n<=500，0.3<alpha<π/2；

题解：

给出一颗树。。求这棵树的阴影面积。。。

什么鬼题！= =

当然这是一道计算几何；

首先我们把树的尖看成半径为0的圆；

三维图形到二维的投影貌似很难啊，不过这题也有很多特殊性；

因为圆对于平行的面的投影都是相同的正圆，所以可以直接投在地面上；

而圆心距则是除了一个tan函数的样子，所以圆已经被扔到待求的平面上了；

曾经的圆台的母线，现在是什么呢？

这些圆的公切线！(这样就从母的变成公的啦233)

公切线的求法就是画图上勾股定理相似乱搞，调两组数据改改就好了；

然后正解似乎是讨论了一堆东西然后分别求面积；

然后我们似乎可以直接上Simpson积分！

具体细节不说了，代码里都有；

这题精度不是特别卡，EPS=1e5就够了？

代码：

#include<math.h>
#include<stdio.h>
#include<string.h>
#include<algorithm>
#define N 550
#define pr pair<double,double>
using namespace std;
const double EPS=1e-5;
const double INF=1e100;
struct Point
{
	double x,y;
	Point(){}
	Point(double _,double __):x(_),y(__){}
	friend Point operator +(Point a,Point b)
	{
		return Point(a.x+b.x,a.y+b.y);
	}
	friend Point operator -(Point a,Point b)
	{
		return Point(a.x-b.x,a.y-b.y);
	}
	friend Point operator *(double a,Point b)
	{
		return Point(a*b.x,a*b.y);
	}
};
struct Line
{
	Point p,v;
	Line(){};
	Line(Point _,Point __){p=_,v=__-_;}
	pr f(double x)
	{
		if(x<p.x||x>(p+v).x)	return pr(0,0);
		double ret=(p+(x-p.x)/v.x*v).y;
		return pr(-ret,ret);
	}
}l[N];
struct Circle
{
	Point p;
	double r;
	pr f(double x)
	{
		if(r-fabs(x-p.x)<=0)	return pr(0,0);
		double ret=sqrt(r*r-(x-p.x)*(x-p.x));
		return pr(-ret,ret);
	}
}O[N];
int n,tot;
pr p[N<<1];
void getL(Circle A,Circle B)
{
	double c=B.p.x-A.p.x;
	if(c+A.r<=B.r||c+B.r<=A.r)	return ;
	double a=B.r-A.r;
	double b=sqrt(c*c-a*a);
	Point lp=A.p-Point(A.r/c*a,-A.r/c*b),rp=B.p-Point(B.r/c*a,-B.r/c*b);//rp++
	l[++tot]=Line(lp,rp);
}
double Cut(double x)
{
	int cnt=0;
	double ret=0,last=-INF;
	for(int i=1;i<=tot;i++)
	{
		p[++cnt]=l[i].f(x);
		if(p[cnt]==pr(0,0))	cnt--;
	}
	for(int i=1;i<=n;i++)
	{
		p[++cnt]=O[i].f(x);
		if(p[cnt]==pr(0,0))	cnt--;
	}
	sort(p+1,p+cnt+1);
	for(int i=1;i<=cnt;i++)
	{
		if(p[i].first>last)
			ret+=p[i].second-p[i].first,last=p[i].second;
		else if(p[i].second>last)
			ret+=p[i].second-last,last=p[i].second;
	}
	return ret;
}
double Simpson(double l,double r,double mid,double Cl,double Cr,double Cm)
{
	double tCl=Cut((l+mid)/2),tCr=Cut((mid+r)/2);
	double ans=(r-l)*(Cl+Cr+4*Cm)/6;
	double lans=(mid-l)*(Cl+Cm+4*tCl)/6,rans=(r-mid)*(Cm+Cr+4*tCr)/6;
	if(fabs(lans+rans-ans)<EPS)
		return ans;
	else
		return Simpson(l,mid,(l+mid)/2,Cl,Cm,tCl)+Simpson(mid,r,(mid+r)/2,Cm,Cr,tCr);
}
int main()
{
	int m,i,j,k;
	double alpha,h,l,r;
	scanf("%d%lf",&n,&alpha);
	n++;
	alpha=1.0/tan(alpha);
	for(i=1;i<=n;i++)
	{
		scanf("%lf",&h);
		O[i].p.x=O[i-1].p.x+h*alpha;
	}
	for(i=1,l=INF,r=-INF;i<n;i++)
	{
		scanf("%lf",&O[i].r);
		l=min(l,O[i].p.x-O[i].r);
		r=max(r,O[i].p.x+O[i].r);
	}
	l=min(l,O[n].p.x-O[n].r);
	r=max(r,O[n].p.x+O[n].r);
	for(i=1;i<n;i++)
	{
		getL(O[i],O[i+1]);
	}
	printf("%.2lf\n",Simpson(l,r,(l+r)/2,0,0,Cut((l+r)/2)));
	return 0;
}