# Untitled

unknown
plain_text
2 months ago
1.6 kB
2
Indexable
Never
#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])