BelfrySCAD / BOSL2

The Belfry OpenScad Library, v2.0. An OpenSCAD library of shapes, masks, and manipulators to make working with OpenSCAD easier. BETA
https://github.com/BelfrySCAD/BOSL2/wiki
BSD 2-Clause "Simplified" License
940 stars 111 forks source link

Support of BSpline and BSpline Surface #1268

Open KeithSloan opened 1 year ago

KeithSloan commented 1 year ago

I would like the library to add support for BSpline and BSpline Surfaces Describe the solution you'd like The library currently supports a Bezier Curve and Surface and I would like to see support added for BSpline Curves and Surfaces Describe alternatives you've considered I have a specific application in mind and it involves reading geometries create with OpenCASCADE and creating the appropriate scad source that calls a new functions in a BOSL2 library to create the equivalent geometries with OpenSCAD.

Documentation on OCCT BSpline structure is available at

adrianVmariano commented 1 year ago

@revarbat Didn't you write some bspline code? This is something we were going to defer to post BOSL2 release, I think?

revarbat commented 1 year ago

I have some code for B-Splines that uses a fixed set of knots whose values I derived from code I saw on the net. I don't know why they work. This is called voodoo programming and is bad.

To be brutally honest, I still don't understand what knots are, and how they affect the curve. I have tried researching this several times, but I've not found a site that explains it clearly to me.

If you want to take a crack at it, @adrianvmariano, go for it.

KeithSloan commented 1 year ago

Knots control the curves.

When I access the OpenCASCADE geometry the info on the object has details of the number of knots and their properties. What I am after is being able to create the same objects with OpenSCAD code.

adrianVmariano commented 1 year ago

The statement "knots control the curves" is so vauge as to be without meaning.

For a proper implementation we need to really understand the math of B splines and whether it's possible to implement them in an efficient matrix representation like I did for Beziers. If B-splines are actually useful (I couldn't ever tell that) then we should do it at some point, but like I said above, it's probably a couple months or more long project, so we should defer it to after the BOSL2 stable release.

KeithSloan commented 1 year ago

"but like I said above, it's probably a couple months or more long project, so we should defer it to after the BOSL2 stable release."

@revarbat Please could you share the code you found and the creation with openscad so that I could create a test library and use in the mean time keith[at]sloan-home.co.uk

KeithSloan commented 1 year ago

Have you looked at B-Spline

A Bezier is just a single bezier line, B-Spline is in effect made of several Bezier curves, where the control points effect just part of the line. If you adjust a control point on a Bezier it effects the whole line.

adrianVmariano commented 1 year ago

If you just want a sequence of beziers that is already supported in beziers.scad as a "bezier path".

A B-spline is a generalization of a bezier, not a sequence of bezier curves. A nth order b-spline is like an nth order bezier, still a degree n polynomial. Revar is the one with the b-spline code.

KeithSloan commented 1 year ago

"A B-spline is a generalization of a bezier, not a sequence of bezier curves"

To quote from the reference

"Hence, the knot points divide a B-spline curve into curve segments, each of which is defined on a knot span. We shall show that these curve segments are all Bézier curve of degree p on the curve subdivision page."

So a B-spline is made up of a number of curve segments that are all Beziers with the same degree.

The point is you can use control points to change sections of the curve, rather than change a control point that effects the whole curve.

What I would like to do is update the FreeCAD OpenSCAD exporter, so that people can design their curves interactively with say the FreeCAD Curves workbench, then be able to export the design to scad code.

adrianVmariano commented 1 year ago

Yes, but for the purpose of rendering curves and surfaces, the decomposition of a B-spline into beziers seems to be of theoretical interest rather than practical utility. That decomposition leads to an explosion in control points. An important consideration in the implementation is making it fast so that you can create surfaces without waiting all day. For beziers the efficient matrix implementation was I think 20x faster than the de casteljau subdivision method. A solid understanding of all of this is part of the reason that it's not a 10 minute project.

revarbat commented 1 year ago

Here's what code I do have, which was more or less a translation of some other code I had found in another programming language. I no longer remember where to find the original code. Some of the variable names were best guesses on my part.

include <BOSL2/std.scad>

function bspline(path, closed=false, splinesteps=8) =
    assert(is_path(path) && len(path)>3)
    assert(is_finite(splinesteps))
    assert(floor(ln(splinesteps)/ln(2)) == (ln(splinesteps)/ln(2)), "splinesteps must be an integer power of 2.")
    let(
        lev = max(0, ceil(ln(splinesteps) / ln(2)) ),
        path = closed? path :
            concat([path[0]], path, [last(path)])
    )
    _bspline_recurse(path, closed=closed, lev=lev);

