您的位置:首页 > 其它

HDU4629扫描线求三角形重叠K次面积

2013-08-01 08:28 411 查看

扫描线求三角形的重合k次的面积

先把所有的三角形的点和他们之间产生的交点求出来,按X轴从小到大排序,然后每相邻的两个点会产生一个区间,在该区间中,从Y轴最上面往下扫描,如果碰到进入三角形的边,当前state++,表示重合数+1,如果碰到一个离开三角形的边,当前state--,表示重合数-1,对于每次碰到一个三角形的边,求出该边在这个区间中的上方的距离,要么是梯形,要么是三角形。



#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cmath>
#include <algorithm>
#include <cstring>
#define eps 1.0e-8
#define zero(x) (((x)>0?(x):-(x))<eps)

using namespace std;
struct node
{
    double x,y;
}point[111111];
struct line{
    node a,b;
};

struct Tri{
    node p[3];
}tri[55];

int parallel(line u,line v){
    return zero((u.a.x-u.b.x)*(v.a.y-v.b.y)-(v.a.x-v.b.x)*(u.a.y-u.b.y));
}

double xmul(node a, node b, node c){
    return fabs((a.x-c.x)*(b.y-c.y)-(a.y-c.y)*(b.x-c.x));
}
node intersection(line u,line v){
    node ret = u.a;
    double t = ((u.a.x-v.a.x)*(v.a.y-v.b.y)-(u.a.y-v.a.y)*(v.a.x-v.b.x))
                /((u.a.x-u.b.x)*(v.a.y-v.b.y)-(u.a.y-u.b.y)*(v.a.x-v.b.x));
    ret.x += (u.b.x-u.a.x)*t;
    ret.y += (u.b.y-u.a.y)*t;
    return ret;
}

bool cmp(node a,node b){
    if(zero(a.x-b.x)){
        return a.y < b.y;
    }
    return a.x < b.x;
}
struct SEG{
    node a;
    node b;
    bool in;
}seg[1111];
bool cmpy(SEG a,SEG b){
    if(zero(max(a.a.y,a.b.y) - max(b.a.y,b.b.y)))
    {
        return min(a.a.y,a.b.y) > min(b.a.y,b.b.y);
    }
    return max(a.a.y,a.b.y) > max(b.a.y,b.b.y);
}

double ans[111];
int main()
{
    int T;
    scanf("%d",&T);
    while(T--){
        int n;
        int cnt = 0;
        scanf("%d",&n);
        memset(ans, 0, sizeof(ans));
        for(int i = 0 ; i < n ; i ++){
            for(int j = 0 ; j < 3 ; j ++){
                scanf("%lf%lf", &tri[i].p[j].x, &tri[i].p[j].y);
                point[cnt++] = tri[i].p[j];
            }
        }
        for(int i = 0 ; i < n ; i ++){
            if(xmul(tri[i].p[0],tri[i].p[1],tri[i].p[2]) < eps)
                continue;
            for(int j = 0 ; j < n ; j ++){
                if(xmul(tri[j].p[0],tri[j].p[1],tri[j].p[2]) < eps)
                    continue;
                for(int k = 0 ; k < 3 ; k ++){
                    for(int l = 0 ; l < 3 ; l ++){
                        line u,v;
                        u.a = tri[i].p[k];
                        u.b = tri[i].p[(k+1)%3];
                        v.a = tri[j].p[l];
                        v.b = tri[j].p[(l+1)%3];
                        if(!parallel(u,v)){
                            point[cnt++] = intersection(u,v);
                        }
                    }
                }
            }
        }
        sort(point,point+cnt,cmp);
        for(int i = 0 ; i < cnt - 1; i ++){
            line left;
            left.a = point[i];
            left.b.x = point[i].x;
            //cout<<point[i].x<<" "<<point[i].y<<endl;
            left.b.y = 200;
            line right;
            right.a = point[i+1];
            right.b.x = point[i+1].x;
            //cout<<point[i+1].x<<" "<<point[i+1].y<<endl;
            right.b.y = 200;
            double d = right.a.x - left.a.x;
            if(d < eps)
                continue;
            int now = 0;
            for(int j = 0; j < n; j++){
                SEG tp[3];
                int num = 0;
                if(xmul(tri[j].p[0],tri[j].p[1],tri[j].p[2]) < eps)
                    continue;
                for(int k = 0; k < 3; k++){
                    line temp;
                    temp.a = tri[j].p[k];
                    temp.b = tri[j].p[(k+1)%3];
                    if(min(temp.a.x,temp.b.x) < left.a.x + eps &&
                       max(temp.a.x,temp.b.x) > left.a.x - eps &&
                       min(temp.a.x,temp.b.x) < right.a.x + eps &&
                       max(temp.a.x,temp.b.x) > right.a.x - eps ){
                        tp[num].in = 0;
                        tp[num].a = intersection(temp,left);
                        tp[num].b = intersection(temp,right);
                        //cout<<tp[num].a.x<<" "<<tp[num].a.y<<endl;
                        //cout<<tp[num].b.x<<" "<<tp[num].b.y<<endl;
                        //cout<<left.a.x<<" "<<right.a.x<<endl;
                        //cout<<temp.a.x<<" "<<temp.a.y<<endl;
                        //cout<<temp.b.x<<" "<<temp.b.y<<endl;
                        //cout<<endl;
                        num++;

                    }
                }
                if(num > 0){
                    if(tp[0].a.y + tp[0].b.y > tp[1].a.y + tp[1].b.y){
                        tp[0].in = 1;
                    }
                    else
                        tp[1].in = 1;
                    seg[now++] = tp[0];
                    seg[now++] = tp[1];

                }
            }
            sort(seg, seg + now, cmpy);
            int state = 0;
            for(int j = 0; j < now - 1; j++){
                if(seg[j].in){
                    state++;
                }
                else
                    state--;
                ans[state] += d*(seg[j].a.y-seg[j+1].a.y + seg[j].b.y-seg[j+1].b.y)/2.0;
            }
        }

        for(int i = 1; i <= n; i ++){
            printf("%.10f\n",ans[i]);
        }

    }
    return 0;
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: