• 多边形与圆相交面积计算


    POJ 3675 裸的题了。直接上模板就行了。注意的是,求出的是有向面积,有可能是负数。

    #include <iostream>
    #include <cstdio>
    #include <cmath>
    #include <vector>
    #include <cstring>
    #include <algorithm>
    #include <string>
    #include <set>
    #include <functional>
    #include <numeric>
    #include <sstream>
    #include <stack>
    
    #define CL(arr, val) memset(arr, val, sizeof(arr))
    #define REP(i, n)for((i) = 0; (i) < (n); ++(i))
    #define FOR(i, l, h)for((i) = (l); (i) <= (h); ++(i))
    #define FORD(i, h, l)for((i) = (h); (i) >= (l); --(i))
    #define L(x)(x) << 1
    #define R(x)(x) << 1 | 1
    #define MID(l, r)(l + r) >> 1
    #define Min(x, y)x < y ? x : y
    #define Max(x, y)x < y ? y : x
    #define E(x)(1 << (x))
    #define iabs(x) (x) < 0 ? -(x) : (x)
    #define OUT(x)printf("%I64d
    ", x)
    #define lowbit(x)(x)&(-x)
    #define Read()freopen("data.in", "r", stdin)
    #define Write()freopen("data.out", "w", stdout);
    
    const double eps = 1e-8;
    typedef long long LL;
    const int inf = ~0u>>2;
    
    using namespace std;
    
    inline double max (double a, double b) { if (a > b) return a; else return b; }
    inline double min (double a, double b) { if (a < b) return a; else return b; }
    inline int fi (double a){
        if (a > eps) return 1;
        else if (a >= -eps) return 0;
        else return -1;
    }
    class Vector{
    public:
        double x, y;
        Vector (void) {}
        Vector (double x0, double y0) : x(x0), y(y0) {}
        double operator * (const Vector& a) const { return x * a.y - y * a.x; }
        double operator % (const Vector& a) const { return x * a.x + y * a.y; }
        Vector verti (void) const { return Vector(-y, x); }
        double length (void) const { return sqrt(x * x + y * y); }
        Vector adjust (double len) {
            double ol = len / length();
            return Vector(x * ol, y * ol);
        }
        Vector oppose (void) { return Vector(-x, -y); }
    };
    class point{
    public:
        double x, y;
        point (void) {}
        point (double x0, double y0) : x(x0), y(y0) {}
        Vector operator - (const point& a) const { return Vector(x - a.x, y - a.y); }
        point operator + (const Vector& a) const { return point(x + a.x, y + a.y); }
    };
    class segment{
    public:
        point a, b;
        segment (void) {}
        segment (point a0, point b0) : a(a0), b(b0) {}
        point intersect (const segment& s) const {
            Vector v1 = s.a - a, v2 = s.b - a, v3 = s.b - b, v4 = s.a - b;
            double s1 = v1 * v2, s2 = v3 * v4;
            double se = s1 + s2;
            s1 /= se, s2 /= se;
            return point(a.x * s2 + b.x * s1, a.y * s2 + b.y * s1);
        }
        point pverti (const point& p) const {
            Vector t = (b - a).verti();
            segment uv(p, p + t);
            return intersect(uv);
        }
        bool on_segment (const point& p) const {
            if (fi(min(a.x, b.x) - p.x) <= 0 && fi(p.x - max(a.x, b.x)) <= 0 &&
                fi(min(a.y, b.y) - p.y) <= 0 && fi(p.y - max(a.y, b.y)) <= 0) return true;
            else return false;
        }
    };
    
    double radius;
    point polygon[110];
    
    double kuras_area (point a, point b, double cx, double cy){
        point ori(cx, cy);
        int sgn = fi((b - ori) * (a - ori));
        double da = (a - ori).length(), db = (b - ori).length();
        int ra = fi(da - radius), rb = fi(db - radius);
        double angle = acos(((b - ori) % (a - ori)) / (da * db));
        segment t(a, b); point h, u; Vector seg;
        double ans, dlt, mov, tangle;
    
        if (fi(da) == 0 || fi(db) == 0) return 0;
        else if (sgn == 0) return 0;
        else if (ra <= 0 && rb <= 0) return fabs((b - ori) * (a - ori)) / 2 * sgn;
        else if (ra >= 0 && rb >= 0){
            h = t.pverti(ori);
            dlt = (h - ori).length();
            if (!t.on_segment(h) || fi(dlt - radius) >= 0)
                return radius * radius * (angle / 2) * sgn;
            else{
                ans = radius * radius * (angle / 2);
                tangle = acos(dlt / radius);
                ans -= radius * radius * tangle;
                ans += radius * sin(tangle) * dlt;
                return ans * sgn;
            }
        }
        else{
            h = t.pverti(ori);
            dlt = (h - ori).length();
            seg = b - a;
            mov = sqrt(radius * radius - dlt * dlt);
            seg = seg.adjust(mov);
            if (t.on_segment(h + seg)) u = h + seg;
            else u = h + seg.oppose();
            if (ra == 1) swap(a, b);
            ans = fabs((a - ori) * (u - ori)) / 2;
            tangle = acos(((u - ori) % (b - ori)) / ((u - ori).length() * (b - ori).length()));
            ans += radius * radius * (tangle / 2);
            return ans * sgn;
        }
    }
    
    const double pi = acos(-1.0);
    
    int main (){
    	int n;
    	while(scanf("%lf",&radius)!=EOF){
    		scanf("%d",&n);
    		for(int i=0;i<n;i++){
    			scanf("%lf%lf",&polygon[i].x,&polygon[i].y);
    		}
    		double ans=0;
    		for(int i=0;i<n;i++){
    			ans+=kuras_area(polygon[i],polygon[(i+1)%n],0,0);
    		}
    		printf("%.2lf
    ",fabs(ans));
    	}
    }
    

      

    HDU 3982。很明显,求出切割的半平面交之后就可以得出一个多边形。然后求多边形与圆面积相交部分即可。

    #include <iostream>
    #include <cstdio>
    #include <cstring>
    #include <algorithm>
    #include <cmath>
    using namespace std;
    
    const double eps = 1e-6;
    typedef long long LL;
    const int inf = ~0u>>2;
    const int MAXN=4222; 
    
    inline double max (double a, double b) { if (a > b) return a; else return b; }
    inline double min (double a, double b) { if (a < b) return a; else return b; }
    inline int fi (double a){
        if (a > eps) return 1;
        else if (a >= -eps) return 0;
        else return -1;
    }
    
    struct Point {
    	double x,y;
    	Point(double xx,double yy){x=xx; y=yy;}
    	Point(){}
    }convex[MAXN],cherry;
    
    class Vector{
    public:
        double x, y;
        Vector (void) {}
        Vector (double x0, double y0) : x(x0), y(y0) {}
        double operator * (const Vector& a) const { return x * a.y - y * a.x; }		//向量叉乘 
        double operator % (const Vector& a) const { return x * a.x + y * a.y; }		//向量点乘
        Vector verti (void) const { return Vector(-y, x); }
        double length (void) const { return sqrt(x * x + y * y); }
        Vector adjust (double len) {
            double ol = len / length();
            return Vector(x * ol, y * ol);
        }
        Vector oppose (void) { return Vector(-x, -y); }
    };
    class point{
    public:
        double x, y;
        point (void) {}
        point (double x0, double y0) : x(x0), y(y0) {}
        Vector operator - (const point& a) const { return Vector(x - a.x, y - a.y); }
        point operator + (const Vector& a) const { return point(x + a.x, y + a.y); }
        void _copy(Point t){
        	x=t.x; y=t.y;
        }
    };
    class segment{
    public:
        point a, b;
        segment (void) {}
        segment (point a0, point b0) : a(a0), b(b0) {}
        point intersect (const segment& s) const {
            Vector v1 = s.a - a, v2 = s.b - a, v3 = s.b - b, v4 = s.a - b;
            double s1 = v1 * v2, s2 = v3 * v4;
            double se = s1 + s2;
            s1 /= se, s2 /= se;
            return point(a.x * s2 + b.x * s1, a.y * s2 + b.y * s1);
        }
        point pverti (const point& p) const {
            Vector t = (b - a).verti();
            segment uv(p, p + t);
            return intersect(uv);
        }
        bool on_segment (const point& p) const {
            if (fi(min(a.x, b.x) - p.x) <= 0 && fi(p.x - max(a.x, b.x)) <= 0 &&
                fi(min(a.y, b.y) - p.y) <= 0 && fi(p.y - max(a.y, b.y)) <= 0) return true;
            else return false;
        }
    };
    
    double radius;
    point polygon[MAXN];//多边形 
    
    double kuras_area (point a, point b, double cx, double cy){
        point ori(cx, cy);
        int sgn = fi((b - ori) * (a - ori));
        double da = (a - ori).length(), db = (b - ori).length();
        int ra = fi(da - radius), rb = fi(db - radius);
        double angle = acos(((b - ori) % (a - ori)) / (da * db));
        segment t(a, b); point h, u; Vector seg;
        double ans, dlt, mov, tangle;
    
        if (fi(da) == 0 || fi(db) == 0) return 0;
        else if (sgn == 0) return 0;
        else if (ra <= 0 && rb <= 0) return fabs((b - ori) * (a - ori)) / 2 * sgn;
        else if (ra >= 0 && rb >= 0){
            h = t.pverti(ori);
            dlt = (h - ori).length();
            if (!t.on_segment(h) || fi(dlt - radius) >= 0)
                return radius * radius * (angle / 2) * sgn;
            else{
                ans = radius * radius * (angle / 2);
                tangle = acos(dlt / radius);
                ans -= radius * radius * tangle;
                ans += radius * sin(tangle) * dlt;
                return ans * sgn;
            }
        }
        else{
            h = t.pverti(ori);
            dlt = (h - ori).length();
            seg = b - a;
            mov = sqrt(radius * radius - dlt * dlt);
            seg = seg.adjust(mov);
            if (t.on_segment(h + seg)) u = h + seg;
            else u = h + seg.oppose();
            if (ra == 1) swap(a, b);
            ans = fabs((a - ori) * (u - ori)) / 2;
            tangle = acos(((u - ori) % (b - ori)) / ((u - ori).length() * (b - ori).length()));
            ans += radius * radius * (tangle / 2);
            return ans * sgn;
        }
    }
    
    const double pi = acos(-1.0);
    
    struct 	Line {
    	Point a,b;
    	double ang;
    }line[MAXN/2],st[MAXN/2];
    int n,cnt;
    
    double operator *(const Point x,const Point y){
    	return x.x*y.y-x.y*y.x;
    }
    Point operator -(Point x,const Point y){
    	x.x-=y.x; x.y-=y.y;
    	return x;
    }
    Point operator * (const Line &x,const Line &y){	//交点 
    	double a1=(y.b-x.a)*(y.a-x.a),a2=(y.a-x.b)*(y.b-x.b);
    	Point r;
    	r.x=(x.a.x*a2+x.b.x*a1)/(a2+a1);
    	r.y=(x.a.y*a2+x.b.y*a1)/(a2+a1);
    	return r;
    }
    
    bool operator ==(const Point &a,const Point &b){
    	return fabs(a.x-b.x)<eps&&fabs(a.y-b.y)<eps;
    }
    
    bool operator <(const Line &x,const Line &y){
    	if(fabs(x.ang-y.ang)<eps)
    	return (y.b-x.a)*(x.b-y.a)>eps;
    	return x.ang <y.ang;
    }
    
    bool JudgeOut(const Line &x,const Point &p){
    	return (p-x.a)*(x.b-x.a)>eps;
    }
    
    bool Parellel(const Line &x,const Line &y){
    	return fabs((x.b-x.a)*(y.b-y.a))<eps;
    }
    
    void HplaneI(){
    	sort(line,line+n);
    	int top=1,bot=0,tmp=1;
    	for(int i=1;i<n;i++){
    		if(line[i].ang-line[i-1].ang>eps) line[tmp++]=line[i];
    	}
    	n=tmp;
    	st[0]=line[0],st[1]=line[1];
    	for(int i=2;i<n;i++){
    		if(Parellel(st[top],st[top-1])||Parellel(st[bot],st[bot+1]))return ;
    		while(bot<top&&JudgeOut(line[i],st[top]*st[top-1])) top--;
    		while(bot<top&&JudgeOut(line[i],st[bot]*st[bot+1])) bot++;
    		st[++top]=line[i];
    	}
    	while(bot<top&&JudgeOut(st[bot],st[top]*st[top-1])) top--;
    	while(bot<top&&JudgeOut(st[top],st[bot]*st[bot+1])) bot++;
    	if(top<=bot+1) return ;
    	st[++top]=st[bot]; cnt=0;
    	for(int i=bot;i<top;i++){
    		polygon[cnt++]._copy(st[i]*st[i+1]);
    	}
    }
    
    int main(){
    	int T,icase=0;
    	scanf("%d",&T);
    	while(T--){
    		scanf("%lf%d",&radius,&n);
    		for(int i=0;i<n;i++){
    			scanf("%lf%lf%lf%lf",&line[i].a.x,&line[i].a.y,&line[i].b.x,&line[i].b.y);
    		}
    		scanf("%lf%lf",&cherry.x,&cherry.y);
    		for(int i=0;i<n;i++){
    			if((line[i].a-cherry)*(line[i].b-cherry)<0){
    				swap(line[i].a,line[i].b);
    			}
    			line[i].ang=atan2(line[i].b.y-line[i].a.y,line[i].b.x-line[i].a.x);
    		}
    		line[n].a=Point(-radius,-radius); line[n].b=Point(radius,-radius); line[n++].ang=0;
    		line[n].a=Point(radius,-radius); line[n].b=Point(radius,radius); line[n++].ang=pi/2;
    		line[n].a=Point(radius,radius); line[n].b=Point(-radius,radius); line[n++].ang=pi;
    		line[n].a=Point(-radius,radius); line[n].b=Point(-radius,-radius); line[n++].ang=-pi/2;
    		cnt=0;
    		HplaneI();
    		double ans=0;
    
    		for(int i=0;i<cnt;i++){
    	//		printf("%lf %lf 
    ",polygon[i].x,polygon[i].y);
    			ans+=kuras_area(polygon[i],polygon[(i+1)%cnt],0,0);
    		}
    		printf("Case %d: ",++icase);
    		printf("%.5lf",fabs(fabs(ans)-pi*radius*radius)<eps?100:fabs(ans)/(pi*radius*radius)*100);
    		puts("%");
    	}
    	return 0;
    }
    

      

  • 相关阅读:
    使用Linq to Sqlite 出现异常Object already attached
    CSS 嵌套DIV布局
    《面试笔记》——MySQL终结篇(30问与答)
    PotPlayer播放器下载
    博客圆的RSS怎么不能用呢
    OPC在自控系统的应用
    TAPI的使用
    刷iPAQ为Linux(zz HiPDA)
    再论软工
    Silverlight的大小自适应中存在的一个问题
  • 原文地址:https://www.cnblogs.com/jie-dcai/p/4509854.html
Copyright © 2020-2023  润新知