function _bspline_recurse(path, closed=false, lev) =
    lev == 0 ? path :
    let(
        endknots = [
            [1.0, 0.0,  0.0,  0.0],
            [0.5, 0.5,  0.0,  0.0],
            [0.0, 0.75, 0.25, 0.0],
            [0.0, 0.1875, 0.6875, 0.125],
        ],
        midknots = [
            [0.5,   0.5,  0.0  ],
            [0.125, 0.75, 0.125],
        ],
        plen = len(path)
    )
    closed
      ? _bspline_recurse([for(i=idx(path), j=[0,1]) midknots[j] * select(path,i,i+2)], closed, lev=lev-1)
      : _bspline_recurse(
            [
                for(i=[0:1:min(3,plen-3)])   endknots[i] * select(path,0,3),
                for(i=[3:1:plen-4], j=[0,1]) midknots[j] * select(path,i-1,i+1),
                midknots[0] * select(path,-4,-2),
                for(i=[min(3,plen-2):-1:0])  endknots[i] * select(path,[-1:-1:-4]),
            ],
            closed,
            lev=lev-1
        );

function bspline_patch(patch, splinesteps=8, col_wrap=false, row_wrap=false) =
    let(
        bswall1 = [for (i = idx(patch[0])) bspline([for (row=patch) row[i]], closed=col_wrap, splinesteps=splinesteps)],
        bswall2 = [for (i = idx(bswall1[0])) bspline([for (row=bswall1) row[i]], closed=row_wrap, splinesteps=splinesteps)]
    )
    vnf_vertex_array(bswall2, col_wrap=col_wrap, row_wrap=row_wrap);

//p = star(r=50,n=5,step=2); closed=true;
//color("#f00") stroke(p, closed=closed, width=0.5, joints="dot");
//color("#00f") stroke(bspline(p, closed=closed), width=0.5, joints="dot");

rsteps = 96;
csteps = 48;
ssteps = 4;

patch = [
    for (u = [0:1/rsteps:1]) [
        for (v=[0:1/csteps:1])
        zrot(u*360, p=[50,0,0] + (20+sin(v*360*6)*sin(u*360*12)*4)*[cos(360*v),0,sin(360*v)])
    ]
];

//color("blue") move_copies(flatten(patch)) cube(1,center=true);
vnf = vnf_vertex_array(patch, style="quincunx", col_wrap=true, row_wrap=true);
vnfb = bspline_patch(patch, splinesteps=ssteps, col_wrap=true, row_wrap=true);
//vnf_polyhedron(vnf);
vnf_polyhedron(vnfb);
revarbat commented 1 year ago

The reason for the specific numbers in the endknots and midknots tables is a mystery to me.

KeithSloan commented 1 year ago

Many Thanks

adrianVmariano commented 12 months ago

I happened to read a bit about B-splines today. I have a question for Keith: what type of B-spline is important for your application? It appears that there are three basic categories, the uniform B-spline, the clamped (also "open") B-spline and the general B-spline. For the uniform version, it appears that it will be possible to decompose it into pieces and compute the pieces separately using the matrix approach. For the general case it appears that no simplification is likely and the De Boors recursive algorithm will be needed. That's probably going to be at least 20x slower than the current bezier code. The clamped case will be intermediate. I think it can be decomposed and the ends treated as general, but then the middle section handled using the fast matrix method.

Revar, it seems that a B-spline of order k can have as many control points as you want, so it acts more like a bezier path than a bezier. The B-spline can be understood as an average of a sequence of bezier curves, and the knots tell the proportionality of the transition from one curve to the next. The cases that appears to be of interest are as noted above, the uniform case first, where the knots are 0,1,2,3,4,...,N. The number of knots is determined by the number of control points and the degree, and scaling and translation of the knots values doesn't matter. One observation: if you want to create a closed curve for a bezier of order (degree?) k then you need to match k control points at the end. For generic B-splines the main utility is doubling up the knots. I understand that this enables you to force the curve to pass through a control point, and also to reduce the continuity at a point, even down to creating a corner, or discontinuity. If you want to create an open curve that starts and end at a control point you need to create a knot pattern like 0,0,0,0,1,2,3,4,4,4,4 where the number of repeats at the ends depends on the degree, so at the ends it's a general case but the middle is uniform. It may be possible to hard code the matrices for the particular end cases that arise up to low order. Something like 4 matrices for order 4, and so on. I don't know what the typical order is for b-splines. If people are doing cubic, or higher or lower degree. (I guess the code above is cubic, since everything seems to be in groups of four.) I have a paper with a recursive formula for computing these matrices. It's hard to imagine it will make sense to use that formula to compute the matrix for each spline point calculation.

I suppose in principle the same scheme might be adapted to the generic case where you have repeated knots. There would be more cases.

