Problem: Find the area of intersection between a circle $c$ and a polygon $g$. The center coordinate of the circle is ($0, 0$).

The polygon $g$ is represented by a sequence of points $p_1$, $p_2$,..., $p_n$ where line segments connecting $p_i$ and $p_{i+1}$ ($1 \leq i \leq n−1$) are sides of the polygon. The line segment connecting $p_n$ and $p_1$ is also a side of the polygon.

Note that the polygon is not necessarily convex.

Constraints

* $3 \leq n \leq 100$
* $1 \leq r \leq 100$
* $-100 \leq x_i, y_i \leq 100$

Input

The input is given in the following format.

$n$ $r$
$x_1$ $y_1$
$x_2$ $y_2$
:
$x_n$ $y_n$


In the first line, an integer n representing the number of points in the polygon is given. The coordinate of a point $p_i$ is given by two integers $x_i$ and $y_i$. The coordinates of the points are given in the order of counter-clockwise visit of them. All input values are given in integers.

Output

Print the area of intersection in a line. The output values should be in a decimal fraction with an error less than 0.00001.

Examples

Input

3 5
1 1
4 1
5 5


Output

4.639858417607


Input

4 5
0 0
-3 -6
1 -3
5 -4


Output

11.787686807576

Solution (in CPP):

Solution:
#include<iostream>
#include<cstdio>
#include<cmath>
#include<cstring>
#include<vector>
#include<set>
#include<queue>
#include<cstdlib>
#include<algorithm>

using namespace std;
#define pb push_back
#define mp make_pair
#define SZ(x) (int)(x).size()

typedef long long LL;
typedef pair<int, int> pii;
typedef pair<double, double> pdd;
typedef pair<LL, LL> pll;

const double eps = 1e-10;
const double PI = acos(-1.0);
//const double PI = 3.14159265358979323846264338327950288419716939937510
int dcmp(const double &x) { if (fabs(x) < eps) return 0; return x < 0 ? -1 : 1; }
const int MOD = 1e9 + 7;
const int INF = 2e9;
const double INF_d = 1e64;

template<class T> 
T qmod(T a, T b){ T ret = 1; while (b){ if (b & 1) ret *= a; b >>= 1; a *= a; } return ret; }
template<class T>
T gcd(T a, T b){ return !b ? a : gcd(b, a % b); }
template<class T>
T ex_gcd(T a, T b, T &x, T &y){
    if (!b){ x = 1, y = 0; return a; }
    T t, ret;
    ret = ex_gcd(b, a % b, x, y);
    t = x, x = y, y = t - a / b * y;
    return ret;
}
template<class T>
T inv(T t, T p){ return t == 1 ? 1 : (p - p / t) * inv(p % t, p) % p; }
// head

typedef double db;

struct point{
    db x, y;
    point(){}
    point(db a, db b): x(a), y(b){}
    point operator + (const point &p){ return point(x + p.x, y + p.y); }
    point operator - (const point &p){ return point(x - p.x, y - p.y); }
    point operator * (const db &k){ return point(x * k, y * k); }
    point operator / (const db &k){ return point(x / k, y / k); }
    db operator ^ (const point &p){ return x * p.y - y * p.x; } // cross product
    db operator * (const point &p){ return x * p.x + y * p.y; } // dot product
    bool operator == (const point &p){ return !dcmp(x - p.x) && !dcmp(y - p.y); }
    bool operator != (const point &p){ return dcmp(x - p.x) || dcmp(y - p.y); }
    bool operator < (const point &p){ return x == p.x ? y < p.y : x < p.x; }
    void read(){
        scanf("%lf%lf", &x, &y);
    }
    void print(){
        printf("%.10f %.10f\n", x, y);
    }
    point rotate(db &ang){ return point(x * cos(ang) - y * sin(ang), y * cos(ang) + x * sin(ang)); }
    point norm() { return point(-y, x); }
    db len() { return sqrt(x * x + y * y); }
    db len2(){ return x * x + y * y; }
    point unit(){ db k = (*this).len(); return point(x / k, y / k); }
};

