+ //================================================================================
+ Intersector::Intersector( SMDS_Mesh* mesh, double tol, const std::vector< gp_XYZ >& normals )
+ {
+ myAlgo = new Algo( mesh, tol, normals );
+ }
+ //================================================================================
+ Intersector::~Intersector()
+ {
+ delete myAlgo;
+ }
+ //================================================================================
+ //! compute cut of two faces of the mesh
+ void Intersector::Cut( const SMDS_MeshElement* face1,
+ const SMDS_MeshElement* face2,
+ const int nbCommonNodes )
+ {
+ myAlgo->Cut( face1, face2, nbCommonNodes );
+ }
+ //================================================================================
+ //! store a face cut by a line given by its ends
+ // accompanied by indices of intersected face edges.
+ // Edge index is <0 if a line end is inside the face.
+ void Intersector::Cut( const SMDS_MeshElement* face,
+ SMESH_NodeXYZ& lineEnd1,
+ int edgeIndex1,
+ SMESH_NodeXYZ& lineEnd2,
+ int edgeIndex2 )
+ {
+ myAlgo->Cut( face, lineEnd1, edgeIndex1, lineEnd2, edgeIndex2 );
+ }
+ //================================================================================
+ //! split all face intersected by Cut() methods
+ void Intersector::MakeNewFaces( SMESH_MeshAlgos::TElemIntPairVec& theNew2OldFaces,
+ SMESH_MeshAlgos::TNodeIntPairVec& theNew2OldNodes,
+ const double theSign,
+ const bool theOptimize )
+ {
+ myAlgo->MakeNewFaces( theNew2OldFaces, theNew2OldNodes, theSign, theOptimize );
+ }
+ //================================================================================
+ //! Cut a face by planes, whose normals point to parts to keep
+ bool Intersector::CutByPlanes(const SMDS_MeshElement* theFace,
+ const std::vector< gp_Ax1 > & thePlanes,
+ const double theTol,
+ std::vector< TFace > & theNewFaceConnectivity )
+ {
+ theNewFaceConnectivity.clear();
+
+ // check if theFace is wholly cut off
+ std::vector< SMESH_NodeXYZ > facePoints( theFace->begin_nodes(), theFace->end_nodes() );
+ facePoints.resize( theFace->NbCornerNodes() );
+ for ( size_t iP = 0; iP < thePlanes.size(); ++iP )
+ {
+ size_t nbOut = 0;
+ const gp_Pnt& O = thePlanes[iP].Location();
+ for ( size_t i = 0; i < facePoints.size(); ++i )
+ {
+ gp_Vec Op( O, facePoints[i] );
+ nbOut += ( Op * thePlanes[iP].Direction() <= 0 );
+ }
+ if ( nbOut == facePoints.size() )
+ return true;
+ }
+
+ // copy theFace into a temporary mesh
+ SMDS_Mesh mesh;
+ Bnd_B3d faceBox;
+ std::vector< const SMDS_MeshNode* > faceNodes;
+ faceNodes.resize( facePoints.size() );
+ for ( size_t i = 0; i < facePoints.size(); ++i )
+ {
+ const SMESH_NodeXYZ& n = facePoints[i];
+ faceNodes[i] = mesh.AddNode( n.X(), n.Y(), n.Z() );
+ faceBox.Add( n );
+ }
+ const SMDS_MeshElement* faceToCut = 0;
+ switch ( theFace->NbCornerNodes() )
+ {
+ case 3:
+ faceToCut = mesh.AddFace( faceNodes[0], faceNodes[1], faceNodes[2] );
+ break;
+ case 4:
+ faceToCut = mesh.AddFace( faceNodes[0], faceNodes[1], faceNodes[2], faceNodes[3] );
+ break;
+ default:
+ faceToCut = mesh.AddPolygonalFace( faceNodes );
+ }
+
+ std::vector< gp_XYZ > normals( 2 + thePlanes.size() );
+ SMESH_MeshAlgos::FaceNormal( faceToCut, normals[ faceToCut->GetID() ]);
+
+ // add faces corresponding to thePlanes
+ std::vector< const SMDS_MeshElement* > planeFaces;
+ double faceSize = Sqrt( faceBox.SquareExtent() );
+ gp_XYZ center = 0.5 * ( faceBox.CornerMin() + faceBox.CornerMax() );
+ for ( size_t i = 0; i < thePlanes.size(); ++i )
+ {
+ gp_Ax2 plnAx( thePlanes[i].Location(), thePlanes[i].Direction() );
+ gp_XYZ O = plnAx.Location().XYZ();
+ gp_XYZ X = plnAx.XDirection().XYZ();
+ gp_XYZ Y = plnAx.YDirection().XYZ();
+ gp_XYZ Z = plnAx.Direction().XYZ();
+
+ double dot = ( O - center ) * Z;
+ gp_XYZ o = center + Z * dot; // center projected to a plane
+
+ gp_XYZ p1 = o + X * faceSize * 2;
+ gp_XYZ p2 = o + Y * faceSize * 2;
+ gp_XYZ p3 = o - (X + Y ) * faceSize * 2;
+
+ const SMDS_MeshNode* n1 = mesh.AddNode( p1.X(), p1.Y(), p1.Z() );
+ const SMDS_MeshNode* n2 = mesh.AddNode( p2.X(), p2.Y(), p2.Z() );
+ const SMDS_MeshNode* n3 = mesh.AddNode( p3.X(), p3.Y(), p3.Z() );
+ planeFaces.push_back( mesh.AddFace( n1, n2, n3 ));
+
+ normals[ planeFaces.back()->GetID() ] = thePlanes[i].Direction().XYZ();
+ }
+
+ // cut theFace
+ Algo algo ( &mesh, theTol, normals );
+ for ( size_t i = 0; i < planeFaces.size(); ++i )
+ {
+ algo.Cut( faceToCut, planeFaces[i], 0 );
+ }
+
+ // retrieve a result
+ SMESH_MeshAlgos::TElemIntPairVec new2OldFaces;
+ SMESH_MeshAlgos::TNodeIntPairVec new2OldNodes;
+ TCutFaceMap::const_iterator cutFacesIt= algo.myCutFaces.cbegin();
+ for ( ; cutFacesIt != algo.myCutFaces.cend(); ++cutFacesIt )
+ {
+ const CutFace& cf = *cutFacesIt;
+ if ( cf.myInitFace != faceToCut )
+ continue;
+
+ if ( !cf.IsCut() )
+ {
+ theNewFaceConnectivity.push_back( facePoints );
+ break;
+ }
+
+ // intersect cut lines
+ algo.IntersectNewEdges( cf );
+
+ // form loops of new faces
+ EdgeLoopSet loopSet;
+ cf.MakeLoops( loopSet, normals[ faceToCut->GetID() ]);
+
+ // erase loops that are cut off by thePlanes
+ const double sign = 1;
+ std::vector< EdgePart > cutOffLinks;
+ TLinkMap cutOffCoplanarLinks;
+ cf.CutOffLoops( loopSet, sign, normals, cutOffLinks, cutOffCoplanarLinks );
+
+ for ( size_t iL = 0; iL < loopSet.myNbLoops; ++iL )
+ {
+ EdgeLoop& loop = loopSet.myLoops[ iL ];
+ if ( loop.myLinks.size() > 0 )
+ {
+ facePoints.clear();
+ for ( SMDS_NodeIteratorPtr nIt = loop.nodeIterator(); nIt->more(); )
+ {
+ const SMDS_MeshNode* n = nIt->next();
+ facePoints.push_back( n );
+ int iN = faceToCut->GetNodeIndex( n );
+ if ( iN < 0 )
+ facePoints.back()._node = 0; // an intersection point
+ else
+ facePoints.back()._node = theFace->GetNode( iN );
+ }
+ theNewFaceConnectivity.push_back( facePoints );
+ }
+ }
+ break;
+ }
+
+ return theNewFaceConnectivity.empty();
+ }
+
+} // namespace SMESH_MeshAlgos
+
+namespace
+{