Curvature of mesh vertex points via Python?

Each vertex of a mesh has an "angle defect" versus a flat mesh also called Gaussian Curvature. How would I obtain this vertex sharpness value (just the sum of angles around a vertex subtracted from 360 degrees) and +/- convex/concave sign for a given vertex of a mesh in a Python module in Grasshopper? I can't find Rhinocommon commands to even approach this for meshes instead of surfaces (which require a UV location using the Gaussian command) and even a list of connected vertices I am not sure if they will return in useful circular order around the vertex or just unpredictably, if I have to manually calculate angles one at a time around the vertex.

In fact, using:

print MESH.Vertices.GetConnectedVertices(45)

shows the returned vertex order is *not* in nice cyclic order around the vertex, so I don't know how to even calculate angles manually one at a time lest some pairs not even represent a face corner.

Daniel Piker suggests doing this calculation on a mesh here:

http://www.grasshopper3d.com/forum/topics/curvature-of-mesh

"What if you triangulated the mesh first, then measured the Gaussian curvature (using the angle defect) per vertex?

If there are non-planar quads then it is true that which of the 2 possible triangulations is used will slightly affect the value, but I imagine this would average out over the surrounding vertices."

I certainly don't want to have to remesh my efficient model with extra triangles though.

I display my face and vertex numbers outside of Python, as a guide:

Intuitively, I think I may be able to derive it from the face normals surrounding each vertex since on a flat mesh, they will all be the same normal versus the vertex normal and on a slightly curved mesh they will be slightly angled and on a sharp point vertex they will be highly angled versus the vertex normal.

This is all especially difficult since the Rhinocommon CHM file still has so few Python examples instead of just C# and Visual Basic. Also, how point data must be translated into Python lists manually using .X, .Y and .Z modifiers and sometimes back for me to operate on them.

Replies to This Discussion

Permalink Reply by djordje on April 26, 2015 at 4:58pm

Hi Nik,

Gianpaolo Savio developed a very nice open source plugin for Rhino, called Rhino open projects.

It has an extensive set of Rhino commands for mesh processing and analysis.

You can find the python code for the Gaussian Curvature analysis in: rhinoopensourcecode\Mesh\MeshCurvature.py.

Permalink Reply by Nik Willmore on April 26, 2015 at 5:01pm

That will be Plan B then. Plan A is to avoid having to convince colleagues to have to install a plugin in possibly different versions of Rhino/Grasshopper.

Permalink Reply by djordje on April 26, 2015 at 5:05pm

You do not have to install the plugin.
As I said, check the python source code I mentioned in the upper reply, and replicate the functionality by yourself.

Permalink Reply by Nik Willmore on April 26, 2015 at 5:15pm

Ah, quite a Python treasure there, indeed, thanks. It's nice just to see more of that instead of having to decypher C#. Those of us new to object oriented dot notation find it hard to even translate the Rhinocommon command definitions into actual strung together lines of references. Trial and error is the rule for now.

Permalink Reply by Nik Willmore on April 26, 2015 at 6:14pm

Ah, the problem with using just the connected surrounding vertex normals is they are not really local to the center vertex but drastically change according to further surrounding geometry:

Only surrounding face normals are local enough to not depend on further surrounding geometry.

Permalink Reply by Nik Willmore on April 26, 2015 at 11:53pm

I managed to grab the face normal vectors around each vertex and display the sum of those angles around each vertex, and it gives fair results:

I'm still lacking +/- signs for concave/convex, which I need to then move the sharp mesh vertices in their negative normal direction to achieve local efficient mesh smoothing.

I'm still lacking +/- signs for concave/convex, which I need to then move the sharp mesh vertices in their negative normal direction to achieve local efficient mesh smoothing. Ah! but I can just average the face centers around each vertex and obtain either an inner or outer point in more or less the right direction to move the vertex to moderate its sharpness by some amount related to the sharpness value. Some vague local geometry idea like that....

Attachments:

Niks_Kink_Killer_E.gh, 24 KB

Permalink Reply by Daniel Piker on April 28, 2015 at 11:13am

It sounds like you are interested in some version of mean curvature flow.

Various approximations of this are commonly used in mesh smoothing. It's a big subject, as a lot of work has been done on this for applications such as improving noisy 3d scan data, but these are some of the classic references:

Taubin: A Signal Processing Approach To Fair Surface Design

Desbrun et al: Implicit Fairing of Irregular Meshes using Diffusion and Curvature ...

Though if you are just after some simple way of reducing sharp points, the basic 'umbrella operator' Laplacian smoothing may suffice:

As you rightly deduced, you can compare this vector to the standard mesh normal at that point to distinguish between convex and concave regions.

Permalink Reply by Nik Willmore on April 28, 2015 at 11:26pm

Great homework assignment, thanks.

Permalink Reply by Nik Willmore on April 29, 2015 at 1:19am

I'm trying a command variation that lack of documentation makes very confusing about the output format which is just a long list:

SHARP_MESH_FACES = []

print "There are: ",MESH.Faces.Count," faces."
print "ANGLE_LIMIT: ",AngleLimitB # 0-180 input slider
SHARP_MESH_FACES = MESH.Faces.GetConnectedFaces(0,math.radians(180-AngleLimitB),False)
print "Size of sharp mesh face array: ",len(SHARP_MESH_FACES)
print SHARP_MESH_FACES

Results in this sort of thing for a large mesh:

There are: 13394 faces.
ANGLE_LIMIT: 160.0
Size of sharp mesh face array: 46
Array[int]((0, 1, 2, 3, 5, 553, 554, 4, 7162, 7163, 549, 556, 7159, 550, 7158, 7165, 551, 552, 7167, 558, 7160, 7161, 557, 7026, 417, 7166, 560, 7016, 407, 559, 7169, 408, 7017, 7168, 7014, 7164, 405, 555, 8, 406, 7020, 7, 411, 7015, 6, 110))

What is that array?! It's not always even length, so it's not pairs. The command page says about GetConnectedFaces: "Find all connected face indices where adjacent face normals meet the criteria of angleRadians and greaterThanAngle." The first parameter is the face to start on, so I used 0 to start at the beginning. As I move the angle limit slider it does grab more and more hits.

Why doesn't rhinocommon.chm actually tell what commands really output and how to use them in Python? UGH.

Permalink Reply by djordje on April 30, 2015 at 5:55am

Hi Nik,

Rhino 5 uses IronPython (2.7).

IronPython is implementation of Python for the .NET framework.

GetConnectedFaces method returns a .NET array or typed array. You can convert it to a regular python list by using list() built in function:

SHARP_MESH_FACES_List = list(SHARP_MESH_FACES)

RhinoCommon.chm and online pages do not have help material related with particular .NET language issues. For this, take a look at some of ironpython documentation pages.

RSS