point err = point(INF_d, INF_d);

db rad(point A, point B){ return atan2((A ^ B), (A * B)); }// the ang(in rad) of two vectors A and B

point proj_SP(point A, point B, point P){ // projection point of P in segment AB
    point AB = B - A, AP = P - A;
    return A + AB * (AP * AB / AB.len2());
}

point refl_SP(point A, point B, point P){ // reflection point of P int segment AB
    return proj_SP(A, B, P) * 2 - P;
}

struct line{
    point p[2], u;
    line(){}
    line(point s, point t){
        p[0] = s, p[1] = t;
        u = p[1] - p[0];
    }
    point get_point(double t){
        return point(p[0] + u * t);
    }
};

bool parallel(line L1, line L2){ // check whether two lines L1 and L2 are parallel
    return !dcmp(L1.u ^ L2.u);
}

bool orthogonal(line L1, line L2){ // check whether two lines L1 and L2 are orthogonal
    return !dcmp(L1.u * L2.u);
}

pair<bool, point> inter_LL(line L1, line L2){ // calculate the intersection of two lines L1 and L2
    if (parallel(L1, L2)) return mp(false, err);
    double t = ((L2.p[0] - L1.p[0]) ^ L2.u) / (L1.u ^ L2.u);
    return mp(true, L1.get_point(t));
}

bool one_dimention_check(db a, db b, db c, db d){ // check whether two one-dimention segments has intersection
    if (dcmp(a - b) > 0) swap(a, b);
    if (dcmp(c - d) > 0) swap(c, d);
    return dcmp(b - c) >= 0 && dcmp(a - d) <= 0; 
}

bool check_SS(point A, point B, point C, point D){ // check whether two segments AB and CD has intersection
    bool flag = one_dimention_check(A.x, B.x, C.x, D.x) && one_dimention_check(A.y, B.y, C.y, D.y);
    double a = dcmp((A - C) ^ (D - C)), b = dcmp((B - C) ^ (D - C)),
            c = dcmp((C - A) ^ (B - A)), d = dcmp((D - A) ^ (B - A));
    return flag && dcmp(a * b) <= 0 && dcmp(c * d) <= 0;
}
// calculate the intersection point of two segments, first check and then call inter_LL to get the cross point

bool inmid(db a, db b, db c){  // check whether a <= c <= b
    return dcmp(dcmp(b - c) * dcmp(a - c)) <= 0; 
}

bool inmid(point A, point B, point P){ // check whether point P is in the rectangle of segment AB
    return inmid(A.x, B.x, P.x) && inmid(A.y, B.y, P.y);
}

db dis_SP(point A, point B, point P){ // the minimum distance from a point P to segment AB
    point p_proj = proj_SP(A, B, P);
    if (inmid(A, B, p_proj)) return (P - p_proj).len();
    else return min((P - A).len(), (P - B).len());
}

db dis_SS(point A, point B, point C, point D){ // the minimum distance from two segments AB and CD
    db res = INF_d;
    if (check_SS(A, B, C, D)) res = 0;
    res = min(res, dis_SP(A, B, C));
    res = min(res, dis_SP(A, B, D));
    res = min(res, dis_SP(C, D, A));
    res = min(res, dis_SP(C, D, B));
    return res;
}

db poly_Area(point *poly, int n){ // calc the area of a simple polygon poly, n is the number of points.
    //the points should be in clockwise or counter_clockwise.
    db res = 0.0;
    for (int i = 1; i < n - 1; i++){
        res += (poly[i] - poly[0]) ^ (poly[i + 1] - poly[0]);
    }
    return fabs(res) / 2;
}

bool is_Convex(point *poly, int n){ // check whether a polygon is a convex polygon. Make sure points are in counter-clockwised order.
    poly[n] = poly[0];
    if (n < 3) return 1;
    for (int i = 1; i < n; i++){
        int tmp_dir = dcmp((poly[i] - poly[i - 1]) ^ (poly[(i + 1) % n] - poly[i]));
        if (tmp_dir < 0) return 0;
    }
    return 1;
}

bool is_PonS(point A, point B, point P){ // check whether a point P is on a segment AB
    return inmid(A, B, P) && !dcmp((P - A) ^ (B - A));
}

int is_PinPoly(point *poly, int n, point P){ // check the position relation of a point and polygon.
    // 0: outside, 1: on a segment, 2: inside
    int pd = 0;
    for (int i = 0; i < n; i++){
        point u = poly[i], v = poly[(i + 1) % n];
        if (is_PonS(u, v, P)) return 1;
        if (dcmp(u.y - v.y) > 0) swap(u, v);
        if (dcmp(u.y - P.y) >= 0 || dcmp(v.y - P.y) < 0) continue;
        if (dcmp((u - v) ^ (P - v)) < 0) pd ^= 1;
    }
    return pd << 1;
}

bool cmpxy(const point &a, const point &b){ // compare x first, then y. from small -> large
    return a.x == b.x ? a.y < b.y : a.x < b.x;
}

bool cmpyx(const point &a, const point &b){ // compare x first, then y. from small -> large
    return a.y == b.y ? a.x < b.x : a.y < b.y;
}

int get_convexhull(point *p, int n, point *ch){ // get the convexhull of given polygon p
    // pay attention to when m = 0, 1, 2 as the return value
    sort(p, p + n, cmpyx); // may be changed by the requirement of problem
    int m = 0;
    for (int i = 0; i < n; i++){
        while (m > 1 && dcmp((ch[m - 1] - ch[m - 2]) ^ (p[i] - ch[m - 2])) < 0) m--;
        ch[m++] = p[i];
    }
    int k = m;
    for (int i = n - 2; i >= 0; i--){
        while (m > k && dcmp((ch[m - 1] - ch[m - 2]) ^ (p[i] - ch[m - 2])) < 0) m--;
        ch[m++] = p[i];
    }
    if (n > 1) m--;
    return m;
}

db convexhull_diameter(point *ch, int n){ // calculate the diameter of a give convexhull
    if (n < 2) return 0;
    int is, js; is = js = 0;
    for (int i = 0; i < n; i++){
        is = ch[i] < ch[is] ? i : is;
        js = ch[js] < ch[i] ? i : js;
    }
    int i = is, j = js;
    db ret = (ch[is] - ch[js]).len();
    do{
        if (dcmp((ch[(i + 1) % n] - ch[i]) ^ (ch[(j + 1) % n] - ch[j])) > 0){
            j = (j + 1) % n;
        }
        else{
            i = (i + 1) % n;
        }
        ret = max(ret, (ch[i] - ch[j]).len());
    } while (i != is || j != js);
    return ret;
}

int convexhull_cut(point *ch, int n, point P, point Q, point *new_ch){
    int m = 0;
    for (int i = 0; i < n; i++){
        int d1 = dcmp((Q - P) ^ (ch[i] - P)), d2 = dcmp((Q - P) ^ (ch[(i + 1) % n] - P));
        if (d1 >= 0) new_ch[m++] = ch[i];
        if (d1 * d2 < 0) new_ch[m++] = inter_LL(line(P, Q), line(ch[i], ch[(i + 1) % n])).second;
    }
    return m;
}

