GeomPlate_BuildPlateSurface Limitations

Hi,

I am trying to fit a surface to an enclosed tubular shape using cross sections and orthogonal guide curves as constraints (black and red respectively in the attached images). The cross sections are roughly elliptical bsplines and the guide curves pass through the y&z extremities of the sections.

Adding the cross sections as constraints to GeomPlate_BuildPlateSurface seems to give what I need if I use only the first few cross sections, but using more than certain number of cross sections causes it to abort calculation - see PlateSurf_ErrorBoundary.png for the point at which things start to break. The algorithm does also breaks if I add the longitudinal guide curves.

Is this a limitation of BuildPlateSurface? Can anyone tell me a better way of fitting this type of enclosed tubular surface from orthogonal curve network in OCC?

For reference, here is some code which may demonstrate what I'm doing ( I took inspiration mostly from http://opencascade.blogspot.co.uk/2010/03/surface-modeling-part6.html​)​:

 

const Standard_Integer aDeg = 3; ​

const Standard_Integer aNbPnts = 15; ​

const Standard_Integer aNbIter = 5; ​

GeomPlate_BuildPlateSurface Builder(aDeg, aNbPnts, aNbIter);  /*I leave the tolerances etc. as default*/

/* Note: if theBoundaries contains longitudinal curves, the builder aborts calculation */

TColStd_ListIteratorOfListOfTransient anIt (theBoundaries);
for (; anIt.More(); anIt.Next(), i++) {

const Handle(Standard_Transient)& aCur = anIt.Value();​
const Handle(GeomAdaptor_HCurve)& aHC =
Handle(GeomAdaptor_HCurve)::DownCast (aCur);
Handle (GeomPlate_CurveConstraint) aConst =
new GeomPlate_CurveConstraint (aHC, 0 /*GeomAbs_G0*/);
Builder.Add (aConst);}

Builder.Perform()

/* I then make an approximate surface with GeomPlateMakeApprox but the error arises before this */

 

Thanks,
Paul

Paul Chambers's picture

Is anyone able to provide some advice on this issue, as it seems I am at the limit of the plate surface functionality and I can't see another way of doing this. I'm looking for something that has the same behaviour as Rhino's "Network Surface".

I would really appreciate any info on this.

Thanks,

Paul

Arjan Schouten's picture

Hi Paul,

Would BRepOffsetAPI_ThruSections help in this case?

Regards,

Arjan

Paul Chambers's picture

Thanks Arjan,

I had tried the ThruSections algorithm before with no joy, however I made some changes to the cross sections and it now gives the desired result. I suppose there is no enforcing of contact with the longitudinal curves, but GeomAbs_C2 continuity seems to do the right thing. 

I notice also that adding vertices to the algorithm (to give a smoothly curved front surface) is order dependent, i.e. I needed to add the front vertex first, then add the cross section wires from front to back. Perhaps this will help someone else who comes across the same issue.

Paul

zejun zuo's picture

Hi ,

I need to receate the surface using the boundary curve and inner points by GeomPlate_BuildPlateSurface.For reference, here is some code which may demonstrate what I'm doing ( I took inspiration mostly from http://opencascade.blogspot.co.uk/2010/03/surface-modeling-part6.html​)​:

However ,when input the constraint points and constraint curve, the     BPSurf.Perform() doesn't work.

I would really appreciate any info on this.

Thanks,

zejun.

   Standard_Real aDeflection = 0.1;
    Standard_Integer i;
   Standard_Integer aIndex = 1, nbNodes = 0;

    //// define two sequence of points
    TColgp_SequenceOfPnt aPoints, aPoints1;

   
    aPoints.Clear();
    vector<D_3DOT> idots;
    LoadData("C:\\Users\\zuo\\Desktop\\3.xyz",idots);
    for(i=0;i<idots.size();i++)
    {
        aPoints.Append(gp_Pnt(idots[i].x,idots[i].y,idots[i].z));
    }

    // remove double points
    nbNodes = aPoints.Length();
    
    for( i = 1; i <= nbNodes; i++)
    {
      gp_Pnt aPi = aPoints(i);
      Standard_Integer j;
      for( j = i + 1; j < nbNodes; j++)
      {
        gp_Pnt aPj = aPoints(j);
        if(!aPi.IsEqual(aPj,0.9))
          aIndex++;
      }
      if(aIndex == j - 1)
        aPoints1.Append(aPi);

      aIndex = i + 1;
    }

    // find max point
    aIndex = 0;
    gp_Pnt aPntMax = aPoints1(1);
    nbNodes = aPoints1.Length();
    for(i = 2; i <= nbNodes; i++)
    {
      if(aPoints1(i).X() > aPntMax.X())
      {
        aIndex = i;
        aPntMax = aPoints1(aIndex);
      } 
    }

    // clear seguence
    aPoints.Clear();

    Standard_Integer nbLeftNodes = nbNodes;

    // ascending sort - fill aPoints with ascending 
    // by X coordinate points from aPoints1
    for(i = 1; i < nbNodes; i++)
    {
      Standard_Real aMin = aPntMax.X();
      aIndex = 1;
      for( Standard_Integer j = 1; j <= nbLeftNodes; j++)
      {
        if(aPoints1(j).X() < aMin)
        {
          aMin = aPoints1(j).X();
          aIndex = j;
        } 
      }
      aPoints.Append(aPoints1(aIndex));
      aPoints1.Remove(aIndex);
      nbLeftNodes = aPoints1.Length();
    }
   
    vector<D_3DOT> dots;
    D_3DOT            dot;

    for( i = 1; i <= aPoints.Length(); i++)
    {
      gp_Pnt aPi = aPoints(i);
      dot.x=aPi.X();
      dot.y=aPi.Y();
      dot.z=aPi.Z();
      dots.push_back(dot);
    }
    SaveData("C:\\Users\\zuo\\Desktop\\bi.xyz",dots);
        //void Add (const Handle(GeomPlate_CurveConstraint)& Cont);
   TColStd_ListOfTransient theBoundaries;
    LoadData("C:\\Users\\zuo\\Desktop\\bian1.xyz",dots);
    RemoveSame(dots);
    gp_Pnt p1,p2;
    for(i=0;i<dots.size()-1;i++)
    {
        p1.SetX(dots[i].x);
        p1.SetY(dots[i].y);
        p1.SetZ(dots[i].z);
        p2.SetX(dots[i+1].x);
        p2.SetY(dots[i+1].y);
        p2.SetZ(dots[i+1].z);

        Handle(Geom_TrimmedCurve)  aSegment1 = GC_MakeSegment(p1,p2);
        Handle(GeomAdaptor_HCurve) aCurve1 = new GeomAdaptor_HCurve(GeomAdaptor_HCurve(aSegment1));
        theBoundaries.Append(aCurve1);
    }

    LoadData("C:\\Users\\zuo\\Desktop\\bian2.xyz",dots);
    RemoveSame(dots);
    for(i=0;i<dots.size()-1;i++)
    {
        Handle(Geom_TrimmedCurve)  aSegment1 = GC_MakeSegment(gp_Pnt(dots[i].x,dots[i].y,dots[i].z),gp_Pnt(dots[i+1].x,dots[i+1].y,dots[i+1].z));
        Handle(GeomAdaptor_HCurve) aCurve1 = new GeomAdaptor_HCurve(GeomAdaptor_HCurve(aSegment1));
        theBoundaries.Append(aCurve1);
    }
    LoadData("C:\\Users\\zuo\\Desktop\\bian3.xyz",dots);
    RemoveSame(dots);
    for(i=0;i<dots.size()-1;i++)
    {
        Handle(Geom_TrimmedCurve)  aSegment1 = GC_MakeSegment(gp_Pnt(dots[i].x,dots[i].y,dots[i].z),gp_Pnt(dots[i+1].x,dots[i+1].y,dots[i+1].z));
        Handle(GeomAdaptor_HCurve) aCurve1 = new GeomAdaptor_HCurve(GeomAdaptor_HCurve(aSegment1));
        theBoundaries.Append(aCurve1);
    }
        LoadData("C:\\Users\\zuo\\Desktop\\bian4.xyz",dots);
        RemoveSame(dots);
    for(i=0;i<dots.size()-1;i++)
    {
        Handle(Geom_TrimmedCurve)  aSegment1 = GC_MakeSegment(gp_Pnt(dots[i].x,dots[i].y,dots[i].z),gp_Pnt(dots[i+1].x,dots[i+1].y,dots[i+1].z));
        Handle(GeomAdaptor_HCurve) aCurve1 = new GeomAdaptor_HCurve(GeomAdaptor_HCurve(aSegment1));
        theBoundaries.Append(aCurve1);
    }

    // define parameters GeomPlate_BuildPlateSurface
    Standard_Integer Degree = 3;
    Standard_Integer NbPtsOnCur = 10;
    Standard_Integer NbIter = 3;
    Standard_Integer Order = 0;
    Standard_Integer MaxSeg = 9;
    Standard_Integer MaxDegree = 5;
    Standard_Real dmax, anApproxTol = 0.001;
    Standard_Real aConstrTol = Precision::Confusion();

    // define object BuildPlateSurface
    GeomPlate_BuildPlateSurface BPSurf(Degree,NbPtsOnCur,NbIter);

    // add point constraints to GeomPlate_BuildPlateSurface object
    nbNodes = aPoints.Length();
 //   for (i = 1; i <= nbNodes; i++)
 //     BPSurf.Add(new GeomPlate_PointConstraint(aPoints(i), Order, aConstrTol));

    const Standard_Real aTol3d        = 1.e-04; 
    const Standard_Real aTol2d        = 1.e-05;  
    const Standard_Real anAngTol    = 1.e-02;  //angular 
    const Standard_Real aCurvTol    = 1.e-01;  //curvature
    TColStd_ListIteratorOfListOfTransient anIt (theBoundaries); 
    if (anIt.More()) 
    { 
        int i = 1;  
        for (; anIt.More(); anIt.Next(), i++)
        {  
            const Handle(Standard_Transient)& aCur = anIt.Value(); 
            if (aCur.IsNull()) 
            {
            } 
            else if (aCur->IsKind (STANDARD_TYPE (Adaptor3d_HCurveOnSurface))) 
            { 
                //G1 constraint  
                const Handle(Adaptor3d_HCurveOnSurface)& aHCOS = Handle(Adaptor3d_HCurveOnSurface)::DownCast (aCur); 
                Handle (GeomPlate_CurveConstraint) aConst =  new GeomPlate_CurveConstraint (aHCOS, 1 /*GeomAbs_G1*/, NbPtsOnCur, aTol3d, anAngTol, aCurvTol); 
                BPSurf.Add (aConst);  
            } 
            else if (aCur->IsKind (STANDARD_TYPE (GeomAdaptor_HCurve))) 
            { 
                //G0 constraint  
                const Handle(GeomAdaptor_HCurve)& aHC =  Handle(GeomAdaptor_HCurve)::DownCast (aCur); 
                Handle (GeomPlate_CurveConstraint) aConst =  new GeomPlate_CurveConstraint (aHC, 0 /*GeomAbs_G0*/, NbPtsOnCur, aTol3d); 
                BPSurf.Add (aConst); 
            } 
            else 
            {  
            } 
        }  
    } 

    BPSurf.Perform();

    // make PlateSurface
    Handle(GeomPlate_Surface) PSurf;
    Handle(Geom_Surface) aSurf;

    if (BPSurf.IsDone())
    {
      PSurf = BPSurf.Surface();

      // define parameter approximation
      dmax = Max(0.01,10*BPSurf.G0Error());

      // make approximation
      GeomPlate_MakeApprox Mapp(PSurf,anApproxTol, MaxSeg,MaxDegree,dmax);
      aSurf = Mapp.Surface();
    }

