显示代码纯文本
#include<cstdio>
#include<iostream>
#include<cmath>
using namespace std;
const int SIZEN=20010;
const double zero=1e-7;
int N;
class miku//点
{
public:
double x,y,z;
void print()
{
printf("%.4lf %.4lf\n",x,y);
}
}P[SIZEN],B[SIZEN];
miku operator - (miku a,miku b)
{
miku c;
c.x=a.x-b.x;c.y=a.y-b.y;
return c;
}
double operator / (miku a,miku b)
{
return a.x*b.x+a.y*b.y;
}
double operator *(miku a,miku b)//叉乘
{
return a.x*b.y-a.y*b.x;
}
double R;
void read()
{
scanf("%d%lf",&N,&R);
for(int i=1;i<=N;i++)
scanf("%lf%lf%lf",&P[i].x,&P[i].y,&P[i].z);
}
miku a[SIZEN],b[SIZEN];
double dis(miku a,miku b)
{
return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));
}
bool check(double now)//这是我写过的最麻烦的check
{
//现在问题已经被简化为:已知一条上折线和一条下折线,判断是否存在一条可以从上下折线中穿过的直线
//只需要把对上下两条折线各求一个凸包,然后判断两个凸包的距离即可
//进一步分析,由于下折线与上折线平行,所以我们只用求半个凸包
int topa,topb;
topa=topb=0;
for(int i=1;i<=N;i++)
{
miku o;
o=B[i];
if(B[i].x>now)
{
o.x=now;
o.y=B[i-1].y+(B[i].y-B[i-1].y)/(B[i].x-B[i-1].x)*(now-B[i-1].x);
}
while(topa>1&&(a[topa]-a[topa-1])*(o-a[topa-1])<zero) topa--;
a[++topa]=o;
o.y--;
while(topb>1&&(b[topb]-b[topb-1])*(o-b[topb-1])>-zero) topb--;
b[++topb]=o;
if(B[i].x>now) break;
}
//判断是否相交
//for(int i=1;i<=topa;i++) a[i].print();
//cout<<endl;
//cout<<"arib"<<endl;
for(int i=1,q=1;i<=topa&&q<=topb;i++)
{
while(q<topb&&b[q].x<a[i].x+zero) q++;
double tem=b[q-1].y+(b[q].y-b[q-1].y)/(b[q].x-b[q-1].x)*(a[i].x-b[q-1].x);
if(tem>a[i].y+zero) return 0;
}
for(int i=1,q=1;i<=topb&&q<=topa;i++)
{
while(q<topa&&a[q].x<b[i].x+zero) q++;
double tem=a[q-1].y+(a[q].y-a[q-1].y)/(a[q].x-a[q-1].x)*(b[i].x-a[q-1].x);
if(tem<b[i].y-zero) return 0;
}
//
//判断距离
for(int i=1,q=1;i<=topa&&q<=topb;i++)
{
while(q<topb&&(b[q+1]-b[q])/(a[i]-b[q])>zero&&(b[q+1]-b[q])*(a[i]-b[q])>zero) q++;
double dist;
if(q==1||(b[q-1]-b[q])/(a[i]-b[q])<zero) dist=dis(a[i],b[q]);
else dist=abs((b[q]-b[q-1])*(a[i]-b[q-1])/dis(b[q],b[q-1]));
if(dist<2*R-zero) return 0;
}
for(int i=1,q=1;i<=topb&&q<=topa;i++)
{
while(q<topa&&(a[q+1]-a[q])/(b[i]-a[q])>-zero&&(a[q+1]-a[q])*(b[i]-a[q])<zero) q++;
double dist;
if(q==1||(a[q-1]-a[q])/(b[i]-a[q])<zero) dist=dis(b[i],a[q]);
else dist=abs((a[q]-a[q-1])*(b[i]-a[q-1])/dis(a[q],a[q-1]));
if(dist<2*R-zero) return 0;
}
return 1;
}
double getdis()//
{
double l=B[1].x,r=B[N].x;
int tot=0;
while(r-l>zero)
{
double mid=(l+r)/2.0;
if(check(mid)) l=mid;
else r=mid;
}
return l;
}
void work()
{
double tem1,tem2;
for(int i=1;i<=N;i++) B[i].x=P[i].y,B[i].y=P[i].x;//先只考虑y-x平面上的点
tem1=getdis();
for(int i=1;i<=N;i++) B[i].x=P[i].y,B[i].y=P[i].z;//只考虑y-z平面上的点
tem2=getdis();
double ans;
ans=min(tem1,tem2);
printf("%.3lf",ans);
}
int main()
{
freopen("nt2011_cs.in","r",stdin);
freopen("nt2011_cs.out","w",stdout);
read();
work();
return 0;
}