Parameter for adaptive computation of minimal 2D meshing precision added in BRepMesh_IncrementalMesh API.
Corresponding option -adaptive added in DRAW command
All meshing parameters are collected in structure, BRepMesh_FastDiscret::Parameters, which is now used to define and manipulate parameters of the algorithm.
theLSNew = theLS;
return iErr;
}
- //
- BOPCol_ListIteratorOfListOfShape aIt;
+ //
TopExp_Explorer aExp;
//
TopAbs_ShapeEnum aType = theLS.First().ShapeType();
BOPTools_AlgoTools::MakeContainer
((aType == TopAbs_FACE) ? TopAbs_SHELL : TopAbs_WIRE, aShape);
//
- aIt.Initialize(theLS);
- for (; aIt.More(); aIt.Next()) {
+ for (BOPCol_ListIteratorOfListOfShape aIt(theLS); aIt.More(); aIt.Next()) {
const TopoDS_Shape& aS = aIt.Value();
aBB.Add(aShape, aS);
}
aNb = aMG.Extent();
for (i = 1; i <= aNb; ++i) {
const TopoDS_Shape& aSS = aMG(i);
- const TopoDS_Shape& aSNew = anUnify.Generated(aSS);
- if (!aSNew.IsNull() && !aSS.IsSame(aSNew)) {
- myMapGenerated.Bind(aSS, aSNew);
+ const TopoDS_Shape& aSGen = anUnify.Generated(aSS);
+ if (!aSGen.IsNull() && !aSS.IsSame(aSGen)) {
+ myMapGenerated.Bind(aSS, aSGen);
}
}
}
// build all solids from the faces
BOPCol_ListOfShape aLSF;
//
- aIt.Initialize(theLS);
- for (; aIt.More(); aIt.Next()) {
+ for (BOPCol_ListIteratorOfListOfShape aIt(theLS); aIt.More(); aIt.Next()) {
const TopoDS_Shape& aS = aIt.Value();
//
aExp.Init(aS, TopAbs_FACE);
Standard_Integer i, aNb;
//
// map faces and solids
- BOPCol_ListIteratorOfListOfShape aIt(theLSNew);
- for (; aIt.More(); aIt.Next()) {
+ for (BOPCol_ListIteratorOfListOfShape aIt(theLSNew); aIt.More(); aIt.Next()) {
const TopoDS_Shape& aS = aIt.Value();
//
aExp.Init(aS, TopAbs_FACE);
//purpose :
//=======================================================================
BRepMesh_DiscretRoot::BRepMesh_DiscretRoot()
-: myDeflection (0.001),
- myAngle (0.1),
- myIsDone (Standard_False)
+: myIsDone (Standard_False)
{
}
//! Destructor
Standard_EXPORT virtual ~BRepMesh_DiscretRoot();
- //! Setup linear deflection.
- inline void SetDeflection(const Standard_Real theDeflection)
- {
- myDeflection = theDeflection;
- }
-
- //! Returns linear deflection.
- inline Standard_Real Deflection() const
- {
- return myDeflection;
- }
-
- //! Setup angular deflection.
- inline void SetAngle(const Standard_Real theAngle)
- {
- myAngle = theAngle;
- }
-
- //! Returns angular deflection.
- inline Standard_Real Angle() const
- {
- return myAngle;
- }
-
//! Set the shape to triangulate.
inline void SetShape(const TopoDS_Shape& theShape)
{
Standard_EXPORT virtual void init();
- Standard_Real myDeflection;
- Standard_Real myAngle;
TopoDS_Shape myShape;
Standard_Boolean myIsDone;
};
#include <BRepMesh_ShapeTool.hxx>
#include <BRepMesh_Classifier.hxx>
#include <BRepTools.hxx>
+#include <TopExp_Explorer.hxx>
+#include <TopoDS_Wire.hxx>
+#include <TopoDS_Edge.hxx>
+#include <TopoDS.hxx>
+#include <TopoDS_Iterator.hxx>
+#include <BRep_Tool.hxx>
//=======================================================================
myVMax (0.),
myDeltaX (1.),
myDeltaY (1.),
- myStatus (BRepMesh_NoError)
+ myMinStep (-1.),
+ myStatus (BRepMesh_NoError),
+ myAdaptiveMin (Standard_False)
{
init();
}
BRepMesh_FaceAttribute::BRepMesh_FaceAttribute(
const TopoDS_Face& theFace,
const BRepMesh::HDMapOfVertexInteger& theBoundaryVertices,
- const BRepMesh::HDMapOfIntegerPnt& theBoundaryPoints)
+ const BRepMesh::HDMapOfIntegerPnt& theBoundaryPoints,
+ const Standard_Boolean theAdaptiveMin)
: myDefFace (0.),
myUMin (0.),
myUMax (0.),
myVMax (0.),
myDeltaX (1.),
myDeltaY (1.),
+ myMinStep (-1.),
myStatus (BRepMesh_NoError),
+ myAdaptiveMin (theAdaptiveMin),
myBoundaryVertices(theBoundaryVertices),
myBoundaryPoints (theBoundaryPoints),
myFace (theFace)
myFace.Orientation(TopAbs_FORWARD);
BRepTools::UVBounds(myFace, myUMin, myUMax, myVMin, myVMax);
+ if (myAdaptiveMin)
+ {
+ // compute minimal UV distance
+ // between vertices
+
+ myMinStep = RealLast();
+ for (TopExp_Explorer anExp(myFace, TopAbs_WIRE); anExp.More(); anExp.Next())
+ {
+ TopoDS_Wire aWire = TopoDS::Wire(anExp.Current());
+
+ for (TopoDS_Iterator aWireExp(aWire); aWireExp.More(); aWireExp.Next())
+ {
+ TopoDS_Edge anEdge = TopoDS::Edge(aWireExp.Value());
+ if (BRep_Tool::IsClosed(anEdge))
+ continue;
+
+ // Get end points on 2d curve
+ gp_Pnt2d aFirst2d, aLast2d;
+ BRep_Tool::UVPoints(anEdge, myFace, aFirst2d, aLast2d);
+ Standard_Real aDist =aFirst2d.Distance(aLast2d);
+ if (aDist < myMinStep)
+ myMinStep = aDist;
+ }
+ }
+ }
+
BRepAdaptor_Surface aSurfAdaptor(myFace, Standard_False);
mySurface = new BRepAdaptor_HSurface(aSurfAdaptor);
}
const Standard_Real theLastParam) const
{
const Standard_Real aDeflectionUV = 1.e-05;
- return Max(Precision::PConfusion(), (theLastParam - theFirstParam) * aDeflectionUV);
+ Standard_Real aPreci = (theLastParam - theFirstParam) * aDeflectionUV;
+ if(myAdaptiveMin && myMinStep < aPreci)
+ aPreci = myMinStep;
+
+ return Max(Precision::PConfusion(), aPreci);
}
//=======================================================================
{
public:
- //! Default constructor.
- Standard_EXPORT BRepMesh_FaceAttribute();
-
//! Constructor.
//! @param theFace face the attribute is created for.
//! Used for default initialization. Attribute keeps reference
//! to the source face with forward orientation.
//! @param theBoundaryVertices shared map of shape vertices.
//! @param theBoundaryPoints shared discretization points of shape boundaries.
+ //! @param theAdaptiveMin switches on adaptive computation of minimal parametric
+ //! tolerance (if true).
Standard_EXPORT BRepMesh_FaceAttribute(
const TopoDS_Face& theFace,
const BRepMesh::HDMapOfVertexInteger& theBoundaryVertices,
- const BRepMesh::HDMapOfIntegerPnt& theBoundaryPoints);
+ const BRepMesh::HDMapOfIntegerPnt& theBoundaryPoints,
+ const Standard_Boolean theAdaptiveMin);
//! Destructor.
Standard_EXPORT virtual ~BRepMesh_FaceAttribute();
return mySurface;
}
- //! Sets reference face.
- inline void SetFace(const TopoDS_Face& theFace)
- {
- myFace = theFace;
- }
//! Returns forward oriented face to be used for calculations.
inline const TopoDS_Face& Face() const
private:
+ //! Default constructor.
+ BRepMesh_FaceAttribute();
+
//! Assignment operator.
void operator =(const BRepMesh_FaceAttribute& /*theOther*/)
{
Standard_Real myVMax; //!< Describes maximal value in V domain
Standard_Real myDeltaX;
Standard_Real myDeltaY;
+ Standard_Real myMinStep;
Standard_Integer myStatus;
+ Standard_Boolean myAdaptiveMin;
BRepMesh::HDMapOfVertexInteger myBoundaryVertices;
BRepMesh::HDMapOfIntegerPnt myBoundaryPoints;
#define UVDEFLECTION 1.e-05
-
//=======================================================================
//function : BRepMesh_FastDiscret
//purpose :
//=======================================================================
-BRepMesh_FastDiscret::BRepMesh_FastDiscret(
- const Standard_Real theDefle,
- const Standard_Real theAngl,
- const Bnd_Box& theBox,
- const Standard_Boolean theWithShare,
- const Standard_Boolean theInshape,
- const Standard_Boolean theRelative,
- const Standard_Boolean theShapetrigu,
- const Standard_Boolean isInParallel,
- const Standard_Real theMinSize,
- const Standard_Boolean isInternalVerticesMode,
- const Standard_Boolean isControlSurfaceDeflection)
-: myAngle (theAngl),
- myDeflection (theDefle),
- myWithShare (theWithShare),
- myInParallel (isInParallel),
- myRelative (theRelative),
- myShapetrigu (theShapetrigu),
- myInshape (theInshape),
+BRepMesh_FastDiscret::BRepMesh_FastDiscret( const Bnd_Box& theBox,
+ const BRepMesh_FastDiscret::Parameters& theParams)
+ :
myBoundaryVertices(new BRepMesh::DMapOfVertexInteger),
myBoundaryPoints(new BRepMesh::DMapOfIntegerPnt),
- myMinSize(theMinSize),
- myInternalVerticesMode(isInternalVerticesMode),
- myIsControlSurfaceDeflection(isControlSurfaceDeflection)
-{
- if ( myRelative )
+ myParameters(theParams),
+ myDtotale(0.)
+{
+ if ( myParameters.Relative )
BRepMesh_ShapeTool::BoxMaxDimension(theBox, myDtotale);
}
-//=======================================================================
-//function : BRepMesh_FastDiscret
-//purpose :
-//=======================================================================
-BRepMesh_FastDiscret::BRepMesh_FastDiscret(
- const TopoDS_Shape& theShape,
- const Standard_Real theDefle,
- const Standard_Real theAngl,
- const Bnd_Box& theBox,
- const Standard_Boolean theWithShare,
- const Standard_Boolean theInshape,
- const Standard_Boolean theRelative,
- const Standard_Boolean theShapetrigu,
- const Standard_Boolean isInParallel,
- const Standard_Real theMinSize,
- const Standard_Boolean isInternalVerticesMode,
- const Standard_Boolean isControlSurfaceDeflection)
-: myAngle (theAngl),
- myDeflection (theDefle),
- myWithShare (theWithShare),
- myInParallel (isInParallel),
- myRelative (theRelative),
- myShapetrigu (theShapetrigu),
- myInshape (theInshape),
- myBoundaryVertices(new BRepMesh::DMapOfVertexInteger),
- myBoundaryPoints(new BRepMesh::DMapOfIntegerPnt),
- myMinSize(theMinSize),
- myInternalVerticesMode(isInternalVerticesMode),
- myIsControlSurfaceDeflection(isControlSurfaceDeflection)
-{
- if ( myRelative )
- BRepMesh_ShapeTool::BoxMaxDimension(theBox, myDtotale);
-
- Perform(theShape);
-}
-
//=======================================================================
//function : InitSharedFaces
//purpose :
aFaces.push_back(aFace);
}
- OSD_Parallel::ForEach(aFaces.begin(), aFaces.end(), *this, !myInParallel);
+ OSD_Parallel::ForEach(aFaces.begin(), aFaces.end(), *this, !myParameters.InParallel);
}
{
OCC_CATCH_SIGNALS
- BRepMesh_FastDiscretFace aTool(GetAngle(), myMinSize,
- myInternalVerticesMode, myIsControlSurfaceDeflection);
+ BRepMesh_FastDiscretFace aTool(myParameters.Angle, myParameters.MinSize,
+ myParameters.InternalVerticesMode, myParameters.ControlSurfaceDeflection);
aTool.Perform(anAttribute);
}
catch (Standard_Failure)
if (myAttribute.IsNull())
{
myAttribute = new BRepMesh_FaceAttribute(theFace,
- myBoundaryVertices, myBoundaryPoints);
+ myBoundaryVertices, myBoundaryPoints,myParameters.AdaptiveMin);
myAttributes.Bind(theFace, myAttribute);
}
-
+
BRepMesh::HIMapOfInteger& aVertexEdgeMap = myAttribute->ChangeVertexEdgeMap();
BRepMesh::HDMapOfShapePairOfPolygon& aInternalEdges = myAttribute->ChangeInternalEdges();
resetDataStructure();
- if (!myWithShare)
- {
- myEdges.Clear();
- myBoundaryVertices->Clear();
- myBoundaryPoints->Clear();
- }
-
Standard_Real defedge;
Standard_Integer nbEdge = 0;
- Standard_Real savangle = myAngle;
+ Standard_Real savangle = myParameters.Angle;
Standard_Real cdef;
Standard_Real maxdef = 2.* BRepMesh_ShapeTool::MaxFaceTolerance(theFace);
Standard_Real defface = 0.;
- if (!myRelative)
- defface = Max(myDeflection, maxdef);
+ if (!myParameters.Relative)
+ defface = Max(myParameters.Deflection, maxdef);
NCollection_Sequence<EdgePCurve> aPCurves;
NCollection_Sequence<TopoDS_Edge> aFaceEdges;
const TopoDS_Edge& aEdge = TopoDS::Edge(aEdgeIt.Current());
if (!myMapdefle.IsBound(aEdge))
{
- if (myRelative)
+ if (myParameters.Relative)
{
if (myEdges.IsBound(aEdge))
{
else
{
defedge = BRepMesh_ShapeTool::RelativeEdgeDeflection(
- aEdge, myDeflection, myDtotale, cdef);
+ aEdge, myParameters.Deflection, myDtotale, cdef);
- myAngle = savangle * cdef;
+ myParameters.Angle = savangle * cdef;
}
defface += defedge;
}
else
{
- defedge = myDeflection;
+ defedge = myParameters.Deflection;
}
defedge = Max(maxdef, defedge);
else
{
defedge = myMapdefle(aEdge);
- if ( myRelative )
+ if ( myParameters.Relative )
{
defface += defedge;
defface = Max(maxdef, defface);
aFaceEdges.Append(aEdge);
add(aEdge, aPCurve, defedge);
- myAngle = savangle;
+ myParameters.Angle = savangle;
}
}
return myAttribute->GetStatus();
}
- if ( myRelative )
+ if ( myParameters.Relative )
{
defface = defface / nbEdge;
}
else
{
- defface = myDeflection;
+ defface = myParameters.Deflection;
}
- if ( myWithShare )
- defface = Max(maxdef, defface);
+ defface = Max(maxdef, defface);
TopLoc_Location aLoc;
Handle(Poly_Triangulation) aTriangulation = BRep_Tool::Triangulation(aFace, aLoc);
- if (!myShapetrigu || aTriangulation.IsNull())
+ if ( aTriangulation.IsNull() )
{
Standard_Real xCur, yCur;
Standard_Real maxX, minX, maxY, minY;
BRepMesh::HClassifier& aClassifier = myAttribute->ChangeClassifier();
BRepMesh_WireChecker aDFaceChecker(aFace, aTol, aInternalEdges,
aVertexEdgeMap, myAttribute->ChangeStructure(),
- myumin, myumax, myvmin, myvmax, myInParallel );
+ myumin, myumax, myvmin, myvmax, myParameters.InParallel );
aDFaceChecker.ReCompute(aClassifier);
BRepMesh_Status aCheckStatus = aDFaceChecker.Status();
Dv = Max(1.0e0 - (defface/r),0.0e0) ;
Standard_Real oldDv = 2.0 * ACos (Dv);
- oldDv = Min(oldDv, myAngle);
+ oldDv = Min(oldDv, myParameters.Angle);
Dv = 0.9*oldDv; //TWOTHIRD * oldDv;
Dv = oldDv;
{
Du = Max(1.0e0 - (defface/ru),0.0e0);
Du = (2.0 * ACos (Du));
- Du = Min(Du, myAngle);
+ Du = Min(Du, myParameters.Angle);
Standard_Real aa = sqrt(Du*Du + oldDv*oldDv);
if (aa < gp::Resolution())
Du = 0.0;
Du = 2.0 * ACos (Du);
- if (Du > GetAngle())
- Du = GetAngle();
+ if (Du > myParameters.Angle)
+ Du = myParameters.Angle;
deltaX = Du / longv;
deltaY = 1.;
if (aEdgeTool.IsNull())
{
aEdgeTool = new BRepMesh_EdgeTessellator(theEdge, myAttribute,
- mySharedFaces, theDefEdge, myAngle, myMinSize);
+ mySharedFaces, theDefEdge, myParameters.Angle, myParameters.MinSize);
}
Standard_Integer ipf, ivf, isvf, ipl, ivl, isvl;
class BRepMesh_FastDiscret : public Standard_Transient
{
public:
+
+
+ //! Structure storing meshing parameters
+ struct Parameters {
+
+ //! Default constructor
+ Parameters()
+ :
+ Angle(0.1),
+ Deflection(0.001),
+ MinSize(Precision::Confusion()),
+ InParallel(Standard_False),
+ Relative(Standard_False),
+ AdaptiveMin(Standard_False),
+ InternalVerticesMode(Standard_True),
+ ControlSurfaceDeflection(Standard_True)
+ {
+ }
+
+ //! Angular deflection
+ Standard_Real Angle;
+
+ //! Deflection
+ Standard_Real Deflection;
+
+ //! Minimal allowed size of mesh element
+ Standard_Real MinSize;
+
+ //! Switches on/off multy thread computation
+ Standard_Boolean InParallel;
+
+ //! Switches on/off relative computation of edge tolerance<br>
+ //! If trur, deflection used for the polygonalisation of each edge will be
+ //! <defle> * Size of Edge. The deflection used for the faces will be the
+ //! maximum deflection of their edges.
+ Standard_Boolean Relative;
+
+ //! Adaptive parametric tolerance flag. <br>
+ //! If this flag is set to true the minimal parametric tolerance
+ //! is computed taking minimal parametric distance between vertices
+ //! into account
+ Standard_Boolean AdaptiveMin;
+
+ //! Mode to take or ont to take internal face vertices into account
+ //! in triangulation process
+ Standard_Boolean InternalVerticesMode;
+
+ //! Prameter to check the deviation of triangulation and interior of
+ //! the face
+ Standard_Boolean ControlSurfaceDeflection;
+ };
+
+public:
+
- Standard_EXPORT BRepMesh_FastDiscret(
- const Standard_Real defle,
- const Standard_Real angle,
- const Bnd_Box& B,
- const Standard_Boolean withShare = Standard_True,
- const Standard_Boolean inshape = Standard_False,
- const Standard_Boolean relative = Standard_False,
- const Standard_Boolean shapetrigu = Standard_False,
- const Standard_Boolean isInParallel = Standard_False,
- const Standard_Real theMinSize = Precision::Confusion(),
- const Standard_Boolean isInternalVerticesMode = Standard_True,
- const Standard_Boolean isControlSurfaceDeflection = Standard_True);
-
- //! if the boolean <relative> is True, the <br>
- //! deflection used for the polygonalisation of <br>
- //! each edge will be <defle> * Size of Edge. <br>
- //! the deflection used for the faces will be the maximum <br>
- //! deflection of their edges. <br>
- //! <br>
- //! if <shapetrigu> is True, the triangulation, if exists <br>
- //! with a correct deflection, can be used to re-triangulate <br>
- //! the shape. <br>
- //! <br>
- //! if <inshape> is True, the calculated <br>
- //! triangulation will be stored in the shape. <br>
- Standard_EXPORT BRepMesh_FastDiscret(
- const TopoDS_Shape& shape,
- const Standard_Real defle,
- const Standard_Real angle,
- const Bnd_Box& B,
- const Standard_Boolean withShare = Standard_True,
- const Standard_Boolean inshape = Standard_False,
- const Standard_Boolean relative = Standard_False,
- const Standard_Boolean shapetrigu = Standard_False,
- const Standard_Boolean isInParallel = Standard_False,
- const Standard_Real theMinSize = Precision::Confusion(),
- const Standard_Boolean isInternalVerticesMode = Standard_True,
- const Standard_Boolean isControlSurfaceDeflection = Standard_True);
+ //! Constructor.
+ //! Sets the meshing parameters and updates
+ //! relative defletion according to bounding box
+ //! @param B - bounding box encompasing shape
+ //! @param theParams - meshing algo parameters
+ Standard_EXPORT BRepMesh_FastDiscret (const Bnd_Box& B,
+ const Parameters& theParams);
//! Build triangulation on the whole shape.
Standard_EXPORT void Perform(const TopoDS_Shape& shape);
//! parallel threads.
Standard_EXPORT void Process(const TopoDS_Face& face) const;
- void operator ()(const TopoDS_Face& face) const
+ void operator () (const TopoDS_Face& face) const
{
Process(face);
}
- //! Request algorithm to launch in multiple threads <br>
- //! to improve performance (should be supported by plugin). <br>
- inline void SetParallel(const Standard_Boolean theInParallel)
+ //! Returns parameters of meshing
+ inline const Parameters& MeshParameters() const
{
- myInParallel = theInParallel;
- }
-
- //! Returns the multi-threading usage flag. <br>
- inline Standard_Boolean IsParallel() const
- {
- return myInParallel;
- }
-
- //! returns the deflection value. <br>
- inline Standard_Real GetDeflection() const
- {
- return myDeflection;
+ return myParameters;
}
- //! returns the deflection value. <br>
- inline Standard_Real GetAngle() const
+ //! Returns modificable mesh parameters
+ inline Parameters& ChangeMeshParameters()
{
- return myAngle;
+ return myParameters;
}
+
- inline Standard_Boolean WithShare() const
- {
- return myWithShare;
- }
-
- inline Standard_Boolean InShape() const
- {
- return myInshape;
- }
-
- inline Standard_Boolean ShapeTrigu() const
- {
- return myShapetrigu;
- }
-
Standard_EXPORT void InitSharedFaces(const TopoDS_Shape& theShape);
inline const TopTools_IndexedDataMapOfShapeListOfShape& SharedFaces() const
private:
TopoDS_Face myFace;
- Standard_Real myAngle;
- Standard_Real myDeflection;
- Standard_Real myDtotale;
- Standard_Boolean myWithShare;
- Standard_Boolean myInParallel;
+
BRepMesh::DMapOfShapePairOfPolygon myEdges;
BRepMesh::DMapOfFaceAttribute myAttributes;
- Standard_Boolean myRelative;
- Standard_Boolean myShapetrigu;
- Standard_Boolean myInshape;
TopTools_DataMapOfShapeReal myMapdefle;
// Data shared for whole shape
Handle(BRepMesh_FaceAttribute) myAttribute;
TopTools_IndexedDataMapOfShapeListOfShape mySharedFaces;
- Standard_Real myMinSize;
- Standard_Boolean myInternalVerticesMode;
- Standard_Boolean myIsControlSurfaceDeflection;
+ Parameters myParameters;
+
+ Standard_Real myDtotale;
};
DEFINE_STANDARD_HANDLE(BRepMesh_FastDiscret, Standard_Transient)
//purpose :
//=======================================================================
BRepMesh_IncrementalMesh::BRepMesh_IncrementalMesh()
-: myRelative (Standard_False),
- myInParallel(Standard_False),
- myMinSize (Precision::Confusion()),
- myInternalVerticesMode(Standard_True),
- myIsControlSurfaceDeflection(Standard_True)
+: myMaxShapeSize(0.),
+ myModified(Standard_False),
+ myStatus(0)
{
}
//function : Constructor
//purpose :
//=======================================================================
-BRepMesh_IncrementalMesh::BRepMesh_IncrementalMesh(
- const TopoDS_Shape& theShape,
- const Standard_Real theLinDeflection,
- const Standard_Boolean isRelative,
- const Standard_Real theAngDeflection,
- const Standard_Boolean isInParallel)
- : myRelative (isRelative),
- myInParallel(isInParallel),
- myMinSize (Precision::Confusion()),
- myInternalVerticesMode(Standard_True),
- myIsControlSurfaceDeflection(Standard_True)
+BRepMesh_IncrementalMesh::BRepMesh_IncrementalMesh( const TopoDS_Shape& theShape,
+ const Standard_Real theLinDeflection,
+ const Standard_Boolean isRelative,
+ const Standard_Real theAngDeflection,
+ const Standard_Boolean isInParallel,
+ const Standard_Boolean adaptiveMin)
+: myMaxShapeSize(0.),
+ myModified(Standard_False),
+ myStatus(0)
+{
+ myParameters.Deflection = theLinDeflection;
+ myParameters.Relative = isRelative;
+ myParameters.Angle = theAngDeflection;
+ myParameters.InParallel = isInParallel;
+ myParameters.AdaptiveMin = adaptiveMin;
+
+ myShape = theShape;
+ Perform();
+}
+
+//=======================================================================
+//function : Constructor
+//purpose :
+//=======================================================================
+BRepMesh_IncrementalMesh::BRepMesh_IncrementalMesh(const TopoDS_Shape& theShape,
+ const BRepMesh_FastDiscret::Parameters& theParameters)
+ : myParameters(theParameters)
{
- myDeflection = theLinDeflection;
- myAngle = theAngDeflection;
myShape = theShape;
Perform();
BRepMesh_ShapeTool::BoxMaxDimension(aBox, myMaxShapeSize);
- myMesh = new BRepMesh_FastDiscret(myDeflection,
- myAngle, aBox, Standard_True, Standard_True,
- myRelative, Standard_True, myInParallel, myMinSize,
- myInternalVerticesMode, myIsControlSurfaceDeflection);
-
+ myMesh = new BRepMesh_FastDiscret (aBox, myParameters);
+
myMesh->InitSharedFaces(myShape);
}
update(aFaceIt.Value());
// Mesh faces
- OSD_Parallel::ForEach(myFaces.begin(), myFaces.end(), *myMesh, !myInParallel);
+ OSD_Parallel::ForEach(myFaces.begin(), myFaces.end(), *myMesh, !myParameters.InParallel);
commit();
clear();
BRepAdaptor_Curve aCurve(aEdge);
GCPnts_TangentialDeflection aDiscret(aCurve, aCurve.FirstParameter(),
- aCurve.LastParameter(), myAngle, aEdgeDeflection, 2,
- Precision::PConfusion(), myMinSize);
+ aCurve.LastParameter(), myParameters.Angle, aEdgeDeflection, 2,
+ Precision::PConfusion(), myParameters.MinSize);
Standard_Integer aNodesNb = aDiscret.NbPoints();
TColgp_Array1OfPnt aNodes (1, aNodesNb);
}
aPoly3D = new Poly_Polygon3D(aNodes, aUVNodes);
- aPoly3D->Deflection(myDeflection);
+ aPoly3D->Deflection(myParameters.Deflection);
BRep_Builder aBuilder;
aBuilder.UpdateEdge(aEdge, aPoly3D);
return myEdgeDeflection(theEdge);
Standard_Real aEdgeDeflection;
- if (myRelative)
+ if ( myParameters.Relative )
{
Standard_Real aScale;
aEdgeDeflection = BRepMesh_ShapeTool::RelativeEdgeDeflection(theEdge,
- myDeflection, myMaxShapeSize, aScale);
+ myParameters.Deflection, myMaxShapeSize, aScale);
}
else
- aEdgeDeflection = myDeflection;
+ aEdgeDeflection = myParameters.Deflection;
myEdgeDeflection.Bind(theEdge, aEdgeDeflection);
return aEdgeDeflection;
Standard_Real BRepMesh_IncrementalMesh::faceDeflection(
const TopoDS_Face& theFace)
{
- if (!myRelative)
- return myDeflection;
+ if ( !myParameters.Relative )
+ return myParameters.Deflection;
Standard_Integer aEdgesNb = 0;
Standard_Real aFaceDeflection = 0.;
aFaceDeflection += edgeDeflection(aEdge);
}
- return (aEdgesNb == 0) ? myDeflection : (aFaceDeflection / aEdgesNb);
+ return (aEdgesNb == 0) ? myParameters.Deflection : (aFaceDeflection / aEdgesNb);
}
//=======================================================================
BRepMesh_DiscretRoot* &theAlgo)
{
BRepMesh_IncrementalMesh* anAlgo = new BRepMesh_IncrementalMesh();
- anAlgo->SetDeflection(theDeflection);
- anAlgo->SetAngle (theAngle);
+ anAlgo->ChangeParameters().Deflection = theDeflection;
+ anAlgo->ChangeParameters().Angle = theAngle;
+ anAlgo->ChangeParameters().InParallel = IS_IN_PARALLEL;
anAlgo->SetShape (theShape);
- anAlgo->SetParallel (IS_IN_PARALLEL);
theAlgo = anAlgo;
return 0; // no error
}
Standard_EXPORT BRepMesh_IncrementalMesh(
const TopoDS_Shape& theShape,
const Standard_Real theLinDeflection,
- const Standard_Boolean isRelative = Standard_False,
+ const Standard_Boolean isRelative = Standard_False,
const Standard_Real theAngDeflection = 0.5,
- const Standard_Boolean isInParallel = Standard_False);
+ const Standard_Boolean isInParallel = Standard_False,
+ const Standard_Boolean adaptiveMin = Standard_False);
+
+ //! Constructor.
+ //! Automatically calls method Perform.
+ //! @param theShape shape to be meshed.
+ //! @param theParameters - parameters of meshing
+ Standard_EXPORT BRepMesh_IncrementalMesh (const TopoDS_Shape& theShape,
+ const BRepMesh_FastDiscret::Parameters& theParameters);
//! Performs meshing ot the shape.
Standard_EXPORT virtual void Perform();
public: //! @name accessing to parameters.
- //! Enables using relative deflection.
- //! @param isRelative if TRUE deflection used for discretization of
- //! each edge will be <theLinDeflection> * <size of edge>. Deflection
- //! used for the faces will be the maximum deflection of their edges.
- inline void SetRelative(const Standard_Boolean isRelative)
+ //! Returns meshing parameters
+ inline const BRepMesh_FastDiscret::Parameters& Parameters() const
{
- myRelative = isRelative;
+ return myParameters;
}
-
- //! Returns relative deflection flag.
- inline Standard_Boolean IsRelative() const
+
+ //! Returns modifiable meshing parameters
+ inline BRepMesh_FastDiscret::Parameters& ChangeParameters()
{
- return myRelative;
+ return myParameters;
}
//! Returns modified flag.
return myStatus;
}
- //! Request algorithm to launch in multiple threads to improve performance.
- inline void SetParallel(const Standard_Boolean isInParallel)
- {
- myInParallel = isInParallel;
- }
-
- //! Returns the multi-threading usage flag.
- inline Standard_Boolean IsParallel() const
- {
- return myInParallel;
- }
-
- //! Sets min size parameter.
- inline void SetMinSize(const Standard_Real theMinSize)
- {
- myMinSize = Max(theMinSize, Precision::Confusion());
- }
-
- //! Returns min size parameter.
- inline Standard_Real GetMinSize() const
- {
- return myMinSize;
- }
-
- //! Enables/disables internal vertices mode (enabled by default).
- inline void SetInternalVerticesMode(const Standard_Boolean isEnabled)
- {
- myInternalVerticesMode = isEnabled;
- }
-
- //! Returns flag indicating is internal vertices mode enabled/disabled.
- inline Standard_Boolean IsInternalVerticesMode() const
- {
- return myInternalVerticesMode;
- }
-
- //! Enables/disables control of deflection of mesh from real surface
- //! (enabled by default).
- inline void SetControlSurfaceDeflection(const Standard_Boolean isEnabled)
- {
- myIsControlSurfaceDeflection = isEnabled;
- }
-
- //! Returns flag indicating is adaptive reconfiguration
- //! of mesh enabled/disabled.
- inline Standard_Boolean IsControlSurfaceDeflection() const
- {
- return myIsControlSurfaceDeflection;
- }
public: //! @name plugin API
protected:
- Standard_Boolean myRelative;
- Standard_Boolean myInParallel;
BRepMesh::DMapOfEdgeListOfTriangulationBool myEdges;
Handle(BRepMesh_FastDiscret) myMesh;
- Standard_Boolean myModified;
TopTools_DataMapOfShapeReal myEdgeDeflection;
+ NCollection_Vector<TopoDS_Face> myFaces;
+
+ BRepMesh_FastDiscret::Parameters myParameters;
+
Standard_Real myMaxShapeSize;
+ Standard_Boolean myModified;
Standard_Integer myStatus;
- NCollection_Vector<TopoDS_Face> myFaces;
- Standard_Real myMinSize;
- Standard_Boolean myInternalVerticesMode;
- Standard_Boolean myIsControlSurfaceDeflection;
};
DEFINE_STANDARD_HANDLE(BRepMesh_IncrementalMesh,BRepMesh_DiscretRoot)
if (!strcmp(dout.GetType(id),"PERS")) focal = dout.Focal(id);
Standard_Real Ang,Def;
HLRBRep::PolyHLRAngleAndDeflection(myAng,Ang,Def);
- BRepMesh_IncrementalMesh MESH(myShape, Def, Standard_True, Ang);
+ BRepMesh_FastDiscret::Parameters aMeshParams;
+ aMeshParams.Relative = Standard_True;
+ aMeshParams.Deflection = Def;
+ aMeshParams.Angle = Ang;
+ BRepMesh_IncrementalMesh MESH(myShape, aMeshParams);
Standard_Boolean recompute = Standard_True;
// find if the view must be recomputed
DBRep_ListIteratorOfListOfHideData it(myHidData);
(enabled by default)\n\
-surf_def_off disables control of deflection of mesh from real\n\
surface (enabled by default)\n\
- -parallel enables parallel execution (switched off by default)\n";
+ -parallel enables parallel execution (switched off by default)\n\
+ -adaptive enables adaptive computation of minimal value in parametric space\n";
return 0;
}
Standard_Boolean isInParallel = Standard_False;
Standard_Boolean isIntVertices = Standard_True;
Standard_Boolean isControlSurDef = Standard_True;
+ Standard_Boolean isAdaptiveMin = Standard_False;
if (nbarg > 3)
{
isIntVertices = Standard_False;
else if (aOpt == "-surf_def_off")
isControlSurDef = Standard_False;
+ else if (aOpt == "-adaptive")
+ isAdaptiveMin = Standard_True;
else if (i < nbarg)
{
Standard_Real aVal = Draw::Atof(argv[i++]);
di << "Incremental Mesh, multi-threading "
<< (isInParallel ? "ON" : "OFF") << "\n";
- BRepMesh_IncrementalMesh aMesher;
- aMesher.SetShape (aShape);
- aMesher.SetDeflection(aLinDeflection);
- aMesher.SetRelative (isRelative);
- aMesher.SetAngle (aAngDeflection);
- aMesher.SetParallel (isInParallel);
- aMesher.SetMinSize (aMinSize);
- aMesher.SetInternalVerticesMode(isIntVertices);
- aMesher.SetControlSurfaceDeflection(isControlSurDef);
- aMesher.Perform();
+ BRepMesh_FastDiscret::Parameters aMeshParams;
+ aMeshParams.Deflection = aLinDeflection;
+ aMeshParams.Angle = aAngDeflection;
+ aMeshParams.Relative = isRelative;
+ aMeshParams.InParallel = isInParallel;
+ aMeshParams.MinSize = aMinSize;
+ aMeshParams.InternalVerticesMode = isIntVertices;
+ aMeshParams.ControlSurfaceDeflection = isControlSurDef;
+ aMeshParams.AdaptiveMin = isAdaptiveMin;
+
+ BRepMesh_IncrementalMesh aMesher (aShape, aMeshParams);
di << "Meshing statuses: ";
Standard_Integer statusFlags = aMesher.GetStatusFlags();
const Standard_Real d = Draw::Atof(argv[2]);
- Standard_Boolean WithShare = Standard_True;
- if (nbarg > 3) WithShare = Draw::Atoi(argv[3]);
-
Bnd_Box B;
BRepBndLib::Add(S,B);
- BRepMesh_FastDiscret MESH(d,0.5,B,WithShare,Standard_True,Standard_False,Standard_True);
+ BRepMesh_FastDiscret::Parameters aParams;
+ aParams.Deflection = d;
+ aParams.Angle = 0.5;
+ BRepMesh_FastDiscret MESH(B,aParams);
//Standard_Integer NbIterations = MESH.NbIterations();
//if (nbarg > 4) NbIterations = Draw::Atoi(argv[4]);
di<<"Starting FastDiscret with :"<<"\n";
di<<" Deflection="<<d<<"\n";
di<<" Angle="<<0.5<<"\n";
- di<<" SharedMode="<< (Standard_Integer) WithShare<<"\n";
- //di<<" NbIterations="<<NbIterations<<"\n";
Handle(Poly_Triangulation) T;
BRep_Builder aBuilder;
theCommands.Add("incmesh","Builds triangular mesh for the shape, run w/o args for help",__FILE__, incrementalmesh, g);
theCommands.Add("MemLeakTest","MemLeakTest",__FILE__, MemLeakTest, g);
- theCommands.Add("fastdiscret","fastdiscret shape deflection [shared [nbiter]]",__FILE__, fastdiscret, g);
+ theCommands.Add("fastdiscret","fastdiscret shape deflection",__FILE__, fastdiscret, g);
theCommands.Add("mesh","mesh result Shape deflection",__FILE__, triangule, g);
theCommands.Add("addshape","addshape meshname Shape [deflection]",__FILE__, addshape, g);
//theCommands.Add("smooth","smooth meshname",__FILE__, smooth, g);
{
myMesher = theMesher;
if (!myMesher.IsNull())
- myDeflection = myMesher->Deflection();
+ myDeflection = myMesher->Parameters().Deflection;
}
//=======================================================================
if (myMesher.IsNull())
{
myMesher = new BRepMesh_IncrementalMesh;
- myMesher->SetDeflection(myDeflection);
- myMesher->SetAngle(0.5);
+ myMesher->ChangeParameters().Deflection = myDeflection;
+ myMesher->ChangeParameters().Angle = 0.5;
}
myMesher->SetShape(theShape);
if(aShape.IsNull()) {di << "OCC369 FAULTY. Entry shape is NULL \n"; return 0;}
// 3. Build mesh
- BRepMesh_IncrementalMesh aMesh(aShape, 0.2, Standard_True, M_PI / 6);
+ BRepMesh_FastDiscret::Parameters aMeshParams;
+ aMeshParams.Relative = Standard_True;
+ aMeshParams.Deflection = 0.2;
+ aMeshParams.Angle = M_PI / 6;
+ BRepMesh_IncrementalMesh aMesh(aShape, aMeshParams);
}
catch (Standard_Failure) {di << "OCC369 Exception \n" ;return 0;}
{
const Standard_Boolean aRel = aDrawer->TypeOfDeflection() == Aspect_TOD_RELATIVE;
Standard_Real aDef = aRel ? aDrawer->HLRDeviationCoefficient() : aDrawer->MaximalChordialDeviation();
- BRepMesh_IncrementalMesh mesh(aShape, aDef, aRel, aDrawer->HLRAngle());
+ BRepMesh_FastDiscret::Parameters aMeshParams;
+ aMeshParams.Relative = aRel;
+ aMeshParams.Angle = aDrawer->HLRAngle();
+ aMeshParams.Deflection = aDef;
+ BRepMesh_IncrementalMesh mesh(aShape, aMeshParams);
}
Handle(HLRBRep_PolyAlgo) hider = new HLRBRep_PolyAlgo(aShape);
Standard_Integer iErr;
//
iErr=0;
- theAlgo=new BRepMesh_IncrementalMesh;
- theAlgo->SetDeflection(theDeflection);
- theAlgo->SetAngle(theAngle);
- theAlgo->SetShape(theShape);
+ BRepMesh_IncrementalMesh* anAlgo = new BRepMesh_IncrementalMesh;
+ anAlgo->ChangeParameters().Deflection = theDeflection;
+ anAlgo->ChangeParameters().Angle = theAngle;
+ anAlgo->SetShape(theShape);
+ theAlgo = anAlgo;
return iErr;
}
--- /dev/null
+puts "========"
+puts "OCC26664"
+puts "========"
+puts ""
+#################################
+# Triangulating a very small polygon fails
+#################################
+
+restore [locate_data_file bug26664_f.brep] a
+
+#1
+vinit
+vsetdispmode 1
+vdisplay a
+vfit
+vdump ${imagedir}/${casename}_1.png
+
+#2
+incmesh a 1e-4
+vdisplay a
+vdump ${imagedir}/${casename}_2.png
+
+#3
+set log2 [incmesh a 1e-4 -adaptive]
+vdisplay a
+vdump ${imagedir}/${casename}_3.png
+
+if { [regexp "NoError" ${log2}] == 1 } {
+ puts "OK : Triangulating a very small polygon is good"
+} else {
+ puts "Error : Triangulating a very small polygon fails"
+}