BZOJ3564 信号增幅仪

http://www.lydsy.com/JudgeOnline/problem.php?id=3564

思路:先旋转坐标系，再缩进x坐标，把椭圆变成圆，然后做最小圆覆盖。

还有，为什么用srand()又错了啊。。。

#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<iostream>
const double Pi=acos(-1);
const double eps=1e-7;
struct Point{
    double x,y;
    Point(){}
    Point(double x0,double y0):x(x0),y(y0){}
}p[200005];
int n;
struct Line{
    Point s,e;
    Line(){}
    Line(Point s0,Point e0):s(s0),e(e0){}
};
int read(){
    int t=0,f=1;char ch=getchar();
    while (ch<'0'||ch>'9'){if (ch=='-') f=-1;ch=getchar();}
    while ('0'<=ch&&ch<='9'){t=t*10+ch-'0';ch=getchar();}
    return t*f;
}
Point operator +(Point p1,Point p2){
    return Point(p1.x+p2.x,p1.y+p2.y);
}
Point operator -(Point p1,Point p2){
    return Point(p1.x-p2.x,p1.y-p2.y);
}
double operator *(Point p1,Point p2){
    return p1.x*p2.y-p1.y*p2.x;
}
Point operator /(Point p1,double x){
    return Point(p1.x/x,p1.y/x);
}
Point operator *(Point p,double x){
    return Point(p.x*x,p.y*x);
}
double sqr(double x){return x*x;}
double dis(Point p){
    return sqrt(sqr(p.x)+sqr(p.y));
}
double dis(Point p1,Point p2){
    return dis(p1-p2);
}
Point turn(Point p,double ang){
    double Cos=cos(ang);
    double Sin=sin(ang);
    double x=p.x*Cos-p.y*Sin;
    double y=p.x*Sin+p.y*Cos;
    return Point(x,y);
}
Point inter(Line p1,Line p2){
    double k1=(p2.e-p1.s)*(p1.e-p1.s);
    double k2=(p1.e-p1.s)*(p2.s-p1.s);
    double t=(k2)/(k1+k2);
    double x=p2.s.x+(p2.e.x-p2.s.x)*t;
    double y=p2.s.y+(p2.e.y-p2.s.y)*t;
    return Point(x,y);
}
Point solve(Point p1,Point p2,Point p3){
    Point a=(p1+p2)/2.0;
    Point b=(p2+p3)/2.0;
    return inter(Line(a,a+turn(p2-a,Pi/2.0)),Line(b,b+turn(p3-b,Pi/2.0)));
}
int main(){
    //srand(233);
    n=read();
    for (int i=1;i<=n;i++)
        p[i].x=read(),p[i].y=read();
    double ang=read();ang/=180.0;ang*=Pi;
    double len=read();
    for (int i=1;i<=n;i++)
       p[i]=turn(p[i],-ang);
    for (int i=1;i<=n;i++)
       p[i].x/=len;
    for (int i=1;i<=n;i++)
        std::swap(p[rand()%n+1],p[rand()%n+1]);
    Point O=p[1];double r=0;
    for (int i=2;i<=n;i++)
      if (dis(p[i],O)+eps>r){
        O=p[i];r=0.0;
        for (int j=1;j<i;j++)
            if (dis(p[j],O)+eps>r){
               O=(p[i]+p[j])/2.0;
               r=dis(O,p[j]);
               for (int k=1;k<j;k++)
                  if (dis(p[k],O)+eps>r){
                    O=solve(p[i],p[j],p[k]);
                    r=dis(O,p[k]);
                  }
            }
      } 
    printf("%.3lf
",r);
}