]> SALOME platform Git repositories - plugins/ghs3dplugin.git/commitdiff
Salome HOME
Merge from PHASE_25_BR 14/11/2010
authorvsr <vsr@opencascade.com>
Mon, 15 Nov 2010 06:54:33 +0000 (06:54 +0000)
committervsr <vsr@opencascade.com>
Mon, 15 Nov 2010 06:54:33 +0000 (06:54 +0000)
src/GHS3DPlugin_GHS3D.cxx
src/Makefile.am

index 16a00f87f1024a2b39c643ca0a2196375e58bbc8..28ee3aa9bff3b2a521c960c69f4d59e92aad9eeb 100644 (file)
@@ -22,7 +22,6 @@
 // Created   : 
 // Author    : Edward AGAPOV, modified by Lioka RAZAFINDRAZAKA (CEA) 09/02/2007
 // Project   : SALOME
-// $Header$
 //=============================================================================
 //
 #include "GHS3DPlugin_GHS3D.hxx"
@@ -42,6 +41,8 @@
 #include "SMDS_FaceOfNodes.hxx"
 #include "SMDS_VolumeOfNodes.hxx"
 
+#include <StdMeshers_QuadToTriaAdaptor.hxx>
+
 #include <BRepAdaptor_Surface.hxx>
 #include <BRepBndLib.hxx>
 #include <BRepClass3d_SolidClassifier.hxx>
@@ -99,6 +100,8 @@ extern "C"
 
 #define HOLE_ID -1
 
