LA 3218 - Find the Border PSLG 平面直线区域划分

最新推荐文章于 2018-10-24 10:30:35 发布

奇怪的猪猪

最新推荐文章于 2018-10-24 10:30:35 发布

阅读量1.3k

点赞数

分类专栏： ACM_计算几何

本文链接：https://blog.csdn.net/jingqi814/article/details/25575569

版权

ACM_计算几何专栏收录该内容

74 篇文章 0 订阅

订阅专栏

题目地址：LA3218

额直接学习LRJ的代码：

// LA3218/UVa1340 Find the Border
// Rujia Liu
// 注意：本题可以直接使用“卷包裹”法求出外轮廓。本程序只是为了演示PSLG的实现
#include<cstdio>
#include<vector>
#include<cmath>
#include<algorithm>
#include<cstring>
#include<cassert>
#include<iostream>

using namespace std;

const double eps = 1e-8;
double dcmp(double x) {
    if(fabs(x) < eps) return 0; else return x < 0 ? -1 : 1;
}

struct Point {
    double x, y;
    Point(double x=0, double y=0):x(x),y(y) { }
};

typedef Point Vector;

Vector operator + (Vector A, Vector B) {
    return Vector(A.x+B.x, A.y+B.y);
}

Vector operator - (Point A, Point B) {
    return Vector(A.x-B.x, A.y-B.y);
}

Vector operator * (Vector A, double p) {
    return Vector(A.x*p, A.y*p);
}

Vector operator / (Vector A, double p) {
    return Vector(A.x/p, A.y/p);
}


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

bool operator == (const Point& a, const Point &b) {
    return dcmp(a.x-b.x) == 0 && dcmp(a.y-b.y) == 0;
}

double Dot(Vector A, Vector B) { return A.x*B.x + A.y*B.y; }
double Cross(Vector A, Vector B) { return A.x*B.y - A.y*B.x; }
double Length(Vector A) { return sqrt(Dot(A, A)); }

typedef vector<Point> Polygon;

Point GetLineIntersection(const Point& P, const Vector& v, const Point& Q, const Vector& w) {
    Vector u = P-Q;
    double t = Cross(w, u) / Cross(v, w);
    return P+v*t;
}

bool SegmentProperIntersection(const Point& a1, const Point& a2, const Point& b1, const Point& b2) {
    double c1 = Cross(a2-a1,b1-a1), c2 = Cross(a2-a1,b2-a1),
    c3 = Cross(b2-b1,a1-b1), c4=Cross(b2-b1,a2-b1);
    return dcmp(c1)*dcmp(c2)<0 && dcmp(c3)*dcmp(c4)<0;
}

bool OnSegment(Point p, Point a1, Point a2) {
    return dcmp(Cross(a1-p, a2-p)) == 0 && dcmp(Dot(a1-p, a2-p)) < 0;
}

typedef   vector<Point>   Polygon;

// 多边形的有向面积

double PolygonArea(Polygon poly)    //   有向面积
{
    double area=0;
    int n=poly.size();
    for(int i=1;i<n-1;i++)
    {
        area+=Cross(poly[i]-poly[0],poly[(i+1)%n]-poly[0]);
        
    }
    
    return area/2;
}



struct Edge
{
    int from;
    int  to;
    double  ang;
};

const int maxn=10000+10;  // 最大边数

struct PSLG
{
    int n,m,face_cnt;
    double  x[maxn],y[maxn];
    vector<Edge> edges;
    vector<int> G[maxn];
    
    int vis[maxn*2];
    int left[maxn*2];
    int prev[maxn*2];
    
    vector<Polygon> faces;
    double area[maxn];
    
    void init(int n)
    {
        this->n=n;
        for(int i=0;i<n;i++)
        {
            G[i].clear();
        }
        
        edges.clear();
        faces.clear();
        
    }
    
    double getAngle(int  from,int to)
    {
        return   atan2(y[to]-y[from],x[to]-x[from]);
    }
    
