Back to index

salome-smesh  6.5.0
StdMeshers_Projection_3D.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D, OPEN CASCADE
00002 //
00003 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
00004 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
00005 //
00006 // This library is free software; you can redistribute it and/or
00007 // modify it under the terms of the GNU Lesser General Public
00008 // License as published by the Free Software Foundation; either
00009 // version 2.1 of the License.
00010 //
00011 // This library is distributed in the hope that it will be useful,
00012 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00013 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00014 // Lesser General Public License for more details.
00015 //
00016 // You should have received a copy of the GNU Lesser General Public
00017 // License along with this library; if not, write to the Free Software
00018 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00019 //
00020 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00021 //
00022 
00023 //  SMESH SMESH : implementaion of SMESH idl descriptions
00024 // File      : StdMeshers_Projection_3D.cxx
00025 // Module    : SMESH
00026 // Created   : Fri Oct 20 11:37:07 2006
00027 // Author    : Edward AGAPOV (eap)
00028 //
00029 #include "StdMeshers_Projection_3D.hxx"
00030 #include "StdMeshers_ProjectionSource3D.hxx"
00031 
00032 #include "StdMeshers_ProjectionUtils.hxx"
00033 
00034 #include "SMDS_PolyhedralVolumeOfNodes.hxx"
00035 #include "SMDS_VolumeTool.hxx"
00036 #include "SMESHDS_Hypothesis.hxx"
00037 #include "SMESHDS_SubMesh.hxx"
00038 #include "SMESH_Block.hxx"
00039 #include "SMESH_Comment.hxx"
00040 #include "SMESH_Gen.hxx"
00041 #include "SMESH_Mesh.hxx"
00042 #include "SMESH_MesherHelper.hxx"
00043 #include "SMESH_Pattern.hxx"
00044 #include "SMESH_subMesh.hxx"
00045 #include "SMESH_subMeshEventListener.hxx"
00046 
00047 #include "utilities.h"
00048 
00049 #include <TopExp.hxx>
00050 #include <TopExp_Explorer.hxx>
00051 #include <TopoDS.hxx>
00052 
00053 #define RETURN_BAD_RESULT(msg) { MESSAGE(")-: Error: " << msg); return false; }
00054 #define gpXYZ(n) gp_XYZ(n->X(),n->Y(),n->Z())
00055 #define SHOWYXZ(msg, xyz) // {\
00056 // gp_Pnt p (xyz); \
00057 // cout << msg << " ("<< p.X() << "; " <<p.Y() << "; " <<p.Z() << ") " <<endl;\
00058 // }
00059 
00060 typedef StdMeshers_ProjectionUtils TAssocTool;
00061 
00062 
00063 //=======================================================================
00064 //function : StdMeshers_Projection_3D
00065 //purpose  : 
00066 //=======================================================================
00067 
00068 StdMeshers_Projection_3D::StdMeshers_Projection_3D(int hypId, int studyId, SMESH_Gen* gen)
00069   :SMESH_3D_Algo(hypId, studyId, gen)
00070 {
00071   _name = "Projection_3D";
00072   _shapeType = (1 << TopAbs_SHELL) | (1 << TopAbs_SOLID);  // 1 bit per shape type
00073 
00074   _compatibleHypothesis.push_back("ProjectionSource3D");
00075   _sourceHypo = 0;
00076 }
00077 
00078 //================================================================================
00082 //================================================================================
00083 
00084 StdMeshers_Projection_3D::~StdMeshers_Projection_3D()
00085 {}
00086 
00087 //=======================================================================
00088 //function : CheckHypothesis
00089 //purpose  : 
00090 //=======================================================================
00091 
00092 bool StdMeshers_Projection_3D::CheckHypothesis(SMESH_Mesh&                          aMesh,
00093                                                const TopoDS_Shape&                  aShape,
00094                                                SMESH_Hypothesis::Hypothesis_Status& aStatus)
00095 {
00096   // check aShape that must be a 6 faces block
00097 /*  PAL16229
00098   if ( TAssocTool::Count( aShape, TopAbs_SHELL, 1 ) != 1 ||
00099        TAssocTool::Count( aShape, TopAbs_FACE , 1 ) != 6 ||
00100        TAssocTool::Count( aShape, TopAbs_EDGE , 1 ) != 12 ||
00101        TAssocTool::Count( aShape, TopAbs_WIRE , 1 ) != 6 )
00102   {
00103     aStatus = HYP_BAD_GEOMETRY;
00104     return false;
00105   }
00106 */
00107   list <const SMESHDS_Hypothesis * >::const_iterator itl;
00108 
00109   const list <const SMESHDS_Hypothesis * >&hyps = GetUsedHypothesis(aMesh, aShape);
00110   if ( hyps.size() == 0 )
00111   {
00112     aStatus = SMESH_Hypothesis::HYP_MISSING;
00113     return false;  // can't work with no hypothesis
00114   }
00115 
00116   if ( hyps.size() > 1 )
00117   {
00118     aStatus = SMESH_Hypothesis::HYP_ALREADY_EXIST;
00119     return false;
00120   }
00121 
00122   const SMESHDS_Hypothesis *theHyp = hyps.front();
00123 
00124   string hypName = theHyp->GetName();
00125 
00126   aStatus = SMESH_Hypothesis::HYP_OK;
00127 
00128   if (hypName == "ProjectionSource3D")
00129   {
00130     _sourceHypo = static_cast<const StdMeshers_ProjectionSource3D *>(theHyp);
00131     // Check hypo parameters
00132 
00133     SMESH_Mesh* srcMesh = _sourceHypo->GetSourceMesh();
00134     SMESH_Mesh* tgtMesh = & aMesh;
00135     if ( !srcMesh )
00136       srcMesh = tgtMesh;
00137 
00138     // check vertices
00139     if ( _sourceHypo->HasVertexAssociation() )
00140     {
00141       // source vertices
00142       TopoDS_Shape edge = TAssocTool::GetEdgeByVertices
00143         ( srcMesh, _sourceHypo->GetSourceVertex(1), _sourceHypo->GetSourceVertex(2) );
00144       if ( edge.IsNull() ||
00145            !SMESH_MesherHelper::IsSubShape( edge, srcMesh ) ||
00146            !SMESH_MesherHelper::IsSubShape( edge, _sourceHypo->GetSource3DShape() ))
00147       {
00148         SCRUTE((edge.IsNull()));
00149         SCRUTE((SMESH_MesherHelper::IsSubShape( edge, srcMesh )));
00150         SCRUTE((SMESH_MesherHelper::IsSubShape( edge, _sourceHypo->GetSource3DShape() )));
00151         aStatus = SMESH_Hypothesis::HYP_BAD_PARAMETER;
00152       }
00153       else
00154       {
00155         // target vertices
00156         edge = TAssocTool::GetEdgeByVertices
00157           ( tgtMesh, _sourceHypo->GetTargetVertex(1), _sourceHypo->GetTargetVertex(2) );
00158         if ( edge.IsNull() ||
00159              !SMESH_MesherHelper::IsSubShape( edge, tgtMesh ) ||
00160              !SMESH_MesherHelper::IsSubShape( edge, aShape ))
00161         {
00162           SCRUTE((edge.IsNull()));
00163           SCRUTE((SMESH_MesherHelper::IsSubShape( edge, tgtMesh )));
00164           SCRUTE((SMESH_MesherHelper::IsSubShape( edge, aShape )));
00165           aStatus = SMESH_Hypothesis::HYP_BAD_PARAMETER;
00166         }
00167       }
00168     }
00169     // check a source shape
00170     if ( !SMESH_MesherHelper::IsSubShape( _sourceHypo->GetSource3DShape(), srcMesh ) ||
00171          ( srcMesh == tgtMesh && aShape == _sourceHypo->GetSource3DShape()))
00172     {
00173       SCRUTE((SMESH_MesherHelper::IsSubShape( _sourceHypo->GetSource3DShape(), srcMesh)));
00174       SCRUTE((srcMesh == tgtMesh));
00175       SCRUTE((aShape == _sourceHypo->GetSource3DShape()));
00176       aStatus = SMESH_Hypothesis::HYP_BAD_PARAMETER;
00177     }
00178   }
00179   else
00180   {
00181     aStatus = SMESH_Hypothesis::HYP_INCOMPATIBLE;
00182   }
00183   return ( aStatus == HYP_OK );
00184 }
00185 
00186 //=======================================================================
00187 //function : Compute
00188 //purpose  : 
00189 //=======================================================================
00190 
00191 bool StdMeshers_Projection_3D::Compute(SMESH_Mesh& aMesh, const TopoDS_Shape& aShape)
00192 {
00193   if ( !_sourceHypo )
00194     return false;
00195 
00196   SMESH_Mesh * srcMesh = _sourceHypo->GetSourceMesh();
00197   SMESH_Mesh * tgtMesh = & aMesh;
00198   if ( !srcMesh )
00199     srcMesh = tgtMesh;
00200 
00201   SMESHDS_Mesh * srcMeshDS = srcMesh->GetMeshDS();
00202   SMESHDS_Mesh * tgtMeshDS = tgtMesh->GetMeshDS();
00203 
00204   // get shell from shape3D
00205   TopoDS_Shell srcShell, tgtShell;
00206   TopExp_Explorer exp( _sourceHypo->GetSource3DShape(), TopAbs_SHELL );
00207   int nbShell;
00208   for ( nbShell = 0; exp.More(); exp.Next(), ++nbShell )
00209     srcShell = TopoDS::Shell( exp.Current() );
00210   if ( nbShell != 1 )
00211     return error(COMPERR_BAD_SHAPE,
00212                  SMESH_Comment("Source shape must have 1 shell but not ") << nbShell);
00213 
00214   exp.Init( aShape, TopAbs_SHELL );
00215   for ( nbShell = 0; exp.More(); exp.Next(), ++nbShell )
00216     tgtShell = TopoDS::Shell( exp.Current() );
00217   if ( nbShell != 1 )
00218     return error(COMPERR_BAD_SHAPE,
00219                  SMESH_Comment("Target shape must have 1 shell but not ") << nbShell);
00220 
00221   // Check that shapes are blocks
00222   if ( TAssocTool::Count( tgtShell, TopAbs_FACE , 1 ) != 6 ||
00223        TAssocTool::Count( tgtShell, TopAbs_EDGE , 1 ) != 12 ||
00224        TAssocTool::Count( tgtShell, TopAbs_WIRE , 1 ) != 6 )
00225     return error(COMPERR_BAD_SHAPE, "Target shape is not a block");
00226   if ( TAssocTool::Count( srcShell, TopAbs_FACE , 1 ) != 6 ||
00227        TAssocTool::Count( srcShell, TopAbs_EDGE , 1 ) != 12 ||
00228        TAssocTool::Count( srcShell, TopAbs_WIRE , 1 ) != 6 )
00229     return error(COMPERR_BAD_SHAPE, "Source shape is not a block");
00230 
00231   // Assure that mesh on a source shape is computed
00232 
00233   SMESH_subMesh* srcSubMesh = srcMesh->GetSubMesh( _sourceHypo->GetSource3DShape() );
00234   //SMESH_subMesh* tgtSubMesh = tgtMesh->GetSubMesh( aShape );
00235 
00236   if ( tgtMesh == srcMesh && !aShape.IsSame( _sourceHypo->GetSource3DShape() )) {
00237     if ( !TAssocTool::MakeComputed( srcSubMesh ))
00238       return error(COMPERR_BAD_INPUT_MESH,"Source mesh not computed");
00239   }
00240   else {
00241     if ( !srcSubMesh->IsMeshComputed() )
00242       return error(COMPERR_BAD_INPUT_MESH,"Source mesh not computed");
00243   }
00244 
00245   // Find 2 pairs of corresponding vertices
00246 
00247   TopoDS_Vertex tgtV000, tgtV100, srcV000, srcV100;
00248   TAssocTool::TShapeShapeMap shape2ShapeMap;
00249 
00250   if ( _sourceHypo->HasVertexAssociation() )
00251   {
00252     tgtV000 = _sourceHypo->GetTargetVertex(1);
00253     tgtV100 = _sourceHypo->GetTargetVertex(2);
00254     srcV000 = _sourceHypo->GetSourceVertex(1);
00255     srcV100 = _sourceHypo->GetSourceVertex(2);
00256   }
00257   else
00258   {
00259     if ( !TAssocTool::FindSubShapeAssociation( tgtShell, tgtMesh, srcShell, srcMesh,
00260                                                shape2ShapeMap) )
00261       return error(COMPERR_BAD_SHAPE,"Topology of source and target shapes seems different" );
00262 
00263     exp.Init( tgtShell, TopAbs_EDGE );
00264     TopExp::Vertices( TopoDS::Edge( exp.Current() ), tgtV000, tgtV100 );
00265 
00266     if ( !shape2ShapeMap.IsBound( tgtV000 ) || !shape2ShapeMap.IsBound( tgtV100 ))
00267       return error("Association of sub-shapes failed" );
00268     srcV000 = TopoDS::Vertex( shape2ShapeMap( tgtV000 ));
00269     srcV100 = TopoDS::Vertex( shape2ShapeMap( tgtV100 ));
00270     if ( !SMESH_MesherHelper::IsSubShape( srcV000, srcShell ) ||
00271          !SMESH_MesherHelper::IsSubShape( srcV100, srcShell ))
00272       return error("Incorrect association of sub-shapes" );
00273   }
00274 
00275   // Load 2 SMESH_Block's with src and tgt shells
00276 
00277   SMESH_Block srcBlock, tgtBlock;
00278   TopTools_IndexedMapOfOrientedShape scrShapes, tgtShapes;
00279   if ( !tgtBlock.LoadBlockShapes( tgtShell, tgtV000, tgtV100, tgtShapes ))
00280     return error(COMPERR_BAD_SHAPE, "Can't detect block sub-shapes. Not a block?");
00281 
00282   if ( !srcBlock.LoadBlockShapes( srcShell, srcV000, srcV100, scrShapes ))
00283     return error(COMPERR_BAD_SHAPE, "Can't detect block sub-shapes. Not a block?");
00284 
00285   // Find matching nodes of src and tgt shells
00286 
00287   TNodeNodeMap src2tgtNodeMap;
00288   for ( int fId = SMESH_Block::ID_FirstF; fId < SMESH_Block::ID_Shell; ++fId )
00289   {
00290     // Corresponding sub-shapes
00291     TopoDS_Face srcFace = TopoDS::Face( scrShapes( fId ));
00292     TopoDS_Face tgtFace = TopoDS::Face( tgtShapes( fId ));
00293     if ( _sourceHypo->HasVertexAssociation() ) { // associate face sub-shapes
00294       shape2ShapeMap.Clear();
00295       vector< int > edgeIdVec;
00296       SMESH_Block::GetFaceEdgesIDs( fId, edgeIdVec );
00297       for ( int i = 0; i < edgeIdVec.size(); ++i ) {
00298         int eID = edgeIdVec[ i ];
00299         shape2ShapeMap.Bind( tgtShapes( eID ), scrShapes( eID ));
00300         if ( i < 2 ) {
00301           vector< int > vertexIdVec;
00302           SMESH_Block::GetEdgeVertexIDs( eID, vertexIdVec );
00303           shape2ShapeMap.Bind( tgtShapes( vertexIdVec[0] ), scrShapes( vertexIdVec[0] ));
00304           shape2ShapeMap.Bind( tgtShapes( vertexIdVec[1] ), scrShapes( vertexIdVec[1] ));
00305         }
00306       }
00307     }
00308     // Find matching nodes of tgt and src faces
00309     TNodeNodeMap faceMatchingNodes;
00310     if ( ! TAssocTool::FindMatchingNodesOnFaces( srcFace, srcMesh, tgtFace, tgtMesh, 
00311                                                  shape2ShapeMap, faceMatchingNodes ))
00312       return error(COMPERR_BAD_INPUT_MESH,SMESH_Comment("Mesh on faces #")
00313                    << srcMeshDS->ShapeToIndex( srcFace ) << " and "
00314                    << tgtMeshDS->ShapeToIndex( tgtFace ) << " seems different" );
00315 
00316     // put found matching nodes of 2 faces to the global map
00317     src2tgtNodeMap.insert( faceMatchingNodes.begin(), faceMatchingNodes.end() );
00318   }
00319 
00320   // ------------------
00321   // Make mesh
00322   // ------------------
00323 
00324   SMDS_VolumeTool volTool;
00325   SMESH_MesherHelper helper( *tgtMesh );
00326   helper.IsQuadraticSubMesh( aShape );
00327 
00328   SMESHDS_SubMesh* srcSMDS = srcSubMesh->GetSubMeshDS();
00329   SMDS_ElemIteratorPtr volIt = srcSMDS->GetElements();
00330   while ( volIt->more() ) // loop on source volumes
00331   {
00332     const SMDS_MeshElement* srcVol = volIt->next();
00333     if ( !srcVol || srcVol->GetType() != SMDSAbs_Volume )
00334         continue;
00335     int nbNodes = srcVol->NbNodes();
00336     SMDS_VolumeTool::VolumeType  volType = volTool.GetType( nbNodes );
00337     if ( srcVol->IsQuadratic() )
00338       nbNodes = volTool.NbCornerNodes( volType );
00339 
00340     // Find or create a new tgt node for each node of a src volume
00341 
00342     vector< const SMDS_MeshNode* > nodes( nbNodes );
00343     for ( int i = 0; i < nbNodes; ++i )
00344     {
00345       const SMDS_MeshNode* srcNode = srcVol->GetNode( i );
00346       const SMDS_MeshNode* tgtNode = 0;
00347       TNodeNodeMap::iterator sN_tN = src2tgtNodeMap.find( srcNode );
00348       if ( sN_tN != src2tgtNodeMap.end() ) // found
00349       {
00350         tgtNode = sN_tN->second;
00351       }
00352       else // Create a new tgt node
00353       {
00354         // compute normalized parameters of source node in srcBlock
00355         gp_Pnt srcCoord = gpXYZ( srcNode );
00356         gp_XYZ srcParam;
00357         if ( !srcBlock.ComputeParameters( srcCoord, srcParam ))
00358           return error(SMESH_Comment("Can't compute normalized parameters ")
00359                        << "for source node " << srcNode->GetID());
00360         // compute coordinates of target node by srcParam
00361         gp_XYZ tgtXYZ;
00362         if ( !tgtBlock.ShellPoint( srcParam, tgtXYZ ))
00363           return error("Can't compute coordinates by normalized parameters");
00364         // add node
00365         SMDS_MeshNode* newNode = tgtMeshDS->AddNode( tgtXYZ.X(), tgtXYZ.Y(), tgtXYZ.Z() );
00366         tgtMeshDS->SetNodeInVolume( newNode, helper.GetSubShapeID() );
00367         tgtNode = newNode;
00368         src2tgtNodeMap.insert( make_pair( srcNode, tgtNode ));
00369       }
00370       nodes[ i ] = tgtNode;
00371     }
00372 
00373     // Create a new volume
00374 
00375     SMDS_MeshVolume * tgtVol = 0;
00376     int id = 0, force3d = false;
00377     switch ( volType ) {
00378     case SMDS_VolumeTool::TETRA     :
00379     case SMDS_VolumeTool::QUAD_TETRA:
00380       tgtVol = helper.AddVolume( nodes[0],
00381                                  nodes[1],
00382                                  nodes[2],
00383                                  nodes[3], id, force3d); break;
00384     case SMDS_VolumeTool::PYRAM     :
00385     case SMDS_VolumeTool::QUAD_PYRAM:
00386       tgtVol = helper.AddVolume( nodes[0],
00387                                  nodes[1],
00388                                  nodes[2],
00389                                  nodes[3],
00390                                  nodes[4], id, force3d); break;
00391     case SMDS_VolumeTool::PENTA     :
00392     case SMDS_VolumeTool::QUAD_PENTA:
00393       tgtVol = helper.AddVolume( nodes[0],
00394                                  nodes[1],
00395                                  nodes[2],
00396                                  nodes[3],
00397                                  nodes[4],
00398                                  nodes[5], id, force3d); break;
00399     case SMDS_VolumeTool::HEXA      :
00400     case SMDS_VolumeTool::QUAD_HEXA :
00401       tgtVol = helper.AddVolume( nodes[0],
00402                                  nodes[1],
00403                                  nodes[2],
00404                                  nodes[3],
00405                                  nodes[4],
00406                                  nodes[5],
00407                                  nodes[6],
00408                                  nodes[7], id, force3d); break;
00409     default: // polyhedron
00410       const SMDS_VtkVolume * poly =
00411         dynamic_cast<const SMDS_VtkVolume*>( srcVol );
00412       if ( !poly )
00413         RETURN_BAD_RESULT("Unexpected volume type");
00414       if ( !poly->IsPoly())
00415         RETURN_BAD_RESULT("Unexpected volume type");
00416       tgtVol = tgtMeshDS->AddPolyhedralVolume( nodes, poly->GetQuantities() );
00417     }
00418     if ( tgtVol ) {
00419       tgtMeshDS->SetMeshElementOnShape( tgtVol, helper.GetSubShapeID() );
00420     }
00421   } // loop on volumes of src shell
00422 
00423   return true;
00424 }
00425 
00426 
00427 //=======================================================================
00428 //function : Evaluate
00429 //purpose  : 
00430 //=======================================================================
00431 
00432 bool StdMeshers_Projection_3D::Evaluate(SMESH_Mesh& aMesh,
00433                                         const TopoDS_Shape& aShape,
00434                                         MapShapeNbElems& aResMap)
00435 {
00436   if ( !_sourceHypo )
00437     return false;
00438 
00439   SMESH_Mesh * srcMesh = _sourceHypo->GetSourceMesh();
00440   SMESH_Mesh * tgtMesh = & aMesh;
00441   if ( !srcMesh )
00442     srcMesh = tgtMesh;
00443 
00444   // get shell from shape3D
00445   TopoDS_Shell srcShell, tgtShell;
00446   TopExp_Explorer exp( _sourceHypo->GetSource3DShape(), TopAbs_SHELL );
00447   int nbShell;
00448   for ( nbShell = 0; exp.More(); exp.Next(), ++nbShell )
00449     srcShell = TopoDS::Shell( exp.Current() );
00450   if ( nbShell != 1 )
00451     return error(COMPERR_BAD_SHAPE,
00452                  SMESH_Comment("Source shape must have 1 shell but not ") << nbShell);
00453 
00454   exp.Init( aShape, TopAbs_SHELL );
00455   for ( nbShell = 0; exp.More(); exp.Next(), ++nbShell )
00456     tgtShell = TopoDS::Shell( exp.Current() );
00457   if ( nbShell != 1 )
00458     return error(COMPERR_BAD_SHAPE,
00459                  SMESH_Comment("Target shape must have 1 shell but not ") << nbShell);
00460 
00461   // Check that shapes are blocks
00462   if ( TAssocTool::Count( tgtShell, TopAbs_FACE , 1 ) != 6 ||
00463        TAssocTool::Count( tgtShell, TopAbs_EDGE , 1 ) != 12 ||
00464        TAssocTool::Count( tgtShell, TopAbs_WIRE , 1 ) != 6 )
00465     return error(COMPERR_BAD_SHAPE, "Target shape is not a block");
00466   if ( TAssocTool::Count( srcShell, TopAbs_FACE , 1 ) != 6 ||
00467        TAssocTool::Count( srcShell, TopAbs_EDGE , 1 ) != 12 ||
00468        TAssocTool::Count( srcShell, TopAbs_WIRE , 1 ) != 6 )
00469     return error(COMPERR_BAD_SHAPE, "Source shape is not a block");
00470 
00471   // Assure that mesh on a source shape is computed
00472 
00473   SMESH_subMesh* srcSubMesh = srcMesh->GetSubMesh( _sourceHypo->GetSource3DShape() );
00474 
00475   if ( !srcSubMesh->IsMeshComputed() )
00476     return error(COMPERR_BAD_INPUT_MESH,"Source mesh not computed");
00477 
00478 
00479   std::vector<int> aVec(SMDSEntity_Last);
00480   for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aVec[i] = 0;
00481 
00482   aVec[SMDSEntity_Node] = srcSubMesh->GetSubMeshDS()->NbNodes();
00483 
00484   //bool quadratic = false;
00485   SMDS_ElemIteratorPtr elemIt = srcSubMesh->GetSubMeshDS()->GetElements();
00486   while ( elemIt->more() ) {
00487     const SMDS_MeshElement* E  = elemIt->next();
00488     if( E->NbNodes()==4 ) {
00489       aVec[SMDSEntity_Tetra]++;
00490     }
00491     else if( E->NbNodes()==5 ) {
00492       aVec[SMDSEntity_Pyramid]++;
00493     }
00494     else if( E->NbNodes()==6 ) {
00495       aVec[SMDSEntity_Penta]++;
00496     }
00497     else if( E->NbNodes()==8 ) {
00498       aVec[SMDSEntity_Hexa]++;
00499     }
00500     else if( E->NbNodes()==10 && E->IsQuadratic() ) {
00501       aVec[SMDSEntity_Quad_Tetra]++;
00502     }
00503     else if( E->NbNodes()==13 && E->IsQuadratic() ) {
00504       aVec[SMDSEntity_Quad_Pyramid]++;
00505     }
00506     else if( E->NbNodes()==15 && E->IsQuadratic() ) {
00507       aVec[SMDSEntity_Quad_Penta]++;
00508     }
00509     else if( E->NbNodes()==20 && E->IsQuadratic() ) {
00510       aVec[SMDSEntity_Quad_Hexa]++;
00511     }
00512     else {
00513       aVec[SMDSEntity_Polyhedra]++;
00514     }
00515   }
00516 
00517   SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
00518   aResMap.insert(std::make_pair(sm,aVec));
00519 
00520   return true;
00521 }
00522 
00523 
00524 //=============================================================================
00534 //=============================================================================
00535 
00536 void StdMeshers_Projection_3D::SetEventListener(SMESH_subMesh* subMesh)
00537 {
00538   TAssocTool::SetEventListener( subMesh,
00539                                 _sourceHypo->GetSource3DShape(),
00540                                 _sourceHypo->GetSourceMesh() );
00541 }
00542