Как найти максимальный радиус окружности которую можно вписать между объектами на плоскости в произвольной точке?

Question

Есть плоскость, а на ней много разных объектов. Если это важно - пусть будут круги с разным радиусом, или квадраты - не принципиально. Скажем, круги.

Задача для произвольной точки быстро определить, какая окрестность не занята. То есть просто для кругов смотрим расстояние от центра до этой точки минус радиус круга, и ищем минимум.

Но кругов много, а времени мало. Как эту задачу решить проще всего? Смотрел на quadtree, на R-дерево. Но или я плохо понимаю их идеи, или они годятся для поиска точек, но не для кругов. Да и для точек - а если ближайшие точки в quadtree находятся в соседних квадрантах - чем это нам может помочь? Все равно, получается, надо пересматривать почти все дерево?

Подскажите, какие алгоритмы в данном случае стоит применять? Чтобы проверять как можно меньше разных кругов?

Для начала можно считать, что круги не пересекаются (если это может помочь). Но вообще говоря, не факт.

Answer 1

Теория

Построение поискового дерева

Для круга и квадрата ориентированного по осям научимся определять наменьшее расстояние от точки в квадрате до круга: d_min. Если они пересекаются, примем это расстояние нулевым. Также определим наибольшее расстояние от точки квадрата до круга: d_max. Если квадрат целиком внутри круга, положим его нулевым.

Пусть у нас есть квадрат и набор кругов. По всем кругам выберем минимум наибольших расстояний: h = min{d_max}. Индексы кругов с условием d_min ≤ h запишем квадрат. Если точка попадёт в квадрат, нужно будет вычислять расстояния только до этого набора кругов, остальные точно дальше.

Накроем квадратом область в которой требуется искать ближайший круг. Вычислим список подходящих кругов. Если список слишком длинный, разобъём квадрат на четыре меньших квадрата и повторим построение списков рекурсивно. Рекурсия останавливается если список кругов не длинее какого-то порога. Рекурсия также должна быть ограничена по глубине.

Локализация

Локализация точки делается так: точка локализуется относительно квадрата, затем выбирается его четверть и так далее, пока не дойдём до листа. Вычислим минимум расстояний от точки до кругов из списка кругов в листе.

Время локализации складывается из спуска по квадратам (не более максимальной глубины рекурсии) и вычисления расстояния до кругов из списка (обычно не более порога на длину списка).

Метод не ограничен кругами. Если вы умеет строить минимальное и максимальное расстояния от квадрата до фигуры, фигура может принять участие в локализации.

Недостаток

Если на плоскости есть точки для которых максимальная пустая окружность касается одновременно k препятствий, список в квадрате будет длиной не меньше k. Другими словами, локализация тратит много памяти и работает медленно если в диаграмме Вороного есть вершины большой степени. Если мы решаем задачу точно, то этот недостаток не обойти. Если можно выдавать ответ с некоторой точностью, список можно сократить до требуемого размера.

На практике, если вы не конструируете сцену специальным образом, то вероятность что, скажем, десять окружностей соберутся вокруг общей точки исчезающе мала. Кроме того у вас есть все средства чтобы обнаружить проблему на этапе построения поискового дерева и решить задачу локализации в одном плохом и очень маленьком квадрате другим способом, если вы придумаете такой способ, конечно.

Пример

const sub = ([x1, y1], [x2, y2]) => [x1 - x2, y1 - y2];
const dot = ([x1, y1], [x2, y2]) => x1 * x2 + y1 * y2;
const mul = ([x, y], f) => [f * x, f * y];
const dist = (p1, p2) => Math.hypot(...sub(p1, p2));

const minDistPointSegment = (q, [p1, p2]) => {
    const dp = sub(p2, p1);
    const dq = sub(q, p1);

    const a = Math.max(0, Math.min(1, dot(dq, dp) / dot(dp, dp)));
    return dist(dq, mul(dp, a));
};

const minDist = ([[x1, y1], [x2, y2]], [x, y, r]) => {
    if (x1 <= x && x <= x2 && y1 <= y && y <= y2) {
        return 0;
    }

    return Math.max(0, Math.min(
        minDistPointSegment([x, y], [[x1, y1], [x2, y1]]),
        minDistPointSegment([x, y], [[x2, y1], [x2, y2]]),
        minDistPointSegment([x, y], [[x2, y2], [x1, y2]]),
        minDistPointSegment([x, y], [[x1, y2], [x1, y1]])
    ) - r);
};

const maxDist = ([[x1, y1], [x2, y2]], [x, y, r]) => Math.max(0, Math.max(
    dist([x, y], [x1, y1]),
    dist([x, y], [x2, y1]),
    dist([x, y], [x2, y2]),
    dist([x, y], [x1, y2])
) - r);

const aMin = a => a.reduce((a, b) => Math.min(a, b), Infinity);

const selectDisks = (rect, disks) => {
    const h = aMin(disks.map(d => maxDist(rect, d)));
    if (h == 0) {
        for (const d of disks) {
            if (maxDist(rect, d) == 0) {
                return [d];
            }
        }
    }
    return disks.filter(d => minDist(rect, d) <= h);
};

