Graham's Scan法求凸包

凸包

凸包(Convex Hull)是一个计算几何(图形学)中的概念。
在一个实数向量空间V中,对于给定集合X,所有包含X的凸集的交集S被称为X的凸包。
X的凸包可以用X内所有点(X1,…Xn)的线性组合来构造.

上面的都没用.
给定二维平面上的点集,凸包就是将最外层的点连接起来构成的凸多边型,它能包含点集中所有的点。
可以想象成一圈橡皮筋收紧成的一圈
这里写图片描述

GrahamsScan G r a h a m ′ s S c a n 算法

用于求解凸包,O(n)
首先我们找到一个最低的(y轴坐标最小的)标准点,然后按照极角排个序,当两点与标准点共线时,距离小的排前面.
极角就是这个点与标准点所成线与x轴的所夹角
就像这样
这里写图片描述

注意,实际上我们不需要求得其夹角(因为涉及到三角函数效率低而且会有精度问题),用叉积判断即可
叉积公式与方向的判断
这里写图片描述
叉积的长度 |a×b| 可以解释成以a和b为邻边的平行四边形的面积。证明这里
这里写图片描述

然后从左边开始进栈,初始先进两个(图中连边点)
这里写图片描述
然后进第三个(d点),我们要判断一下如果直接取BD,那么还需不需要BC
这里写图片描述

怎么判断呢?
其实就是判断BD在不在C的外侧1,如果BD就可以包括C了,那凸包加上C也没有意义了,在栈中把C退掉,然后在对A,B重复上述过程直到不满足上述条件.
判断外侧的话还是要用到叉积,当其叉积为正时我们知道p在q的顺时针
一直进完其他点,最后留在栈中的就是这个点集的凸包了.

CODE

type
    dcm=real;
    point=record
        x,y:dcm;
    end;
var
    pnts,tmp:array[0..10000] of point;
    stack:array[0..10000] of longint;
    len,size,i,j,n,m,k:longint;
    ans:dcm;
    std:point;
function getDis(a,b:point):dcm; //平面两点距离
begin
    exit(sqrt(sqr(a.x-b.x)+sqr(a.y-b.y)));
end;

function DPT(a,b,mid:point):dcm; //向量叉积,mid是原点
begin
    exit((a.x-mid.x)*(b.y-mid.y)-(b.x-mid.x)*(a.y-mid.y));
end;

function inOrder(left,right:point):boolean; //如left在right的逆时针则返回true
var result:dcm;
begin
    result:=DPT(left,right,std);
    if (result<0)or(result=0)and(getDis(left,std)<getDis(right,std)) then exit(true)  //left-->right
        else exit(false); //q-->p
end;

procedure qsort(l,r: longint);
var
    i,j,y: longint;
    x:point;
begin
    i:=l;
    j:=r;
    x:=pnts[(l+r) shr 1]; //这里不能直接写下标,因为被排序时值会交换,这样中间值就会变.
    repeat
        while inOrder(pnts[i],x) do inc(i);
        while inOrder(x,pnts[j]) do dec(j);
        if not(i>j) then
        begin
            pnts[0]:=pnts[i];
            pnts[i]:=pnts[j];
            pnts[j]:=pnts[0];
            inc(i);dec(j);
        end;
    until i>j;
    if l<j then qsort(l,j);
    if i<r then qsort(i,r);
end;

procedure init;
begin
    readln(n);
    std.x:=maxlongint;
    std.y:=maxlongint; //标准点STD
    for i:=1 to n do
        with tmp[i] do
        begin
            readln(x,y);
            if (y<std.y)or((y=std.y)and(x<std.x)) then
            begin
                std:=tmp[i];
            end;
        end; 

    for i:=1 to n do
        if (tmp[i].x<>std.x)or(tmp[i].y<>std.y) then
        begin
            inc(len);
            with pnts[len] do
            begin
                x:=tmp[i].x;
                y:=tmp[i].y;
            end; //其实可以将std放到0下标直接排序就好了  这里写的略麻烦
        end;
    dec(n);
    qsort(1,n);
    pnts[0]:=std;

end;
procedure main;
begin
    stack[1]:=0;
    stack[2]:=1;
    stack[0]:=2;
    for i:=2 to n do
    begin
        while (DPT(pnts[stack[stack[0]]],pnts[i],pnts[stack[stack[0]-1]])>=0) do
            dec(stack[0]);
        inc(stack[0]);
        stack[stack[0]]:=i;
    end;
end;
procedure print;
begin
    ans:=0;
    for i:=2 to stack[0] do
        ans:=ans+getDis(pnts[stack[i]],pnts[stack[i-1]]);

    writeln(ans+getDis(pnts[stack[stack[0]]],pnts[stack[1]]):0:2);
end;

begin
    init();
    main();
    print();
end.

  1. 这个要看顺时针还是逆时针
  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 2
    评论
评论 2
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值