Gym 101917 E 简单计算几何，I 最大流

  1 #include <bits/stdc++.h>
  2 #define Vector Point
  3 using namespace std;
  4 typedef long long ll;
  5 const double eps = 1e-9;
  6 const double PI=acos(-1);
  7 const int maxn = 1e5+10,inf=0x3f3f3f3f;
  8 //head-------------------------------------------------------------------
  9 int dcmp(double x) { return fabs(x) < eps ? 0 : (x < 0 ? -1 : 1); }
 10 struct Point {
 11     double x, y;
 12 
 13     Point(const Point& rhs): x(rhs.x), y(rhs.y) { } //拷贝构造函数
 14     Point(double x = 0.0, double y = 0.0): x(x), y(y) { }   //构造函数
 15 
 16     friend istream& operator >> (istream& in, Point& P) { return in >> P.x >> P.y; }
 17     friend ostream& operator << (ostream& out, const Point& P) { return out << P.x << ‘ ‘ << P.y; }
 18 
 19     friend Vector operator + (const Vector& A, const Vector& B) { return Vector(A.x+B.x, A.y+B.y); }
 20     friend Vector operator - (const Point& A, const Point& B) { return Vector(A.x-B.x, A.y-B.y); }
 21     friend Vector operator * (const Vector& A, const double& p) { return Vector(A.x*p, A.y*p); }
 22     friend Vector operator / (const Vector& A, const double& p) { return Vector(A.x/p, A.y/p); }
 23     friend bool operator == (const Point& A, const Point& B) { return dcmp(A.x-B.x) == 0 && dcmp(A.y-B.y) == 0; }
 24     friend bool operator < (const Point& A, const Point& B) { return dcmp(A.x - B.x)<0 || (dcmp(A.x - B.x)==0 && dcmp(A.y - B.y)<0); }
 25     void in() { scanf("%lf%lf", &x, &y); }
 26     void out() { printf("%.10f %.10f\n", x, y); }
 27 };
 28 
 29 template <class T> T sqr(T x) { return x * x;}
 30 double Dot(const Vector& A, const Vector& B) { return A.x*B.x + A.y*B.y; }  //点积
 31 double Length(const Vector& A){ return sqrt(Dot(A, A)); } //向量长度
 32 double Angle(const Vector& A, const Vector& B) { return acos(Dot(A, B)/Length(A)/Length(B)); }  //AB向量夹角
 33 double Cross(const Vector& A, const Vector& B) { return A.x*B.y - A.y*B.x; }    //AB叉积 有向面积
 34 double Area(const Point& A, const Point& B, const Point& C) { return fabs(Cross(B-A, C-A)); }
 35 //向量旋转 rad正表示逆时针 反之顺时针
 36 Vector Rotate(Vector A,double rad){ return Vector(A.x*cos(rad)-A.y*sin(rad),A.x*sin(rad)+A.y*cos(rad));}
 37 //A不能是0向量
 38 Vector normal(Vector A) { return Point(-A.y, A.x) / Length(A);} //向量A的单位法向量，即A左转90度 以后把长度归归一化
 39 double angle(Vector A) { return atan2(A.y, A.x);} //向量极角 向量A与x轴的夹角
 40 Vector vecunit(Vector A){ return A / Length(A);} //单位向量
 41 
 42 struct Line {
 43     Point P;    //直线上一点
 44     Vector dir; //方向向量(半平面交中该向量左侧表示相应的半平面)
 45     double ang; //极角，即从x正半轴旋转到向量dir所需要的角（弧度）
 46 
 47     Line() { }  //构造函数
 48     Line(const Line& L): P(L.P), dir(L.dir), ang(L.ang) { }
 49     Line(const Point& P, const Vector& dir): P(P), dir(dir) { ang = atan2(dir.y, dir.x); }
 50 
 51     bool operator < (const Line& L) const { //极角排序
 52         return ang < L.ang;
 53     }
 54     Point point(double t) { return P + dir*t; }
 55 };
 56 
 57 
 58 //直线交点1 P+tv Q+tw
 59 Point GetLineIntersection(Point P,Vector v,Point Q,Vector w)
 60 {
 61     Vector u=P-Q;
 62     double t = Cross(w,u)/Cross(v,w);
 63     return P+v*t;
 64 }
 65 //直线交点2
 66 Point GetIntersection(Line a, Line b)
 67 {
 68     Vector u = a.P-b.P;
 69     double t = Cross(b.dir, u) / Cross(a.dir, b.dir);
 70     return a.P + a.dir*t;
 71 }
 72 bool SegmentProperInntersection(Point a1,Point a2,Point b1,Point b2)//线段相交判定
 73 {
 74     double c1=Cross(a2-a1,b1-a1),c2=Cross(a2-a1,b2-a1),
 75            c3=Cross(b2-b1,a1-b1),c4=Cross(b2-b1,a2-b1);
 76     return dcmp(c1)*dcmp(c2)<0&&dcmp(c3)*dcmp(c4)<0;
 77 }
 78 double DistanceToSegment(Point P,Point A,Point B) //点P到线段AB的距离
 79 {
 80     if(A==B) return Length(P-A);
 81     Vector v1 = B - A,v2 = P - A,v3 = P - B;
 82     if(dcmp(Dot(v1,v2))<0) return Length(v2);
 83     else if(dcmp(Dot(v1,v3))>0) return Length(v3);
 84     else return fabs(Cross(v1,v2))/Length(v1);
 85 }
 86 Point GetLineProjection(Point P,Point A,Point B) //点P在直线AB上的投影
 87 {
 88     Vector v = B-A;
 89     return A+v*(Dot(v,P-A)/Dot(v,v));
 90 }
 91 bool OnSegment(Point p, Point a1, Point a2) //判断点是否在线段a1a2上 不含端点
 92 {
 93     return dcmp(Cross(a1-p, a2-p)) == 0 && dcmp(Dot(a1-p, a2-p)) < 0;
 94 }
 95 
 96 struct Circle {
 97     Point c;    //圆心
 98     double r;   //半径
 99 
100     Circle() { }
101     Circle(const Circle& rhs): c(rhs.c), r(rhs.r) { }
102     Circle(const Point& c, const double& r): c(c), r(r) { }
103 
104     Point point(double ang) const { return Point(c.x + cos(ang)*r, c.y + sin(ang)*r); } //圆心角所对应的点
105     double area(void) const { return PI * r * r; }
106 };
107 bool InCircle(Point x, Circle c)
108 {
109     return dcmp(c.r*c.r - Length(c.c - x)*Length(c.c - x)) >= 0;
110 }
111 //直线与圆的交点
112 int getLineCircleIntersection(Line L, Circle C, Point* sol) //函数返回交点个数 sol数组存放交点
113 {
114     Vector nor = normal(L.dir);
115     Line pl = Line(C.c, nor);
116     Point ip = GetIntersection(pl, L);
117     double dis = Length(ip - C.c);
118     if (dcmp(dis - C.r) > 0) return 0;
119     Point dxy = vecunit(L.dir) * sqrt(C.r*C.r - dis*dis);
120     int ret = 0;
121     sol[ret] = ip + dxy;
122     if (OnSegment(sol[ret], L. P, L.point(1))) ret++;
123     sol[ret] = ip - dxy;
124     if (OnSegment(sol[ret], L.P, L.point(1))) ret++;
125     return ret;
126 }
127 //圆与圆的交点
128 int getCircleCircleIntersection(Circle C1,Circle C2,vector<Point>& sol)//函数返回交点个数 sol数组存放交点
129 {
130     double d=Length(C1.c-C2.c);
131     if(dcmp(d)==0){
132         if(dcmp(C1.r-C2.r)==0) return -1;//两圆重合
133         return 0;
134     }
135     if(dcmp(C1.r+C2.r-d)<0) return 0;
136     if(dcmp(fabs(C1.r-C2.r)-d)>0) return 0;
137 
138     double a=angle(C2.c-C1.c);
139     double da=acos((C1.r*C1.r+d*d-C2.r*C2.r)/(2*C1.r*d));
140     Point p1=C1.point(a-da),p2=C1.point(a+da);
141     sol.push_back(p1);
142     if(p1==p2) return 1;
143     sol.push_back(p2);
144     return 2;
145 }
146 //过点p到圆c的切线
147 int getTangents(Point p,Circle C,Vector* v)//函数返回条数,v[i]是第i条切线的向量
148 {
149     Vector u=C.c-p;
150     double dist=Length(u);
151     if(dist<C.r) return 0;  //点在圆内
152     else if(dcmp(dist-C.r)==0){//p在圆c上，只有一条切线
153         v[0]=Rotate(u,PI/2);
154         return 1;
155     }else{
156         double ang = asin(C.r/dist);
157         v[0]=Rotate(u,-ang);
158         v[1]=Rotate(u,ang);
159         return 2;
160     }
161 }
162 //两圆的公切线
163 int getTangents(Circle A,Circle B,Point* a,Point* b)
164 {
165     int cnt=0;
166     if(A.r<B.r){swap(A,B);swap(a,b);}
167     double d2=(A.c.x-B.c.x)*(A.c.x-B.c.x)+(A.c.y-B.c.y)*(A.c.y-B.c.y);
168     double rdiff=A.r - B.r;
169     double rsum=A.r + B.r;
170     if(d2<rdiff*rdiff) return 0;//内含
171 
172     double base = atan2(B.c.y-A.c.y,B.c.x-A.c.x);
173     if(dcmp(d2)==0&&dcmp(A.r-B.r)==0) return -1;//两圆重合,无数条切线
174     if(dcmp(d2-rdiff*rdiff)==0){                 //内切,1条切线
175         a[cnt]=A.point(base);b[cnt]=B.point(base);cnt++;
176         return 1;
177     }
178     //有外公切线
179     double ang=acos((A.r-B.r)/sqrt(d2));
180     a[cnt] = A.point(base+ang);b[cnt] = B.point(base+ang);cnt++;
181     a[cnt] = A.point(base-ang);b[cnt] = B.point(base-ang);cnt++;
182     if(dcmp(d2-rsum*rsum)==0){       //一条内公切线
183         a[cnt]=b[cnt]=A.point(base);cnt++;
184     }
185     else if(dcmp(d2-rsum*rsum)>0){   //两条内公切线
186         ang = acos((A.r+B.r)/sqrt(d2));
187         a[cnt] = A.point(base+ang);b[cnt] = B.point(PI+base+ang);cnt++;
188         a[cnt] = A.point(base-ang);b[cnt] = B.point(PI+base-ang);cnt++;
189     }
190     return cnt;
191 }
192 double SegCircleArea(Circle C, Point a, Point b) //线段切割圆
193 {
194     double a1 = angle(a - C.c);
195     double a2 = angle(b - C.c);
196     double da = fabs(a1 - a2);
197     if (da > PI) da = PI * 2.0 - da;
198     return dcmp(Cross(b - C.c, a - C.c)) * da * sqr(C.r) / 2.0;
199 }
200 
201 double PolyCiclrArea(Circle C, Point *p, int n)//多边形与圆相交面积
202 {
203     double ret = 0.0;
204     Point sol[2];
205     p[n] = p[0];
206     for(int i=0;i<n;i++)
207     {
208         double t1, t2;
209         int cnt = getLineCircleIntersection(Line(p[i], p[i+1]-p[i]), C, sol);
210         if (cnt == 0)
211         {
212             if (!InCircle(p[i], C) || !InCircle(p[i+1], C)) ret += SegCircleArea(C, p[i], p[i+1]);
213             else ret += Cross(p[i+1] - C.c, p[i] - C.c) / 2.0;
214         }
215         if (cnt == 1)
216         {
217             if (InCircle(p[i], C) && !InCircle(p[i+1], C)) ret += Cross(sol[0] - C.c, p[i] - C.c) / 2.0, ret += SegCircleArea(C, sol[0], p[i+1]);
218             else ret += SegCircleArea(C, p[i], sol[0]), ret += Cross(p[i+1] - C.c, sol[0] - C.c) / 2.0;
219         }
220         if (cnt == 2)
221         {
222             if ((p[i] < p[i + 1]) ^ (sol[0] < sol[1])) swap(sol[0], sol[1]);
223             ret += SegCircleArea(C, p[i], sol[0]);
224             ret += Cross(sol[1] - C.c, sol[0] - C.c) / 2.0;
225             ret += SegCircleArea(C, sol[1], p[i+1]);
226         }
227     }
228     return fabs(ret);
229 }
230 double PolygonArea(Point *po, int n) //多边形面积
231 {
232     double area = 0.0;
233     for(int i = 1; i < n-1; i++) {
234         area += Cross(po[i]-po[0], po[i+1]-po[0]);
235     }
236     return area * 0.5;
237 }
238 //-----------------------------------------------------------------------------------
239 Point p[maxn];
240 int main()
241 {
242     double a,w,h;
243     int n;
244     cin>>a>>w>>h>>n;
245     double ang=a/180*PI;
246     double minx=inf,maxx=-1.0;
247     double miny=inf,maxy=-1.0;
248     for(int i=0;i<n;i++)
249     {
250         p[i].in();
251         p[i]=p[0]+Rotate(p[i]-p[0],ang);
252         minx=min(minx,p[i].x);
253         miny=min(miny,p[i].y);
254         maxx=max(maxx,p[i].x);
255         maxy=max(maxy,p[i].y);
256     }
257     double dux=w/(maxx-minx);
258     double duy=h/(maxy-miny);
259     for(int i=0;i<n;i++)
260     {
261         p[i].x-=minx;
262         p[i].y-=miny;
263         p[i].x*=dux;
264         p[i].y*=duy;
265         p[i].out();
266     }
267 }