Arbitrary law-based curve and surface modeling. Part 2

(Continuation of part1)
 Same approach works for surfaces. In this case GeomConvert_ApproxSurface and Adaptor3d_Surface classes should be used in exact same manner.

Below are a few examples of creating law surfaces.

The first case is an example of creating a variable offset surface, where a surface is defined as follows:
S(u,v) = B(u,v) + Offset(u,v) * N (u,v), where
B(u,v) is a basis surface,
N (u,v) is a unit normal to the basis surface,
Offset (u,v) is a function C + u ^ 2 + v ^ 2, where C is constant.

The two below examples apply such offset laws to plane and sphere respectively. The basis surface B is in red, and the resulting surface S – in green.

The other example below demonstrates surface warping, where a planar face is twisted along one of its directions:

Enjoy! :-)

Arbitrary law-based curve and surface modeling

A recent case with enhancing CAD Exchanger ACIS importer to broaden support of ACIS primitives inspired me to write this post.

Like other modeling kernels, Open CASCADE comes with a finite set of supported types of curves and surfaces. For instance, for curves this includes lines and conic curves (circles, ellipses, parabolas and hyperbolas), B-Splines, Bezier curves, offset curves plus explicitly trimmed curves. OCC supports parametrics definition where 3D coordinate (x,y,z) is evaluated via parameter t, for surface it is calculated from a pair (u,v).

For instance, for line it is:
C(t) = O + Dt, where O is an origin and D is a unit vector.

Thus, a line in OCC is parametrized by its length.

ACIS and Parasolid additionally support so called procedural geometries (e.g. intersection curve or rolling ball surface) when there is no explicit formula but each point is still unambiguously calculated from parameters t or (u,v).

However, a limited set of explicitly supported types does not allow to express other possible types of curves and surfaces which could be easily represented in parametric definition. For instance, there were a few questions on the OCC forum how to model a helix curve with the help of OCC.

Let us consider how you could that indeed.

Say, a helix is parametrized as follows:

X(t) = O + Rcos(u) X
Y(t) = O + Rsin(u) Y
Z(t) = O + (pitch/(2*PI) * v Z

Where {O ,X, Y, Z} is a local axis system (origin and 3 unit vectors). Pitch – is a helix step along the Z axis. See the below screenshot (taken from ACIS documentation):

Now given this explicit definition (or law) how could you map this to Open CASCADE ?

One option would be to subclass Geom_Curve and implement respective methods (D0(), D1(), ..., Continuity(), etc). That would be fine if you only need a limited time span of such an object and won’t feed it into various OCC modeling algorithms. The OCC classes ShapeExtend_ComplexCurve and _CompositeSurface follow this approach.

However, although you could create an edge on such a curve, you would not be able to save it in a .brep file right away (you would need to enable a vehicle for saving/retrieving user-defined types). Some algorithms may throw an exception on an unrecognized type and so on.

Another option is to do one time approximation with B-Spline and keep the B-Spline representation after that. You would certainly lose original definition and evaluation of a point and derivative would be times more expensive. However you could keep this representation in persistent representation and confidently use it anywhere. By the way ACIS does support helix type and combines both the original definition and its B-Spline approximation.

To approximate with B-Spline, the GeomConvert_ApproxCurve class will do the job.

GeomConvert_ApproxCurve accepts a curve adaptor which essentially implements an Adaptor design pattern and produces a B-Spline curve. You would essentially need to create a subclass of Adaptor3d_Curve and implement respective methods (D0(), D1(), Intervals(), etc). You could possibly mix this approach with the former and create GeomAdaptor_Curve which would accept your Geom_Curve subclass.

Below is an excerpt from CAD Exchanger that approximates a helix with OCC B-Spline:
/*! Uses adaptor classes and invokes GeomConvert_ApproxCurve to approximate with a B-Spline.
    Created B-Spline is polynomial and is of C2-continuity.
    Returns true if the B-Spline has been successfully created and false otherwise.
bool ACISAlgo_Helix::MakeHelix (const ACISGeom_HelixData& theSource,
    Handle_Geom_BSplineCurve& theTarget)
    Handle_ACISAlgo_HHelixCurveAdaptor anAdaptor = new ACISAlgo_HHelixCurveAdaptor (theSource);
    Standard_Real aTol = Precision::Confusion();
    GeomAbs_Shape aContinuity = GeomAbs_C2 /*highest supported continuity*/;
    Standard_Integer aMaxSeg = 10000, /*max number of spans*/
                     aMaxDeg = 9; /*max degree, consistent with settings in Algo*/
    GeomConvert_ApproxCurve anApprox (anAdaptor, aTol,
    if (anApprox.HasResult()) {
        theTarget = anApprox.Curve();
        Base_Debugger& aDebugger = Base_Debugger::GlobalInstance();
        aDebugger.Save (theTarget, "curve");
    return !theTarget.IsNull();
//! Defines data elements that can be reused for both common ACIS 'helix' and ASM 'helix_int_cur'.
struct ACISGeom_HelixData
    ACISGeom_HelixData() :
        myRangeMin (0.),
        myRangeMax (2 * M_PI),
        myScaleFactor (1.)

    __CADEX_DEFINE_PROPERTY(gp_Ax3,Position) //can be right- or left-handed
    __CADEX_DEFINE_PRIMITIVE_PROPERTY(double,Pitch)  //must be >= 0, if = 0 then is planar
    __CADEX_DEFINE_PRIMITIVE_PROPERTY(double,Taper) //if > 0, helix widens along the Z-axis, if 0 - then lies on a cylinder

/*! A few methods in OCC 6.9.0 have been made const.*/
#if OCC_VERSION_HEX < 0x060900

/* \class ACISAlgo_HelixCurveAdaptor
   \brief Defines an adaptor to represent a helix curve.

   Helix data is defined by ACISGeom_HelixData.
   Evaluation is performed in the Evaluator subclass which can either represent a helix lying on a
   cylinder (if taper is 0) or on a cone (if taper is not 0).

   Helix can have distinct radii along X and Y axes, i.e. to have an elliptical section.
class ACISAlgo_HelixCurveAdaptor : public Adaptor3d_Curve

    class Evaluator;

    //! Constructor.
    ACISAlgo_HelixCurveAdaptor (const ACISGeom_HelixData& theData);

    //! Constructor
    ACISAlgo_HelixCurveAdaptor (const std::shared_ptr& theEvaluator,
        Standard_Real theMin,
        Standard_Real theMax);

    virtual Standard_Real FirstParameter()const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual Standard_Real LastParameter() const __CADEX_OVERRIDE_ATTRIBUTE;

    virtual GeomAbs_Shape Continuity() const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual Standard_Integer NbIntervals (const GeomAbs_Shape S) __CADEX_ADAPTOR3D_CURVE_CONST
    virtual void Intervals (TColStd_Array1OfReal& T, const GeomAbs_Shape S) __CADEX_ADAPTOR3D_CURVE_CONST
    virtual Handle(Adaptor3d_HCurve) Trim (const Standard_Real First,
        const Standard_Real Last,
        const Standard_Real Tol) const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual Standard_Boolean IsClosed() const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual Standard_Boolean IsPeriodic() const __CADEX_OVERRIDE_ATTRIBUTE;

    virtual gp_Pnt Value (const Standard_Real U) const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual void D0 (const Standard_Real U, gp_Pnt& P) const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual void D1 (const Standard_Real U, gp_Pnt& P, gp_Vec& V) const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual void D2 (const Standard_Real U, gp_Pnt& P, gp_Vec& V1, gp_Vec& V2) const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual void D3 (const Standard_Real U, gp_Pnt& P, gp_Vec& V1, gp_Vec& V2, gp_Vec& V3) const
    virtual gp_Vec DN (const Standard_Real U, const Standard_Integer N) const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual Standard_Real Resolution (const Standard_Real R3d) const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual GeomAbs_CurveType GetType() const __CADEX_OVERRIDE_ATTRIBUTE;

    std::shared_ptr  myEvaluator;
    Standard_Real               myMin;
    Standard_Real               myMax;
class ACISAlgo_HHelixCurveAdaptor : public Adaptor3d_HCurve

    //! Constructor.
    ACISAlgo_HHelixCurveAdaptor (const ACISGeom_HelixData& theData) : myAdaptor (theData) {}

    //! Constructor.
    ACISAlgo_HHelixCurveAdaptor (const std::shared_ptr& theEvaluator,
        Standard_Real theMin,
        Standard_Real theMax) : myAdaptor (theEvaluator, theMin, theMax) {}

    //! Returns the adaptor as Adaptor3d_Curve.
    /*! Return the internal ACISAlgo_HelixCurveAdaptor object.*/
    virtual const Adaptor3d_Curve& Curve() const __CADEX_OVERRIDE_ATTRIBUTE { return myAdaptor; }

    //! Returns the adaptor as Adaptor3d_Curve.
    /*! Return the internal ACISAlgo_HelixCurveAdaptor object.*/
    virtual Adaptor3d_Curve& GetCurve() __CADEX_OVERRIDE_ATTRIBUTE  { return myAdaptor; }

    ACISAlgo_HelixCurveAdaptor  myAdaptor;


/*! \class ACISAlgo_HelixCurveAdaptor::Evaluator
    \brief Base abstract class to evaluate helix.
class ACISAlgo_HelixCurveAdaptor::Evaluator

    Evaluator (const ACISGeom_HelixData& theData) : myData (theData), myVCoef (1.) {}
    virtual ~Evaluator() {}

    const ACISGeom_HelixData& Data() const { return myData; }

    double VParameter (Standard_Real U) const { return U * myVCoef; }
    virtual void D0 (Standard_Real U, gp_Pnt& P) const = 0;
    virtual void D1 (Standard_Real U, gp_Pnt& P, gp_Vec& V) const = 0;
    virtual void D2 (Standard_Real U, gp_Pnt& P, gp_Vec& V1, gp_Vec& V2) const = 0;
    virtual void D3 (Standard_Real U, gp_Pnt& P, gp_Vec& V1, gp_Vec& V2, gp_Vec& V3) const = 0;
    virtual gp_Vec DN (Standard_Real U, Standard_Integer N) const = 0;

    const gp_XYZ& Loc()  const { return myData.Position().Location().XYZ(); }
    const gp_XYZ& XDir() const { return myData.Position().XDirection().XYZ(); }
    const gp_XYZ& YDir() const { return myData.Position().YDirection().XYZ(); }
    const gp_XYZ& ZDir() const { return myData.Position().Direction().XYZ(); }

    ACISGeom_HelixData  myData;
    double              myVCoef; //coefficient to multiply U to get a V parameter on a respective surface
/*! \class ACISAlgo_HelixCurveAdaptor_CylinderEvaluator
    \brief Evaluates a helix lying on a cylinder.
class ACISAlgo_HelixCurveAdaptor_CylinderEvaluator : public ACISAlgo_HelixCurveAdaptor::Evaluator
    ACISAlgo_HelixCurveAdaptor_CylinderEvaluator (const ACISGeom_HelixData& theData);
    virtual void D0 (Standard_Real U, gp_Pnt& P) const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual void D1 (Standard_Real U, gp_Pnt& P, gp_Vec& V) const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual void D2 (Standard_Real U, gp_Pnt& P, gp_Vec& V1, gp_Vec& V2) const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual void D3 (Standard_Real U, gp_Pnt& P, gp_Vec& V1, gp_Vec& V2, gp_Vec& V3) const __CADEX_OVERRIDE_ATTRIBUTE;
    virtual gp_Vec DN (Standard_Real U, Standard_Integer N) const __CADEX_OVERRIDE_ATTRIBUTE;

ACISAlgo_HelixCurveAdaptor_CylinderEvaluator::ACISAlgo_HelixCurveAdaptor_CylinderEvaluator (
    const ACISGeom_HelixData& theData) :
    ACISAlgo_HelixCurveAdaptor::Evaluator (theData)
    myVCoef = theData.Pitch() * theData.ScaleFactor() / (2 * M_PI);

void ACISAlgo_HelixCurveAdaptor_CylinderEvaluator::D0 (Standard_Real U,
    gp_Pnt& P) const
    Standard_Real v    = VParameter (U);
    Standard_Real Rx   = myData.XRadius();
    Standard_Real Ry   = myData.YRadius();
    Standard_Real sinU = sin (U);
    Standard_Real cosU = cos (U);
    P = Rx * cosU * XDir() + Ry * sinU * YDir() + v * ZDir() + Loc();

void ACISAlgo_HelixCurveAdaptor_CylinderEvaluator::D1 (
    Standard_Real U,
    gp_Pnt& P,
    gp_Vec& V) const
    Standard_Real v    = VParameter (U);
    Standard_Real Rx   = myData.XRadius();
    Standard_Real Ry   = myData.YRadius();
    Standard_Real sinU = sin (U);
    Standard_Real cosU = cos (U);
    P = Rx * cosU * XDir() + Ry * sinU * YDir() + v * ZDir() + Loc();

    Standard_Real k    = myVCoef;
    V = -Rx * sinU * XDir() + Ry * cosU * YDir() + k * ZDir();

void ACISAlgo_HelixCurveAdaptor_CylinderEvaluator::D2 (
    Standard_Real U,
    gp_Pnt& P,
    gp_Vec& V1,
    gp_Vec& V2) const
    Standard_Real v    = VParameter (U);
    Standard_Real Rx   = myData.XRadius();
    Standard_Real Ry   = myData.YRadius();
    Standard_Real sinU = sin (U);
    Standard_Real cosU = cos (U);
    P = Rx * cosU * XDir() + Ry * sinU * YDir() + v * ZDir() + Loc();

    Standard_Real k    = myVCoef;
    V1 = -Rx * sinU * XDir() + Ry * cosU * YDir() + k * ZDir();
    V2 = -Rx * cosU * XDir() - Ry * sinU * YDir();

ACISAlgo_HelixCurveAdaptor::ACISAlgo_HelixCurveAdaptor (const ACISGeom_HelixData& theData) :
    myMin (theData.RangeMin()),
    myMax (theData.RangeMax())
    if (std::fabs (theData.Taper()) < Precision::Confusion()) { //cylinder
        myEvaluator.reset (new ACISAlgo_HelixCurveAdaptor_CylinderEvaluator (theData));

    } else { //cone
        myEvaluator.reset (new ACISAlgo_HelixCurveAdaptor_ConeEvaluator (theData));

/*! Used when trimming*/
ACISAlgo_HelixCurveAdaptor::ACISAlgo_HelixCurveAdaptor (const std::shared_ptr& theEvaluator,
    Standard_Real theMin,
    Standard_Real theMax) :
    myEvaluator (theEvaluator),
    myMin (theMin),
    myMax (theMax)
    __CADEX_ASSERT_INVALID_VALUE(myEvaluator->Data().RangeMin() <= theMin);
    __CADEX_ASSERT_INVALID_VALUE(myEvaluator->Data().RangeMax() >= theMax);

Standard_Real ACISAlgo_HelixCurveAdaptor::FirstParameter() const
    return myMin;

Standard_Real ACISAlgo_HelixCurveAdaptor::LastParameter() const
    return myMax;

GeomAbs_Shape ACISAlgo_HelixCurveAdaptor::Continuity() const
    return GeomAbs_CN;

Standard_Integer ACISAlgo_HelixCurveAdaptor::NbIntervals (const GeomAbs_Shape /*S*/) __CADEX_ADAPTOR3D_CURVE_CONST
    return 1;

void ACISAlgo_HelixCurveAdaptor::Intervals (TColStd_Array1OfReal& T, const GeomAbs_Shape /*S*/) __CADEX_ADAPTOR3D_CURVE_CONST
    T (T.Lower()) = FirstParameter();
    T (T.Upper()) = LastParameter();

Handle(Adaptor3d_HCurve) ACISAlgo_HelixCurveAdaptor::Trim (const Standard_Real First,
    const Standard_Real Last,
    const Standard_Real /*Tol*/) const
    return new ACISAlgo_HHelixCurveAdaptor (myEvaluator, First, Last);

Standard_Boolean ACISAlgo_HelixCurveAdaptor::IsClosed() const
    return Standard_False;

Standard_Boolean ACISAlgo_HelixCurveAdaptor::IsPeriodic() const
    return Standard_False;

gp_Pnt ACISAlgo_HelixCurveAdaptor::Value (const Standard_Real U) const
    gp_Pnt aP;
    D0 (U, aP);
    return aP;

void ACISAlgo_HelixCurveAdaptor::D0 (const Standard_Real U, gp_Pnt& P) const
    myEvaluator->D0 (U, P);

void ACISAlgo_HelixCurveAdaptor::D1 (const Standard_Real U, gp_Pnt& P, gp_Vec& V) const
    myEvaluator->D1 (U, P, V);

void ACISAlgo_HelixCurveAdaptor::D2 (const Standard_Real U, gp_Pnt& P, gp_Vec& V1, gp_Vec& V2) const
    myEvaluator->D2 (U, P, V1, V2);

Standard_Real ACISAlgo_HelixCurveAdaptor::Resolution (const Standard_Real R3d) const
    //see GeomAdaptor_Curve::Resolution()
    const auto& aData = myEvaluator->Data();
    Standard_Real R = std::max (aData.XRadius(), aData.YRadius());
    if (R3d < 2 * R)
        return 2 * ASin (R3d / (2 * R));
        return 2 * M_PI;

GeomAbs_CurveType ACISAlgo_HelixCurveAdaptor::GetType() const
    return GeomAbs_OtherCurve;

Below are two screenshot of approximated helices in DRAW – the first one has different Rx and Ry radii and is lying on a cylindrical surface, the second is of equal Rx and Ry radii lying on a conical surface.

Hope this post will give you some hints on how to approximate arbitrary curves and surfaces using B-Spline approximation techniques.

If you have any interesting examples to share that would certainly be good to know.




Applying Vectorization Techniques for B-Spline Surface Evaluation

Last year I mentored an Intel Summer school project which was dedicated to demonstrating vectorization parallelism with the help of Intel C/C++ compiler. We chose Open CASCADE and its NURBS surface evaluation algorithms as a target.

The outcome was quite nice - up to 16x speedup of the kernel computational functions. If you are interested, feel free to have a look at the paper published at Intel Developer Zone here - https://software.intel.com/en-us/articles/applying-vectorization-techniques-for-b-spline-surface-evaluation.


AIS. Connecting objects.

(This post was written 2+ months ago but got stuck in my drafts folder. Sorry for that)

The Nokia's once famous motto 'Connecting people' has probably gone forever (especially given recent acquisition by Microsoft), but let us grab its idea. This post will be about a nice 'connecting' concept offered by the AIS (Application Interactive Services) package, part of the visualization mechanism of Open CASCADE. Unfortunately, like many other gems scattered across the product, this one sits virtually undocumented and therefore its powerful features are significantly underutilized in OCC-based apps.

This is about AIS_ConnectedInteractive and AIS_MultipleConnected, which allows you to construct derived visual representations from already computed (or to be computed) ones. Imagine an assembly consisting of various sub-assemblies, which in their turn may consist of further sub-assemblies, and so on. Eventually each sub-assembly is broken down into part(s). Each part or sub-assembly are 'instantiated' in a parent assembly, where an instance is a reference to a referred part or assembly plus an attached transformation. Here is an example of an assembly (from STEP test suite as1*.stp):
The next image is an assembly tree where you can see instances of various sub-assemblies.
Then also goes a sub-assembly l-bracket-assembly consisting of instance of a part named l_bracket and 3 instances of nut-bolt-assembly (each with different transformation), where each is a combination of two instantiated parts (bolt and nut).
AIS_ConnectedInteractive and AIS_MultipleConnected offer efficient way to compute visual representations of such complex assemblies.

AIS_ConnectedInteractive essentially implements the proxy pattern (or alike) and contains a reference to another AIS_InteractiveObject and transformation matrix.

Handle_AIS_InteractiveObject aRefObject =
 CreateRepresentation (...); TopLoc_Location aLoc = ...; Handle_AIS_ConnectedInteractive anInstance =  new AIS_ConnectedInteractive; anInstance->Connect (theObject, aLoc);

Note that you don't need to know details of aRefObject representation, you just attach it to anInstance.
Thus, AIS_ConnectedInteractive is well suited to create visual representation of an instance in the assembly hierarchy.

AIS_MultipleConnected follows the composite pattern and allows to combine multiple representations into one. Therefore it is an efficient way to construct visual representation of the (sub-)assembly.

Handle_AIS_MultipleConnectedInteractive anAssebly =
 new AIS_MultipleConnectedInteractive; anAssebly->Connect (aChild1); anAssebly->Connect (aChild2);

Again, this allows to abstract from internal details of each child object.

I use the above approach to construct visual representations of objects in a new data model designed in CAD Exchanger (currently planned to be made part of public API in version 3.0), and so far it works great.

To further abstract the creation of a final representation of the root object (or any interim one selected in a tree browser), there is a factory method that returns a handle to an object of the base type AIS_InteractiveObject.

Handle_AIS_InteractiveObject anObject =
 ModelPrs_InteractiveObjectFactory::Create (...);

Thus, the callers do not know the real type of the object that will be created – it can be either AIS_MultipleConnected, or AIS_ ConnectedInteractive, or AIS_Shape, or any other subclass of AIS_InteractiveObject.

The key benefits of *Connected* are in reusing representations of referred objects, instead of recomputing them. This allows to:
  • Reduce computation time (as creation of part representations takes the greatest time)
  • Reduce memory footprint (you do not have to create extra objects supporting the representations)
  • Reduce code size (you do not have to design extra classes for composite objects, only parts representations are really necessary)
  • Achieve greater flexibility and abstract approach (internal implementation details of referred objects are hidden and can change independently from instances and assemblies) 

To have a quick test, you might want to try DRAWEXE:
pload ALL
wedge w 1 2 3 0.5
vdisplay w
vconnect iw -5 0 0 1 0 0 0 0 1 w #creates AIS_ConnectedInteractive of AIS_Shape
box b 2 0 0 3 1 2
vconnect a -5 0 0 1 0 0 0 0 1 w b #creates AIS_ConnectedInteractive of AIS_MultipleConnected

At last, some limitations to be aware of:
  • The referred objects must belong to the same AIS_InteractiveContext. This is unfortunate but is not specific to *Connected*.
  • Setting attributes to a referred object (e.g. part) affects the referring representation (which is sort of expected). Setting an attribute to a referring object seems to have no effect.
  • Some glitches with selection (as reproduced in DRAW). 

Currently these are not affecting CAD Exchanger development, so I did not investigate these in greater details.

Anyway, hopefully these hints will be helpful for those dealing with complex data structures and GUI. If you already worked with the *Connected* and have some experience to share, it would be great to hear your comments.