I do wonder if there's some input scheme for denoting repeated knots (which maybe also require repeated control points? Not sure). The structured cases of the clamped curve and the closed uniform curve could be handled with flags rather than the user explicitly repeating data. But internal repeated knots at arbitrary places seem harder to manage.

Also there's NURBS which appear to be such a trivial generalization that we may as well include it: to do NURBS you add weights. And this is done by just adding an extra dimension to the data, running the regular b-spline computation, and then dividing by the extra component when you're done.

In any case, it's apparent that there complexity here both in implementation and in interface. And it's possible I overlooked something or misunderstood something.

KeithSloan commented 12 months ago

What I would like to implement is support for export of surfaces and lines from FreeCAD to OpenSCAD scad.

The Alternate_OpenSCAD importer https://github.com/KeithSloan/OpenSCAD_Alt_Import, has despite the name support for export. It handles objects created via the Part Workbench to scad statements and if the Object type is not supported it Meshes the Shape and calls a Mesh to Polyhedron function. The user could achieve the same by exporting an STL and using that in OpenSCAD, see Image 13-10-2023 at 13 01 Which is a screen shot from a FreeCAD issue

Now what I would like to do is process the FreeCAD Shape and where it finds a Spline or Surface export the SCAD code to recreate the shape, so Spline, Surface and other graphic items.

FreeCAD has a number of workbenches for creating curves and surfaces plus a rudimentary importer for OpenNURBS files. ie. export from Rhino etc.

So I like the idea of enabling a facility for people to design their surfaces interactively etc in FreeCAD and then have the ability to export to SCAD code and use however they like in OpenSCAD.

FreeCAD is built on OpenCASCADE and there are some difference between using FreeCAD functions and direct OCC but FreeCAD Shapes are a subset of OCC.

The OCC documentation for Splines and Surfaces are as follows

adrianVmariano commented 12 months ago

From my perspective, the larger context is a distraction. It sounds like you need an implementation for Geom_BSplineCurve? You don't need all the other stuff in that API, do you? Like the stuff for modifying splines by inserting knots, for example?

That implementation is fully general, meaning it includes rational bsplines without any constraings (NURBS). I have a hard time understanding exactly what happens when you specify periodic. Does this automatically repeat the "poles" as well as the knots? The way the document uses bullet lists is sometimes rather confusing (information in one bullet item actually applies to the next one).

There is guidance about the interface from that document, namely that they separate knot values from knot multiplicities. In retrospect it's clear that this is the right way to do things. Then you can simply not give the knots if they are uniform.

I can't tell from that document what matters for your use case: will your splines be uniform? Will they be rational or non-rational? What degrees will occur? Generally speaking, what are the normal b-spline use cases? The question basically is whether the common use cases can be made fast.

KeithSloan commented 12 months ago

Trouble is there are a number of curve, surface creating work benches and it would take a lot of effort to extract the essence of all.

Here are some examples of where I have writen related code, don't know if this helps

        bs = Part.BSplineCurve()
        bs.buildFromPolesMultsKnots(pts, mu, ku, periodic, nc.Degree, weights)
        if mu[0] < (nc.Degree + 1):
            bs.setPeriodic()
        return bs

Where pts = list of control points mu and mv come from knots

        ku, mu = self.getFCKnots(nu.KnotsU)
        kv, mv = self.getFCKnots(nu.KnotsV)
  bs = Part.BSplineSurface()
  bs.buildFromPolesMultsKnots(
        pts,
        mu,
        mv,
        ku,
        kv,
        uperiodic,
        vperiodic,
        nu.Degree(0),
        nu.Degree(1),
        weights,
    )
    if mu[0] < (nu.Degree(0) + 1):
        bs.setUPeriodic()
    if mv[0] < (nu.Degree(1) + 1):
        bs.setVPeriodic()
    return bs
adrianVmariano commented 12 months ago

Those examples look like they are using a (possibly) generic possibly rational B-spline. I assume that a generic bspline calculator would be able to do that. But it provides no guidance on what types of spline are common or important to optimize.

revarbat commented 11 months ago

I get mu and ku, what what is nu?

adrianVmariano commented 11 months ago

I think nu must be a 2-vector holding the degree of the spline in the two coordinate directions. Unlike with a bezier, this cannot be inferred from the control points.

KeithSloan commented 11 months ago

As Adrian stated its degrees.

Source (https://github.com/FreeCAD/FreeCAD/blob/main/src/Mod/Part/App/BSplineCurvePyImp.cpp) has

PyObject* BSplineCurvePy::buildFromPolesMultsKnots(PyObject *args, PyObject *keywds)
{
    static const std::array<const char *, 8> kwlist{"poles", "mults", "knots", "periodic", "degree", "weights",
                                                    "CheckRational", nullptr};