C++ code for cotangent weights over a triangular mesh

$$w_{ij} = \cot(\alpha) + \cot(\beta)$$

Here is some code to compute the cotangent weights described in the article Discrete Differential-Geometry Operators for Triangulated 2-Manifolds. Cotangent weights have many application such as mesh smoothing, computing harmonic weights over the mesh, or diffusing values on its surface. I'm leaving this here for other tutorial to reference.

A cotangent weight can be thought as barycentric weights, for every vertex \( i \) we associate the weights to each of its direct neighbor \( j \). The code below does not compile but hoppefuly it can serve as a reference to adapt it to your own project.

struct Mesh {

    /// For every vertex, a list of direct (1st ring) neighbors.
    /// first_ring[index_vertex].size() == number of direct neighbors to 'index_vertex'
    /// first_ring[index_vertex][ith_neighbor] == neighbor_index
    std::vector< std::vector<Vert_idx> > first_ring();

    Vec3 vertex(Vert_idx idx) const { return _vertices[idx]; }
    int nb_vertices() const { return (int)_vertices.size(); }

    std::vector<Vec3> _vertices;
};

// -----------------------------------------------------------------------------

void laplacian_cotan_weights(const Mesh& m,
                             std::vector< std::vector >& per_edge_cotan_weights)
{
    int nb_verts = m.geom().nb_vertices();
    per_edge_cotan_weights.resize( nb_verts );
    for(int vertex_index = 0; vertex_index < nb_verts; ++vertex_index)
    {
        int nb_neighbors = m.first_ring()[vertex_index].size();
        per_edge_cotan_weights[i].resize( nb_neighs );
        for(int n = 0; n < nb_neighbors; ++n)
        {
            float weight = laplacian_cotan_weight(m, vertex_index, n);
            //std::cout << "vert: " << i << " cotan wij " << w << std::endl;
            per_edge_cotan_weights[vertex_index][n] = weight;
        }
    }
}

// -----------------------------------------------------------------------------

/// compute the cotangent weight for the edge (v_i, v_j)
float laplacian_cotan_weight(const Mesh& m,
                             Vert_idx vertex_i,
                             Edge_idx edge_j )
{
    int nb_neighs = m.first_ring()[vertex_i].size();

    Vert_idx vertex_j      = m.first_ring()[vertex_i][  edge_j                                       ];
    Vert_idx vertex_j_next = m.first_ring()[vertex_i][ (edge_j+1) % nb_neighs                        ];
    Vert_idx vertex_j_prev = m.first_ring()[vertex_i][ (edge_j-1) >= 0 ? (edge_j-1) : nb_neighs-1 ];

    Vec3 v1 = m.vertex(vertex_i) - m.vertex(vertex_j_prev);
    Vec3 v2 = m.vertex(vertex_j) - m.vertex(vertex_j_prev);
    Vec3 v3 = m.vertex(vertex_i) - m.vertex(vertex_j_next);
    Vec3 v4 = m.vertex(vertex_j) - m.vertex(vertex_j_next);

    // wij = (cot(alpha) + cot(beta)),
    // for boundary edge, there is only one such edge
    float cotan1 = 0.0f;
    float cotan2 = 0.0f;

    // If the mesh is not water-tight (e.g. a piece of cloth without thickness)
    // we must check for edges (v_i, v_j) that lies on the sides/boundaries of
    // the mesh
    if(!m.is_vert_on_side(vertex_i) || !m.is_vert_on_side(vertex_j))
    {
        // general case: not a boundary
        cotan1 = (v1.dot(v2)) / (v1.cross(v2)).norm();
        cotan2 = (v3.dot(v4)) / (v3.cross(v4)).norm();
    }
    else // boundary edge, only have one such angle
    {
        if( vertex_j_next == vertex_j_prev )
        {
            // two angles are the same, e.g. corner of a square
            cotan1 = (v1.dot(v2)) / (v1.cross(v2)).norm();
        }
        else
        {
            // find the angle not on the boundary
            if( !m.is_vert_on_side(vertex_i) || !m.is_vert_on_side(vertex_j_next) )
                cotan2 = (v3.dot(v4)) / (v3.cross(v4)).norm();
            else
                cotan1 = (v1.dot(v2)) / (v1.cross(v2)).norm();
        }
    }

    return (cotan1 + cotan2);
}

Note: recall that
$$
\cot \theta = \frac{\cos \theta}{\sin \theta}
$$
$$
\frac{\vec u . \vec v}{ \vec u \times \vec v} = \frac{\||\vec u\||.\||\vec v\|| cos(\theta) }{\||\vec u\||.\||\vec v\|| sin(\theta)} = \cot(\theta)
$$

Some other references making use of cotangent weights:

One comment

this is cool

yes - 16-04-’18 21:18
All html tags except <b> and <i> will be removed from your comment. You can make links by just typing the url or mail-address.
Spam bot question: