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

У меня есть 2D-плоскость, разделенная на n-сторонние выпуклые многоугольники. Я использую алгоритм WRF PNPOLY для включения полигонов в убедитесь, что точка принадлежит одному и только одному многоугольнику.

Есть ли алгоритм, который я могу использовать для обрезки сегмента линии PO до заданного многоугольника на плоскости, предполагая, что pnpoly(O) == true, так что pnpoly(P') всегда будет истинным ?

Моя текущая реализация clipToPoly выполняет тест пересечения линии с сегментом и каждым краем многоугольника, а затем использует точку пересечения (как подробно описано в это ТАК ответ), но это не всегда дает точку, которая удовлетворяет PNPOLY.

function clipPointToPoly(p, o, poly) {
    var i, j, n = poly.length,
        q, r = {}, s = {}, pq = {},
        rxs, t, u;

    function cross2(v, w) {
        return v.x * w.y - v.y * w.x;
    }

    for (i = 0, j = n - 1; i < n; j = i++) {
        q = poly[i];
        s.x = poly[j].x - q.x;
        s.y = poly[j].y - q.y;

        r.x = o.x - p.x;
        r.y = o.y - p.y;

        rxs = cross2(r, s);
        if (rxs !== 0) {
            pq.x = q.x - p.x;
            pq.y = q.y - p.y;

            t = cross2(pq, s) / rxs;
            u = cross2(pq, r) / rxs;
            if (0 < u && u < 1 && 0 < t && t < 1) {
                p.x = p.x + t * r.x;
                p.y = p.y + t * r.y;
                return true;
            }
        }
    }
    return false;
};

Вот моя реализация PNPOLY:

function pnpoly(p, poly) {
    var i, j, e0, e1,
        n = polygon.length,
        inside = false;
    for (i = 0, j = n - 1; i < n; j = i++) {
        e0 = poly[i];
        e1 = poly[j];
        if ( ((p0.y > p.y) !== (p1.y > p.y)) &&
             ((p.x < (p1.x - p0.x) * (p.y - p0.y) / (p1.y - p0.y) + p0.x)) ) {
            inside = !inside;
        }
    }
    return inside;
};

Я не думаю, что понимаю достаточно о симуляции простоты или о том, как работать с полуоткрытыми наборами PNPOLY при использовании чисел с плавающей запятой для правильной обработки крайних случаев.

Например:

poly: [(1,1), (-1,1), (-1,-1), (1,-1)]
p: (5,5)
o: (0,0)
p' = (1,1)

Это не удается, потому что (1,1) не включено в соответствии с PNPOLY (оно находится на открытой стороне набора), но clipToPoly не принимает это во внимание. Я полагаю, что мог бы подтолкнуть его на эпсилон, если бы знал, что он находится на открытом конце набора, но я бы предпочел более стабильное решение.

Еще один пример:

poly: [-995.9592341908675, -88.48705014724577
       -1040.5031753180106, -176.53192722405026
       -549.9211095905894, -330.8462151682281
       -653.7143990581328, -211.59193148034612]
p: -1032.3773586525654, -208.3586379393678
o: -957.4172402148379, -202.6668958854324

В этом случае clipToPoly дает сбой, потому что O находится так близко к краю многоугольника, что даже не обнаруживает пересечение из-за неточности с плавающей запятой.

t: 1.0000000000000002 u: 0.8306380503739466

Есть ли способ получить неточность с плавающей запятой clipToPoly, чтобы она соответствовала PNPOLY, чтобы оба были согласованы?