+typedef const list<const SMDS_MeshFace*> TTriaList;
+
 static void removeFile( const TCollection_AsciiString& fileName )
 {
   try {
@@ -242,9 +245,11 @@ static int countShape( Mesh* mesh, Shape shape ) {
 //purpose  : 
 //=======================================================================
 
-static bool writeFaces (ofstream &            theFile,
-                        SMESHDS_Mesh *        theMesh,
-                        const map <int,int> & theSmdsToGhs3dIdMap)
+static bool writeFaces (ofstream &                            theFile,
+                        SMESHDS_Mesh*                         theMesh,
+                        const TopoDS_Shape&                   theShape,
+                        vector<StdMeshers_QuadToTriaAdaptor>& theQuad2Trias,
+                        const map <int,int> &                 theSmdsToGhs3dIdMap)
 {
   // record structure:
   //
@@ -252,10 +257,6 @@ static bool writeFaces (ofstream &            theFile,
   // Loop from 1 to NB_ELEMS
   // NB_NODES NODE_NB_1 NODE_NB_2 ... (NB_NODES + 1) times: DUMMY_INT
 
-  int nbShape = countShape( theMesh, TopAbs_FACE );
-
-  int *tabID;             tabID    = new int[nbShape];
-  TopoDS_Shape *tabShape; tabShape = new TopoDS_Shape[nbShape];
   TopoDS_Shape aShape;
   SMESHDS_SubMesh* theSubMesh;
   const SMDS_MeshElement* aFace;
@@ -263,46 +264,97 @@ static bool writeFaces (ofstream &            theFile,
   const int   dummyint = 0;
   map<int,int>::const_iterator itOnMap;
   SMDS_ElemIteratorPtr itOnSubMesh, itOnSubFace;
-  int shapeID, nbNodes, aSmdsID;
-  bool idFound;
+  int nbNodes, aSmdsID;
 
-  // count faces bound to geometry
-  int nbFaces = 0;
+  // count triangles bound to geometry
+  int nbTriangles = 0;
 
-  TopExp_Explorer expface( theMesh->ShapeToMesh(), TopAbs_FACE );
-  for ( int i = 0; expface.More(); expface.Next(), i++ ) {
-    tabID[i] = 0;
-    aShape   = expface.Current();
-    shapeID  = theMesh->ShapeToIndex( aShape );
-    idFound  = false;
-    for ( int j=0; j<=i; j++) {
-      if ( shapeID == tabID[j] ) {
-        idFound = true;
-        break;
-      }
-    }
-    if ( ! idFound ) {
-      tabID[i]    = shapeID;
-      tabShape[i] = aShape;
-      theSubMesh  = theMesh->MeshElements( aShape );
-      if ( theSubMesh )
-        nbFaces += theSubMesh->NbElements();
+  TopTools_IndexedMapOfShape facesMap;
+  TopExp::MapShapes( theShape, TopAbs_FACE, facesMap );
+
+  // 2 adaptors for each face in facesMap, as a face can belong to 2 solids
+  typedef vector< StdMeshers_QuadToTriaAdaptor* > TwoAdaptors;
+  vector< TwoAdaptors > qttaByFace;
+  if ( theQuad2Trias.empty() )
+  {
+    // case w/o quadrangles
+    for ( int i = 1; i <= facesMap.Extent(); ++i )
+      if (( theSubMesh  = theMesh->MeshElements( facesMap(i))))
+        nbTriangles += theSubMesh->NbElements();
+  }
+  else
+  {
+    // case with quadrangles
+    qttaByFace.resize( facesMap.Extent() );
+    for ( unsigned i = 0; i < theQuad2Trias.size(); ++i )
+    {
+      TopoDS_Shape solid = theQuad2Trias[i].GetShape();
+      TopExp_Explorer expface( solid, TopAbs_FACE );
+      for ( ; expface.More(); expface.Next() )
+        if (( theSubMesh = theMesh->MeshElements( expface.Current()) ))
+        {
+          const int faceIndex = facesMap.Add( expface.Current() );
+          TwoAdaptors& aTwoAdaptors = qttaByFace[ faceIndex-1 ];
+          const bool newFaceEncounters = aTwoAdaptors.empty();
+          aTwoAdaptors.push_back( & theQuad2Trias[i] );
+
+          // on a shared face encountered for the second time
+          // we count only triangles of pyramids
+          const int countTrias = int( newFaceEncounters );
+          itOnSubMesh = theSubMesh->GetElements();
+          while ( itOnSubMesh->more() )
+          {
+            aFace = itOnSubMesh->next();
+            if ( aFace->NbCornerNodes() != 4 )
+              nbTriangles += countTrias;
+            else if ( TTriaList* trias = theQuad2Trias[i].GetTriangles( aFace ))
+              nbTriangles += trias->size();
+            else
+              nbTriangles += countTrias;
+          }
+        }
     }
   }
-  std::cout << "    " << nbFaces << " shapes of 2D dimension" << std::endl;
+
+  std::cout << "    " << facesMap.Extent() << " shapes of 2D dimension" << std::endl;
   std::cout << std::endl;
 
-  theFile << space << nbFaces << space << dummyint << std::endl;
-
-  for ( int i =0; i < nbShape; i++ ) {
-    if ( tabID[i] != 0 ) {
-      aShape      = tabShape[i];
-      shapeID     = tabID[i];
-      theSubMesh  = theMesh->MeshElements( aShape );
-      if ( !theSubMesh ) continue;
-      itOnSubMesh = theSubMesh->GetElements();
-      while ( itOnSubMesh->more() ) {
-        aFace   = itOnSubMesh->next();
+  theFile << space << nbTriangles << space << dummyint << std::endl;
+
+  vector< const SMDS_MeshElement* > trias;
+  trias.resize( 8 ); // 4 triangles from 2 pyramids basing on one quadranle
+
+  for ( int i = 1; i <= facesMap.Extent(); i++ )
+  {
+    aShape = facesMap(i);
+    theSubMesh = theMesh->MeshElements( aShape );
+    if ( !theSubMesh ) continue;
+    TwoAdaptors& aTwoAdaptors = qttaByFace[ i-1 ];
+    itOnSubMesh = theSubMesh->GetElements();
+    while ( itOnSubMesh->more() )
+    {
+      aFace = itOnSubMesh->next();
+      if ( aFace->NbCornerNodes() == 4 )
+      {
+        nbTriangles = 0;
+        for ( unsigned j = 0; j < aTwoAdaptors.size(); ++j )
+          if ( TTriaList* t = aTwoAdaptors[j]->GetTriangles( aFace ))
+            for ( TTriaList::const_iterator tIt = t->begin(); tIt != t->end(); ++tIt)
+              trias[nbTriangles++] = *tIt;
+        if ( nbTriangles == 0 )
+        {
+          nbTriangles = 1;
+          trias[0] = aFace;
+        }
+      }
+      else
+      {
+        nbTriangles = 1;
+        trias[0] = aFace;
+      }
+      for ( int j = 0; j < nbTriangles; ++j )
+      {
+        aFace = trias[j];
         nbNodes = aFace->NbNodes();
 
         theFile << space << nbNodes;
@@ -330,10 +382,6 @@ static bool writeFaces (ofstream &            theFile,
     }
   }
   
-  
-  delete [] tabID;
-  delete [] tabShape;
-
   return true;
 }
 
@@ -342,8 +390,9 @@ static bool writeFaces (ofstream &            theFile,
 //purpose  : Write Faces in case if generate 3D mesh w/o geometry
 //=======================================================================
 
-static bool writeFaces (ofstream &            theFile,
-                        SMESHDS_Mesh *        theMesh,
+static bool writeFaces (ofstream &                      theFile,
+                        SMESHDS_Mesh *                  theMesh,
+                        StdMeshers_QuadToTriaAdaptor&   theQuad2Trias,
                         vector <const SMDS_MeshNode*> & theNodeByGhs3dId)
 {
   // record structure:
@@ -352,63 +401,83 @@ static bool writeFaces (ofstream &            theFile,
   // Loop from 1 to NB_ELEMS
   //   NB_NODES NODE_NB_1 NODE_NB_2 ... (NB_NODES + 1) times: DUMMY_INT
 
-
-  int nbFaces = 0;
-  list< const SMDS_MeshElement* > faces;
-  list< const SMDS_MeshElement* >::iterator f;
+  int nbNodes, nbTriangles = 0;
   map< const SMDS_MeshNode*,int >::iterator it;
   SMDS_ElemIteratorPtr nodeIt;
   const SMDS_MeshElement* elem;
-  int nbNodes;
 
   const char* space    = "  ";
   const int   dummyint = 0;
 
-  //get all faces from mesh
-  SMDS_FaceIteratorPtr eIt = theMesh->facesIterator();
-  while ( eIt->more() ) {
-    const SMDS_MeshElement* elem = eIt->next();
-    if ( !elem )
-      return false;
-    faces.push_back( elem );
-    nbFaces++;
-  }
+  // count faces
+
+  nbTriangles = theQuad2Trias.TotalNbOfTriangles();
+  if ( nbTriangles == 0 )
+    // theQuad2Trias not computed as there are no quadrangles in mesh
+    nbTriangles = theMesh->NbFaces();
 
-  if ( nbFaces == 0 )
+  if ( nbTriangles == 0 )
     return false;
-  
-  std::cout << "The initial 2D mesh contains " << nbFaces << " faces and ";
+
+  std::cout << "The initial 2D mesh contains " << nbTriangles << " faces and ";
 
   // NB_ELEMS DUMMY_INT
-  theFile << space << nbFaces << space << dummyint << std::endl;
+  theFile << space << nbTriangles << space << dummyint << std::endl;
 
-  // Loop from 1 to NB_ELEMS
 
   map<const SMDS_MeshNode*,int> aNodeToGhs3dIdMap;
-  f = faces.begin();
-  for ( ; f != faces.end(); ++f )
+
+  vector< const SMDS_MeshElement* > trias;
+  trias.resize( 8 ); // 8 == 4 triangles from 2 pyramids basing on one quadranle
+
+  // Loop from 1 to NB_ELEMS
+
+  SMDS_FaceIteratorPtr eIt = theMesh->facesIterator();
+  while ( eIt->more() )
   {
-    // NB_NODES PER FACE
-    elem = *f;
-    nbNodes = elem->NbNodes();
-    theFile << space << nbNodes;
-
-    // NODE_NB_1 NODE_NB_2 ...
-    nodeIt = elem->nodesIterator();
-    while ( nodeIt->more() )
+    elem = eIt->next();
+    // get triangles
+    if ( elem->NbCornerNodes() == 4 )
     {
-      // find GHS3D ID
-      const SMDS_MeshNode* node = castToNode( nodeIt->next() );
-      int newId = aNodeToGhs3dIdMap.size() + 1; // ghs3d ids count from 1
-      it = aNodeToGhs3dIdMap.insert( make_pair( node, newId )).first;
-      theFile << space << it->second;
+      nbTriangles = 0;
+      if ( TTriaList* t = theQuad2Trias.GetTriangles( elem ))
+        for ( TTriaList::const_iterator tIt = t->begin(); tIt != t->end(); ++tIt)
+          trias[nbTriangles++] = *tIt;
+      if ( nbTriangles == 0 )
+      {
+        nbTriangles = 1;
+        trias[0] = elem;
+      }
     }
+    else
+    {
+      nbTriangles = 1;
+      trias[0] = elem;
+    }
+    // write triangles
+    for ( int j = 0; j < nbTriangles; ++j )
+    {
+      // NB_NODES PER FACE
+      elem = trias[j];
+      nbNodes = elem->NbNodes();
+      theFile << space << nbNodes;
+
+      // NODE_NB_1 NODE_NB_2 ...
+      nodeIt = elem->nodesIterator();
+      while ( nodeIt->more() )
+      {
+        // find GHS3D ID
+        const SMDS_MeshNode* node = castToNode( nodeIt->next() );
+        int newId = aNodeToGhs3dIdMap.size() + 1; // ghs3d ids count from 1
+        it = aNodeToGhs3dIdMap.insert( make_pair( node, newId )).first;
+        theFile << space << it->second;
+      }
 
-    // (NB_NODES + 1) times: DUMMY_INT
-    for ( int i=0; i<=nbNodes; i++)
-      theFile << space << dummyint;
-
-    theFile << std::endl;
+      // (NB_NODES + 1) times: DUMMY_INT
+      for ( int i=0; i<=nbNodes; i++)
+        theFile << space << dummyint;
+      theFile << std::endl;
+    }
   }
 
   // put nodes to theNodeByGhs3dId vector
@@ -1212,7 +1281,6 @@ bool GHS3DPlugin_GHS3D::Compute(SMESH_Mesh&         theMesh,
     tabBox[i] = new double[6];
   Standard_Real Xmin, Ymin, Zmin, Xmax, Ymax, Zmax;
 
-  // TopExp_Explorer expBox (meshDS->ShapeToMesh(), TopAbs_SOLID); -- 0020330:...ghs3d as a submesh
   for (expBox.ReInit(); expBox.More(); expBox.Next()) {
     tabShape[iShape] = expBox.Current();
     Bnd_Box BoundingBox;
@@ -1265,8 +1333,17 @@ bool GHS3DPlugin_GHS3D::Compute(SMESH_Mesh&         theMesh,
   SMESH_MesherHelper helper( theMesh );
   helper.SetSubShape( theShape );
 
-  Ok = writePoints( aPointsFile, helper, aSmdsToGhs3dIdMap, aGhs3dIdToNodeMap, enforcedVertices) &&
-       writeFaces ( aFacesFile,  meshDS, aSmdsToGhs3dIdMap );
+  // make prisms on quadrangles
+  vector<StdMeshers_QuadToTriaAdaptor> aQuad2Trias;
+  if ( theMesh.NbQuadrangles() > 0 )
+  {
+    aQuad2Trias.resize( _nbShape );
+    for (_iShape = 0, expBox.ReInit(); expBox.More(); expBox.Next())
+      aQuad2Trias[ _iShape++ ].Compute( theMesh, expBox.Current() );
+  }
+
+  Ok = (writePoints( aPointsFile, helper, aSmdsToGhs3dIdMap, aGhs3dIdToNodeMap, enforcedVertices) &&
+        writeFaces ( aFacesFile,  meshDS, theShape, aQuad2Trias, aSmdsToGhs3dIdMap ));
 
   // Write aSmdsToGhs3dIdMap to temp file
   TCollection_AsciiString aSmdsToGhs3dIdMapFileName;
@@ -1434,7 +1511,11 @@ bool GHS3DPlugin_GHS3D::Compute(SMESH_Mesh&         theMesh,
 
   vector <const SMDS_MeshNode*> aNodeByGhs3dId;
 
-  Ok = (writeFaces ( aFacesFile, meshDS, aNodeByGhs3dId ) &&
+  StdMeshers_QuadToTriaAdaptor aQuad2Trias;
+  if ( theMesh.NbQuadrangles() > 0 )
+    aQuad2Trias.Compute( theMesh );
+
+  Ok = (writeFaces ( aFacesFile, meshDS, aQuad2Trias, aNodeByGhs3dId ) &&
         writePoints( aPointsFile, &theMesh, aNodeByGhs3dId,enforcedVertices));
   
   aFacesFile.close();
index d4487d4e87278ab2aa964615e80f375d44d8f188..44ae6d1ac6f05f14832b47e79f7ca0dd9f229461 100644 (file)
@@ -64,5 +64,5 @@ libGHS3DEngine_la_LDFLAGS  = \
        ../idl/libSalomeIDLGHS3DPLUGIN.la \
        $(CAS_KERNEL) -lTKBRep -lTKG2d -lTKG3d -lTKTopAlgo -lTKGeomBase -lTKGeomAlgo \
        $(MED_LDFLAGS) -lSalomeIDLMED \
-       $(SMESH_LDFLAGS) -lSMESHimpl -lSMESHEngine -lSMESHDS -lSMDS \
+       $(SMESH_LDFLAGS) -lSMESHimpl -lSMESHEngine -lSMESHDS -lSMDS -lStdMeshers \
        $(KERNEL_LDFLAGS) -lSalomeGenericObj -lSALOMELocalTrace -lSALOMEBasics