+//=======================================================================
+//function : GetCentralNode
+//purpose : Return existing or create a new central node for a quardilateral
+// quadratic face given its 8 nodes.
+//@param : force3d - true means node creation in between the given nodes,
+// else node position is found on a geometrical face if any.
+//=======================================================================
+
+const SMDS_MeshNode* SMESH_MesherHelper::GetCentralNode(const SMDS_MeshNode* n1,
+ const SMDS_MeshNode* n2,
+ const SMDS_MeshNode* n3,
+ const SMDS_MeshNode* n4,
+ const SMDS_MeshNode* n12,
+ const SMDS_MeshNode* n23,
+ const SMDS_MeshNode* n34,
+ const SMDS_MeshNode* n41,
+ bool force3d)
+{
+ SMDS_MeshNode *centralNode = 0; // central node to return
+
+ // Find an existing central node
+
+ TBiQuad keyOfMap(n1,n2,n3,n4);
+ std::map<TBiQuad, SMDS_MeshNode* >::iterator itMapCentralNode;
+ itMapCentralNode = myMapWithCentralNode.find( keyOfMap );
+ if ( itMapCentralNode != myMapWithCentralNode.end() )
+ {
+ return (*itMapCentralNode).second;
+ }
+
+ // Get type of shape for the new central node
+
+ TopAbs_ShapeEnum shapeType = TopAbs_SHAPE;
+ int shapeID = -1;
+ int faceID = -1;
+ TopoDS_Shape shape;
+ TopTools_ListIteratorOfListOfShape it;
+
+ std::map< int, int > faceId2nbNodes;
+ std::map< int, int > ::iterator itMapWithIdFace;
+
+ SMESHDS_Mesh* meshDS = GetMeshDS();
+
+ // check if a face lie on a FACE, i.e. its all corner nodes lie either on the FACE or
+ // on sub-shapes of the FACE
+ if ( GetMesh()->HasShapeToMesh() )
+ {
+ const SMDS_MeshNode* nodes[] = { n1, n2, n3, n4 };
+ for(int i = 0; i < 4; i++)
+ {
+ shape = GetSubShapeByNode( nodes[i], meshDS );
+ if ( shape.IsNull() ) break;
+ if ( shape.ShapeType() == TopAbs_SOLID )
+ {
+ shapeID = nodes[i]->getshapeId();
+ shapeType = TopAbs_SOLID;
+ break;
+ }
+ if ( shape.ShapeType() == TopAbs_FACE )
+ {
+ faceID = nodes[i]->getshapeId();
+ itMapWithIdFace = faceId2nbNodes.insert( std::make_pair( faceID, 0 ) ).first;
+ itMapWithIdFace->second++;
+ }
+ else
+ {
+ PShapeIteratorPtr it = GetAncestors(shape, *GetMesh(), TopAbs_FACE );
+ while ( const TopoDS_Shape* face = it->next() )
+ {
+ faceID = meshDS->ShapeToIndex( *face );
+ itMapWithIdFace = faceId2nbNodes.insert( std::make_pair( faceID, 0 ) ).first;
+ itMapWithIdFace->second++;
+ }
+ }
+ }
+ }
+ if ( shapeID < 1 && !faceId2nbNodes.empty() ) // SOLID not found
+ {
+ // find ID of the FACE the four corner nodes belong to
+ itMapWithIdFace = faceId2nbNodes.begin();
+ for ( ; itMapWithIdFace != faceId2nbNodes.end(); ++itMapWithIdFace)
+ {
+ if ( itMapWithIdFace->second == 4 )
+ {
+ shapeType = TopAbs_FACE;
+ faceID = (*itMapWithIdFace).first;
+ break;
+ }
+ }
+ }
+
+ TopoDS_Face F;
+ if ( shapeType == TopAbs_FACE )
+ {
+ F = TopoDS::Face( meshDS->IndexToShape( faceID ));
+ }
+
+ // Create a node
+
+ gp_XY uvAvg;
+ gp_Pnt P;
+ if ( !F.IsNull() )
+ {
+ if ( !force3d )
+ {
+ uvAvg = calcTFI (0.5, 0.5,
+ GetNodeUV(F,n1,n3), GetNodeUV(F,n2,n4),
+ GetNodeUV(F,n3,n1), GetNodeUV(F,n4,n2),
+ GetNodeUV(F,n12,n3), GetNodeUV(F,n23,n4),
+ GetNodeUV(F,n34,n2), GetNodeUV(F,n41,n2));
+ TopLoc_Location loc;
+ Handle( Geom_Surface ) S = BRep_Tool::Surface( F, loc );
+ P = S->Value( uvAvg.X(), uvAvg.Y() ).Transformed( loc );
+ centralNode = meshDS->AddNode( P.X(), P.Y(), P.Z() );
+ if ( mySetElemOnShape )
+ meshDS->SetNodeOnFace( centralNode, faceID, uvAvg.X(), uvAvg.Y() );
+ myMapWithCentralNode.insert( std::make_pair( keyOfMap, centralNode ) );
+ return centralNode;
+ }
+ }
+
+ P = ( SMESH_TNodeXYZ( n1 ) +
+ SMESH_TNodeXYZ( n2 ) +
+ SMESH_TNodeXYZ( n3 ) +
+ SMESH_TNodeXYZ( n4 ) ) / 4;
+ centralNode = meshDS->AddNode( P.X(), P.Y(), P.Z() );
+
+ if ( mySetElemOnShape )
+ {
+ if ( !F.IsNull() )
+ {
+ uvAvg = (GetNodeUV(F,n1,n3) +
+ GetNodeUV(F,n2,n4) +
+ GetNodeUV(F,n3,n1) +
+ GetNodeUV(F,n4,n2)) / 4;
+ CheckNodeUV( F, centralNode, uvAvg, 2*BRep_Tool::Tolerance( F ), /*force=*/true);
+ meshDS->SetNodeOnFace( centralNode, faceID, uvAvg.X(), uvAvg.Y() );
+ }
+ else if ( shapeID > 0 )
+ meshDS->SetNodeInVolume( centralNode, shapeID );
+ else if ( myShapeID > 0 )
+ meshDS->SetMeshElementOnShape( centralNode, myShapeID );
+ }
+
+ myMapWithCentralNode.insert( std::make_pair( keyOfMap, centralNode ) );
+ return centralNode;
+}
+
