#include<bits/stdc++.h>
#define vct point
using namespace std;
const double pi=atan2(0,-1);
const double eps=1e-8;
int sgn(double d){
if(fabs(d)<eps)return 0;
return d<0?-1:1;
}
struct point{
double x,y;
bool operator<(point b)const{
if(sgn(x-b.x))return sgn(x-b.x)<0;
return sgn(y-b.y)<0;
}
bool operator==(point b)const{return sgn(x-b.x)==0&&sgn(y-b.y)==0;}
//基本运算
vct operator+(vct b){return {x+b.x,y+b.y};}
vct operator-(vct b){return {x-b.x,y-b.y};}
vct operator*(double d){return {x*d,y*d};}
vct operator/(double d){return {x/d,y/d};}
double dot(vct b){return x*b.x+y*b.y;}
double cross(vct b){return x*b.y-y*b.x;}
double dis(){return hypot(x,y);}
double dis(point b){return hypot(x-b.x,y-b.y);}
double agl(){return atan2(y,x);}
double agl(vct b){return acos(dot(b)/dis()/b.dis());}
vct rot(double d){return {x*cos(d)-y*sin(d),x*sin(d)+y*cos(d)};}
//线段
int cmp(point a,point b){//点在线段上1,端点上0,外-1
if(*this==a||*this==b)return 0;
vct pa=a-*this,pb=b-*this;
return sgn(pa.cross(pb))==0&&pa.dot(pb)<0?1:-1;
}
double dis(point a,point b){//点到线段ab的距离
if(a==b)return dis(a);
vct ab=b-a,ap=*this-a,bp=*this-b;
if(sgn(ab.dot(ap))<0)return ap.dis();
if(sgn(ab.dot(bp))>0)return bp.dis();
return fabs(ab.cross(ap)/ab.dis());
}
//多边形
template<typename T>
int cmp(T first,T last){//点在多边形内1,上0,外-1
int wn=0;
for(T i=first;i<last;i++){
point a=*i,b=(i+1)==last?*first:*(i+1);
if(cmp(a,b)==0)return 0;
int k=sgn((b-a).cross(*this-a));
int d1=sgn(a.y-y);
int d2=sgn(b.y-y);
if(k>0&&d1<=0&&d2>0)wn++;
if(k<0&&d2<=0&&d1>0)wn--;
}
return wn!=0?1:-1;
}
};
struct line{
point p,v;
line(point p,vct v):p(p),v(v/v.dis()){}//方向向量必须单位化,小心除0错误
point getp(double t){return p+v*t;}
point intersection(line l){//直线交点
double t=l.v.cross(p-l.p)/v.cross(l.v);//平行分母为0,小心除0错误
return getp(t);
}
double dis(point q){return fabs(v.cross(q-p));}//不取绝对值就是有向距离
point projection(point q){//点在直线上的投影
double t=v.dot(q-p);//向量pq在单位向量上的投影长度
return getp(t);
}
};
//线段
bool intersect(point a,point b,point c,point d){//线段相交
return sgn((b-a).cross(c-a))*sgn((b-a).cross(d-a))<0
&&sgn((d-c).cross(a-c))*sgn((d-c).cross(b-c))<0;
}
double agl(double a,double b,double c){return acos((a*a+b*b-c*c)/2/a/b);}//角c
struct circle{
point o;
double r;
point getp(double d){return o+point{cos(d),sin(d)}*r;}
int intersection(line l,vector<point>&ans){//直线与圆交点
double d=l.dis(o);
if(sgn(d-r)>0)return 0;//相离
point q=l.projection(o);
if(sgn(d-r)==0){
ans.push_back(q);
return 1;//相切
}
double t=sqrt(r*r-d*d);
ans.push_back(q+l.v*t);
ans.push_back(q-l.v*t);
return 2;//相交
}
int intersection(circle c,vector<point>&ans){
double d=o.dis(c.o);
if(sgn(d)==0){
if(sgn(r-c.r)==0)return -1;//重合
return 0;//内含
}
if(sgn(r+c.r-d)<0)return 0;//外离
if(sgn(fabs(r-c.r)-d)>0)return 0;//内含
double a=(o-c.o).agl();//o1o2极角
double b=agl(r,c.r,d);
point p1=getp(a+b),p2=getp(a-b);
ans.push_back(p1);
if(p1==p2)return 1;//相切
ans.push_back(p2);
return 2;//相交
}
int tangent(point p,vector<line>&ans){//过点p做圆的切线
vct po=o-p;
double d=p.dis(o);
if(sgn(d-r)<0)return 0;
if(sgn(d-r)==0){//点在圆上
ans.push_back(line(p,po.rot(pi/2)));
return 1;
}
double a=asin(d/r);
ans.push_back(line(p,po.rot(a)));
ans.push_back(line(p,po.rot(-a)));
return 2;
}
};
int tangent(circle a,circle b,vector<point>&va,vector<point>&vb){//va,vb是切点
if(a.r<b.r)swap(a,b);//使a半径大
double d=a.o.dis(b.o);
double rsub=a.r-b.r;
double rsum=a.r+b.r;
if(sgn(rsub-d)>0)return 0;//内含
if(sgn(d)==0&&sgn(rsub)==0)return -1;//重合
double base=(b.o-a.o).agl();
if(sgn(rsub-d)==0){//内切
va.push_back(a.getp(base));
vb.push_back(b.getp(base));
return 1;
}
//有外公切线
double ang=acos((rsub)/d);
va.push_back(a.getp(base+ang));
vb.push_back(b.getp(base+ang));
va.push_back(a.getp(base-ang));
vb.push_back(b.getp(base-ang));
if(sgn(d-rsum)==0){//外切,一条内公切线
va.push_back(a.getp(base));
vb.push_back(b.getp(base));
return 3;
}
if(sgn(d-rsum)>0){//外离,两条外公切线
double ang=acos(rsum/d);
va.push_back(a.getp(base+ang));
vb.push_back(b.getp(base+ang));
va.push_back(a.getp(base-ang));
vb.push_back(b.getp(base-ang));
return 4;
}
return 2;//相交
}
计算几何模板
最新推荐文章于 2024-03-31 15:51:24 发布