Qr Qr's picture

OpenCascade comes with MFC sample on gap filling where Plate tool is used. Here is a code excerpt (slightly cleaned) which demonstrates the approach:

//! Constructs TPS (Thin Plate Spline) approximation for the passed points.
//! \param points [in]  point set to build the approximation surface for.
//! \param result [out] approximation surface.
//! \return true in case of success, false -- otherwise.
bool BuildPlateOnPoints(const std::vector<gp_Pnt>& points, Handle(Geom_BSplineSurface)& result)
{
  const int    Degree      = 3;
  const int    NbPtsOnCur  = 10;
  const int    NbIter      = 3;
  const double Tol2d       = 0.00001;
  const double Tol3d       = 0.0001;
  const double TolAng      = 0.01;
  const double TolCurv     = 0.1;
  const bool   Anisotropie = false;

  /* =================================
   *  STAGE 1: build an average plane
   * ================================= */

  gp_Pln plane;
  PlaneOnPoints planeAlgo;
  //
  if ( !planeAlgo.Build(points, plane) )
  {
#if defined COUT_DEBUG
    std::cout << "Error: cannot build average plane" << std::endl;
#endif
    SendLogMessage(LogErr(Normal) << "Cannot build average plane" );
    return false;
  }
  //
  Handle(Geom_Plane) planeSurf = new Geom_Plane(plane);

  /* ======================
   *  STAGE 2: build plate
   * ====================== */

  GeomPlate_BuildPlateSurface plateAlgo(planeSurf,
                                        Degree,
                                        NbPtsOnCur,
                                        NbIter,
                                        Tol2d,
                                        Tol3d,
                                        TolAng,
                                        TolCurv,
                                        Anisotropie);

  // Add pinpoint constraints
  for ( size_t k = 0; k < points.size(); ++k )
  {
    const gp_Pnt& P = points[k];

    Handle(GeomPlate_PointConstraint)
      PC = new GeomPlate_PointConstraint( P, GeomAbs_C0, Precision::Confusion() );
    //
    plateAlgo.Add(PC);
  }

  // Run plate algorithm
  plateAlgo.Perform();
  //
  if ( !plateAlgo.IsDone() )
  {
    SendLogMessage(LogErr(Normal) << "Plating failed");
    return false;
  }
  //
  Handle(GeomPlate_Surface) plateSurf = plateAlgo.Surface();

  /* ====================================
   *  STAGE 3: approximate plate surface
   * ==================================== */

  int    maxNumOfConstraints = 1000;
  double tol                 = Precision::Confusion();
  int    nbP                 = 4;
  int    maxNumOfSegments    = maxNumOfConstraints / 10;
  //
  TColgp_SequenceOfXY  seqXY;
  TColgp_SequenceOfXYZ seqXYZ;

  plateAlgo.Disc2dContour(nbP, seqXY);
  plateAlgo.Disc3dContour(nbP, GeomAbs_C0, seqXYZ);
  //
  GeomPlate_PlateG0Criterion aPlateCriterion(seqXY, seqXYZ, tol);
  GeomPlate_MakeApprox plateApprox(plateSurf, aPlateCriterion, tol, maxNumOfSegments, Degree);
  result = plateApprox.Surface();

  return true;
}

The code starts with average plane construction. However, note that Plate algorithm can construct this average plane internally, so this piece of code can be skipped.