Computational Geometry Templates




//author: birdstorm
#include<bits/stdc++.h>
using namespace std;
/*************useful TIPS in complex*********//**
*1. use eps carefully
*2. remember useful functions in complex:
*   abs() for distance,
*   arg() for angle,
*   conj() for conjunction,
*   norm() for |v|^2,
*   polar(double mag, double ang=0.0) for creating complex
*   real(), imag(), ...
*3. useful functions 2:
*   cos(), sin(), acos(), asin(), tan(), atan(), cosh(), sinh(), tanh(),
*   pow(), exp(), sqrt(), log(), ...
*   can also be used in complex directly!
*4. debug with cout<<v<<endl;
*5. we can also use "type _Complex" in definition(C99/C++11)
*   e.g. double _Complex p; float _Complex q;
*6. To Be Continued...
**************************************************/
/*********CONSTANTS**************/
const double eps=1.0e-8;
const double pi=acos(-1.0);
/***************DEFINATION*********************/
typedef complex<double> point;
typedef vector<point> polygon;
typedef point Vec;
struct line : public vector<point> {
    line() {}
    line(const point& a, const point& b) { push_back(a); push_back(b); }
};
typedef vector<line> Obstacle;
struct Circle{
    point center;
    double radius;
    Circle(const point& c, const double& r):center(c),radius(r) {}
};
struct polygon_convex{
    polygon p;
    polygon_convex(int sz=0){
        p.resize(sz);
    }
};
/******************COMMONLY-USED FUNCTIONS*************/
point unit(const point& v) { return v/abs(v);}
point ortho(const point& v) { return v*point(0,1);}
point spin(const point& v, const double &angle) {return v*point(cos(angle),sin(angle));}
inline point vec(const line& l) { return l[1]-l[0];}
bool equal(const double a, const double b) { return abs(a-b)<eps;}
bool equal(const point& a, const point& b) { return abs(a-b)<eps;}
inline double dot  (const point& a, const point& b) { return (a*conj(b)).real();}
inline double cross(const point& a, const point& b) { return (conj(a)*b).imag();}
inline int dlcmp(const double& x){
    return x<-eps?-1:x>eps;
}
/**************COMMONLY-USED FUNCTIONS 2***************/
inline int ccw(const point& a, const point& b, const point& c) {
    point u=b-a, v=c-a;
    if(cross(u,v)>0) return 1;
    if(cross(u,v)<0) return -1;
    if(dot(u,v)<0) return 2;
    if(abs(u)<abs(v)) return -2;
    return 0;
}
inline int ccw(const line& s, const point& p) {
    return ccw(s[0], s[1], p);
}
/**************COMPARISON********************/
namespace std{
    bool operator<(const point& a, const point& b) {
        return a.real()!=b.real()?a.real()<b.real():a.imag()<b.imag();
    }
    bool operator==(const point& a, const point& b) {
        return abs(a-b)<eps;
    }
}
bool cmp_less(const point& a, const point& b){
    return a.real()!=b.real()?a.real()<b.real():a.imag()<b.imag();
}
/**************INTERSECTIONS Judging****************/
bool intersectLP(const line& l, const point& p) {
    return abs(cross(l[1]-p, l[0]-p)) < eps;
}
bool intersectSP(const line& s, const point& p) {
    return abs(s[0]-p)+abs(s[1]-p) < abs(s[1]-s[0])+eps;
}
bool intersectLL(const line& l, const line& m) {
    return cross(vec(l), vec(m)) > eps || cross(vec(l), m[0]-l[0]) < eps;
}
bool intersectLS(const line& l, const line& s) {
    return cross(vec(l),s[0]-l[0]) * cross(vec(l),s[1]-l[0]) <= 0;
}
bool intersectSS(const line& s, const line& t) {
    return ccw(s,t[0])*ccw(s,t[1]) <= 0 && ccw(t,s[0])*ccw(t,s[1]) <= 0;
}
/************INTERSECTIONS Calculating***************/
point crosspoint(const line& l, const line& m) {    ///use Intersection Judging first
    double A = cross(vec(l), vec(m));
    double B = cross(vec(l), l[1]-m[0]);
    if(abs(A)<eps) {    // parallel
        return m[0];    // sameline
    }
    return m[0] + B/A*vec(m);
}
point intersection_point(const line& A, const Circle& C){
    double dx=real(A[1]-A[0]), dy=imag(A[1]-A[0]);
    double rx=real(A[0]-C.center), ry=imag(A[0]-C.center);
    double a=dx*dx+dy*dy, b=2*dx*rx+2*dy*ry, c=rx*rx+ry*ry-C.radius*C.radius;
    double delta=b*b-4*a*c;
    if(delta<-eps) return point(-10000,-10000);
    double ans=(-b+sqrt(delta))/(2.0*a);
    return point(A[0].real()+ans*dx,A[0].imag()+ans*dy);
}
point projection(const line &l, const point &p) {
    double t = dot(p-l[0], l[0]-l[1]) / norm(l[0]-l[1]);
    return l[0] + t*(l[0]-l[1]);
}
double distanceSP(const line &s, const point &p) {
    const point r = projection(s, p);
    if (intersectSP(s, r)) return abs(r - p);
    return min(abs(s[0] - p), abs(s[1] - p));
}
line generateSP(const line &s, const point &p) {
    const point r = projection(s, p);
    if(intersectSP(s, r)) return line(r, p);
    else{
        if(abs(s[0]-p) < abs(s[1] - p)) return line(s[0], p);
        else return line(s[1], p);
    }
    return s;
}
bool isSameLine(const point& a, const point& b, const line& l){
    double x = distanceSP(l, a);
    double y = distanceSP(l, b);
    return x < eps && y < eps;
}
/**************more complicated APPLICATION...***********/
bool comp(const point& a, const point& b){
    return abs(a-b) < eps;
}
vector<point> crosspointSO(const line& sight, const vector<line>& obstacle){
    vector<point> res;
    for(int i=0; i<obstacle.size(); i++) {
        if(intersectSS(sight, obstacle[i]))
            res.push_back(crosspoint(sight, obstacle[i]));
    }
    sort(res.begin(), res.end());
    res.erase(unique(res.begin(), res.end(), comp), res.end());
    if(res.size() <= 1) return vector<point>();
    for(int i=0; i<obstacle.size(); i++) {
        if(isSameLine(res[0], res[1], obstacle[i])) return vector<point>();
    }
    return res;
}
line calcSight(const line& sight, const vector<Obstacle>& obstacles){
    vector<pair<double,point> > res;
    for(int i=0; i<obstacles.size(); i++) {
        vector<point> p = crosspointSO(sight, obstacles[i]);
        for(int j=0; j<p.size(); j++) {
            res.push_back(make_pair(abs(p[j]-sight[0]), p[j]));
        }
    }
    sort(res.begin(), res.end());
    if(res.size() == 0) return sight;
    return line(sight[0], res[0].second);
}
bool reachable(const line& sight, const vector<Obstacle>& obstacles){
    for(int i=0; i<obstacles.size(); i++) {
        vector<point> p = crosspointSO(sight, obstacles[i]);
        if(p.size() != 0) return false;
    }
    return true;
}
double TriAngleCircleInsection(Circle C, point A, point B){
    Vec OX = A-C.center, OY = B-C.center;
    Vec YX = A-B, YO = C.center-B;
    Vec XY = B-A, XO = C.center-A;
    double da = abs(OX), db = abs(OY),dc = abs(XY), r = C.radius;
    if(dlcmp(cross(OX,OY)) == 0) return 0;
    if(dlcmp(da-r) < 0 && dlcmp(db-r) < 0) return cross(OX,OY)*0.5;
    else if(db < r && da >= r) {
        double x = (dot(YX,YO) + sqrt(r*r*dc*dc-cross(YX,YO)*cross(YX,YO)))/dc;
        double TS = cross(OX,OY)*0.5;
        return asin(TS*(1-x/dc)*2/r/da)*r*r*0.5+TS*x/dc;
    }
    else if(db >= r && da < r) {
        double y = (dot(XY,XO)+sqrt(r*r*dc*dc-cross(XY,XO)*cross(XY,XO)))/dc;
        double TS = cross(OX,OY)*0.5;
        return asin(TS*(1-y/dc)*2/r/db)*r*r*0.5+TS*y/dc;
    }
    else if(fabs(cross(OX,OY)) >= r*dc || dot(XY,XO) <= 0 || dot(YX,YO) <= 0) {
        if(dot(OX,OY) < 0) {
            if(cross(OX,OY) < 0) return (-acos(-1.0)-asin(cross(OX,OY)/da/db))*r*r*0.5;
            else                 return ( acos(-1.0)-asin(cross(OX,OY)/da/db))*r*r*0.5;
        }
        else                     return asin(cross(OX,OY)/da/db)*r*r*0.5;
    }
    else {
        double x = (dot(YX,YO)+sqrt(r*r*dc*dc-cross(YX,YO)*cross(YX,YO)))/dc;
        double y = (dot(XY,XO)+sqrt(r*r*dc*dc-cross(XY,XO)*cross(XY,XO)))/dc;
        double TS = cross(OX,OY)*0.5;
        return (asin(TS*(1-x/dc)*2/r/da)+asin(TS*(1-y/dc)*2/r/db))*r*r*0.5 + TS*((x+y)/dc-1);
    }
}
/********************POLYGONS**********************/
polygon_convex convex_hull(polygon a){
    polygon_convex res(2*a.size()+5);
    sort(a.begin(),a.end(),cmp_less);
    int m=0, n=a.size();
    for(int i=0; i<n; i++){
        while(m>1&&dlcmp(cross(res.p[m-1]-res.p[m-2],a[i]-res.p[m-2]))<=0) --m;
        res.p[m++]=a[i];
    }
    int k=m;
    for(int i=n-2;i>=0;--i){
        while(m>k&&dlcmp(cross(res.p[m-1]-res.p[m-2],a[i]-res.p[m-2]))<=0) --m;
        res.p[m++]=a[i];
    }
    res.p.resize(m);
    if(a.size()>1) res.p.resize(m-1);
    return res;
}
bool contain(const polygon_convex& a, const point& b){
    int n=a.p.size();
    #define next(i) ((i+1)%n)
    int sign=0;
    for(int i=0; i<n; i++){
        int x=dlcmp(cross(a.p[i]-b,a.p[next(i)]-b));
        if(x){
            if(sign){
                if(sign!=x) return false;
            }
            else sign=x;
        }
    }
    return true;
}
/******************END********************/
int main(){
    return 0;
}





