%--------------------------
 % @Author: Jingqiao Hu
 % @Date: 2022-01-24 13:48:58
 % @LastEditTime: 2022-01-24 14:32:52

 % Compute the barycentric weights for a point p in an n-gon Q
 % Assumes p is strictly within Q and the vertices qj of Q are ordered.
 % assumes weight w = 1

 % PARAS:
 % polygon: ordered vertices qj of Q
%--------------------------
function w = compute_barycentric(polygon, p)

    n = size(polygon, 1);
    w = zeros(n, 1);
    wsum = 0;

    for j = 0:n-1
        prev = mod(j + n - 1, n) + 1;
        next = mod(j + 1, n) + 1;
        jj = j+1;

        qj = polygon(jj, :);
        qp = polygon(prev, :);
        qn = polygon(next, :);

        ct = cotangent(p, qj, qp) + cotangent(p, qj, qn);
        w(jj) = ct / (norm(p - qj)^2 + 1e-9);
        wsum = wsum + w(jj);
    end
    w = w ./ wsum;
end

function ct = cotangent(a, b, c)

    ba = a-b;
    bc = c-b;

    ct = dot(bc, ba) / (norm(cross(bc, ba)) + 1e-9);
end