OCCT-计算二维图形的中轴线

方案一：三角剖分、重心的连接

// 获取TopoDS_Shape中所有三角形的重心并连接成中轴线
TopoDS_Shape GetTriangleCentroids(const TopoDS_Shape& shape)
{std::vector<gp_Pnt> points;TopExp_Explorer explorer(shape, TopAbs_FACE);while (explorer.More()){const TopoDS_Face& face = TopoDS::Face(explorer.Current());CalculateTriangleCentroid(points,face);explorer.Next();}// 使用BRepBuilderAPI_MakePolygon类构造PolylineBRepBuilderAPI_MakePolygon makePolygon;TColgp_Array1OfPnt PointsArray(1, points.size());int i = 0;for (auto & onept:points){makePolygon.Add(onept);PointsArray.SetValue(++i, onept);}// 获取构造的Polyline对象if (makePolygon.IsDone()){return makePolygon.Shape();}{Standard_Integer Dmin = 3;Standard_Integer Dmax = 8;Standard_Real Tol3d = 1.e-3;Handle(Geom_BSplineCurve) TheCurve;GeomAPI_PointsToBSpline aPointToBSpline(PointsArray);TheCurve = aPointToBSpline.Curve();// 创建一个AIS_Shape对象TopoDS_Shape TheCurveShape = BRepBuilderAPI_MakeEdge(TheCurve);return TheCurveShape;}return makePolygon.Edge();}

// 计算三角形的重心
void CalculateTriangleCentroid(std::vector<gp_Pnt> &points,const TopoDS_Face& face)
{TopLoc_Location loc;Handle(Poly_Triangulation) triangulation = BRep_Tool::Triangulation(face, loc);const Poly_Array1OfTriangle& triangles = triangulation->Triangles();for (size_t i = 1; i <= triangles.Size(); i++){const Poly_Triangle& triangle = triangles(i); // 假设只有一个三角形const gp_Pnt& p1 = triangulation->Node(triangle.Value(1));const gp_Pnt& p2 = triangulation->Node(triangle.Value(2));const gp_Pnt& p3 = triangulation->Node(triangle.Value(3));// 计算三个顶点的坐标平均值gp_Pnt centroid((p1.X() + p2.X() + p3.X()) / 3.0,(p1.Y() + p2.Y() + p3.Y()) / 3.0,(p1.Z() + p2.Z() + p3.Z()) / 3.0);points.emplace_back(centroid);}
}

方案二：曲线连续性

// 定义一个函数，用于提取曲线的中轴线
TopoDS_Shape ExtractCurveAxis(const Handle(Geom_Curve)& curve)
{int numPoints = 3;// 创建一个GeomLProp_CLProps对象，并将待提取中轴线的曲线作为参数传入GeomLProp_CLProps props(curve,1,1e-6);// 调用GeomLProp_CLProps类的SetParameter函数，设置参数，如精度等props.SetParameter(1e-6);// 调用GeomLProp_CLProps类的Curvature函数，获取曲线上每个点的曲率信息std::vector<Standard_Real> curvatures;for (Standard_Integer i = 1; i <= numPoints; i++) {Standard_Real u = curve->FirstParameter() + (i - 1) * (curve->LastParameter() - curve->FirstParameter()) / (numPoints - 1);props.SetParameter(u);Standard_Real curvature = props.Curvature();curvatures.push_back(curvature);}// 调用GeomLProp_CLProps类的D1函数，获取曲线上每个点的一阶导数信息std::vector<gp_Vec> derivatives;for (Standard_Integer i = 1; i <= numPoints; i++){Standard_Real u = curve->FirstParameter() + (i - 1) * (curve->LastParameter() - curve->FirstParameter()) / (numPoints - 1);props.SetParameter(u);gp_Vec derivative=props.D1();derivatives.push_back(derivative);}// 根据曲率和一阶导数信息，计算中轴线的方向和位置std::vector<gp_Pnt> axisPoints;for (Standard_Integer i = 0; i < numPoints; i++){Standard_Real u = curve->FirstParameter() + i * (curve->LastParameter() - curve->FirstParameter()) / (numPoints - 1);gp_Pnt point = curve->Value(u);gp_Vec derivative = derivatives[i];Standard_Real curvature = curvatures[i];gp_Vec axisDirection = derivative.Normalized().Crossed(gp_Vec(0, 0, 1));gp_Pnt axisPoint = point.Translated(axisDirection.Multiplied(-curvature));axisPoints.push_back(axisPoint);}// 使用BRepBuilderAPI_MakePolygon类构造PolylineBRepBuilderAPI_MakePolygon makePolygon;TColgp_Array1OfPnt PointsArray(1, axisPoints.size());int i = 0;for (auto& onept : axisPoints){makePolygon.Add(onept);PointsArray.SetValue(++i, onept);}// 获取构造的Polyline对象if (makePolygon.IsDone()){//return makePolygon.Shape();}// 根据计算得到的中轴线信息，创建一个Handle(Geom_Curve)对象//  Handle(Geom_BSplineCurve) axisCurve = new Geom_BSplineCurve(axisPoints.data(), numPoints, degree, GeomAbs_C2, knots, multiplicities, weights);// Make a BSpline curve from the points array.TColgp_Array1OfPnt aPoles(1, axisPoints.size());int nNo = 0;for (auto& pt : axisPoints){aPoles.SetValue(++nNo, pt);}Handle(Geom_BSplineCurve) aBSplineCurve = GeomAPI_PointsToBSpline(aPoles).Curve();{// 将Geom_BSplineCurve对象转换为可视化对象TopoDS_Shape BSplineShape = BRepBuilderAPI_MakeEdge(aBSplineCurve);return BSplineShape;}
}