-// Copyright (C) 2007-2013 CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2014 CEA/DEN, EDF R&D, OPEN CASCADE
//
// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
-// version 2.1 of the License.
+// version 2.1 of the License, or (at your option) any later version.
//
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
#include "SMESH_MeshAlgos.hxx"
+#include "SMDS_FaceOfNodes.hxx"
#include "SMDS_LinearEdge.hxx"
-#include "SMDS_VolumeTool.hxx"
#include "SMDS_Mesh.hxx"
+#include "SMDS_PolygonalFaceOfNodes.hxx"
+#include "SMDS_VolumeTool.hxx"
#include "SMESH_OctreeNode.hxx"
#include <GC_MakeSegment.hxx>
void getElementsInSphere ( const gp_XYZ& center,
const double radius, TIDSortedElemSet& foundElems);
size_t getSize() { return std::max( _size, _elements.size() ); }
- ~ElementBndBoxTree();
+ virtual ~ElementBndBoxTree();
protected:
ElementBndBoxTree():_size(0) {}
*/
//=======================================================================
+SMESH_ElementSearcher::~SMESH_ElementSearcher()
+{
+}
+
struct SMESH_ElementSearcherImpl: public SMESH_ElementSearcher
{
SMDS_Mesh* _mesh;
SMESH_ElementSearcherImpl( SMDS_Mesh& mesh, SMDS_ElemIteratorPtr elemIt=SMDS_ElemIteratorPtr())
: _mesh(&mesh),_meshPartIt(elemIt),_ebbTree(0),_nodeSearcher(0),_tolerance(-1),_outerFacesFound(false) {}
- ~SMESH_ElementSearcherImpl()
+ virtual ~SMESH_ElementSearcherImpl()
{
if ( _ebbTree ) delete _ebbTree; _ebbTree = 0;
if ( _nodeSearcher ) delete _nodeSearcher; _nodeSearcher = 0;
{
const SMDS_MeshElement* closestElem = 0;
- if ( type == SMDSAbs_Face )
+ if ( type == SMDSAbs_Face || type == SMDSAbs_Volume )
{
if ( !_ebbTree || _elementType != type )
{
{
gp_Pnt boxCenter = 0.5 * ( _ebbTree->getBox()->CornerMin() +
_ebbTree->getBox()->CornerMax() );
- double radius;
+ double radius = -1;
if ( _ebbTree->getBox()->IsOut( point.XYZ() ))
radius = point.Distance( boxCenter ) - 0.5 * _ebbTree->maxSize();
- else
+ if ( radius < 0 )
radius = _ebbTree->maxSize() / pow( 2., _ebbTree->getHeight()) / 2;
while ( suspectElems.empty() )
{
TIDSortedElemSet::iterator elem = suspectElems.begin();
for ( ; elem != suspectElems.end(); ++elem )
{
- double dist = SMESH_MeshAlgos::GetDistance( dynamic_cast<const SMDS_MeshFace*>(*elem),
- point );
+ double dist = SMESH_MeshAlgos::GetDistance( *elem, point );
if ( dist < minDist + 1e-10)
{
minDist = dist;
}
}
+//=======================================================================
+/*!
+ * \brief Return minimal distance from a point to an element
+ *
+ * Currently we ignore non-planarity and 2nd order of face
+ */
+//=======================================================================
+
+double SMESH_MeshAlgos::GetDistance( const SMDS_MeshElement* elem,
+ const gp_Pnt& point )
+{
+ switch ( elem->GetType() )
+ {
+ case SMDSAbs_Volume:
+ return GetDistance( dynamic_cast<const SMDS_MeshVolume*>( elem ), point);
+ case SMDSAbs_Face:
+ return GetDistance( dynamic_cast<const SMDS_MeshFace*>( elem ), point);
+ case SMDSAbs_Edge:
+ return GetDistance( dynamic_cast<const SMDS_MeshEdge*>( elem ), point);
+ case SMDSAbs_Node:
+ return point.Distance( SMESH_TNodeXYZ( elem ));
+ }
+ return -1;
+}
+
//=======================================================================
/*!
* \brief Return minimal distance from a point to a face
return badDistance;
}
+//=======================================================================
+/*!
+ * \brief Return minimal distance from a point to an edge
+ */
+//=======================================================================
+
+double SMESH_MeshAlgos::GetDistance( const SMDS_MeshEdge* edge, const gp_Pnt& point )
+{
+ throw SALOME_Exception(LOCALIZED("not implemented so far"));
+}
+
+//=======================================================================
+/*!
+ * \brief Return minimal distance from a point to a volume
+ *
+ * Currently we ignore non-planarity and 2nd order
+ */
+//=======================================================================
+
+double SMESH_MeshAlgos::GetDistance( const SMDS_MeshVolume* volume, const gp_Pnt& point )
+{
+ SMDS_VolumeTool vTool( volume );
+ vTool.SetExternalNormal();
+ const int iQ = volume->IsQuadratic() ? 2 : 1;
+
+ double n[3], bc[3];
+ double minDist = 1e100, dist;
+ for ( int iF = 0; iF < vTool.NbFaces(); ++iF )
+ {
+ // skip a facet with normal not "looking at" the point
+ if ( !vTool.GetFaceNormal( iF, n[0], n[1], n[2] ) ||
+ !vTool.GetFaceBaryCenter( iF, bc[0], bc[1], bc[2] ))
+ continue;
+ gp_XYZ bcp = point.XYZ() - gp_XYZ( bc[0], bc[1], bc[2] );
+ if ( gp_XYZ( n[0], n[1], n[2] ) * bcp < 1e-6 )
+ continue;
+
+ // find distance to a facet
+ const SMDS_MeshNode** nodes = vTool.GetFaceNodes( iF );
+ switch ( vTool.NbFaceNodes( iF ) / iQ ) {
+ case 3:
+ {
+ SMDS_FaceOfNodes tmpFace( nodes[0], nodes[ 1*iQ ], nodes[ 2*iQ ] );
+ dist = GetDistance( &tmpFace, point );
+ break;
+ }
+ case 4:
+ {
+ SMDS_FaceOfNodes tmpFace( nodes[0], nodes[ 1*iQ ], nodes[ 2*iQ ], nodes[ 3*iQ ]);
+ dist = GetDistance( &tmpFace, point );
+ break;
+ }
+ default:
+ vector<const SMDS_MeshNode *> nvec( nodes, nodes + vTool.NbFaceNodes( iF ));
+ SMDS_PolygonalFaceOfNodes tmpFace( nvec );
+ dist = GetDistance( &tmpFace, point );
+ }
+ minDist = Min( minDist, dist );
+ }
+ return minDist;
+}
+
+//================================================================================
+/*!
+ * \brief Returns barycentric coordinates of a point within a triangle.
+ * A not returned bc2 = 1. - bc0 - bc1.
+ * The point lies within the triangle if ( bc0 >= 0 && bc1 >= 0 && bc0+bc1 <= 1 )
+ */
+//================================================================================
+
+void SMESH_MeshAlgos::GetBarycentricCoords( const gp_XY& p,
+ const gp_XY& t0,
+ const gp_XY& t1,
+ const gp_XY& t2,
+ double & bc0,
+ double & bc1)
+{
+ const double // matrix 2x2
+ T11 = t0.X()-t2.X(), T12 = t1.X()-t2.X(),
+ T21 = t0.Y()-t2.Y(), T22 = t1.Y()-t2.Y();
+ const double Tdet = T11*T22 - T12*T21; // matrix determinant
+ if ( Abs( Tdet ) < std::numeric_limits<double>::min() )
+ {
+ bc0 = bc1 = 2.;
+ return;
+ }
+ // matrix inverse
+ const double t11 = T22, t12 = -T12, t21 = -T21, t22 = T11;
+ // vector
+ const double r11 = p.X()-t2.X(), r12 = p.Y()-t2.Y();
+ // barycentric coordinates: mutiply matrix by vector
+ bc0 = (t11 * r11 + t12 * r12)/Tdet;
+ bc1 = (t21 * r11 + t22 * r12)/Tdet;
+}
+
//=======================================================================
//function : FindFaceInSet
//purpose : Return a face having linked nodes n1 and n2 and which is