db cloest_pair(point *p, int l, int r){ // find the distance of cloest pair of points in p
//please make sure that p is sorted by x coordinate | sort(p, p + n, cmpxy);
//this is nlognlogn, if need nlogn, replace the sort with a merge sort.
    if (r - l <= 5){ db cloest = INF_d;
        for (int i = l; i < r; i++) for (int j = i + 1; j < r; j++) cloest = min(cloest, (p[i] - p[j]).len());
        return cloest;
    }
    else{
        int mid = (l + r) >> 1; db cloest = min(cloest_pair(p, l, mid), cloest_pair(p, mid, r));
        vector<point> v; v.clear();
        for (int i = l; i < r; i++) if (fabs(p[i].x - p[mid].x) <= cloest) v.pb(p[i]); sort(v.begin(), v.end(), cmpyx);
        for (int i = 0; i < SZ(v); i++) for (int j = i + 1; j < SZ(v) && dcmp(v[j].y - v[i].y - cloest) < 0; j++)
                cloest = min(cloest, (v[i] - v[j]).len());
        return cloest;
    }
}

struct circle{
    point o; double r;
    circle(){}
    circle(point a, double b): o(a), r(b){}
    void read(){
        o.read(), scanf("%lf", &r);
    }
    int inside(point p){ // 1: inside, 2: on the circle, 0: outside
        int sgn = dcmp((o - p).len() - r); if (sgn < 0) return 1; if (!sgn) return 2; return 0;
    }
};

int check_CC(circle C1, circle C2){ // return the number of common tangents of two circles C1 and C2
//4: xiang_li, 3: wai_qie, 2: xiang_jiao, 1: nei_qie, 0: nei_han
    db d =(C1.o - C2.o).len();
    if (dcmp(d - C1.r - C2.r) > 0) return 4;
    if (!dcmp(d - C1.r - C2.r)) return 3;
    if (!dcmp(d - fabs(C1.r - C2.r))) return 1;
    if (dcmp(d - fabs(C1.r - C2.r)) > 0 && dcmp(d - C1.r - C2.r) < 0) return 2;
    return 0;
}

vector<point> inter_CL(circle C, point A, point B){ // return the intersection of a circle and a line
// the result give in the direction from A -> B
    point p_proj = proj_SP(A, B, C.o);
    db d = C.r * C.r - (C.o - p_proj).len2();
    if (dcmp(d) < 0) return {};
    point delta = (B - A).unit() * sqrt(max((db)0, d)); return {p_proj - delta, p_proj + delta};
}

vector<point> inter_CL(circle C, line L){ // the same as the former one, but argument is different
    return inter_CL(C, L.p[0], L.p[1]);
}

vector<point> inter_CC(circle C1, circle C2){ // return the cross points in counter-clockwised order
    int sgn = check_CC(C1, C2);
    if (!sgn || sgn == 4) return {};
    point dir = C2.o - C1.o;
    db d = dir.len(); db cosA = (d * d + C1.r * C1.r - C2.r * C2.r) / (2 * C1.r * d);
    db sinA = sqrt(1 - cosA * cosA); db h = C1.r * sinA; db base = C1.r * cosA;
    point k = dir.unit(), b = k.norm() * h;
    return {C1.o + k * base - b, C1.o + k * base + b};
}

vector<point> tange_CP(circle C, point P){ // get then tangent point of a point P for a circle C, return in counter-clockwise order
    point OP = P - C.o;
    if (dcmp(OP.len() - C.r) < 0) return {};
    db cosA = C.r / OP.len(); db base = C.r * cosA; db h = sqrt(max((db)0, C.r * C.r - base * base));
    point dir = OP.unit().norm();
    return {C.o + OP.unit() * base - dir * h, C.o + OP.unit() * base + dir * h}; 
}

vector<line> tange_inCC(circle C1, circle C2){ // find the inner tangent lines of two circles
    int sgn = check_CC(C1, C2);
    if (sgn <= 2) return {}; 
    if (sgn == 3){
        point itp = inter_CC(C1, C2)[0];
        return { line(itp, itp) };
    }
    else if (sgn == 4){
        point itp = (C2.o * C1.r + C1.o * C2.r) / (C1.r + C2.r);
        vector<point> p1 = tange_CP(C1, itp);
        vector<point> p2 = tange_CP(C2, itp);
        vector<line> res; res.clear();
        for (int i = 0; i < SZ(p1); i++) res.pb(line(p1[i], p2[i]));
    }
}

