OpenCascade中读取STL文件的相同节点合并器(源码分析)
一、概述
STL文件的格式是没有节点索引的，是以数值来表示节点的，这就导致如果读入的时候不做去重的话，读取后节点会出现重复的情况，不仅浪费内存，而且可能在运行某些算法的时候出现惊人的错误。
STL的面片格式如下
facet normal 5.223066e-02 -1.736872e-01 9.834148e-01
outer loop
vertex 3.002048e+02 9.027633e+01 1.639474e+02
vertex 2.440949e+02 8.656459e+01 1.662719e+02
vertex 2.513544e+02 6.242381e+01 1.616227e+02
endloop
endfacet
二、代码分析
OpenCascade中的RWStl_Reader.cxx实现了对STL文件的读取，并且做到了去掉重复的点，让我们看一下源码
Class MergeNodeTool用于合并相同的点
首先分析一下MergeNodeTool的父类Poly_MergeNodesTool
构造函数：
public:
//! @param[in] theSmoothAngle smooth angle in radians or 0.0 to disable merging by angle
//! @param[in] theMergeTolerance node merging maximum distance
//! @param[in] theNbFacets estimated number of facets for map preallocation
Standard_EXPORT Poly_MergeNodesTool (const double theSmoothAngle,
const double theMergeTolerance = 0.0,
const int theNbFacets = -1);
类成员定义：
private:
Handle(Poly_Triangulation) myPolyData; //!< output triangulation
MergedNodesMap myNodeIndexMap; //!< map of merged nodes
NCollection_Map<NCollection_Vec4<int>, MergedElemHasher>
myElemMap; //!< map of elements
NCollection_Vec4<int> myNodeInds; //!< current element indexes
NCollection_Vec3<float> myTriNormal; //!< current triangle normal
gp_XYZ myPlaces[4]; //!< current triangle/quad coordinates to push
Standard_Real myUnitFactor; //!< scale factor to apply
Standard_Integer myNbNodes; //!< number of output nodes
Standard_Integer myNbElems; //!< number of output elements
Standard_Integer myNbDegenElems; //!< number of degenerated elements
Standard_Integer myNbMergedElems; //!< number of merged elements
Standard_Boolean myToDropDegenerative; //!< flag to filter our degenerate elements
Standard_Boolean myToMergeElems; //!< flag to merge elements
1、如何判断为相同的节点
在vec3AreEqual中使用两种评判标准，如果myInvTol <= 0 则不启用可容忍的误差来判断两点是否相等，仅当点完全一致时为相同顶点，另一种就是加入了可以容忍的误差，当点足够近的时候，也判断为相同的节点（注：代码中的myInvTol = 1.0/可容忍的误差值，避免因为太小而被判断为 == 0.0 导致没有加入可容忍的误差值）
// =======================================================================
inline bool Poly_MergeNodesTool::MergedNodesMap::vec3AreEqual (const NCollection_Vec3<float>& theKey1,
const NCollection_Vec3<float>& theKey2) const
{
if (myInvTol <= 0.0f)
{
//! 如果myInvTol <= 0 则不启用可容忍的误差来判断两点是否相等，仅当点完全一致时为相同.
//bool IsEqual (const NCollection_Vec3& theOther) const
//{
//return v[0] == theOther.v[0]
//&& v[1] == theOther.v[1]
//&& v[2] == theOther.v[2];
//}
//调用了NCollection中的IsEqual方法 对比每一个值来判断是否相等
return theKey1.IsEqual (theKey2);
}
/// tolerance should be smaller than triangle size to avoid artifacts
//const CellVec3i anIndex1 = vec3ToCell (theKey1);
//const CellVec3i anIndex2 = vec3ToCell (theKey2);
//return anIndex1.IsEqual (anIndex2);
float aVal = theKey1.x() - theKey2.x();
if (aVal < 0) { aVal = -aVal; }
if (aVal > myTolerance) { return false; }
aVal = theKey1.y() - theKey2.y();
if (aVal < 0) { aVal = -aVal; }
if (aVal > myTolerance) { return false; }
aVal = theKey1.z() - theKey2.z();
if (aVal < 0) { aVal = -aVal; }
if (aVal > myTolerance) { return false; }
return true;
}
isEqual方法先是调用了之前点是否相等，其次是对两个点的法向量进行了点乘，由于都是单位向量，点乘的结果等于两个法向量夹角的余弦值，如果结果大于我们设定的cos值时，则说明时同一个顶点，应该合并。
当aCosinus <= -myAngleCos时，如果设定了myToMergeOpposite = true(允许合并具有相反方向法线的节点) 则合并。
（注：这玩意在等会的插入三角形的时候会用到，如果节点判断已经出现过了，则会直接返回id，而不会新增一个节点的id）
// =======================================================================
// function : MergedNodesMap::isEqual
// purpose :
// =======================================================================
inline bool Poly_MergeNodesTool::MergedNodesMap::isEqual (const Vec3AndNormal& theKey1,
const NCollection_Vec3<float>& thePos2,
const NCollection_Vec3<float>& theNorm2,
bool& theIsOpposite) const
{
if (!vec3AreEqual (theKey1.Pos, thePos2))
{
return false;
}
const float aCosinus = theKey1.Norm.Dot (theNorm2);
if (aCosinus >= myAngleCos)
{
//theIsOpposite = false;
return true;
}
else if (myToMergeOpposite
&& aCosinus <= -myAngleCos)
{
theIsOpposite = true;
return true;
}
return false;
}
2、Poly_MergeNodesTool中的MergeNodes方法
Handle(Poly_Triangulation) Poly_MergeNodesTool::MergeNodes (const Handle(Poly_Triangulation)& theTris,
const gp_Trsf& theTrsf,
const Standard_Boolean theToReverse,
const double theSmoothAngle,
const double theMergeTolerance,
const bool theToForce)
{
//判断Poly_Triangulation是否合法
if (theTris.IsNull()
|| theTris->NbNodes() < 3
|| theTris->NbTriangles() < 1)
{
//不合法，直接返回一个空的handle
return Handle(Poly_Triangulation)();
}
//初始化一个合并工具
Poly_MergeNodesTool aMergeTool (theSmoothAngle, theMergeTolerance, theTris->NbTriangles());
aMergeTool.AddTriangulation (theTris, theTrsf, theToReverse);
//theToForce为true的时候代表强制放回新的Poly_Triangulation，即便点的数量和面的数量没有改变 (没有发生合并)
if (!theToForce
&& aMergeTool.NbNodes() == theTris->NbNodes()
&& aMergeTool.NbElements() == theTris->NbTriangles())
{
return Handle(Poly_Triangulation)();
}
//返回结果
return aMergeTool.Result();
}
aMergeTool.AddTriangulation的代码如下
void Poly_MergeNodesTool::AddTriangulation (const Handle(Poly_Triangulation)& theTris,
const gp_Trsf& theTrsf,
const Standard_Boolean theToReverse)
{
if (theTris.IsNull())
{
return;
}
if (!myPolyData.IsNull()
&& myPolyData->NbNodes() == 0)
{
// 提前分配好空间
//设置精度
myPolyData->SetDoublePrecision (theTris->IsDoublePrecision());
//设置点的数量
myPolyData->ResizeNodes (theTris->NbNodes(), false);
//设置三角面的数量
myPolyData->ResizeTriangles (theTris->NbTriangles(), false);
}
//开始遍历所有的三角面
for (int anElemIter = 1; anElemIter <= theTris->NbTriangles(); ++anElemIter)
{
//Poly_Triangle中含有一个int mynodes[3]用来存放三个点的下标
//获得当前遍历的三角形
Poly_Triangle anElem = theTris->Triangle (anElemIter);
//如果需要翻转三角形，则翻转三角形
if (theToReverse)
{
anElem = Poly_Triangle (anElem.Value (1), anElem.Value (3), anElem.Value (2));
}
//遍历这个面的每个顶点，并进行空间的变换到世界坐标下
for (int aTriNodeIter = 0; aTriNodeIter < 3; ++aTriNodeIter)
{
const gp_Pnt aNode = theTris->Node (anElem.Value (aTriNodeIter + 1)).Transformed (theTrsf);
//myPlaces是类成员,可以用来存放当前三角面的三个顶点或者四边形网格的四个顶点
myPlaces[aTriNodeIter] = aNode.XYZ();
}
//开始把三角形数据插入到我们类成员的myPolyData
//这一步非常的关键
PushLastTriangle();
}
}
PushLastTriangle();这一步非常的重要，是合并代码的核心！具体流程如下
- 调用PushLastElement (3);
- 判断点的数量是否合法
- 判断是否设置了可容忍的误差(Angle和Tolerance)
3.1. 如果开启了，判断是否可以跳过角度检查（angle接近90度 ），如果不能跳过，则计算三角形的法向量并保存。调用pushNodeCheck方法对每一个点进行检查，查找在哈希表中是否存在相同或者相近的顶点，如果存在则直接放回index值，如果不存在则需要新增顶点，保存在哈希表和myPolyData中。
3.2. 如果没有开启，则直接调用pushNodeNoMerge对每个点进行插入，不进行去重操作。 - 其他一些选项的检查
- 在myPolyData中新增三角面
代码如下：
void Poly_MergeNodesTool::PushLastElement (int theNbNodes)
{
if (theNbNodes != 3
&& theNbNodes != 4)
{
throw Standard_ProgramError ("Poly_MergeNodesTool::PushLastElement() - Internal error");
}
bool isOpposite = false;
myNodeInds[3] = -1;
//判断是否具有合并操作
if (myNodeIndexMap.HasMergeAngle()
|| myNodeIndexMap.HasMergeTolerance())
{
if (!myNodeIndexMap.ToMergeAnyAngle())
{
myTriNormal = computeTriNormal();
}
//合并节点操作
pushNodeCheck (isOpposite, 0);
pushNodeCheck (isOpposite, 1);
pushNodeCheck (isOpposite, 2);
if (theNbNodes == 4)
{
pushNodeCheck (isOpposite, 3);
}
}
else
{
//直接插入操作
pushNodeNoMerge (0);
pushNodeNoMerge (1);
pushNodeNoMerge (2);
if (theNbNodes == 4)
{
pushNodeNoMerge (3);
}
}
//判断是否有面退化成边了
if (myToDropDegenerative)
{
// warning - removing degenerate elements may produce unused nodes
if (myNodeInds[0] == myNodeInds[1]
|| myNodeInds[0] == myNodeInds[2]
|| myNodeInds[1] == myNodeInds[2])
{
if (theNbNodes == 4)
{
//
}
else
{
++myNbDegenElems;
return;
}
}
}
//如果开启了合并三角面
if (myToMergeElems)
{
NCollection_Vec4<int> aSorted = myNodeInds;
std::sort (aSorted.ChangeData(), aSorted.ChangeData() + theNbNodes);
if (!myElemMap.Add (aSorted))
{
++myNbMergedElems;
return;
}
}
//在myPolyData中新增新的三角形
++myNbElems;
if (!myPolyData.IsNull())
{
if (myPolyData->NbTriangles() < myNbElems)
{
myPolyData->ResizeTriangles (myNbElems * 2, true);
}
myPolyData->SetTriangle (myNbElems, Poly_Triangle (myNodeInds[0] + 1, myNodeInds[1] + 1, myNodeInds[2] + 1));
if (theNbNodes == 4)
{
++myNbElems;
if (myPolyData->NbTriangles() < myNbElems)
{
myPolyData->ResizeTriangles (myNbElems * 2, true);
}
myPolyData->SetTriangle (myNbElems, Poly_Triangle (myNodeInds[0] + 1, myNodeInds[2] + 1, myNodeInds[3] + 1));
}
}
}
检查是否需要合并顶点操作并新增顶点
void pushNodeCheck (bool& theIsOpposite,
const int theTriNode)
{
int aNodeIndex = myNbNodes;
const gp_XYZ& aPlace = myPlaces[theTriNode];
const NCollection_Vec3<float> aVec3 ((float )aPlace.X(), (float )aPlace.Y(), (float )aPlace.Z());
if (myNodeIndexMap.Bind (aNodeIndex, theIsOpposite, aVec3, myTriNormal))
{
++myNbNodes;
if (!myPolyData.IsNull())
{
if (myPolyData->NbNodes() < myNbNodes)
{
myPolyData->ResizeNodes (myNbNodes * 2, true);
}
myPolyData->SetNode (myNbNodes, aPlace * myUnitFactor);
}
}
myNodeInds[theTriNode] = aNodeIndex;
}
// =======================================================================
// function : MergedNodesMap::Bind
// purpose :
// =======================================================================
inline bool Poly_MergeNodesTool::MergedNodesMap::Bind (int& theIndex,
bool& theIsOpposite,
const NCollection_Vec3<float>& thePos,
const NCollection_Vec3<float>& theNorm)
{
if (Resizable())
{
ReSize (Extent());
}
//NCollection_ListNode** NCollection_BaseMap::myData1
//是基础map中的一个二级指针
DataMapNode** aData = (DataMapNode** )myData1;
const int aHash = hashCode (thePos, theNorm, NbBuckets());
for (DataMapNode* aNodeIter = aData[aHash]; aNodeIter != NULL;
aNodeIter = (DataMapNode* )aNodeIter->Next())
{
//比较两个点是否相同
if (isEqual (aNodeIter->Key(), thePos, theNorm, theIsOpposite))
{
//访问可变的value值
//TheItemType& ChangeValue () { return myValue; }
theIndex = aNodeIter->ChangeValue();
//已经存在节点，无需绑定
return false;
}
}
//没有找到完全相同的点，查看是否允许误差的存在，如果存在利用误差进行判断
if (myInvTol > 0.0f)
{
//非常妙的判断，利用可容忍度的倒数，将点的坐标转化成int，并且在xyz方向上增减1来判断是否相近
//myInvTol = 1.0/可容忍度
//例：可容忍度为0.001 则myInvTol=1000，点的坐标全部乘1000
// 当点的坐标为0.45333时，会被放大到(int)453,此时454或453或452都会被认为是相同
// 的点，下面的THE_NEIGHBRS就是用来帮助判断的。
static const CellVec3i THE_NEIGHBRS[26] =
{
CellVec3i(-1, 0, 0),CellVec3i( 1, 0, 0),CellVec3i( 0,-1, 0),CellVec3i( 0, 1, 0),CellVec3i( 0, 0,-1),CellVec3i( 0, 0, 1),
CellVec3i(-1,-1, 0),CellVec3i( 1,-1, 0),CellVec3i( 1, 1, 0),CellVec3i(-1, 1, 0),
CellVec3i( 0,-1,-1),CellVec3i( 0, 1,-1),CellVec3i( 0, 1, 1),CellVec3i( 0,-1, 1),
CellVec3i(-1, 0,-1),CellVec3i( 1, 0,-1),CellVec3i( 1, 0, 1),CellVec3i(-1, 0, 1),
CellVec3i(-1,-1,-1),CellVec3i( 1,-1,-1),CellVec3i(-1, 1,-1),CellVec3i( 1, 1,-1),CellVec3i(-1,-1, 1),CellVec3i( 1,-1, 1),CellVec3i(-1, 1, 1),CellVec3i(1, 1, 1)
};
const CellVec3i anIndexCnt = vec3ToCell (thePos);
for (int aNeigIter = 0; aNeigIter < 26; ++aNeigIter)
{
const CellVec3i anIndex = anIndexCnt + THE_NEIGHBRS[aNeigIter];
const int aHashEx = vec3iHashCode (anIndex, NbBuckets());
for (DataMapNode* aNodeIter = aData[aHashEx]; aNodeIter != NULL;
aNodeIter = (DataMapNode* )aNodeIter->Next())
{
if (isEqual (aNodeIter->Key(), thePos, theNorm, theIsOpposite))
{
theIndex = aNodeIter->ChangeValue();
return false;
}
}
}
}
//theIsOpposite = false;
//如果上述的情况都不存在，则新增一个节点
aData[aHash] = new (this->myAllocator) DataMapNode (thePos, theNorm, theIndex, aData[aHash]);
//map中的mysize++
Increment();
return true;
}
至此，Poly_MergeNodesTool的大致运行机制都已经了解了，就是利用map对点来进行合并，对每个面都进行遍历，将认为是相同的点给设置相同的index，这样就保证了不会出现重合的顶点，接下来就是对Poly_MergeNodesTool的继承实现。
3、MergeNodeTool（继承于Poly_MergeNodesTool）
MergeNodeTool中重写了AddTriangle方法，主要是在原来的方法上新增了一些对于点索引的维护，同时新增了两个成员变量
private:
//一个读的指针
RWStl_Reader* myReader;
//维护顶点索引的map
NCollection_DataMap<Standard_Integer, Standard_Integer> myNodeIndexMap;
class MergeNodeTool : public Poly_MergeNodesTool
{
public:
//! Constructor
MergeNodeTool (RWStl_Reader* theReader,
const Standard_Integer theNbFacets = -1)
: Poly_MergeNodesTool (theReader->MergeAngle(), 0.0, theNbFacets),
myReader (theReader),
myNodeIndexMap (1024, new NCollection_IncAllocator (1024 * 1024))
{
// avoid redundant allocations as final triangulation is managed by RWStl_Reader subclass
ChangeOutput().Nullify();
}
//! Add new triangle
void AddTriangle (const gp_XYZ theElemNodes[3])
{
//调用父类的增加三角形方法
Poly_MergeNodesTool::AddTriangle (theElemNodes);
// remap node indices returned by RWStl_Reader::AddNode();
// 对于顺序增加是一个很耗费时间的方法, 但是保留了，因为这是STLAPI的接口
//ElementNodeIndex返回的是当前处理的三角形的顶点下标
int aNodesSrc[3] = { ElementNodeIndex (0), ElementNodeIndex (1), ElementNodeIndex (2) };
int aNodesRes[3] = { -1, -1, -1 };
for (int aNodeIter = 0; aNodeIter < 3; ++aNodeIter)
{
// 查找map中是否有对应下标已经存储过了，如果有则直接填入aNodesRes中
if (!myNodeIndexMap.Find (aNodesSrc[aNodeIter], aNodesRes[aNodeIter]))
{
//此时没有相同的点，则增加新的节点到当前模型并返回下标
aNodesRes[aNodeIter] = myReader->AddNode (theElemNodes[aNodeIter]);
//绑定原来的index和当前模型下的index
myNodeIndexMap.Bind (aNodesSrc[aNodeIter], aNodesRes[aNodeIter]);
}
}
//增加面
if (aNodesRes[0] != aNodesRes[1]
&& aNodesRes[1] != aNodesRes[2]
&& aNodesRes[2] != aNodesRes[0])
{
myReader->AddTriangle (aNodesRes[0], aNodesRes[1], aNodesRes[2]);
}
}
private:
RWStl_Reader* myReader;
NCollection_DataMap<Standard_Integer, Standard_Integer> myNodeIndexMap;
};
三、总结
整体来说这一块的设计还是非常妙的，但是代码有点难读，如果不理解的话可以自己再看看源码，基本上用到下面的文件
- RWStl_Reader.hxx ｜ RWStl_Reader.cxx
- Poly_MergeNodesTool.hxx ｜ Poly_MergeNodesTool.hxx
- NCollection_DataMap.hxx ｜ NCollection_DataMap.cxx
- Poly_Triangulation.hxx ｜ Poly_Triangulation.cxx
980
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
