Untitled

#Reset positions to non-noisy values
pos[:] = pos_orig

# Compute the symmetric Laplacian
def create_laplacian(m, symmetric=False):
    # symmetric is a flag that you can choose in order to make the matrix L symmetric or not
    no_verts = m.no_allocated_vertices()
    L = zeros((no_verts,no_verts)) # The cot Laplacian that we end up outputting 
    A = zeros((no_verts)) # Use this vector to store 1-ring areas

    for i in m.vertices():
        for f in m.circulate_face(i,mode = 'f'):
            A[i] += m.area(f)
    
    for i in m.vertices():
        neighbour_edges = m.circulate_vertex(i, 'h')
        for h in neighbour_edges:
            j = m.incident_vertex(h)
            next = m.incident_vertex(m.next_halfedge(m.opposite_halfedge(h)))
            prev = m.incident_vertex(m.next_halfedge(h))
            cot_alpha = cotan_angle(m, i, j, prev)
            cot_beta = cotan_angle(m, i, j, next)
            if symmetric==True:
                L[i,j] = (cot_alpha+cot_beta)/(2*np.sqrt(A[i]*A[j]))
            else:
                L[i,j] = (cot_alpha+cot_beta)/(2*A[i])
        L[i,i] = -(L[i,:].sum())
    return L

L = create_laplacian(m, symmetric=True)

# Next, find the eigenvectors and eigenvalues of the symmetric Laplacian
# In fact, we need to use -L because that way the eigenvectors are positive.


(W,V) = eigh(-1*L)
print(V)

# Now display some of the eigenvectors. The code below uses the OpenGL viewer
# uncomment to use.

# viewer = gl.Viewer()
# for i in range(20):
#    print(W[i])
#    viewer.display(m,mode='s',smooth=True,data=V[:,i])
# del viewer
jd.display(m,data=V[:,1])