    void AddEdge(int from,int to)
    {
        edges.push_back((Edge) {from,to,getAngle(from, to)});
        edges.push_back((Edge) {to,from,getAngle(to, from)});
        m=edges.size();
        G[from].push_back(m-2);
        G[to].push_back(m-1);
        
    }
    // 图已经建好  AddEdge() 已经完成
     void Build()
    {
        for(int u=0;u<n;u++)
        {
            int d=G[u].size();
            for(int i=0;i<d;i++)
                for(int j=i+1;j<d;j++)
                {
                    if(edges[G[u][i]].ang>edges[G[u][j]].ang)
                    {
                        swap(G[u][i],G[u][j]);
                    }
                }
            
            for(int i=0;i<d;i++)
            {
                prev[G[u][(i+1)%d]]=G[u][i];
            }
            
        }
        
        memset(vis,0, sizeof(vis));
        face_cnt=0;
        // 从0号顶点开始
        for(int u=0;u<n;u++)
        {
            for(int i=0;i<G[u].size();i++)
            {
                int e=G[u][i];
                if(!vis[e])
                {
                    face_cnt++;  // 从1开始计数
                    Polygon poly;
                    for(;;)
                    {
                        vis[e]=1;
                        left[e]=face_cnt;
                        int from=edges[e].from;
                        Point P(x[from],y[from]);
                        poly.push_back(P);
                        e=prev[e^1];  // 反向边顺时针旋转的第一条
                        if(e==G[u][i])  break;
                        
                    }
                    faces.push_back(poly);
                }
            }
        }
        
        for(int i=0;i<face_cnt;i++)
        {
            area[i]=PolygonArea(faces[i]);
        }
    
   }
    
    
};

PSLG g;
const int maxp=100+5;

Point P[maxp];
Point V[maxp*(maxp-1)/2+maxp];

int n,c;
int ID(Point P)
{
    return lower_bound(V,V+c,P)-V;
    
}

Polygon  simplify(const Polygon & poly)   //  去共线
{
    Polygon ans;
    int n=poly.size();
    for(int i=0;i<n;i++)
    {
        Point a=poly[i];
        Point b=poly[(i+1)%n];
        Point c=poly[(i+2)%n];
        if(dcmp(Cross(b-a,c-b))!=0)
        {
            ans.push_back(b);
            
        }
           
    }
    
    return ans;
}

void  build_gragh()
{
    for(int i=0;i<n;i++)
    {
        V[i]=P[i];
    }
    
    vector<double> dist[maxp];
    
    c=n;
    for(int i=0;i<n;i++)
        for(int j=i+1;j<n;j++)
        {
             if(SegmentProperIntersection(P[i], P[(i+1)%n], P[j], P[(j+1)%n]))
             {
                 Point ip=GetLineIntersection(P[i], P[(i+1)%n]-P[i], P[j], P[(j+1)%n]-P[j]);
                 V[c++]=ip;
                 dist[i].push_back(Length(ip-P[i]));
                 dist[j].push_back(Length(ip-P[j]));
                 
             }
        }
    
    sort(V,V+c);
    c=unique(V,V+c)-V;
    
    g.init(c);
    
    for(int i=0;i<c;i++)
    {
        g.x[i]=V[i].x;
        g.y[i]=V[i].y;
    }
    
    
    for(int i=0;i<n;i++)   // 按在那条边上来统计
    {
        Vector v=P[(i+1)%n]-P[i];
        double len=Length(v);
        v=v/len;
        
        dist[i].push_back(0);
        dist[i].push_back(len);
        
        sort(dist[i].begin(),dist[i].end());
        
        int sz=dist[i].size();
        
        for(int j=1;j<sz;j++)
        {
            Point a=P[i]+v*dist[i][j-1];
            Point b=P[i]+v*dist[i][j];
            if(a==b) continue;
            g.AddEdge(ID(a),ID(b));
            
        }
        
    }
    
    g.Build();
    
    
    
}

int main()
{
    while(cin>>n&&n)
    {
        
        for(int i=0;i<n;i++)
        {
            scanf("%lf%lf",&P[i].x,&P[i].y);
        }
        build_gragh();
        
        
//        for(int i=0;i<g.faces.size();i++)
//        {
//            cout<<g.faces[i].size()<<endl;
//            for(int j=0;j<g.faces[i].size();j++)
//            {
//                cout<<g.faces[i][j].x<<' '<<g.faces[i][j].y<<' ';
//            }
//            cout<<endl;
//        }
        Polygon poly;
        
        for(int i=0;i<g.faces.size();i++)
        {
             if(g.area[i]<0)
             {
                 
                 poly=g.faces[i];
                 reverse(poly.begin(),poly.end());
                 poly=simplify(poly);
                 break;
             }
        }
        
     
        
        int start=0;
        int m=poly.size();
        cout<<m<<endl;
        for(int i=0;i<m;i++)
        {
            if(poly[i]<poly[start])
            {
                start=i;
            }
        }
        
        for(int i=start;i<m;i++)
        {
            printf("%.4f %.4f\n",poly[i].x,poly[i].y);
        }
        
        for(int i=0;i<start;i++)
        {
            printf("%.4f %.4f\n",poly[i].x,poly[i].y);
        }
        
    }
}