圆和多边形面积交模板

Posted 2020-11-10 walfy

hdu5130

//#pragma comment(linker, "/stack:200000000")
//#pragma GCC optimize("Ofast,no-stack-protector")
//#pragma GCC target("sse,sse2,sse3,ssse3,sse4,popcnt,abm,mmx,avx,tune=native")
//#pragma GCC optimize("unroll-loops")
#include<bits/stdc++.h>
#define fi first
#define se second
#define mp make_pair
#define pb push_back
//#define pi acos(-1.0)
#define ll long long
#define vi vector<int>
#define mod 1000000007
#define ld long double
//#define C 0.5772156649
#define ls l,m,rt<<1
#define rs m+1,r,rt<<1|1
#define pil pair<int,ll>
#define pli pair<ll,int>
#define pii pair<int,int>
#define cd complex<double>
#define ull unsigned long long
#define base 1000000000000000000
#define fio ios::sync_with_stdio(false);cin.tie(0)

using namespace std;

const double eps=1e-6,PI = acos( -1.0 ) ;
const int N=500000+10,maxn=20000+10,inf=0x3f3f3f3f,INF=0x3f3f3f3f3f3f3f3f;

inline double sqr( double x ){ return x * x ; }
inline int sgn( double x ){
    if ( fabs(x) < eps ) return 0 ;
    return x > 0? 1 : -1 ;
}

struct Point{
    double x , y ;
    Point(){}
    Point( double _x , double _y ): x(_x) , y(_y) {}
    void input() { scanf( "%lf%lf" ,&x ,&y ); }
    double norm() { return sqrt( sqr(x) + sqr(y) ); }

    friend Point operator + ( const Point &a , const Point &b ) { return Point( a.x + b.x , a.y + b.y ) ; }
    friend Point operator - ( const Point &a , const Point &b ) { return Point( a.x - b.x , a.y - b.y ) ; }
    friend Point operator * ( const Point &a , const double &b ) { return Point( a.x * b , a.y * b ) ; }
    friend Point operator * ( const double &a , const Point &b ) { return Point( b.x * a , b.y * a ) ; }
    friend Point operator / ( const Point &a , const double &b ) { return Point( a.x / b , a.y / b ) ; }
    friend bool operator == ( const Point &a , const Point &b ) { return sgn( a.x - b.x ) == 0 && sgn( a.y - b.y ) == 0 ; }

    bool operator < ( const Point &a )const{
        return ( sgn( x - a.x ) < 0 ) || ( sgn( x - a.x ) == 0 && sgn( y - a.y ) < 0 ) ;
    }
};

double dot( Point a , Point b ) { return a.x * b.x + a.y * b.y ; }
double det( Point a , Point b ) { return a.x * b.y - a.y * b.x ; }
double dist( Point a , Point b ) { return ( a - b ).norm() ; }

int n ;
double k ;
Point A,B ;
Point p[505] ;
Point o ;
double r ;

int CircleInterLine( Point a, Point b, Point o, double r, Point *p )
{
    Point p1 = a - o ;
    Point d = b - a ;
    double A = dot( d, d ) ;
    double B = 2 * dot( d, p1 ) ;
    double C = dot( p1, p1 ) - sqr(r) ;

    double delta = sqr(B) - 4*A*C ;
    if ( sgn(delta) < 0 ) return 0 ;//相离
    if ( sgn(delta) == 0 ) { //相切
        double t = -B / (2*A) ; // 0 <= t <= 1说明交点在线段上
        if ( sgn( t - 1 ) <= 0 && sgn( t ) >= 0 ) {
            p[0] = a + t * d ;
            return 1 ;
        }
    }
    if ( sgn(delta) > 0 ) { //相交
        double t1 = ( -B - sqrt(delta) ) / (2*A) ;
        double t2 = ( -B + sqrt(delta) ) / (2*A) ; //0 <= t1, t2 <= 1说明交点在线段上
        int k = 0 ;
        if ( sgn( t1 - 1 ) <= 0 && sgn( t1 ) >= 0 )
            p[k++] = a + t1 * d ;
        if ( sgn( t2 - 1 ) <= 0 && sgn( t2 ) >= 0 )
            p[k++] = a + t2 * d ;
        return k ;
    }
}
double Triangle_area( Point a, Point b )
{
    return fabs( det( a , b ) ) / 2.0  ;
}
double Sector_area( Point a, Point b )
{
    double ang = atan2( a.y , a.x ) - atan2( b.y, b.x  ) ;
    while ( ang <= 0 ) ang += 2 * PI ;
    while ( ang > 2 * PI ) ang -= 2 * PI ;
    ang = min( ang, 2*PI - ang ) ;
    return sqr(r) * ang/2 ;
}
double calc( Point a , Point b , double r )
{
    Point pi[2] ;
    if ( sgn( a.norm() - r ) < 0 ) {
        if ( sgn( b.norm() - r ) < 0 ) {
            return Triangle_area( a, b ) ;
        }
        else {
            CircleInterLine( a, b, Point(0,0), r, pi) ;
            return Sector_area( b, pi[0] ) + Triangle_area( a, pi[0] ) ;
        }
    }
    else {
        int cnt = CircleInterLine( a, b, Point(0,0), r, pi ) ;
        if ( sgn( b.norm() - r ) < 0 ) {
            return Sector_area( a, pi[0] ) + Triangle_area( b, pi[0] ) ;
        }
        else {
            if ( cnt == 2 )
                return Sector_area( a, pi[0] ) + Sector_area( b, pi[1] ) + Triangle_area( pi[0], pi[1] ) ;
            else
                return Sector_area( a, b ) ;
        }
    }
}
double area_CircleAndPolygon( Point *p , int n , Point o , double r )
{
    double res = 0 ;
    p[n] = p[0] ;
    for ( int i = 0 ; i < n ; i++ ) {
        int tmp = sgn( det( p[i] - o , p[i+1] - o ) ) ;
        if ( tmp )
            res += tmp * calc( p[i] - o , p[i+1] - o , r ) ;
    }
    return fabs( res ) ;
}

void gao()
{
    double a1=1.0-k*k,b1=2.0*(sqr(k)*A.x-B.x),c=sqr(B.x)+sqr(B.y)-sqr(k)*(sqr(A.x)+sqr(A.y));
    double a2=1.0-k*k,b2=2.0*(sqr(k)*A.y-B.y);
    o.x=-b1/2.0/a1,o.y=-b2/2.0/a2;
    r=sqrt(sqr(o.x)+sqr(o.y)+(sqr(k*A.x)-sqr(B.x)+sqr(k*A.y)-sqr(B.y))/a1);
//    printf("%.12f %.12f %.12f\n",o.x,o.y,r);
}
int main()
{
    int cnt=1;
    while(~scanf("%d%lf",&n,&k))
    {
        for(int i=0;i<n;i++)
            p[i].input();
        A.input(),B.input();
        gao();
        printf("Case %d: %.12f\n",cnt++,area_CircleAndPolygon(p,n,o,r));
    }
    return 0;
}
/********************

********************/