mikedh / trimesh

Python library for loading and using triangular meshes.
https://trimesh.org
MIT License
3k stars 580 forks source link

Is there an issue with the discretize_bspline or am i defining my bspline incorrectly #2143

Closed LyndonAlcock closed 8 months ago

LyndonAlcock commented 8 months ago

Context

I am currently writing a function that is intended to take faces of a shape in freeCAD and convert it into a 2D path in trimesh so that i can obtain its medial axis tranform (MAT) , i will include the code in the appendix of this issue for completion however i have developed a minimal reproducible example for clarity:

Minimum Reproducible Example

I am using Python311 and all librarys are installed correctly I am trying to create a quadratic bspline with a single knot as follows:

p2d = Path2D(
    entities=[
    BSpline(points=np.arange(2),knots=[(38.0, 25.0)],closed=False)
    ], 
    vertices=[(0.0, 30.0), (52.0, 0.0)])

When i try to do p2d.show() it throws this error:

---> [24] p2d.show()
     [25] #np.arcsin
     [26] #circle = entities.Arc(,closed=True)
     [28] degree

File [<path-to-site-packages>/site-packages/trimesh/path/path.py:928], in Path2D.show(self, annotations)
    [926]     self.plot_discrete(show=True, annotations=annotations)
    [927] else:
--> [928]     self.plot_entities(show=True, annotations=annotations)

File [<path-to-site-packages>/site-packages/trimesh/path/path.py:1373), in Path2D.plot_entities(self, show, annotations, color)
   [1371]     continue
   [1372] # otherwise plot the discrete curve
-> [1373] discrete = entity.discrete(self.vertices)
   [1374] # a unique key for entities
   [1375] e_key = entity.__class__.__name__ + str(int(entity.closed))

File [<path-to-site-packages>/site-packages/trimesh/path/entities.py:779], in BSpline.discrete(self, vertices, count, scale)
    [761] def discrete(self, vertices, count=None, scale=1.0):
    [762]     """
...
--> [313]         raise ValueError("0<=der=%d<=k=%d must hold" % (der, k))
    [314]     if ext not in (0, 1, 2, 3):
    [315]         raise ValueError("ext = %s not in (0, 1, 2, 3) " % ext)

ValueError: 0<=der=0<=k=-2 must hold

Stepping through the code I can see that it executes something analagous to the following :

vertices = np.array([(0.0, 30.0), (52.0, 0.0)])
knots=[(38.0, 25.0)]
knots = np.asanyarray(knots,dtype=np.float64)
points = np.arange(2)
points = np.asanyarray(points, dtype=np.int64)
print(points)
control = vertices[points]
control = np.asanyarray(control,dtype=np.float64)

degree = len(knots) - len(control) - 1

So in this specific example the degree is -2 instead of 2, which feels like it is incorrect but i cant tell what would be correct formatting of this code

thanks for any help :) From Lyndon

Appendix

This is the code that I am currently working on that throws an issue


def faceCenterNormal(face):
    center_of_mass = face.CenterOfMass
    uv_center = face.Surface.parameter(center_of_mass)
    return face.normalAt(*uv_center)

upward_face = next(filter(
    lambda face: faceCenterNormal(face) == FreeCAD.Vector(0, 0, 1),
    FreeCAD.activeDocument().Body.Shape.Faces))

for i, edge in enumerate(upward_face.Edges):
    curve = edge.Curve
    # code for GeomLine and Curve
    if curve.TypeId == 'Part::GeomBSplineCurve':
        print("index: ", i)
        knot_number = curve.Degree - 1
        knots = [tuple(curve.getPole(knot_index)[0:2]) for knot_index in np.arange(knot_number)+2]
        print("knots: ", knots)
        points = [tuple(vertex.Point[0:2]) for vertex in edge.Vertexes]
        print("points: ", points)
        closed = len(points) == 1
        print("closed: ", closed)
        for point in points:
                if index_dict.get(point,None) is None:
                    index_dict[point] = current_index
                    vertices.append(point)
                    current_index+=1

        print("point indices: ", [index_dict[point] for point in points])
        entities.append(BSpline(
            points=[index_dict[point] for point in points],
            knots=knots,
            closed=closed
        ))
mikedh commented 8 months ago

Hey, I'd look at what comes in from a DXF file as that is the only format that supports B-splines currently:

In [1]: p = trimesh.load('models/2D/spline.DXF')
DEBUG              load.py:69    loaded <trimesh.Path2D(vertices.shape=(34, 2), len(entities)=10)> in 0.0065s

In [3]: p.entities[4]
Out[3]: <trimesh.path.entities.BSpline at 0x7fcd57f8aec0>

In [4]: p.entities[4].knots
Out[4]: array([0., 0., 0., 1., 1., 1.])

In [5]: p.entities[4].points
Out[5]: array([18, 29, 21])

I have never had any luck manually defining knot vectors, but your mileage may vary haha. This runs and fixes the data type errors (knots are a flat vector). It is kind of suspicious, knots probably need to be normalized to 0.0-1.0, why are they not increasing in order, etc:

    p = trimesh.path.Path2D(
        entities=[trimesh.path.entities.BSpline(
            points=[0,1],knots=[38.0, 38.0, 25.0, 25.0],closed=False)],
    vertices=[(0.0, 30.0), (52.0, 0.0)])
LyndonAlcock commented 8 months ago

Hello, just an update, it turns out this:

In [4]: p.entities[4].knots
Out[4]: array([0., 0., 0., 1., 1., 1.])

is the output of when you ask for knots = curve.KnotSequence in FreeCAD however, i am unsure what this means in relation to the knot weights as this does not change with the varying knot weights, but i believe that the problem goes a lot deeper than I expected as varying the knot weights in FreeCAD does not change the 3D geometry on my laptop,

long story short, I've kinda got it to work but I am practically treating the BSplines as Bezier Curves but I am behind on my dissertation so this will have to do :expressionless: I have managed to get the format transfer to work with this code:

if curve.TypeId == 'Part::GeomBSplineCurve':
        knots = curve.KnotSequence
        print(knots)
        points = [pole[0:2] for pole in curve.getPoles()]

        closed = len(edge.Vertexes) == 1

        for point in points:
                if index_dict.get(point,None) is None:
                    index_dict[point] = current_index
                    vertices.append(point)
                    current_index+=1

        entities.append(BSpline(
            points=[index_dict[point] for point in points],
            knots=knots,
            closed=closed
        ))