const buildSearchTree = (depth, length, rect, disks) => {
    const selected = selectDisks(rect, disks);
    if (selected.length <= length || depth <= 0) {
        const maxRadius = (x, y) => [
            Math.max(0, aMin(selected.map(
                ([cx, cy, r]) => Math.hypot(x - cx, y - cy) - r
            ))),
            depth,
            selected.length
        ];
        return [maxRadius, rect, []];
    }
    const [[x1, y1], [x2, y2]] = rect;
    const xm = (x1 + x2) / 2;
    const ym = (y1 + y2) / 2;

    const subtrees = [
        buildSearchTree(depth - 1, length, [[x1, y1], [xm, ym]], selected),
        buildSearchTree(depth - 1, length, [[xm, y1], [x2, ym]], selected),
        buildSearchTree(depth - 1, length, [[x1, ym], [xm, y2]], selected),
        buildSearchTree(depth - 1, length, [[xm, ym], [x2, y2]], selected)
    ];

    const maxRadius = (x, y) => subtrees[
        (x < xm) ? ((y < ym)? 0 : 2) : ((y < ym)? 1 : 3)
    ][0](x, y);

    return [maxRadius, rect, subtrees];
};

const makeRect = (svg, [[x1, y1], [x2, y2]]) => {
    const rect = document.createElementNS(
        'http://www.w3.org/2000/svg',
        'rect'
    );
    rect.setAttribute('x', x1);
    rect.setAttribute('y', y1);
    rect.setAttribute('width', x2 - x1);
    rect.setAttribute('height', y2 - y1);
    rect.setAttribute('stroke', 'grey');
    rect.setAttribute('fill', 'none');
    svg.appendChild(rect);
};

const makeCircle = (svg, [x, y, r]) => {
    const circle = document.createElementNS(
        'http://www.w3.org/2000/svg',
        'circle'
    );
    circle.setAttribute('cx', x);
    circle.setAttribute('cy', y);
    circle.setAttribute('r', r);
    circle.setAttribute('stroke', 'red');
    circle.setAttribute('fill', 'none');
    svg.appendChild(circle);
};

const makeDisk = svg => {

    const group = document.createElementNS('http://www.w3.org/2000/svg', 'g');
    svg.appendChild(group);

    const moveGroup = (x, y) => {
        const ctm = group.getCTM();
        ctm.e = x;
        ctm.f = y;

        const t = svg.createSVGTransform();
        t.setMatrix(ctm);

        group.transform.baseVal.initialize(t);
    };

    const center = document.createElementNS(
        'http://www.w3.org/2000/svg',
        'circle'
    );
    center.setAttribute('r', 2);
    center.setAttribute('stroke', 'black');
    center.setAttribute('fill', 'black');
    group.appendChild(center);

    const circle = document.createElementNS(
        'http://www.w3.org/2000/svg',
        'circle'
    );
    circle.setAttribute('stroke', 'black');
    circle.setAttribute('fill', 'none');
    group.appendChild(circle);

    const text = document.createElementNS(
        'http://www.w3.org/2000/svg',
        'text'
    );
    group.appendChild(text);

    return {
        'move': moveGroup,
        'setRadius': r => circle.setAttribute('r', r),
        'setText': t => text.textContent = t
    };
};

const x3 = x => (-2 * x + 3) * x * x;
const random2 = () => x3(x3(Math.random()));

const main = () => {
    const body = document.getElementsByTagName('body')[0];

    const svg = document.createElementNS(
        'http://www.w3.org/2000/svg',
        'svg'
    );
    const width = 600;
    const height = 200;
    svg.setAttribute('width', width);
    svg.setAttribute('height', height);
    body.appendChild(svg);

    // const n = 100000;
    const n = 100;
    const maxR = 40;
    const disks = [];
    for (let i = 0; i < n; ++i) {
        disks.push([
            width * random2(),
            height * random2(),
            maxR * Math.random()
        ]);
    }

    const maxDepth = 10;
    const size = Math.max(width, height);
    const tree = buildSearchTree(
        maxDepth,
        5,
        [[0, 0], [size, size]],
        disks
    );

    const drawTree = tree => {
        makeRect(svg, tree[1]);
        tree[2].forEach(drawTree);
    };
    drawTree(tree);

    disks.forEach(d => makeCircle(svg, d));

    const disk = makeDisk(svg);

    svg.onmousemove = e => {
        const p = svg.createSVGPoint();
        p.x = e.clientX;
        p.y = e.clientY;
        const {x, y} = p.matrixTransform(svg.getScreenCTM().inverse());
        disk.move(x, y);

        const [r, depth, length] = tree[0](x, y);
        disk.setRadius(r);
        disk.setText(`${maxDepth - depth}+${length}`);
    };
};

document.addEventListener('DOMContentLoaded', main);

Красные круги - препятствия. Чёрный круг следует за курсором мыши и показывает максимальный радиус в данной точке. Внутри чёрного круга два числа: глубина узла в дереве и количество дисков-препятствий в этом узле. Запустите код и подвигайте курсор мыши по картинке. В примере сто случайных препятствий, для любого положения курсора вычисляется растояние до пяти препятствий максимум.

Картинка

Чтобы нарисовать взвешенную диаграмму Вороного, порог для длины списка ставится в единицу. Тогда рекурсия останавливается только если квадрат внутри ячейки Вороного одного круга. Рядом с рёбрами диаграммы деление продолжается до предела глубины рекурсии. Мелкие серые квадраты собираются вдоль рёбер и вершин диаграммы и они становятся видны.

Границы ячеек между кругами разного радиуса искривлены - это куски гипербол.

БЛОГ НА HUSL