vector<line> tange_outCC(circle C1, circle C2){ // find the outter tangent lines of two circles
    int sgn = check_CC(C1, C2);
    if (!sgn) return {};
    if (sgn == 1){ point itp = inter_CC(C1, C2)[0]; return { line(itp, itp) }; }
    if (!dcmp(C1.r - C2.r)){
        point dir = (C2.o - C1.o).unit().norm();
        return { line(C1.o - dir * C1.r, C2.o - dir * C2.r), line(C1.o + dir * C1.r, C2.o + dir * C2.r)};
    }
    else{
        point itp = (C2.o * C1.r - C1.o * C2.r) / (C1.r - C2.r);
        vector<point> p1 = tange_CP(C1, itp); vector<point> p2 = tange_CP(C2, itp);
        vector<line> res; res.clear();
        for (int i = 0; i < SZ(p1); i++){ res.pb(line(p1[i], p2[i])); }
    }
}

vector<line> tange_CC(circle C1, circle C2){ // find the tanget lines of two circles
    bool flag = 0; if (dcmp(C1.r - C2.r) < 0) { swap(C1, C2); flag = 1; }
    int sgn = check_CC(C1, C2);
    if (!sgn) return {};       
    vector<line> sol2 = tange_outCC(C1, C2);
    vector<line> sol1 = tange_inCC(C1, C2);
    for (int i = 0; i < SZ(sol2); i++){
        sol1.push_back(sol2[i]);
    }
    if (flag) for (line &k : sol1) swap(k.p[0], k.p[1]);
    return sol1;
}

vector<point> unique_point(vector<point> &v){ // return the unique_point of the points v
    vector<point> ret; ret.clear();
    if (!SZ(v)) return {};
    ret.pb(v[0]);
    for (int i = 1; i < SZ(v); i++){
        if (v[i] == v[i - 1]) continue;
        else ret.pb(v[i]);
    }
    return ret;
}

db CT_area(circle C, point A, point B){ // calculate the intersection area of a triangle and a circle
// one point of the triangle is the center of the circle
    point t = C.o; C.o = C.o - t, A = A - t, B = B - t;
    int sgn1 = C.inside(A), sgn2 = C.inside(B);
    vector<point> itp = inter_CL(C, line(A, B));
    if (sgn1 >= 1){
        if (sgn2 >= 1) return ((A - C.o) ^ (B - C.o)) / 2;
        else return ((A - C.o) ^ (itp[1] - C.o)) / 2 + rad(itp[1], B) * C.r * C.r / 2;
    }
    else{
        if (sgn2 >= 1) return ((itp[0] - C.o) ^ (B - C.o)) / 2 + rad(A, itp[0]) * C.r * C.r / 2;
        else{
            int sgn = dcmp(C.r - dis_SP(A, B, C.o));
            if (sgn <= 0) return rad(A, B) * C.r * C.r / 2; 
            else return (rad(A, itp[0]) + rad(itp[1], B)) * C.r * C.r / 2 + ((itp[0] - C.o) ^ (itp[1] - C.o)) / 2;
        }
    }
}

db PolyC_area(point *p, int n, circle C){
    db res = 0;
    for (int i = 0; i < n; i++){
        point u = p[i], v = p[(i + 1) % n];
        res += CT_area(C, u, v);
    }
    return fabs(res);
}

const int N = 1e5 + 7;
int n;
point p[N], x;
db R;

int main(){
    scanf("%d%lf", &n, &R);
    for (int i = 0; i < n; i++){
        p[i].read();
    }
    circle C = circle(point(0, 0), R);
    printf("%.12f\n", PolyC_area(p, n, C));
    return 0;
}