N圆面积并

#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cctype>
#include <cmath>
#include <algorithm>
#include <vector>
#include <queue>
#include <stack>
#include <map>

using namespace std;

const int inf = 0x3f;
const int INF = 0x3f3f3f3f;
const int maxn = 1005;
const double eps = 1e-6;   

typedef pair<double, double> PDD;
int n, L, R;
int x[maxn], y[maxn], r[maxn];
int nx[maxn], ny[maxn], nr[maxn];
int xl[maxn], xr[maxn];
int s[maxn];
PDD se[maxn];

inline bool cmp(int a, int b)
{
    if(x[a] - r[a] == x[b] - r[b])
        return x[a] + r[a] < x[b] + r[b];
    return x[a] - r[a] < x[b] - r[b];
}

inline bool cmp0(int a, int b)
{
    return r[a] > r[b];
}

inline double f(double v)
{
    int sz = 0, i, j;
    double ret = 0.0;
    for(int i = L; i < R; i++)
    {
        if(v <= xl[i] || v >= xr[i]) continue;
        j = s[i];
        double d = sqrt(r[j] - (v - x[j]) * (v - x[j]));
        se[sz].first = y[j] - d;
        se[sz].second = y[j] + d;
        ++sz;
    }
    sort(se, se + sz);
    for(i = 0; i < sz; i++)
    {
        double nowl, nowr;
        nowl = se[i].first;
        nowr = se[i].second;
        for(j = i + 1; j < sz; ++j)
        {
            if(se[j].first > nowr) break;
            else nowr = max(nowr, se[j].second);
        }
        ret += nowr - nowl;
        i = j - 1;          
    }
    return ret;
}

inline double rsimp(double l, double m, double r, double sl, 
        double sm, double sr, double tot)
{
    double m1 = (l + m) * 0.5, m2 = (m + r) * 0.5;
    double s0 = f(m1), s2 = f(m2);
    double gl = (sl + sm + s0 * 4.0)*(m-l), 
           gr = (sm + sr + s2 * 4.0)*(r-m);
    if(fabs(gl + gr - tot) < eps) return gl + gr;
    return rsimp( l, m1, m, sl, s0, sm, gl) + rsimp( m, m2,r, sm, s2, sr, gr);         
}
void work()
{
    sort(s, s+n, cmp);
    double lo, hi, ans = 0.0;
    int i, j;
    for(i = 0; i < n; i++)
    {
        xl[i] = x[s[i]] - r[s[i]];
        xr[i] = x[s[i]] + r[s[i]];
        r[s[i]] *= r[s[i]];
    }
    for(i = 0; i < n; i++)
    {
        int ilo, ihi;
        ilo = xl[i];
        ihi = xr[i];
        for(j = i+1; j < n; j++)
        {
            if(xl[j] > ihi) break;
            ihi = max(ihi, xr[j]);
        }
        lo = ilo, hi = ihi;
        L = i; R = j;
        double mid = (lo + hi) * 0.5;
        double sl = f(lo), sm= f(mid), sr = f(hi);
        double tot = sl + sr + sm * 4.0;
        ans += rsimp(lo, mid, hi, sl, sm, sr, tot);
        i = j - 1;
    }
    printf("%.3f\n", ans / 6.0);
}
int main()
{
    int i, j, k; 
    while(~scanf("%d", &n))
    {
        i = 0; j = 0; k = 0;
        for(i = 0; i < n; i++)
        {
            scanf("%d %d %d", &x[i], &y[i], &r[i]);
            s[i] = i;
        }
        sort(s, s+n, cmp0);
        for(i = 0; i < n; i++)
        {
            for(k = 0; k < j; k++)
            {
                if((nx[k]-x[s[i]])*(nx[k]-x[s[i]])+
                   (ny[k]-y[s[i]])*(ny[k]-y[s[i]])
                   <= (nr[k]-r[s[i]])*(nr[k]-r[s[i]]))
                    break;
            }
            if(k == j)
            {
                nx[j] = x[s[i]];
                ny[j] = y[s[i]];
                nr[j] = r[s[i]];
                s[j] = j;
                j++;
            }
        }       
        n = j;
        for(i = 0; i < n; i++)
        {
            x[i] = nx[i];
            y[i] = ny[i];
            r[i] = nr[i];
        }
        work();
    }
    return 0;
}