poj3525 二分+半平面交

该博客介绍了如何解决POJ3525问题,即计算给定凸多边形中距离边界最远的点到边界的最大距离。通过使用二分搜索结合半平面交的方法,逐步将多边形的每条边向内平移,判断是否存在内部点,从而找到最远距离。
摘要由CSDN通过智能技术生成

题目链接:http://poj.org/problem?id=3525

题意:给定一个凸多边形,求多边形中距离边界最远的点到边界的距离。

思路 : 每次将凸多边形每条边往里平移d,判断是否存在核;二分d即可。

#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cmath>
using namespace std;

const double eps = 1e-10;
const int maxn = 105;

int dq[maxn], top, bot, pn, order[maxn], ln;
struct Point {
    double x, y;
} p[maxn];

struct Line {
    Point a, b;
    double angle;
} l[maxn], tmp[maxn];

int dblcmp(double k) {
    if (fabs(k) < eps) return 0;
    return k > 0 ? 1 : -1;
}

double multi(Point p0, Point p1, Point p2) {
    return (p1.x-p0.x)*(p2.y-p0.y)-(p1.y-p0.y)*(p2.x-p0.x);
}

bool cmp(int u, int v) {
    int d = dblcmp(l[u].angle-l[v].angle);
    if (!d) return dblcmp(multi(l[u].a, l[v].a, l[v].b)) > 0;
    return d < 0;
}

void getIntersect(Line l1, Line l2, Point& p) {
    double dot1,dot2;
    dot1 = multi(l2.a, l1.b, l1.a);
    dot2 = multi(l1.b, l2.b, l1.a);
    p.x = (l2.a.x * dot2 + l2.b.x * dot1) / (dot2 + dot1);
    p.y = (l2.a.y * dot2 + l2.b.y * dot1) / (dot2 + dot1);
}


bool judge(Line l0, Line l1, Line l2) {
    Point p;
    getIntersect(l1, l2, p);
    return dblcmp(multi(p, l0.a, l0.b)) < 0;
}

void addLine(double x1, double y1, double x2, double y2) {
    l[ln].a.x = x1; l[ln].a.y = y1;
    l[ln].b.x = x2; l[ln].b.y = y2;
    l[ln].angle = atan2(y2-y1, x2-x1);
    ln++;
}

bool halfPlaneIntersection(Line l[], int n) {
    int i, j;
    for (i = 0; i < n; i++) order[i] = i;
    sort(order, order+n, cmp);
    for (i = 1, j = 0; i < n; i++)
        if (dblcmp(l[order[i]].angle-l[order[j]].angle) > 0)
            order[++j] = order[i];
    n = j + 1;
    dq[0] = order[0];
    dq[1] = order[1];
    bot = 0;
    top = 1;
    for (i = 2; i < n; i++) {
        while (bot < top && judge(l[order[i]], l[dq[top-1]], l[dq[top]])) top--;
        while (bot < top && judge(l[order[i]], l[dq[bot+1]], l[dq[bot]])) bot++;
        dq[++top] = order[i];
    }
    while (bot < top && judge(l[dq[bot]], l[dq[top-1]], l[dq[top]])) top--;
    while (bot < top && judge(l[dq[top]], l[dq[bot+1]], l[dq[bot]])) bot++;
    if (bot + 1 >= top) return false; //当dq中少于等于两条边时,说明半平面无交集
    return true;
}

double getDis(Point a, Point b) {
    return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));
}

void changePolygon(double h) {
    double len, dx, dy;
    for (int i = 0; i < ln; i++) {
        len = getDis(l[i].a, l[i].b);
        dx = (l[i].a.y - l[i].b.y) / len * h;
        dy = (l[i].b.x - l[i].a.x) / len * h;
        tmp[i].a.x = l[i].a.x + dx;
        tmp[i].a.y = l[i].a.y + dy;
        tmp[i].b.x = l[i].b.x + dx;
        tmp[i].b.y = l[i].b.y + dy;
        tmp[i].angle = l[i].angle;
    }
}

double BSearch() {
    double l = 0, r = 20000, mid;
    while (l + eps < r) {
        mid = (l + r) / 2;
        changePolygon(mid);
        if (halfPlaneIntersection(tmp, ln))
            l = mid;
        else r = mid;
    }
    return l;
}

int main()
{
    int i;

    while (scanf ("%d", &pn) && pn) {
        for (i = 0; i < pn; i++)
            scanf ("%lf%lf", &p[i].x, &p[i].y);
        for (i = ln = 0; i < pn-1; i++)
            addLine(p[i].x, p[i].y, p[i+1].x, p[i+1].y);
        addLine(p[i].x, p[i].y, p[0].x, p[0].y);

        printf ("%.6lf\n", BSearch());
    }
    return 0;
}




评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值