Uniform Knot Spacing BSpline Interpolation

Hello all,

I'm trying to interpolate a degree 3 periodic B-spline given 4 points on the curve. The points define the extremities of the roughly circular curve, and the resulting bspline should be symmetric in the vertical central line connecting the top and bottom points , i.e. the XZ plane.

I can achieve the desired shape in Rhino using a 'periodic uniform knot spacing' interpolation, however I can't find how to set this interpolation strategy in OCC and the result is skewed/assymetric (see attached images). I also cannot find how to interpolate a spline to have uniform multiplicity i.e. no repeating knots such that control points are not the start/end points of curve.
 
Can anyone please advise how I can achieve the desired b-spline interpolation from OpenCascade?

​Thanks

Paul

 

Attachments: 
Roman Lygin's picture

Hello Paul,

Can you give a snippet of the code you use to interpolate the curve ?
I recently came across the issue in GeomAPI_Interpolate (see #27112) so this might be related.

Regards,

Roman

Paul Chambers's picture

Thanks Roman,

I am using Python-OCC (though the user forum for that seems quiet), so this snippet is translated roughly but hopefully you get the idea:

TColgp_Array1OfPnt Pnts(1,5);
Pnts.SetValue(1,gp_Pnt(24.838275, 2.81719693646, 1.3169089144));
Pnts.SetValue(2,gp_Pnt(24.838275, 0.0, 4.163813147));
Pnts.SetValue(3,gp_Pnt(24.838275, -2.81719693646, 1.3169089144));
Pnts.SetValue(4,gp_Pnt(24.838275, 0.0, -1.5299419524));
Pnts.SetValue(5,gp_Pnt(24.838275, 2.81719693646, 1.3169089144));     //Repeated Start point

Handle(Geom_BSplineCurve) SPL = GeomAPI_PointsToBSpline(Pnts, 3, 3, GeomAbs_C1).Curve();
SPL.GetObject().SetPeriodic()

I think I need C2 2nd derivative continuity rather than C1, however this also gives a diagonally squashed spline. I also tried GeomAPI_Interpolate:

TColgp_HArray1OfPnt harray(1, 4);
harray.SetValue(1,gp_Pnt(24.838275, 2.81719693646, 1.3169089144));
harray.SetValue(2,gp_Pnt(24.838275, 0.0, 4.163813147));
harray.SetValue(3,gp_Pnt(24.838275, -2.81719693646, 1.3169089144));
harray.SetValue(4,gp_Pnt(24.838275, 0.0, -1.5299419524));

GeomAPI_Interpolate interp(harray.GetHandle(), Standard_True, 0.01);
interp.Perform();
Handle(Geom_BSplineCurve) SPL = interp.Curve();

I'm not sure loading tangents will help, as I think the main issue is that the multiplicity is default so that the first and last pole are the start and end point of the curve, which I can't seem to prevent. Any ideas?

Paul

Roman Lygin's picture

Paul,

I made some further experiments with GeomAPI_Interpolate trying to interpolate a circle using 4 distinct points (+ periodicity flag). 4 points are sampled at PI/2, PI, 3/2PI, 2PI angle.

Results are attached (the original circle is in red):

- no tangents,

- with first tangent provided (as (-1,0,0))

- with all tangents provided.

It is seen that when no tangents are provided the result is totally weird. So I bet the default tangent computed as described #27112 has the greatest impact.

Perhaps the only meaningful work-around is to provide a greater number of sample points (i.e. to reduce chordal and angular deviation from expected result) and to precompute boundary tangents.

GeomAPI_Interpolate interpolates using a cubic spline and provides C2-continuous result when only points or point + boundary tangents are provided. If interim tangents are computed, then respective knots will have multiplicity 2 (thus leading to C1 continuity only).

Hope this helps.

Roman

 

Paul Chambers's picture

Hi Roman,

That looks promising, however I can't reproduce either of the images you show there. In the case of all tangents, are you referring to the GeomAPI_Interpolate Load function using a tangents vector array containing a gp_Vec at each point? Here is my attempt and an image of the output (control points in yellow, sample points in blue, black curve is the fitted curve):

// harray is the Handle(TColgp_HArray1OfPnt) containing the sampled circle points in the YZ plane at PI/4, PI/2, 3/2PI and 2PI

GeomAPI_Interpolate interp(harray, Standard_True, 0.001);

TColgp_Array1OfVec tangents(1, 4);

tangets.SetValue(1, gp_Vec(0, -1, 0));

tangents.SetValue(2, gp_Vec(0, 0, -1));

tangents.SetValue(3, gp_Vec(0, 1, 0));

tangents.SetValue(4, gp_Vec(0, 0, 1));

TColstd_HArray1OfBoolean flags(1,4);

flags.SetValue(1, Standard_True);

flags.SetValue(2, Standard_True);

flags.SetValue(3, Standard_True);

flags.SetValue(4, Standard_True);

interp.Load(tangents, flags.GetHandle());

interp.Perform();

Have I missed something? If you have the sample code, could you please share it here?

Thanks

Paul

Roman Lygin's picture

Hi Paul,

Please see below - for all normals. You can easily tweak it to only use required tangents. Pay attention to set theScale option to false in Load().

int n = 4;
Handle_TColgp_HArray1OfPnt harray = new TColgp_HArray1OfPnt (1, n);
Handle_Geom_Circle aC = new Geom_Circle (gp_Ax2(), 1.);
for (int i = 1; i <= n; ++i) {
    double p = 2. * M_PI / n * i;
    harray->SetValue(i, aC->Value (p));
}

Handle_TColStd_HArray1OfBoolean farray = new TColStd_HArray1OfBoolean (1,n);
farray->Init (Standard_True);

TColgp_Array1OfVec tarray (1, n);
tarray.SetValue (1, gp_Vec(-1,0,0)); //or can rather use aC->D1() to populate tarray in the loop above
tarray.SetValue (2, gp_Vec(0,-1,0));
tarray.SetValue (3, gp_Vec(1,0,0));
tarray.SetValue (4, gp_Vec(0,1,0));

GeomAPI_Interpolate interp(harray, Standard_True /*periodic*/, 0.01);
interp.Load (tarray, farray, Standard_False /*theScale*/);
interp.Perform();

Good luck!

Roman

Paul Chambers's picture

Thats great, it was the Scale=False that I had missed in the GeomAPI_Interpolate Load tangents. It was not clear from the documentation that this was to be used if specifying a unit vector tangent - my curve is now as I had intended.

Thanks Roman!

P.

Roman Lygin's picture

Yep, glad this helped ;-).

Beware - theScaling option does not exist in Geom2dAPI (yet). So we work-around this by creating 3D data (vector and points) with Z=0, interpolate with GeomAPI, project back using GeomAPI::To2d().

Arjan Schouten's picture

Hello Roman,

Maybe it is related to this? http://wiki.mcneel.com/developer/onsuperfluousknot