圆交多边形:

#include <bits/stdc++.h>
using namespace std;
typedef complex<double> point;
const double pi=acos(-1.0);
const double eps=1.0e-8;
int dcmp(double k){
	return k<-eps?-1:k>eps;
}
const int MAXN=105;
double cross(const point& a, const point& b){ return (conj(a)*b).imag(); }
double dot  (const point& a, const point& b){ return (a*conj(b)).real(); }
point crosspt(const point& a, const point& b, const point& p, const point& q){
	double a1=cross(b-a,p-a);
	double a2=cross(b-a,q-a);
	return (p*a2-q*a1)/(a2-a1);
}
double mysqrt(double x){
	return sqrt(max(0.0,x));
}
double sqr(double x){
	return x*x;
}
double sector_area(const point& a, const point& b, const double& r){
	double theta=arg(a)-arg(b);
	while(theta<=0) theta+=2*pi;
	while(theta>2*pi) theta-=2*pi;
	theta=min(theta,2*pi-theta);
	return r*r*theta/2;
}
void circle_cross_line(point a, point b, point o, double r, point p[], int& num){
	double x0=o.real(), y0=o.imag();
	double x1=a.real(), y1=a.imag();
	double x2=b.real(), y2=b.imag();
	double dx=x2-x1, dy=y2-y1;
	double A=dx*dx+dy*dy;
	double B=2*dx*(x1-x0)+2*dy*(y1-y0);
	double C=sqr(x1-x0)+sqr(y1-y0)-sqr(r);
	double delta=B*B-4*A*C;
	num=0;
	if(dcmp(delta)>=0){
		double t1=(-B-mysqrt(delta))/(2*A);
		double t2=(-B+mysqrt(delta))/(2*A);
		if(dcmp(t1-1)<=0&&dcmp(t1)>=0){
			p[num++]=point(x1+t1*dx,y1+t1*dy);
		}
		if(dcmp(t2-1)<=0&&dcmp(t2)>=0){
			p[num++]=point(x1+t2*dx,y1+t2*dy);
		}
	}
}
double calc(const point& a, const point& b, const double& r){
	point p[2];
	int num=0;
	int ina=dcmp(abs(a)-r)<0;
	int inb=dcmp(abs(b)-r)<0;
	if(ina){
		if(inb){
			return abs(cross(a,b))/2.0;
		}
		else{
			circle_cross_line(a,b,point(0,0),r,p,num);
			return sector_area(b,p[0],r)+abs(cross(a,p[0]))/2.0;
		}
	}
	else{
		if(inb){
			circle_cross_line(a,b,point(0,0),r,p,num);
			return sector_area(p[0],a,r)+abs(cross(p[0],b))/2.0;
		}
		else{
			circle_cross_line(a,b,point(0,0),r,p,num);
			if(num==2){
				return sector_area(a,p[0],r)+sector_area(p[1],b,r)+abs(cross(p[0],p[1]))/2.0;
			}
			else{
				sector_area(a,b,r);
			}
		}
	}
}
point p[MAXN];
int n;
double area(const double& r){
	double ret=0;
	for(int i=0; i<n; i++){
		int sgn=dcmp(cross(p[i],p[i+1]));
		if(sgn){
			ret+=sgn*calc(p[i],p[i+1],r);
		}
	}
	return ret;
}
int main(){
	double x0, y0, v0, theta, t, g, R;
	while(~scanf("%lf%lf%lf%lf%lf%lf%lf",&x0,&y0,&v0,&theta,&t,&g,&R)){
		if(dcmp(x0)==0&&dcmp(y0)==0&&dcmp(R)==0&&dcmp(g)==0&&dcmp(t)==0&&dcmp(theta)==0&&dcmp(v0)==0) break;
		theta=theta/180.0*pi;
		double dx=v0*cos(theta), dy=v0*sin(theta);
		x0+=dx*t;
		y0+=dy*t-0.5*g*t*t;
		cin>>n;
		double x, y;
		for(int i=0; i<n; i++){
			scanf("%lf%lf",&x,&y);
			x-=x0; y-=y0;
            p[i]=point(x,y);
		}
		p[n]=p[0];
		printf("%.2f\n",abs(area(R)));
	}
	return 0;
}





评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值