mirror of
https://github.com/Open-Cascade-SAS/OCCT.git
synced 2026-05-10 09:30:48 +08:00
6f6be66842
To input parameters of ShapeConstruct_ProjectCurveOnSurface::insertAdditionalPointOrAdjust() new arrays were assigned. The size of new arrays are guaranteed to be 1 element larger that original array. NCollection_Array1::Assign() is guaranteed to throw exception when arrays size mismatch. So, given that assignment code is reached, function will always throw. Fixed by calling move assignment instead - size mismatch is allowed in this case. Also it is more optimal.
2786 lines
82 KiB
C++
2786 lines
82 KiB
C++
// Copyright (c) 1999-2014 OPEN CASCADE SAS
|
|
//
|
|
// This file is part of Open CASCADE Technology software library.
|
|
//
|
|
// This library is free software; you can redistribute it and/or modify it under
|
|
// the terms of the GNU Lesser General Public License version 2.1 as published
|
|
// by the Free Software Foundation, with special exception defined in the file
|
|
// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
|
|
// distribution for complete text of the license and disclaimer of any warranty.
|
|
//
|
|
// Alternatively, this file may be used under the terms of Open CASCADE
|
|
// commercial license or contractual agreement.
|
|
|
|
#include <ShapeConstruct_ProjectCurveOnSurface.hxx>
|
|
|
|
#include <Approx_CurveOnSurface.hxx>
|
|
#include <ElCLib.hxx>
|
|
#include <Geom2d_BSplineCurve.hxx>
|
|
#include <Geom2d_Circle.hxx>
|
|
#include <Geom2d_Curve.hxx>
|
|
#include <Geom2d_Ellipse.hxx>
|
|
#include <Geom2d_Hyperbola.hxx>
|
|
#include <Geom2d_Line.hxx>
|
|
#include <Geom2d_Parabola.hxx>
|
|
#include <Geom2d_TrimmedCurve.hxx>
|
|
#include <Geom2dAdaptor.hxx>
|
|
#include <Geom2dAPI_Interpolate.hxx>
|
|
#include <Geom_BezierSurface.hxx>
|
|
#include <Geom_BSplineCurve.hxx>
|
|
#include <Geom_BSplineSurface.hxx>
|
|
#include <Geom_Curve.hxx>
|
|
#include <Geom_OffsetSurface.hxx>
|
|
#include <Geom_Plane.hxx>
|
|
#include <Geom_RectangularTrimmedSurface.hxx>
|
|
#include <Geom_SphericalSurface.hxx>
|
|
#include <Geom_Surface.hxx>
|
|
#include <Geom_SurfaceOfLinearExtrusion.hxx>
|
|
#include <Geom_TrimmedCurve.hxx>
|
|
#include <GeomAdaptor_Curve.hxx>
|
|
#include <GeomAdaptor_Surface.hxx>
|
|
#include <GeomAPI_Interpolate.hxx>
|
|
#include <GeomAPI_PointsToBSpline.hxx>
|
|
#include <GeomProjLib.hxx>
|
|
#include <gp_Pnt2d.hxx>
|
|
#include <IntCurve_IntConicConic.hxx>
|
|
#include <IntRes2d_Domain.hxx>
|
|
#include <NCollection_IncAllocator.hxx>
|
|
#include <NCollection_Sequence.hxx>
|
|
#include <NCollection_DynamicArray.hxx>
|
|
#include <Precision.hxx>
|
|
#include <ProjLib_ProjectedCurve.hxx>
|
|
#include <ShapeAnalysis.hxx>
|
|
#include <ShapeAnalysis_Curve.hxx>
|
|
#include <ShapeAnalysis_Surface.hxx>
|
|
#include <ShapeExtend.hxx>
|
|
#include <Standard_ErrorHandler.hxx>
|
|
#include <Standard_Failure.hxx>
|
|
#include <Standard_Type.hxx>
|
|
#include <gp_Pnt.hxx>
|
|
#include <NCollection_Array1.hxx>
|
|
#include <NCollection_HArray1.hxx>
|
|
#include <Standard_Integer.hxx>
|
|
|
|
#include <algorithm>
|
|
|
|
IMPLEMENT_STANDARD_RTTIEXT(ShapeConstruct_ProjectCurveOnSurface, Standard_Transient)
|
|
|
|
namespace
|
|
{
|
|
//! Default number of control points for discretization.
|
|
constexpr int THE_NCONTROL = 23;
|
|
|
|
//! Utility class for projecting points onto a surface with B-spline corner cache optimization.
|
|
//! For clamped B-spline surfaces, caches corner pole positions and their exact UV parameters
|
|
//! to avoid expensive ValueOfUV calls when projected points coincide with surface corners.
|
|
//! Only corner poles are cached because B-spline surfaces pass through corners only
|
|
//! (when end multiplicities equal degree + 1).
|
|
class SurfaceProjectorWithCache
|
|
{
|
|
public:
|
|
//! Constructor - initializes the projector with a surface.
|
|
//! For B-spline surfaces, builds the corner cache automatically.
|
|
//! @param[in] theSurf the surface analysis object
|
|
SurfaceProjectorWithCache(const occ::handle<ShapeAnalysis_Surface>& theSurf)
|
|
: mySurf(theSurf)
|
|
{
|
|
if (theSurf.IsNull())
|
|
{
|
|
return;
|
|
}
|
|
|
|
occ::handle<Geom_BSplineSurface> aBSplineSurf =
|
|
occ::down_cast<Geom_BSplineSurface>(theSurf->Surface());
|
|
if (!aBSplineSurf.IsNull())
|
|
{
|
|
buildCornerCache(aBSplineSurf);
|
|
}
|
|
}
|
|
|
|
//! Projects a 3D point onto the surface.
|
|
//! First checks B-spline corner cache, then falls back to ValueOfUV.
|
|
//! @param[in] thePoint the 3D point to project
|
|
//! @param[in] theTol the tolerance for projection
|
|
//! @param[in] theTolSq squared tolerance for corner matching
|
|
//! @return the UV coordinates on the surface
|
|
gp_Pnt2d ValueOfUV(const gp_Pnt& thePoint, const double theTol, const double theTolSq) const
|
|
{
|
|
gp_Pnt2d aResult;
|
|
if (findInCornerCache(thePoint, theTolSq, aResult))
|
|
{
|
|
// Corner UV is exact, but refine for numerical stability
|
|
return mySurf->NextValueOfUV(aResult, thePoint, theTol, theTol);
|
|
}
|
|
return mySurf->ValueOfUV(thePoint, theTol);
|
|
}
|
|
|
|
//! Projects a 3D point onto the surface using a hint from previous projection.
|
|
//! First checks B-spline corner cache, then falls back to NextValueOfUV.
|
|
//! @param[in] theHint the UV hint from previous projection
|
|
//! @param[in] thePoint the 3D point to project
|
|
//! @param[in] theTol the tolerance for projection
|
|
//! @param[in] theTolSq squared tolerance for corner matching
|
|
//! @param[in] theStep the step tolerance
|
|
//! @return the UV coordinates on the surface
|
|
gp_Pnt2d NextValueOfUV(const gp_Pnt2d& theHint,
|
|
const gp_Pnt& thePoint,
|
|
const double theTol,
|
|
const double theTolSq,
|
|
const double theStep) const
|
|
{
|
|
gp_Pnt2d aResult;
|
|
if (findInCornerCache(thePoint, theTolSq, aResult))
|
|
{
|
|
// Corner UV is exact, but refine for numerical stability
|
|
return mySurf->NextValueOfUV(aResult, thePoint, theTol, theTol);
|
|
}
|
|
return mySurf->NextValueOfUV(theHint, thePoint, theTol, theStep);
|
|
}
|
|
|
|
//! Returns the gap from the last projection
|
|
double Gap() const { return mySurf->Gap(); }
|
|
|
|
private:
|
|
//! Builds cache from B-spline surface corner poles.
|
|
//! For clamped B-splines (multiplicity = degree + 1 at ends), the surface passes
|
|
//! exactly through corner poles, so we can use their UV parameters directly.
|
|
//! @param[in] theSurface the B-spline surface
|
|
void buildCornerCache(const occ::handle<Geom_BSplineSurface>& theSurface)
|
|
{
|
|
const int aNbPolesU = theSurface->NbUPoles();
|
|
const int aNbPolesV = theSurface->NbVPoles();
|
|
const int aDegreeU = theSurface->UDegree();
|
|
const int aDegreeV = theSurface->VDegree();
|
|
|
|
// Check if surface is clamped (end multiplicities = degree + 1)
|
|
const int aFirstUMult = theSurface->UMultiplicity(1);
|
|
const int aLastUMult = theSurface->UMultiplicity(theSurface->NbUKnots());
|
|
const int aFirstVMult = theSurface->VMultiplicity(1);
|
|
const int aLastVMult = theSurface->VMultiplicity(theSurface->NbVKnots());
|
|
|
|
const bool isUFirstClamped = (aFirstUMult >= aDegreeU + 1);
|
|
const bool isULastClamped = (aLastUMult >= aDegreeU + 1);
|
|
const bool isVFirstClamped = (aFirstVMult >= aDegreeV + 1);
|
|
const bool isVLastClamped = (aLastVMult >= aDegreeV + 1);
|
|
|
|
// Get parameter bounds
|
|
const double aUFirst = theSurface->UKnot(1);
|
|
const double aULast = theSurface->UKnot(theSurface->NbUKnots());
|
|
const double aVFirst = theSurface->VKnot(1);
|
|
const double aVLast = theSurface->VKnot(theSurface->NbVKnots());
|
|
|
|
// Cache corner poles where surface passes through them
|
|
// Corner (1, 1) - UFirst, VFirst
|
|
if (isUFirstClamped && isVFirstClamped)
|
|
{
|
|
myCorners3d.Append(theSurface->Pole(1, 1));
|
|
myCorners2d.Append(gp_Pnt2d(aUFirst, aVFirst));
|
|
}
|
|
|
|
// Corner (NbPolesU, 1) - ULast, VFirst
|
|
if (isULastClamped && isVFirstClamped)
|
|
{
|
|
myCorners3d.Append(theSurface->Pole(aNbPolesU, 1));
|
|
myCorners2d.Append(gp_Pnt2d(aULast, aVFirst));
|
|
}
|
|
|
|
// Corner (1, NbPolesV) - UFirst, VLast
|
|
if (isUFirstClamped && isVLastClamped)
|
|
{
|
|
myCorners3d.Append(theSurface->Pole(1, aNbPolesV));
|
|
myCorners2d.Append(gp_Pnt2d(aUFirst, aVLast));
|
|
}
|
|
|
|
// Corner (NbPolesU, NbPolesV) - ULast, VLast
|
|
if (isULastClamped && isVLastClamped)
|
|
{
|
|
myCorners3d.Append(theSurface->Pole(aNbPolesU, aNbPolesV));
|
|
myCorners2d.Append(gp_Pnt2d(aULast, aVLast));
|
|
}
|
|
}
|
|
|
|
//! Finds the closest corner to a 3D point within tolerance.
|
|
//! @param[in] thePoint the 3D point to find
|
|
//! @param[in] theTolSq squared tolerance for matching
|
|
//! @param[out] theUV the UV parameter if found
|
|
//! @return true if a matching corner was found
|
|
bool findInCornerCache(const gp_Pnt& thePoint, const double theTolSq, gp_Pnt2d& theUV) const
|
|
{
|
|
for (int i = 0; i < myCorners3d.Length(); ++i)
|
|
{
|
|
if (myCorners3d(i).SquareDistance(thePoint) < theTolSq)
|
|
{
|
|
theUV = myCorners2d(i);
|
|
return true;
|
|
}
|
|
}
|
|
return false;
|
|
}
|
|
|
|
private:
|
|
occ::handle<ShapeAnalysis_Surface> mySurf; //!< Surface to project on
|
|
NCollection_DynamicArray<gp_Pnt> myCorners3d; //!< 3D positions of B-spline surface corners
|
|
NCollection_DynamicArray<gp_Pnt2d> myCorners2d; //!< UV parameters of B-spline surface corners
|
|
};
|
|
|
|
//=================================================================================================
|
|
|
|
//! Extracts B-spline curve from a possibly nested trimmed curve.
|
|
//! Recursively unwraps trimmed curves to find B-spline basis.
|
|
//! @param[in] theCurve the curve to extract from
|
|
//! @return the extracted B-spline curve, or null if not a B-spline
|
|
occ::handle<Geom_BSplineCurve> extractBSplineCurve(const occ::handle<Geom_Curve>& theCurve)
|
|
{
|
|
occ::handle<Geom_Curve> aCurve = theCurve;
|
|
|
|
// Recursively unwrap trimmed curves
|
|
while (!aCurve.IsNull() && aCurve->IsKind(STANDARD_TYPE(Geom_TrimmedCurve)))
|
|
{
|
|
aCurve = occ::down_cast<Geom_TrimmedCurve>(aCurve)->BasisCurve();
|
|
}
|
|
|
|
return occ::down_cast<Geom_BSplineCurve>(aCurve);
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
//! Adjusts second point to first point considering surface periodicity.
|
|
//! For periodic surfaces, adjusts the second point coordinates to be within
|
|
//! half-period of the first point.
|
|
//! @param[in] theFirstPoint the reference point
|
|
//! @param[in,out] theSecondPoint the point to adjust
|
|
//! @param[in] theSurf the surface to check for periodicity
|
|
void adjustSecondToFirstPoint(const gp_Pnt2d& theFirstPoint,
|
|
gp_Pnt2d& theSecondPoint,
|
|
const occ::handle<Geom_Surface>& theSurf)
|
|
{
|
|
if (theSurf->IsUPeriodic())
|
|
{
|
|
const double anUPeriod = theSurf->UPeriod();
|
|
const double aNewU = ElCLib::InPeriod(theSecondPoint.X(),
|
|
theFirstPoint.X() - anUPeriod / 2,
|
|
theFirstPoint.X() + anUPeriod / 2);
|
|
theSecondPoint.SetX(aNewU);
|
|
}
|
|
if (theSurf->IsVPeriodic())
|
|
{
|
|
const double aVPeriod = theSurf->VPeriod();
|
|
const double aNewV = ElCLib::InPeriod(theSecondPoint.Y(),
|
|
theFirstPoint.Y() - aVPeriod / 2,
|
|
theFirstPoint.Y() + aVPeriod / 2);
|
|
theSecondPoint.SetY(aNewV);
|
|
}
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
//! Fixes possible period jumps in an array of 4 points.
|
|
//! Handles the walking period parameter to ensure smooth transition
|
|
//! across periodic boundaries.
|
|
//! @param[in,out] thePnt array of 4 points to fix
|
|
//! @param[in] theIdx coordinate index (1 for U, 2 for V)
|
|
//! @param[in] thePeriod the period value
|
|
//! @param[in] theSavedPoint index of the reference point (-1 for none)
|
|
//! @param[in] theSavedParam the saved parameter value
|
|
//! @return true if a period jump was detected and fixed
|
|
bool fixPeriodicityTroubles(gp_Pnt2d* thePnt,
|
|
const int theIdx,
|
|
const double thePeriod,
|
|
const int theSavedPoint,
|
|
const double theSavedParam)
|
|
{
|
|
double aSavedParam;
|
|
int aSavedPoint;
|
|
double aMinParam = 0.0, aMaxParam = thePeriod;
|
|
|
|
if (theSavedPoint < 0)
|
|
{
|
|
aSavedParam = 0.5 * thePeriod;
|
|
aSavedPoint = 0;
|
|
}
|
|
else
|
|
{
|
|
aSavedParam = theSavedParam;
|
|
aSavedPoint = theSavedPoint;
|
|
while (aMinParam > aSavedParam)
|
|
{
|
|
aMinParam -= thePeriod;
|
|
aMaxParam -= thePeriod;
|
|
}
|
|
while (aMaxParam < aSavedParam)
|
|
{
|
|
aMinParam += thePeriod;
|
|
aMaxParam += thePeriod;
|
|
}
|
|
}
|
|
|
|
double aFixIsoParam = aMinParam;
|
|
bool isIsoLine = false;
|
|
if (aMaxParam - aSavedParam < Precision::PConfusion()
|
|
|| aSavedParam - aMinParam < Precision::PConfusion())
|
|
{
|
|
aFixIsoParam = aSavedParam;
|
|
isIsoLine = true;
|
|
}
|
|
|
|
// Normalize all coordinates to [aMinParam, aMaxParam)
|
|
for (int i = 0; i < 4; i++)
|
|
{
|
|
double aParam = thePnt[i].Coord(theIdx);
|
|
double aShift = ShapeAnalysis::AdjustToPeriod(aParam, aMinParam, aMaxParam);
|
|
aParam += aShift;
|
|
|
|
if (isIsoLine)
|
|
{
|
|
if (aMaxParam - aParam < Precision::PConfusion()
|
|
|| aParam - aMinParam < Precision::PConfusion())
|
|
{
|
|
aParam = aFixIsoParam;
|
|
}
|
|
}
|
|
else
|
|
{
|
|
if (aMaxParam - aParam < Precision::PConfusion())
|
|
{
|
|
aParam = aMaxParam;
|
|
}
|
|
if (aParam - aMinParam < Precision::PConfusion())
|
|
{
|
|
aParam = aMinParam;
|
|
}
|
|
}
|
|
|
|
thePnt[i].SetCoord(theIdx, aParam);
|
|
}
|
|
|
|
// Find possible period jump and increasing or decreasing coordinates vector
|
|
bool isJump = false;
|
|
double aPrevDiff = 0.0;
|
|
double aSumDiff = 1.0;
|
|
for (int i = 0; i < 3; i++)
|
|
{
|
|
double aDiff = thePnt[i + 1].Coord(theIdx) - thePnt[i].Coord(theIdx);
|
|
if (aDiff < -Precision::PConfusion())
|
|
{
|
|
aSumDiff *= -1.0;
|
|
}
|
|
if (aDiff * aPrevDiff < -Precision::PConfusion())
|
|
{
|
|
isJump = true;
|
|
}
|
|
aPrevDiff = aDiff;
|
|
}
|
|
|
|
if (!isJump)
|
|
{
|
|
return false;
|
|
}
|
|
|
|
if (aSumDiff > 0)
|
|
{
|
|
for (int i = aSavedPoint; i > 0; i--)
|
|
{
|
|
if (thePnt[i].Coord(theIdx) > thePnt[i - 1].Coord(theIdx))
|
|
{
|
|
thePnt[i - 1].SetCoord(theIdx, thePnt[i - 1].Coord(theIdx) + thePeriod);
|
|
}
|
|
}
|
|
for (int i = aSavedPoint; i < 3; i++)
|
|
{
|
|
if (thePnt[i].Coord(theIdx) < thePnt[i + 1].Coord(theIdx))
|
|
{
|
|
thePnt[i + 1].SetCoord(theIdx, thePnt[i + 1].Coord(theIdx) - thePeriod);
|
|
}
|
|
}
|
|
}
|
|
else
|
|
{
|
|
for (int i = aSavedPoint; i > 0; i--)
|
|
{
|
|
if (thePnt[i].Coord(theIdx) < thePnt[i - 1].Coord(theIdx))
|
|
{
|
|
thePnt[i - 1].SetCoord(theIdx, thePnt[i - 1].Coord(theIdx) - thePeriod);
|
|
}
|
|
}
|
|
for (int i = aSavedPoint; i < 3; i++)
|
|
{
|
|
if (thePnt[i].Coord(theIdx) > thePnt[i + 1].Coord(theIdx))
|
|
{
|
|
thePnt[i + 1].SetCoord(theIdx, thePnt[i + 1].Coord(theIdx) + thePeriod);
|
|
}
|
|
}
|
|
}
|
|
|
|
return true;
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
//! Checks if B-spline curve has uneven parameterization requiring special handling.
|
|
//! Computes the ratio of maximum to minimum parameterization speed across knot intervals.
|
|
//! If this ratio exceeds a threshold, the curve should be projected using ProjLib instead
|
|
//! of the standard approximation approach.
|
|
//! @param[in] theCurve the curve to check (may be trimmed)
|
|
//! @param[in] theFirst the first parameter of the working range
|
|
//! @param[in] theLast the last parameter of the working range
|
|
//! @param[out] theBSpline the extracted B-spline curve (if any)
|
|
//! @return true if the curve has uneven parameterization requiring ProjLib
|
|
bool isBSplineCurveInvalid(const occ::handle<Geom_Curve>& theCurve,
|
|
const double theFirst,
|
|
const double theLast,
|
|
occ::handle<Geom_BSplineCurve>& theBSpline)
|
|
{
|
|
theBSpline = extractBSplineCurve(theCurve);
|
|
if (theBSpline.IsNull())
|
|
{
|
|
return false;
|
|
}
|
|
|
|
// Compute parametrization speed on each knot interval inside [theFirst, theLast].
|
|
// If quotient = (MaxSpeed / MinSpeed) >= aMaxQuotientCoeff then use PerformByProjLib.
|
|
double aFirstParam = theFirst;
|
|
double aLastParam = theLast;
|
|
|
|
// Find first knot index
|
|
int anIdx = 1;
|
|
for (; anIdx <= theBSpline->NbKnots() && aFirstParam < theLast; anIdx++)
|
|
{
|
|
if (theBSpline->Knot(anIdx) > theFirst)
|
|
{
|
|
break;
|
|
}
|
|
}
|
|
|
|
GeomAdaptor_Curve aC3DAdaptor(theCurve);
|
|
double aMinParSpeed = Precision::Infinite();
|
|
NCollection_Sequence<double> aKnotCoeffs;
|
|
|
|
for (; anIdx <= theBSpline->NbKnots() && aFirstParam < theLast; anIdx++)
|
|
{
|
|
// Fill current knot interval
|
|
aLastParam = std::min(theLast, theBSpline->Knot(anIdx));
|
|
int aNbIntPnts = THE_NCONTROL;
|
|
|
|
// Adapt number of inner points according to the length of the interval
|
|
// to avoid a lot of calculations on small range of parameters.
|
|
if (anIdx > 1)
|
|
{
|
|
const double aLenThres = 1.e-2;
|
|
const double aLenRatio =
|
|
(aLastParam - aFirstParam) / (theBSpline->Knot(anIdx) - theBSpline->Knot(anIdx - 1));
|
|
if (aLenRatio < aLenThres)
|
|
{
|
|
aNbIntPnts = int(aLenRatio / aLenThres * aNbIntPnts);
|
|
if (aNbIntPnts < 2)
|
|
{
|
|
aNbIntPnts = 2;
|
|
}
|
|
}
|
|
}
|
|
|
|
double aStep = (aLastParam - aFirstParam) / (aNbIntPnts - 1);
|
|
gp_Pnt p3d1, p3d2;
|
|
|
|
// Start filling from first point
|
|
aC3DAdaptor.D0(aFirstParam, p3d1);
|
|
|
|
double aLength3d = 0.0;
|
|
for (int anIntIdx = 1; anIntIdx < aNbIntPnts; anIntIdx++)
|
|
{
|
|
double aParam = aFirstParam + aStep * anIntIdx;
|
|
aC3DAdaptor.D0(aParam, p3d2);
|
|
const double aDist = p3d2.Distance(p3d1);
|
|
|
|
aLength3d += aDist;
|
|
p3d1 = p3d2;
|
|
|
|
aMinParSpeed = std::min(aMinParSpeed, aDist / aStep);
|
|
}
|
|
|
|
const double aCoeff = aLength3d / (aLastParam - aFirstParam);
|
|
if (std::abs(aCoeff) > gp::Resolution())
|
|
{
|
|
aKnotCoeffs.Append(aCoeff);
|
|
}
|
|
aFirstParam = aLastParam;
|
|
}
|
|
|
|
double anEvenlyCoeff = 0;
|
|
if (aKnotCoeffs.Size() > 0)
|
|
{
|
|
anEvenlyCoeff = *std::max_element(aKnotCoeffs.begin(), aKnotCoeffs.end())
|
|
/ *std::min_element(aKnotCoeffs.begin(), aKnotCoeffs.end());
|
|
}
|
|
|
|
const double aMaxQuotientCoeff = 1500.0;
|
|
return (anEvenlyCoeff > aMaxQuotientCoeff && aMinParSpeed > Precision::Confusion());
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
//! Generates discretization points for a curve.
|
|
//! Uses uniform distribution for general curves, and adjusts the number of points
|
|
//! for B-splines based on their knot structure to ensure adequate sampling.
|
|
//! @param[in] theCurve the curve to discretize
|
|
//! @param[in] theFirst the first parameter
|
|
//! @param[in] theLast the last parameter
|
|
//! @param[in] theNbControlPoints the minimum number of points
|
|
//! @param[out] thePoints the generated 3D points (resized as needed)
|
|
//! @param[out] theParams the corresponding parameters (resized as needed)
|
|
//! @return the actual number of generated points
|
|
int generateCurvePoints(const occ::handle<Geom_Curve>& theCurve,
|
|
const double theFirst,
|
|
const double theLast,
|
|
const int theNbControlPoints,
|
|
ShapeConstruct_ProjectCurveOnSurface::ArrayOfPnt& thePoints,
|
|
ShapeConstruct_ProjectCurveOnSurface::ArrayOfReal& theParams)
|
|
{
|
|
occ::handle<Geom_BSplineCurve> aBSpline = extractBSplineCurve(theCurve);
|
|
int aNbPini = theNbControlPoints;
|
|
|
|
if (!aBSpline.IsNull())
|
|
{
|
|
int aUsedKnots = 0;
|
|
for (int i = 1; i < aBSpline->NbKnots(); i++)
|
|
{
|
|
if (aBSpline->Knot(i + 1) > theFirst && aBSpline->Knot(i) < theLast)
|
|
{
|
|
aUsedKnots++;
|
|
}
|
|
}
|
|
int aMinPnt = aUsedKnots * (aBSpline->Degree() + 1);
|
|
while (aNbPini < aMinPnt)
|
|
{
|
|
aNbPini += THE_NCONTROL - 1;
|
|
}
|
|
}
|
|
|
|
thePoints.Resize(1, aNbPini, false);
|
|
theParams.Resize(1, aNbPini, false);
|
|
|
|
const double aDeltaParam = (theLast - theFirst) / (aNbPini - 1);
|
|
|
|
for (int i = 1; i <= aNbPini; ++i)
|
|
{
|
|
double aParam;
|
|
if (i == 1)
|
|
{
|
|
aParam = theFirst;
|
|
}
|
|
else if (i == aNbPini)
|
|
{
|
|
aParam = theLast;
|
|
}
|
|
else
|
|
{
|
|
aParam = theFirst + (i - 1) * aDeltaParam;
|
|
}
|
|
|
|
gp_Pnt aPoint;
|
|
theCurve->D0(aParam, aPoint);
|
|
thePoints.SetValue(i, aPoint);
|
|
theParams.SetValue(i, aParam);
|
|
}
|
|
|
|
return aNbPini;
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
//! Wrapper for ShapeAnalysis_Surface::ProjectDegenerated.
|
|
//! Converts NCollection_Array1 containers to sequences, performs projection,
|
|
//! then copies results back. Required because ShapeAnalysis_Surface uses
|
|
//! sequence containers.
|
|
//! @param[in] theSurf the analysis surface
|
|
//! @param[in] theNbPnt the number of points
|
|
//! @param[in] thePoints the 3D points array
|
|
//! @param[in,out] thePoints2d the 2D points array (adjusted for degenerate regions)
|
|
//! @param[in] thePreci the precision
|
|
//! @param[in] theDirect true to check from start, false to check from end
|
|
void projectDegeneratedPoints(const occ::handle<ShapeAnalysis_Surface>& theSurf,
|
|
const int theNbPnt,
|
|
const NCollection_Array1<gp_Pnt>& thePoints,
|
|
NCollection_Array1<gp_Pnt2d>& thePoints2d,
|
|
const double thePreci,
|
|
const bool theDirect)
|
|
{
|
|
// Use incremental allocator for sequences when there are enough elements
|
|
// to benefit from pooled memory allocation.
|
|
constexpr int THE_ALLOC_THRESHOLD = 100;
|
|
|
|
occ::handle<NCollection_BaseAllocator> anAlloc;
|
|
if (theNbPnt > THE_ALLOC_THRESHOLD)
|
|
{
|
|
anAlloc = new NCollection_IncAllocator(NCollection_IncAllocator::THE_MINIMUM_BLOCK_SIZE);
|
|
}
|
|
|
|
// Convert arrays to sequences for ShapeAnalysis_Surface::ProjectDegenerated
|
|
NCollection_Sequence<gp_Pnt> aPoints3d(anAlloc);
|
|
NCollection_Sequence<gp_Pnt2d> aPoints2d(anAlloc);
|
|
|
|
for (int i = 1; i <= theNbPnt; ++i)
|
|
{
|
|
aPoints3d.Append(thePoints(i));
|
|
aPoints2d.Append(thePoints2d(i));
|
|
}
|
|
|
|
// Call the method that expects sequences
|
|
theSurf->ProjectDegenerated(theNbPnt, aPoints3d, aPoints2d, thePreci, theDirect);
|
|
|
|
// Copy results back to array
|
|
for (int i = 1; i <= theNbPnt; ++i)
|
|
{
|
|
thePoints2d.SetValue(i, aPoints2d.Value(i));
|
|
}
|
|
}
|
|
|
|
} // namespace
|
|
|
|
//=================================================================================================
|
|
|
|
ShapeConstruct_ProjectCurveOnSurface::ShapeConstruct_ProjectCurveOnSurface()
|
|
: myPreci(Precision::Confusion()),
|
|
myStatus(ShapeExtend::EncodeStatus(ShapeExtend_OK)),
|
|
myAdjustOverDegen(1)
|
|
{
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
void ShapeConstruct_ProjectCurveOnSurface::Init(const occ::handle<Geom_Surface>& theSurf,
|
|
const double thePreci)
|
|
{
|
|
Init(new ShapeAnalysis_Surface(theSurf), thePreci);
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
void ShapeConstruct_ProjectCurveOnSurface::Init(const occ::handle<ShapeAnalysis_Surface>& theSurf,
|
|
const double thePreci)
|
|
{
|
|
SetSurface(theSurf);
|
|
SetPrecision(thePreci);
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
void ShapeConstruct_ProjectCurveOnSurface::SetSurface(const occ::handle<Geom_Surface>& theSurf)
|
|
{
|
|
SetSurface(new ShapeAnalysis_Surface(theSurf));
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
void ShapeConstruct_ProjectCurveOnSurface::SetSurface(
|
|
const occ::handle<ShapeAnalysis_Surface>& theSurf)
|
|
{
|
|
if (mySurf == theSurf)
|
|
{
|
|
return;
|
|
}
|
|
mySurf = theSurf;
|
|
myCache = CacheArray();
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
void ShapeConstruct_ProjectCurveOnSurface::SetPrecision(const double thePreci)
|
|
{
|
|
myPreci = thePreci;
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
int& ShapeConstruct_ProjectCurveOnSurface::AdjustOverDegenMode()
|
|
{
|
|
return myAdjustOverDegen;
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
bool ShapeConstruct_ProjectCurveOnSurface::Status(const ShapeExtend_Status theStatus) const
|
|
{
|
|
return ShapeExtend::DecodeStatus(myStatus, theStatus);
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
bool ShapeConstruct_ProjectCurveOnSurface::Perform(const occ::handle<Geom_Curve>& theC3D,
|
|
const double theFirst,
|
|
const double theLast,
|
|
occ::handle<Geom2d_Curve>& theC2D,
|
|
const double theTolFirst,
|
|
const double theTolLast)
|
|
{
|
|
myStatus = ShapeExtend::EncodeStatus(ShapeExtend_OK);
|
|
|
|
if (mySurf.IsNull() || theC3D.IsNull())
|
|
{
|
|
theC2D.Nullify();
|
|
myStatus |= ShapeExtend::EncodeStatus(ShapeExtend_FAIL1);
|
|
return false;
|
|
}
|
|
|
|
// Try analytical projection first
|
|
occ::handle<Geom_Curve> aCurve3DTrim = theC3D;
|
|
if (!theC3D->IsKind(STANDARD_TYPE(Geom_BoundedCurve)))
|
|
{
|
|
aCurve3DTrim = new Geom_TrimmedCurve(theC3D, theFirst, theLast);
|
|
}
|
|
|
|
theC2D = projectAnalytic(aCurve3DTrim);
|
|
if (!theC2D.IsNull())
|
|
{
|
|
myStatus |= ShapeExtend::EncodeStatus(ShapeExtend_DONE1);
|
|
return true;
|
|
}
|
|
|
|
// Handle problematic B-spline curves with uneven parameterization
|
|
occ::handle<Geom_BSplineCurve> aBSpline;
|
|
if (isBSplineCurveInvalid(theC3D, theFirst, theLast, aBSpline))
|
|
{
|
|
// Use ProjLib for curves with highly uneven parameterization speed
|
|
if (PerformByProjLib(theC3D, theFirst, theLast, theC2D))
|
|
{
|
|
return Status(ShapeExtend_DONE);
|
|
}
|
|
}
|
|
|
|
// Generate curve points
|
|
ArrayOfPnt aPoints;
|
|
ArrayOfReal aParams;
|
|
const int aNbPnt = generateCurvePoints(theC3D, theFirst, theLast, THE_NCONTROL, aPoints, aParams);
|
|
|
|
// Approximate pcurve
|
|
ArrayOfPnt2d aPoints2d(1, aNbPnt);
|
|
for (int i = 1; i <= aNbPnt; i++)
|
|
{
|
|
aPoints2d.SetValue(i, gp_Pnt2d(0., 0.));
|
|
}
|
|
|
|
approxPCurve(aNbPnt, theC3D, theTolFirst, theTolLast, aPoints, aParams, aPoints2d, theC2D);
|
|
|
|
const int aNbPini = aPoints.Length();
|
|
if (!theC2D.IsNull())
|
|
{
|
|
myStatus |= ShapeExtend::EncodeStatus(ShapeExtend_DONE2);
|
|
return true;
|
|
}
|
|
|
|
// Interpolate the result
|
|
theC2D = interpolatePCurve(aNbPini, aPoints2d, aParams);
|
|
|
|
myStatus |= ShapeExtend::EncodeStatus(theC2D.IsNull() ? ShapeExtend_FAIL1 : ShapeExtend_DONE2);
|
|
return Status(ShapeExtend_DONE);
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
bool ShapeConstruct_ProjectCurveOnSurface::PerformByProjLib(const occ::handle<Geom_Curve>& theC3D,
|
|
const double theFirst,
|
|
const double theLast,
|
|
occ::handle<Geom2d_Curve>& theC2D)
|
|
{
|
|
theC2D.Nullify();
|
|
if (mySurf.IsNull())
|
|
{
|
|
myStatus = ShapeExtend::EncodeStatus(ShapeExtend_FAIL1);
|
|
return false;
|
|
}
|
|
|
|
try
|
|
{
|
|
OCC_CATCH_SIGNALS
|
|
occ::handle<GeomAdaptor_Surface> aGAS = mySurf->Adaptor3d();
|
|
occ::handle<GeomAdaptor_Curve> aGAC = new GeomAdaptor_Curve(theC3D, theFirst, theLast);
|
|
ProjLib_ProjectedCurve aProjector(aGAS, aGAC);
|
|
|
|
switch (aProjector.GetType())
|
|
{
|
|
case GeomAbs_Line:
|
|
theC2D = new Geom2d_Line(aProjector.Line());
|
|
break;
|
|
case GeomAbs_Circle:
|
|
theC2D = new Geom2d_Circle(aProjector.Circle());
|
|
break;
|
|
case GeomAbs_Ellipse:
|
|
theC2D = new Geom2d_Ellipse(aProjector.Ellipse());
|
|
break;
|
|
case GeomAbs_Parabola:
|
|
theC2D = new Geom2d_Parabola(aProjector.Parabola());
|
|
break;
|
|
case GeomAbs_Hyperbola:
|
|
theC2D = new Geom2d_Hyperbola(aProjector.Hyperbola());
|
|
break;
|
|
case GeomAbs_BSplineCurve:
|
|
theC2D = aProjector.BSpline();
|
|
break;
|
|
default:
|
|
break;
|
|
}
|
|
|
|
if (theC2D.IsNull())
|
|
{
|
|
myStatus = ShapeExtend::EncodeStatus(ShapeExtend_FAIL2);
|
|
return false;
|
|
}
|
|
else
|
|
{
|
|
myStatus = ShapeExtend::EncodeStatus(ShapeExtend_DONE1);
|
|
return true;
|
|
}
|
|
}
|
|
catch (Standard_Failure const& anException)
|
|
{
|
|
#ifdef OCCT_DEBUG
|
|
std::cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::PerformByProjLib(): Exception: ";
|
|
anException.Print(std::cout);
|
|
std::cout << std::endl;
|
|
#endif
|
|
(void)anException;
|
|
myStatus = ShapeExtend::EncodeStatus(ShapeExtend_FAIL3);
|
|
theC2D.Nullify();
|
|
}
|
|
return false;
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
occ::handle<Geom2d_Curve> ShapeConstruct_ProjectCurveOnSurface::projectAnalytic(
|
|
const occ::handle<Geom_Curve>& theC3D) const
|
|
{
|
|
occ::handle<Geom2d_Curve> aResult;
|
|
|
|
occ::handle<Geom_Surface> aSurf = mySurf->Surface();
|
|
occ::handle<Geom_Plane> aPlane = occ::down_cast<Geom_Plane>(aSurf);
|
|
|
|
if (aPlane.IsNull())
|
|
{
|
|
occ::handle<Geom_RectangularTrimmedSurface> aRTS =
|
|
occ::down_cast<Geom_RectangularTrimmedSurface>(aSurf);
|
|
if (!aRTS.IsNull())
|
|
{
|
|
aPlane = occ::down_cast<Geom_Plane>(aRTS->BasisSurface());
|
|
}
|
|
else
|
|
{
|
|
occ::handle<Geom_OffsetSurface> anOS = occ::down_cast<Geom_OffsetSurface>(aSurf);
|
|
if (!anOS.IsNull())
|
|
{
|
|
aPlane = occ::down_cast<Geom_Plane>(anOS->BasisSurface());
|
|
}
|
|
}
|
|
}
|
|
|
|
if (!aPlane.IsNull())
|
|
{
|
|
occ::handle<Geom_Curve> aProjOnPlane =
|
|
GeomProjLib::ProjectOnPlane(theC3D, aPlane, aPlane->Position().Direction(), true);
|
|
occ::handle<GeomAdaptor_Curve> aHC = new GeomAdaptor_Curve(aProjOnPlane);
|
|
ProjLib_ProjectedCurve aProj(mySurf->Adaptor3d(), aHC);
|
|
|
|
aResult = Geom2dAdaptor::MakeCurve(aProj);
|
|
if (aResult.IsNull())
|
|
{
|
|
return aResult;
|
|
}
|
|
|
|
if (aResult->IsKind(STANDARD_TYPE(Geom2d_TrimmedCurve)))
|
|
{
|
|
occ::handle<Geom2d_TrimmedCurve> aTC = occ::down_cast<Geom2d_TrimmedCurve>(aResult);
|
|
aResult = aTC->BasisCurve();
|
|
}
|
|
|
|
return aResult;
|
|
}
|
|
|
|
return aResult;
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
occ::handle<Geom2d_Curve> ShapeConstruct_ProjectCurveOnSurface::getLine(
|
|
const ArrayOfPnt& thePoints,
|
|
const ArrayOfReal& theParams,
|
|
ArrayOfPnt2d& thePoints2d,
|
|
const double theTol,
|
|
bool& theIsRecompute,
|
|
bool& theIsFromCache) const
|
|
{
|
|
const int aNb = thePoints.Length();
|
|
gp_Pnt aP[4];
|
|
aP[0] = thePoints(1);
|
|
aP[1] = thePoints(2);
|
|
aP[2] = thePoints(aNb - 1);
|
|
aP[3] = thePoints(aNb);
|
|
|
|
gp_Pnt2d aP2d[4];
|
|
int i = 0;
|
|
|
|
double aTol2 = theTol * theTol;
|
|
double aTolWorking = theTol;
|
|
const bool isPeriodicU = mySurf->Surface()->IsUPeriodic();
|
|
const bool isPeriodicV = mySurf->Surface()->IsVPeriodic();
|
|
|
|
// Protection against bad tolerance shapes
|
|
if (aTol2 > 1.0)
|
|
{
|
|
aTolWorking = Precision::Confusion();
|
|
aTol2 = aTolWorking * aTolWorking;
|
|
}
|
|
if (aTol2 < Precision::SquareConfusion())
|
|
{
|
|
aTol2 = Precision::SquareConfusion();
|
|
}
|
|
const double anOldTol2 = aTol2;
|
|
|
|
int aSavedPointNum = -1;
|
|
gp_Pnt2d aSavedPoint;
|
|
|
|
// Create projector with B-spline surface pole cache optimization
|
|
SurfaceProjectorWithCache aProjector(mySurf);
|
|
|
|
// Helper lambda to project a point with cache lookup
|
|
auto projectPoint = [&](const gp_Pnt& thePoint, gp_Pnt2d& theResult, const int theIndex) -> void {
|
|
// Try existing endpoint cache first
|
|
const int aNbCache = myCache.Length();
|
|
for (int j = 0; j < aNbCache; ++j)
|
|
{
|
|
const CachePoint& aCachePnt = myCache(j);
|
|
if (aCachePnt.first.SquareDistance(thePoint) < aTol2)
|
|
{
|
|
theResult =
|
|
aProjector.NextValueOfUV(aCachePnt.second, thePoint, aTolWorking, aTol2, aTolWorking);
|
|
aSavedPointNum = theIndex;
|
|
aSavedPoint = aCachePnt.second;
|
|
if (theIndex == 0)
|
|
{
|
|
theIsFromCache = true;
|
|
}
|
|
return;
|
|
}
|
|
}
|
|
|
|
// Fall back to full projection (with B-spline pole cache optimization inside)
|
|
theResult = aProjector.ValueOfUV(thePoint, aTolWorking, aTol2);
|
|
};
|
|
|
|
// Project first and last points
|
|
for (; i < 4; i += 3)
|
|
{
|
|
projectPoint(aP[i], aP2d[i], i);
|
|
|
|
const double aDist = aProjector.Gap();
|
|
const double aCurDist = aDist * aDist;
|
|
if (aTol2 < aCurDist)
|
|
{
|
|
aTol2 = aCurDist;
|
|
}
|
|
}
|
|
|
|
if (isPeriodicU || isPeriodicV)
|
|
{
|
|
// Compute second and last but one c2d points
|
|
for (i = 1; i < 3; i++)
|
|
{
|
|
projectPoint(aP[i], aP2d[i], i);
|
|
|
|
const double aDist = aProjector.Gap();
|
|
const double aCurDist = aDist * aDist;
|
|
if (aTol2 < aCurDist)
|
|
{
|
|
aTol2 = aCurDist;
|
|
}
|
|
}
|
|
|
|
if (isPeriodicU)
|
|
{
|
|
theIsRecompute = fixPeriodicityTroubles(&aP2d[0],
|
|
1,
|
|
mySurf->Surface()->UPeriod(),
|
|
aSavedPointNum,
|
|
aSavedPoint.X());
|
|
}
|
|
|
|
if (isPeriodicV)
|
|
{
|
|
theIsRecompute = fixPeriodicityTroubles(&aP2d[0],
|
|
2,
|
|
mySurf->Surface()->VPeriod(),
|
|
aSavedPointNum,
|
|
aSavedPoint.Y());
|
|
}
|
|
}
|
|
|
|
if (theIsRecompute && mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_SphericalSurface)))
|
|
{
|
|
return nullptr;
|
|
}
|
|
|
|
thePoints2d.SetValue(1, aP2d[0]);
|
|
thePoints2d.SetValue(aNb, aP2d[3]);
|
|
|
|
// Restore old tolerance
|
|
aTol2 = anOldTol2;
|
|
|
|
const double dPar = theParams(aNb) - theParams(1);
|
|
if (std::abs(dPar) < Precision::PConfusion())
|
|
{
|
|
return nullptr;
|
|
}
|
|
|
|
const gp_Vec2d aVec0(aP2d[0], aP2d[3]);
|
|
const gp_Vec2d aVec = aVec0 / dPar;
|
|
occ::handle<Geom_Surface> aSurf = mySurf->Surface();
|
|
bool isNormalCheck = aSurf->IsCNu(1) && aSurf->IsCNv(1);
|
|
|
|
if (isNormalCheck)
|
|
{
|
|
for (i = 1; i <= aNb; i++)
|
|
{
|
|
const gp_XY aCurPoint = aP2d[0].XY() + aVec.XY() * (theParams(i) - theParams(1));
|
|
gp_Pnt aCurP;
|
|
gp_Vec aNormalVec, aDu, aDv;
|
|
aSurf->D1(aCurPoint.X(), aCurPoint.Y(), aCurP, aDu, aDv);
|
|
aNormalVec = aDu ^ aDv;
|
|
if (aNormalVec.SquareMagnitude() < Precision::SquareConfusion())
|
|
{
|
|
isNormalCheck = false;
|
|
break;
|
|
}
|
|
const gp_Lin aNormalLine(aCurP, gp_Dir(aNormalVec));
|
|
const double aDist = aNormalLine.Distance(thePoints(i));
|
|
if (aDist > aTolWorking)
|
|
{
|
|
return nullptr;
|
|
}
|
|
}
|
|
}
|
|
|
|
if (!isNormalCheck)
|
|
{
|
|
const double aFirstPointDist =
|
|
mySurf->Surface()->Value(aP2d[0].X(), aP2d[0].Y()).SquareDistance(thePoints(1));
|
|
aTol2 = std::max(aTol2, aTol2 * 2 * aFirstPointDist);
|
|
for (i = 2; i < aNb; i++)
|
|
{
|
|
const gp_XY aCurPoint = aP2d[0].XY() + aVec.XY() * (theParams(i) - theParams(1));
|
|
gp_Pnt aCurP;
|
|
aSurf->D0(aCurPoint.X(), aCurPoint.Y(), aCurP);
|
|
const double aDist1 = aCurP.SquareDistance(thePoints(i));
|
|
|
|
if (std::abs(aFirstPointDist - aDist1) > aTol2)
|
|
{
|
|
return nullptr;
|
|
}
|
|
}
|
|
}
|
|
|
|
// Check if pcurve can be represented by Geom2d_Line
|
|
const double aLLength = aVec0.Magnitude();
|
|
if (std::abs(aLLength - dPar) <= Precision::PConfusion())
|
|
{
|
|
const gp_XY aDirL = aVec0.XY() / aLLength;
|
|
const gp_Pnt2d aPL(aP2d[0].XY() - theParams(1) * aDirL);
|
|
return new Geom2d_Line(aPL, gp_Dir2d(aDirL));
|
|
}
|
|
|
|
// Create straight bspline
|
|
NCollection_Array1<gp_Pnt2d> aPoles(1, 2);
|
|
aPoles(1) = aP2d[0];
|
|
aPoles(2) = aP2d[3];
|
|
|
|
NCollection_Array1<double> aKnots(1, 2);
|
|
aKnots(1) = theParams(1);
|
|
aKnots(2) = theParams(theParams.Length());
|
|
|
|
NCollection_Array1<int> aMults(1, 2);
|
|
aMults(1) = 2;
|
|
aMults(2) = 2;
|
|
|
|
return new Geom2d_BSplineCurve(aPoles, aKnots, aMults, 1);
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
bool ShapeConstruct_ProjectCurveOnSurface::approxPCurve(const int theNbPnt,
|
|
const occ::handle<Geom_Curve>& theC3D,
|
|
const double theTolFirst,
|
|
const double theTolLast,
|
|
ArrayOfPnt& thePoints,
|
|
ArrayOfReal& theParams,
|
|
ArrayOfPnt2d& thePoints2d,
|
|
occ::handle<Geom2d_Curve>& theC2D)
|
|
{
|
|
// For performance, first try to handle typical case when pcurve is straight
|
|
bool isRecompute = false;
|
|
bool isFromCacheLine = false;
|
|
|
|
theC2D = getLine(thePoints, theParams, thePoints2d, myPreci, isRecompute, isFromCacheLine);
|
|
if (!theC2D.IsNull())
|
|
{
|
|
// Fill cache
|
|
bool aChangeCycle = false;
|
|
if (!myCache.IsEmpty()
|
|
&& myCache(0).first.Distance(thePoints(1)) > myCache(0).first.Distance(thePoints(theNbPnt))
|
|
&& myCache(0).first.Distance(thePoints(theNbPnt)) < Precision::Confusion())
|
|
{
|
|
aChangeCycle = true;
|
|
}
|
|
|
|
myCache.Resize(0, 1, false);
|
|
if (aChangeCycle)
|
|
{
|
|
myCache(0) = CachePoint(thePoints(1), thePoints2d(1));
|
|
myCache(1) = CachePoint(thePoints(theNbPnt), thePoints2d(theNbPnt));
|
|
}
|
|
else
|
|
{
|
|
myCache(0) = CachePoint(thePoints(theNbPnt), thePoints2d(theNbPnt));
|
|
myCache(1) = CachePoint(thePoints(1), thePoints2d(1));
|
|
}
|
|
return true;
|
|
}
|
|
|
|
bool isDone = true;
|
|
|
|
// Test if the curve 3d is a boundary of the surface
|
|
bool isoParam, isoPar2d3d, isoTypeU, p1OnIso, p2OnIso, isoClosed;
|
|
gp_Pnt2d valueP1, valueP2;
|
|
occ::handle<Geom_Curve> cIso;
|
|
double t1, t2;
|
|
|
|
occ::handle<Standard_Type> sType = mySurf->Surface()->DynamicType();
|
|
bool isAnalytic = true;
|
|
if (sType == STANDARD_TYPE(Geom_BezierSurface) || sType == STANDARD_TYPE(Geom_BSplineSurface))
|
|
{
|
|
isAnalytic = false;
|
|
}
|
|
|
|
double uf, ul, vf, vl;
|
|
mySurf->Surface()->Bounds(uf, ul, vf, vl);
|
|
isoClosed = false;
|
|
|
|
ArrayOfReal pout(1, theNbPnt);
|
|
for (int i = 1; i <= theNbPnt; i++)
|
|
{
|
|
pout.SetValue(i, 0.0);
|
|
}
|
|
|
|
isoParam = isAnIsoparametric(theNbPnt,
|
|
thePoints,
|
|
theParams,
|
|
isoTypeU,
|
|
p1OnIso,
|
|
valueP1,
|
|
p2OnIso,
|
|
valueP2,
|
|
isoPar2d3d,
|
|
cIso,
|
|
t1,
|
|
t2,
|
|
pout);
|
|
|
|
// Projection of the points on surfaces
|
|
gp_Pnt p3d;
|
|
gp_Pnt2d p2d;
|
|
double isoValue = 0., isoPar1 = 0., isoPar2 = 0., tPar = 0., tdeb, tfin;
|
|
double Cf, Cl, parf, parl;
|
|
|
|
if (isoParam)
|
|
{
|
|
if (isoTypeU)
|
|
{
|
|
isoValue = valueP1.X();
|
|
isoPar1 = valueP1.Y();
|
|
isoPar2 = valueP2.Y();
|
|
isoClosed = mySurf->IsVClosed(myPreci);
|
|
parf = vf;
|
|
parl = vl;
|
|
}
|
|
else
|
|
{
|
|
isoValue = valueP1.Y();
|
|
isoPar1 = valueP1.X();
|
|
isoPar2 = valueP2.X();
|
|
isoClosed = mySurf->IsUClosed(myPreci);
|
|
parf = uf;
|
|
parl = ul;
|
|
}
|
|
if (!isoPar2d3d && !isAnalytic)
|
|
{
|
|
Cf = cIso->FirstParameter();
|
|
Cl = cIso->LastParameter();
|
|
if (Precision::IsInfinite(Cf))
|
|
{
|
|
Cf = -1000;
|
|
}
|
|
if (Precision::IsInfinite(Cl))
|
|
{
|
|
Cl = +1000;
|
|
}
|
|
|
|
tdeb = pout(2);
|
|
|
|
if (isoClosed && (isoPar1 == parf || isoPar1 == parl))
|
|
{
|
|
if (std::abs(tdeb - parf) < std::abs(tdeb - parl))
|
|
{
|
|
isoPar1 = parf;
|
|
}
|
|
else
|
|
{
|
|
isoPar1 = parl;
|
|
}
|
|
if (isoTypeU)
|
|
{
|
|
valueP1.SetY(isoPar1);
|
|
}
|
|
else
|
|
{
|
|
valueP1.SetX(isoPar1);
|
|
}
|
|
}
|
|
if (isoClosed && (isoPar2 == parf || isoPar2 == parl))
|
|
{
|
|
tfin = pout(theNbPnt - 1);
|
|
if (std::abs(tfin - parf) < std::abs(tfin - parl))
|
|
{
|
|
isoPar2 = parf;
|
|
}
|
|
else
|
|
{
|
|
isoPar2 = parl;
|
|
}
|
|
if (isoTypeU)
|
|
{
|
|
valueP2.SetY(isoPar2);
|
|
}
|
|
else
|
|
{
|
|
valueP2.SetX(isoPar2);
|
|
}
|
|
}
|
|
|
|
if (!isoClosed)
|
|
{
|
|
if ((std::abs(tdeb - isoPar1) > std::abs(tdeb - isoPar2))
|
|
&& (std::abs(pout(theNbPnt - 1) - isoPar2) > std::abs(pout(theNbPnt - 1) - isoPar1)))
|
|
{
|
|
gp_Pnt2d valueTmp = valueP1;
|
|
valueP1 = valueP2;
|
|
valueP2 = valueTmp;
|
|
if (isoTypeU)
|
|
{
|
|
isoValue = valueP1.X();
|
|
isoPar1 = valueP1.Y();
|
|
isoPar2 = valueP2.Y();
|
|
}
|
|
else
|
|
{
|
|
isoValue = valueP1.Y();
|
|
isoPar1 = valueP1.X();
|
|
isoPar2 = valueP2.X();
|
|
}
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
const double Up = ul - uf;
|
|
const double Vp = vl - vf;
|
|
double gap = myPreci;
|
|
bool aChangeCycle = false;
|
|
|
|
bool isFromCache = false;
|
|
gp_Pnt2d aSavedPoint;
|
|
|
|
if (!myCache.IsEmpty()
|
|
&& myCache(0).first.Distance(thePoints(1)) > myCache(0).first.Distance(thePoints(theNbPnt)))
|
|
{
|
|
if (myCache(0).first.Distance(thePoints(theNbPnt)) < Precision::Confusion())
|
|
{
|
|
aChangeCycle = true;
|
|
}
|
|
}
|
|
|
|
bool needResolveUJump = false;
|
|
bool needResolveVJump = false;
|
|
gp_Pnt prevP3d;
|
|
gp_Pnt2d prevP2d;
|
|
|
|
for (int ii = 1; ii <= theNbPnt; ii++)
|
|
{
|
|
const int aPntIndex = aChangeCycle ? (theNbPnt - ii + 1) : ii;
|
|
p3d = thePoints(aPntIndex);
|
|
|
|
if (isoParam)
|
|
{
|
|
if (isoPar2d3d)
|
|
{
|
|
if (isoPar2 > isoPar1)
|
|
{
|
|
tPar = theParams(aPntIndex);
|
|
}
|
|
else
|
|
{
|
|
tPar = t1 + t2 - theParams(aPntIndex);
|
|
}
|
|
}
|
|
else if (!isAnalytic)
|
|
{
|
|
if (aPntIndex == 1)
|
|
{
|
|
tPar = isoPar1;
|
|
}
|
|
else if (aPntIndex == theNbPnt)
|
|
{
|
|
tPar = isoPar2;
|
|
}
|
|
else
|
|
{
|
|
tPar = pout(aPntIndex);
|
|
}
|
|
}
|
|
|
|
if (!isoPar2d3d && isAnalytic)
|
|
{
|
|
if (aPntIndex == 1)
|
|
{
|
|
p2d = valueP1;
|
|
}
|
|
else if (aPntIndex == theNbPnt)
|
|
{
|
|
p2d = valueP2;
|
|
}
|
|
else
|
|
{
|
|
p2d = mySurf->NextValueOfUV(p2d, p3d, myPreci, Precision::Confusion() + 1000 * gap);
|
|
gap = mySurf->Gap();
|
|
}
|
|
}
|
|
else
|
|
{
|
|
if (isoTypeU)
|
|
{
|
|
p2d.SetX(isoValue);
|
|
p2d.SetY(tPar);
|
|
}
|
|
else
|
|
{
|
|
p2d.SetX(tPar);
|
|
p2d.SetY(isoValue);
|
|
}
|
|
}
|
|
}
|
|
else
|
|
{
|
|
if ((aPntIndex == 1) && p1OnIso)
|
|
{
|
|
p2d = valueP1;
|
|
}
|
|
else if ((aPntIndex == theNbPnt) && p2OnIso)
|
|
{
|
|
p2d = valueP2;
|
|
}
|
|
else
|
|
{
|
|
if (aPntIndex == 1 || aPntIndex == theNbPnt)
|
|
{
|
|
if (!isRecompute)
|
|
{
|
|
p2d = thePoints2d(aPntIndex);
|
|
gap = mySurf->Gap();
|
|
if (aPntIndex == 1)
|
|
{
|
|
isFromCache = isFromCacheLine;
|
|
aSavedPoint = p2d;
|
|
}
|
|
continue;
|
|
}
|
|
else
|
|
{
|
|
int j = 0;
|
|
const int aNbCache = myCache.Length();
|
|
for (; j < aNbCache; ++j)
|
|
{
|
|
const CachePoint& aCachePnt = myCache(j);
|
|
if (aCachePnt.first.SquareDistance(p3d) < myPreci * myPreci)
|
|
{
|
|
p2d = mySurf->NextValueOfUV(aCachePnt.second,
|
|
p3d,
|
|
myPreci,
|
|
Precision::Confusion() + gap);
|
|
if (aPntIndex == 1)
|
|
{
|
|
isFromCache = true;
|
|
aSavedPoint = aCachePnt.second;
|
|
}
|
|
break;
|
|
}
|
|
}
|
|
if (j >= aNbCache)
|
|
{
|
|
p2d = mySurf->ValueOfUV(p3d, myPreci);
|
|
}
|
|
}
|
|
}
|
|
else
|
|
{
|
|
p2d = mySurf->NextValueOfUV(p2d, p3d, myPreci, Precision::Confusion() + 1000 * gap);
|
|
}
|
|
gap = mySurf->Gap();
|
|
}
|
|
}
|
|
thePoints2d.SetValue(aPntIndex, p2d);
|
|
|
|
if (theNbPnt > 23 && ii > 2 && ii < theNbPnt)
|
|
{
|
|
if (fabs(p2d.X() - prevP2d.X()) > 0.95 * Up && prevP3d.Distance(p3d) < myPreci
|
|
&& !mySurf->IsUClosed(myPreci) && mySurf->NbSingularities(myPreci) > 0
|
|
&& mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_BSplineSurface)))
|
|
{
|
|
needResolveUJump = true;
|
|
}
|
|
if (fabs(p2d.Y() - prevP2d.Y()) > 0.95 * Vp && prevP3d.Distance(p3d) < myPreci
|
|
&& !mySurf->IsVClosed(myPreci) && mySurf->NbSingularities(myPreci) > 0
|
|
&& mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_BSplineSurface)))
|
|
{
|
|
needResolveVJump = true;
|
|
}
|
|
}
|
|
prevP3d = p3d;
|
|
prevP2d = p2d;
|
|
if (ii > 1)
|
|
{
|
|
if (aChangeCycle)
|
|
{
|
|
p2d.SetXY(2. * p2d.XY() - thePoints2d(aPntIndex + 1).XY());
|
|
}
|
|
else
|
|
{
|
|
p2d.SetXY(2. * p2d.XY() - thePoints2d(aPntIndex - 1).XY());
|
|
}
|
|
}
|
|
}
|
|
|
|
if (!isoPar2d3d)
|
|
{
|
|
projectDegeneratedPoints(mySurf, theNbPnt, thePoints, thePoints2d, myPreci, true);
|
|
projectDegeneratedPoints(mySurf, theNbPnt, thePoints, thePoints2d, myPreci, false);
|
|
}
|
|
|
|
// Check extremities for singularities
|
|
const gp_Pnt aPointFirst = thePoints.First();
|
|
const gp_Pnt aPointLast = thePoints.Last();
|
|
const double aTolFirst = (theTolFirst < 0) ? Precision::Confusion() : theTolFirst;
|
|
const double aTolLast = (theTolLast < 0) ? Precision::Confusion() : theTolLast;
|
|
|
|
for (int i = 1;; i++)
|
|
{
|
|
double aPreci, aFirstPar, aLastPar;
|
|
gp_Pnt aP3d;
|
|
gp_Pnt2d aFirstP2d, aLastP2d;
|
|
bool IsUiso;
|
|
if (!mySurf->Singularity(i, aPreci, aP3d, aFirstP2d, aLastP2d, aFirstPar, aLastPar, IsUiso))
|
|
{
|
|
break;
|
|
}
|
|
if (aPreci <= Precision::Confusion() && aPointFirst.Distance(aP3d) <= aTolFirst)
|
|
{
|
|
correctExtremity(theC3D, theParams, thePoints2d, true, aFirstP2d, IsUiso);
|
|
}
|
|
if (aPreci <= Precision::Confusion() && aPointLast.Distance(aP3d) <= aTolLast)
|
|
{
|
|
correctExtremity(theC3D, theParams, thePoints2d, false, aFirstP2d, IsUiso);
|
|
}
|
|
}
|
|
|
|
// Handle U-closed surfaces
|
|
double dist2d;
|
|
const double TolOnUPeriod = Precision::Confusion() * Up;
|
|
const double TolOnVPeriod = Precision::Confusion() * Vp;
|
|
|
|
if (mySurf->IsUClosed(myPreci) || needResolveUJump)
|
|
{
|
|
double prevX, firstX = thePoints2d(1).X();
|
|
if (!isFromCache)
|
|
{
|
|
while (firstX < uf)
|
|
{
|
|
firstX += Up;
|
|
thePoints2d.ChangeValue(1).SetX(firstX);
|
|
}
|
|
while (firstX > ul)
|
|
{
|
|
firstX -= Up;
|
|
thePoints2d.ChangeValue(1).SetX(firstX);
|
|
}
|
|
}
|
|
|
|
if (mySurf->Surface()->IsUPeriodic() && isFromCache)
|
|
{
|
|
double aMinParam = uf, aMaxParam = ul;
|
|
while (aMinParam > aSavedPoint.X())
|
|
{
|
|
aMinParam -= Up;
|
|
aMaxParam -= Up;
|
|
}
|
|
while (aMaxParam < aSavedPoint.X())
|
|
{
|
|
aMinParam += Up;
|
|
aMaxParam += Up;
|
|
}
|
|
double aShift = ShapeAnalysis::AdjustToPeriod(firstX, aMinParam, aMaxParam);
|
|
firstX += aShift;
|
|
thePoints2d.ChangeValue(1).SetX(firstX);
|
|
}
|
|
prevX = firstX;
|
|
|
|
double minX = firstX, maxX = firstX;
|
|
bool ToAdjust = false;
|
|
|
|
for (int aPntIter = 2; aPntIter <= thePoints2d.Length(); ++aPntIter)
|
|
{
|
|
double CurX = thePoints2d(aPntIter).X();
|
|
dist2d = std::abs(CurX - prevX);
|
|
if (dist2d > (Up / 2))
|
|
{
|
|
insertAdditionalPointOrAdjust(ToAdjust,
|
|
1,
|
|
Up,
|
|
TolOnUPeriod,
|
|
CurX,
|
|
prevX,
|
|
theC3D,
|
|
aPntIter,
|
|
thePoints,
|
|
theParams,
|
|
thePoints2d);
|
|
}
|
|
prevX = CurX;
|
|
if (minX > CurX)
|
|
{
|
|
minX = CurX;
|
|
}
|
|
else if (maxX < CurX)
|
|
{
|
|
maxX = CurX;
|
|
}
|
|
}
|
|
|
|
if (!isFromCache)
|
|
{
|
|
double midX = 0.5 * (minX + maxX);
|
|
double shiftX = 0.;
|
|
if (midX > ul)
|
|
{
|
|
shiftX = -Up;
|
|
}
|
|
else if (midX < uf)
|
|
{
|
|
shiftX = Up;
|
|
}
|
|
if (shiftX != 0.)
|
|
{
|
|
for (int aPntIter = 1; aPntIter <= thePoints2d.Length(); ++aPntIter)
|
|
{
|
|
thePoints2d.ChangeValue(aPntIter).SetX(thePoints2d(aPntIter).X() + shiftX);
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
// Handle V-closed surfaces
|
|
if (mySurf->IsVClosed(myPreci) || needResolveVJump
|
|
|| mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_SphericalSurface)))
|
|
{
|
|
double prevY, firstY = thePoints2d(1).Y();
|
|
if (!isFromCache)
|
|
{
|
|
while (firstY < vf)
|
|
{
|
|
firstY += Vp;
|
|
thePoints2d.ChangeValue(1).SetY(firstY);
|
|
}
|
|
while (firstY > vl)
|
|
{
|
|
firstY -= Vp;
|
|
thePoints2d.ChangeValue(1).SetY(firstY);
|
|
}
|
|
}
|
|
|
|
if (mySurf->Surface()->IsVPeriodic() && isFromCache)
|
|
{
|
|
double aMinParam = vf, aMaxParam = vl;
|
|
while (aMinParam > aSavedPoint.Y())
|
|
{
|
|
aMinParam -= Vp;
|
|
aMaxParam -= Vp;
|
|
}
|
|
while (aMaxParam < aSavedPoint.Y())
|
|
{
|
|
aMinParam += Vp;
|
|
aMaxParam += Vp;
|
|
}
|
|
double aShift = ShapeAnalysis::AdjustToPeriod(firstY, aMinParam, aMaxParam);
|
|
firstY += aShift;
|
|
thePoints2d.ChangeValue(1).SetY(firstY);
|
|
}
|
|
prevY = firstY;
|
|
|
|
double minY = firstY, maxY = firstY;
|
|
bool ToAdjust = false;
|
|
|
|
for (int aPntIter = 2; aPntIter <= thePoints2d.Length(); ++aPntIter)
|
|
{
|
|
double CurY = thePoints2d(aPntIter).Y();
|
|
dist2d = std::abs(CurY - prevY);
|
|
if (dist2d > (Vp / 2))
|
|
{
|
|
insertAdditionalPointOrAdjust(ToAdjust,
|
|
2,
|
|
Vp,
|
|
TolOnVPeriod,
|
|
CurY,
|
|
prevY,
|
|
theC3D,
|
|
aPntIter,
|
|
thePoints,
|
|
theParams,
|
|
thePoints2d);
|
|
}
|
|
prevY = CurY;
|
|
if (minY > CurY)
|
|
{
|
|
minY = CurY;
|
|
}
|
|
else if (maxY < CurY)
|
|
{
|
|
maxY = CurY;
|
|
}
|
|
}
|
|
|
|
if (!isFromCache)
|
|
{
|
|
double midY = 0.5 * (minY + maxY);
|
|
double shiftY = 0.;
|
|
if (midY > vl)
|
|
{
|
|
shiftY = -Vp;
|
|
}
|
|
else if (midY < vf)
|
|
{
|
|
shiftY = Vp;
|
|
}
|
|
if (shiftY != 0.)
|
|
{
|
|
for (int aPntIter = 1; aPntIter <= thePoints2d.Length(); ++aPntIter)
|
|
{
|
|
thePoints2d.ChangeValue(aPntIter).SetY(thePoints2d(aPntIter).Y() + shiftY);
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
// Handle V-closed seam adjustment
|
|
if (mySurf->IsVClosed(myPreci) || mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_SphericalSurface)))
|
|
{
|
|
for (int aPntIter = 2; aPntIter <= thePoints2d.Length(); ++aPntIter)
|
|
{
|
|
dist2d = std::abs(thePoints2d(aPntIter).Y() - thePoints2d(aPntIter - 1).Y());
|
|
if (dist2d > (Vp / 2))
|
|
{
|
|
bool prevOnFirst = false;
|
|
bool prevOnLast = false;
|
|
bool currOnFirst = false;
|
|
bool currOnLast = false;
|
|
|
|
double distPrevVF = std::abs(thePoints2d(aPntIter - 1).Y() - vf);
|
|
double distPrevVL = std::abs(thePoints2d(aPntIter - 1).Y() - vl);
|
|
double distCurrVF = std::abs(thePoints2d(aPntIter).Y() - vf);
|
|
double distCurrVL = std::abs(thePoints2d(aPntIter).Y() - vl);
|
|
|
|
double theMin = distPrevVF;
|
|
prevOnFirst = true;
|
|
if (distPrevVL < theMin)
|
|
{
|
|
theMin = distPrevVL;
|
|
prevOnFirst = false;
|
|
prevOnLast = true;
|
|
}
|
|
if (distCurrVF < theMin)
|
|
{
|
|
theMin = distCurrVF;
|
|
prevOnFirst = false;
|
|
prevOnLast = false;
|
|
currOnFirst = true;
|
|
}
|
|
if (distCurrVL < theMin)
|
|
{
|
|
prevOnFirst = false;
|
|
prevOnLast = false;
|
|
currOnFirst = false;
|
|
currOnLast = true;
|
|
}
|
|
|
|
if (prevOnFirst)
|
|
{
|
|
thePoints2d.ChangeValue(aPntIter - 1).SetY(vf);
|
|
}
|
|
else if (prevOnLast)
|
|
{
|
|
thePoints2d.ChangeValue(aPntIter - 1).SetY(vl);
|
|
}
|
|
else if (currOnFirst)
|
|
{
|
|
thePoints2d.ChangeValue(aPntIter).SetY(vf);
|
|
}
|
|
else if (currOnLast)
|
|
{
|
|
thePoints2d.ChangeValue(aPntIter).SetY(vl);
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
// Handle AdjustOverDegen
|
|
if (myAdjustOverDegen != -1)
|
|
{
|
|
if (mySurf->IsUClosed(myPreci))
|
|
{
|
|
mySurf->IsDegenerated(gp_Pnt(0, 0, 0), myPreci);
|
|
if (mySurf->NbSingularities(myPreci) > 0)
|
|
{
|
|
double PrevX = 0.;
|
|
int OnBound = 0, PrevOnBound = 0;
|
|
int ind;
|
|
bool start = true;
|
|
for (ind = 1; ind <= thePoints2d.Length(); ind++)
|
|
{
|
|
double CurX = thePoints2d(ind).X();
|
|
if (mySurf->IsDegenerated(thePoints(ind), Precision::Confusion()))
|
|
{
|
|
continue;
|
|
}
|
|
OnBound =
|
|
(std::abs(std::abs(CurX - 0.5 * (ul + uf)) - Up / 2) <= Precision::PConfusion());
|
|
if (!start && std::abs(std::abs(CurX - PrevX) - Up / 2) <= 0.01 * Up)
|
|
{
|
|
break;
|
|
}
|
|
start = false;
|
|
PrevX = CurX;
|
|
PrevOnBound = OnBound;
|
|
}
|
|
if (ind <= thePoints2d.Length())
|
|
{
|
|
PrevX = (myAdjustOverDegen ? uf : ul);
|
|
double dU = Up / 2 + Precision::PConfusion();
|
|
if (PrevOnBound)
|
|
{
|
|
thePoints2d.ChangeValue(ind - 1).SetX(PrevX);
|
|
for (int j = ind - 2; j > 0; j--)
|
|
{
|
|
double CurX = thePoints2d(j).X();
|
|
while (CurX < PrevX - dU)
|
|
{
|
|
CurX += Up;
|
|
thePoints2d.ChangeValue(j).SetX(CurX);
|
|
}
|
|
while (CurX > PrevX + dU)
|
|
{
|
|
CurX -= Up;
|
|
thePoints2d.ChangeValue(j).SetX(CurX);
|
|
}
|
|
}
|
|
}
|
|
else if (OnBound)
|
|
{
|
|
thePoints2d.ChangeValue(ind).SetX(PrevX);
|
|
for (int j = ind + 1; j <= thePoints2d.Length(); j++)
|
|
{
|
|
double CurX = thePoints2d(j).X();
|
|
while (CurX < PrevX - dU)
|
|
{
|
|
CurX += Up;
|
|
thePoints2d.ChangeValue(j).SetX(CurX);
|
|
}
|
|
while (CurX > PrevX + dU)
|
|
{
|
|
CurX -= Up;
|
|
thePoints2d.ChangeValue(j).SetX(CurX);
|
|
}
|
|
}
|
|
}
|
|
myStatus |= ShapeExtend::EncodeStatus(ShapeExtend_DONE4);
|
|
}
|
|
}
|
|
}
|
|
else if (mySurf->IsVClosed(myPreci))
|
|
{
|
|
mySurf->IsDegenerated(gp_Pnt(0, 0, 0), myPreci);
|
|
if (mySurf->NbSingularities(myPreci) > 0)
|
|
{
|
|
double PrevY = 0.;
|
|
int OnBound = 0, PrevOnBound = 0;
|
|
int ind;
|
|
bool start = true;
|
|
for (ind = 1; ind <= thePoints2d.Length(); ind++)
|
|
{
|
|
double CurY = thePoints2d(ind).Y();
|
|
if (mySurf->IsDegenerated(thePoints(ind), Precision::Confusion()))
|
|
{
|
|
continue;
|
|
}
|
|
OnBound =
|
|
(std::abs(std::abs(CurY - 0.5 * (vl + vf)) - Vp / 2) <= Precision::PConfusion());
|
|
if (!start && std::abs(std::abs(CurY - PrevY) - Vp / 2) <= 0.01 * Vp)
|
|
{
|
|
break;
|
|
}
|
|
start = false;
|
|
PrevY = CurY;
|
|
PrevOnBound = OnBound;
|
|
}
|
|
if (ind <= thePoints2d.Length())
|
|
{
|
|
PrevY = (myAdjustOverDegen ? vf : vl);
|
|
double dV = Vp / 2 + Precision::PConfusion();
|
|
if (PrevOnBound)
|
|
{
|
|
thePoints2d.ChangeValue(ind - 1).SetY(PrevY);
|
|
for (int j = ind - 2; j > 0; j--)
|
|
{
|
|
double CurY = thePoints2d(j).Y();
|
|
while (CurY < PrevY - dV)
|
|
{
|
|
CurY += Vp;
|
|
thePoints2d.ChangeValue(j).SetY(CurY);
|
|
}
|
|
while (CurY > PrevY + dV)
|
|
{
|
|
CurY -= Vp;
|
|
thePoints2d.ChangeValue(j).SetY(CurY);
|
|
}
|
|
}
|
|
}
|
|
else if (OnBound)
|
|
{
|
|
thePoints2d.ChangeValue(ind).SetY(PrevY);
|
|
for (int j = ind + 1; j <= thePoints2d.Length(); j++)
|
|
{
|
|
double CurY = thePoints2d(j).Y();
|
|
while (CurY < PrevY - dV)
|
|
{
|
|
CurY += Vp;
|
|
thePoints2d.ChangeValue(j).SetY(CurY);
|
|
}
|
|
while (CurY > PrevY + dV)
|
|
{
|
|
CurY -= Vp;
|
|
thePoints2d.ChangeValue(j).SetY(CurY);
|
|
}
|
|
}
|
|
}
|
|
myStatus |= ShapeExtend::EncodeStatus(ShapeExtend_DONE4);
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
// Fill cache
|
|
myCache.Resize(0, 1, false);
|
|
if (aChangeCycle)
|
|
{
|
|
myCache(0) = CachePoint(thePoints(1), thePoints2d(1));
|
|
myCache(1) = CachePoint(thePoints.Last(), thePoints2d.Last());
|
|
}
|
|
else
|
|
{
|
|
myCache(0) = CachePoint(thePoints.Last(), thePoints2d.Last());
|
|
myCache(1) = CachePoint(thePoints(1), thePoints2d(1));
|
|
}
|
|
|
|
return isDone;
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
void ShapeConstruct_ProjectCurveOnSurface::correctExtremity(const occ::handle<Geom_Curve>& theC3D,
|
|
const ArrayOfReal& theParams,
|
|
ArrayOfPnt2d& thePoints2d,
|
|
const bool theIsFirstPoint,
|
|
const gp_Pnt2d& thePointOnIsoLine,
|
|
const bool theIsUIso)
|
|
{
|
|
const int NbPnt = thePoints2d.Length();
|
|
const int IndCoord = (theIsUIso) ? 2 : 1;
|
|
const double SingularityCoord = thePointOnIsoLine.Coord(3 - IndCoord);
|
|
gp_Pnt2d EndPoint = (theIsFirstPoint) ? thePoints2d(1) : thePoints2d(NbPnt);
|
|
const double FinishCoord = EndPoint.Coord(3 - IndCoord);
|
|
|
|
gp_Dir2d aDir = (theIsUIso) ? gp::DY2d() : gp::DX2d();
|
|
gp_Lin2d anIsoLine(EndPoint, aDir);
|
|
IntRes2d_Domain Dom1, Dom2;
|
|
|
|
const bool IsPeriodic =
|
|
(theIsUIso) ? mySurf->Surface()->IsVPeriodic() : mySurf->Surface()->IsUPeriodic();
|
|
|
|
gp_Pnt2d FirstPointOfLine, SecondPointOfLine;
|
|
double FinishParam, FirstParam, SecondParam;
|
|
|
|
if (theIsFirstPoint)
|
|
{
|
|
FirstPointOfLine = thePoints2d(3);
|
|
SecondPointOfLine = thePoints2d(2);
|
|
FinishParam = theParams(1);
|
|
FirstParam = theParams(3);
|
|
SecondParam = theParams(2);
|
|
}
|
|
else
|
|
{
|
|
FirstPointOfLine = thePoints2d(NbPnt - 2);
|
|
SecondPointOfLine = thePoints2d(NbPnt - 1);
|
|
FinishParam = theParams(NbPnt);
|
|
FirstParam = theParams(NbPnt - 2);
|
|
SecondParam = theParams(NbPnt - 1);
|
|
}
|
|
|
|
if (SingularityCoord > FinishCoord && SecondPointOfLine.Coord(3 - IndCoord) > FinishCoord)
|
|
{
|
|
return;
|
|
}
|
|
if (SingularityCoord < FinishCoord && SecondPointOfLine.Coord(3 - IndCoord) < FinishCoord)
|
|
{
|
|
return;
|
|
}
|
|
|
|
{
|
|
const double aPrevDist =
|
|
std::abs(SecondPointOfLine.Coord(IndCoord) - FirstPointOfLine.Coord(IndCoord));
|
|
const double aCurDist = std::abs(EndPoint.Coord(IndCoord) - SecondPointOfLine.Coord(IndCoord));
|
|
if (aCurDist <= 2 * aPrevDist)
|
|
{
|
|
return;
|
|
}
|
|
}
|
|
|
|
gp_Pnt2d FinishPoint = (theIsUIso) ? gp_Pnt2d(FinishCoord, SecondPointOfLine.Y())
|
|
: gp_Pnt2d(SecondPointOfLine.X(), FinishCoord);
|
|
|
|
for (;;)
|
|
{
|
|
if (std::abs(SecondPointOfLine.Coord(3 - IndCoord) - FinishCoord)
|
|
<= 2 * Precision::PConfusion())
|
|
{
|
|
break;
|
|
}
|
|
|
|
gp_Vec2d aVec(FirstPointOfLine, SecondPointOfLine);
|
|
const double aSqMagnitude = aVec.SquareMagnitude();
|
|
if (aSqMagnitude <= 1.e-32)
|
|
{
|
|
break;
|
|
}
|
|
aDir.SetCoord(aVec.X(), aVec.Y());
|
|
|
|
gp_Lin2d aLine(FirstPointOfLine, aDir);
|
|
IntCurve_IntConicConic Intersector(anIsoLine, Dom1, aLine, Dom2, 1.e-10, 1.e-10);
|
|
if (Intersector.IsDone() && !Intersector.IsEmpty())
|
|
{
|
|
IntRes2d_IntersectionPoint IntPoint = Intersector.Point(1);
|
|
FinishPoint = IntPoint.Value();
|
|
}
|
|
else
|
|
{
|
|
FinishPoint = (theIsUIso) ? gp_Pnt2d(FinishCoord, SecondPointOfLine.Y())
|
|
: gp_Pnt2d(SecondPointOfLine.X(), FinishCoord);
|
|
}
|
|
|
|
gp_Pnt2d PrevPoint = FirstPointOfLine;
|
|
FirstPointOfLine = SecondPointOfLine;
|
|
FirstParam = SecondParam;
|
|
SecondParam = (FirstParam + FinishParam) / 2;
|
|
if (std::abs(SecondParam - FirstParam) <= 2 * Precision::PConfusion())
|
|
{
|
|
break;
|
|
}
|
|
gp_Pnt aP3d;
|
|
theC3D->D0(SecondParam, aP3d);
|
|
SecondPointOfLine =
|
|
mySurf->NextValueOfUV(FirstPointOfLine, aP3d, myPreci, Precision::Confusion());
|
|
if (IsPeriodic)
|
|
{
|
|
adjustSecondToFirstPoint(FirstPointOfLine, SecondPointOfLine, mySurf->Surface());
|
|
}
|
|
|
|
const double aPrevDist = std::abs(FirstPointOfLine.Coord(IndCoord) - PrevPoint.Coord(IndCoord));
|
|
const double aCurDist =
|
|
std::abs(SecondPointOfLine.Coord(IndCoord) - FirstPointOfLine.Coord(IndCoord));
|
|
if (aCurDist > 2 * aPrevDist)
|
|
{
|
|
break;
|
|
}
|
|
}
|
|
|
|
if (theIsFirstPoint)
|
|
{
|
|
thePoints2d.ChangeValue(1) = FinishPoint;
|
|
}
|
|
else
|
|
{
|
|
thePoints2d.ChangeValue(NbPnt) = FinishPoint;
|
|
}
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
void ShapeConstruct_ProjectCurveOnSurface::insertAdditionalPointOrAdjust(
|
|
bool& theToAdjust,
|
|
const int theIndCoord,
|
|
const double thePeriod,
|
|
const double theTolOnPeriod,
|
|
double& theCurCoord,
|
|
const double thePrevCoord,
|
|
const occ::handle<Geom_Curve>& theC3D,
|
|
int& theIndex,
|
|
ArrayOfPnt& thePoints,
|
|
ArrayOfReal& theParams,
|
|
ArrayOfPnt2d& thePoints2d)
|
|
{
|
|
const double CorrectedCurCoord =
|
|
ElCLib::InPeriod(theCurCoord, thePrevCoord - thePeriod / 2, thePrevCoord + thePeriod / 2);
|
|
|
|
if (!theToAdjust)
|
|
{
|
|
const double CurPar = theParams(theIndex);
|
|
const double PrevPar = theParams(theIndex - 1);
|
|
double MidPar = (PrevPar + CurPar) / 2;
|
|
gp_Pnt MidP3d;
|
|
theC3D->D0(MidPar, MidP3d);
|
|
gp_Pnt2d MidP2d = mySurf->ValueOfUV(MidP3d, myPreci);
|
|
double MidCoord = MidP2d.Coord(theIndCoord);
|
|
MidCoord =
|
|
ElCLib::InPeriod(MidCoord, thePrevCoord - thePeriod / 2, thePrevCoord + thePeriod / 2);
|
|
double FirstCoord = thePrevCoord, LastCoord = CorrectedCurCoord;
|
|
if (LastCoord < FirstCoord)
|
|
{
|
|
double tmp = FirstCoord;
|
|
FirstCoord = LastCoord;
|
|
LastCoord = tmp;
|
|
}
|
|
if (LastCoord - FirstCoord <= theTolOnPeriod)
|
|
{
|
|
theToAdjust = true;
|
|
}
|
|
else if (FirstCoord <= MidCoord && MidCoord <= LastCoord)
|
|
{
|
|
theToAdjust = true;
|
|
}
|
|
else
|
|
{
|
|
bool Success = true;
|
|
double FirstT = PrevPar;
|
|
double LastT = CurPar;
|
|
MidCoord = MidP2d.Coord(theIndCoord);
|
|
while (std::abs(MidCoord - thePrevCoord) >= thePeriod / 2 - theTolOnPeriod
|
|
|| std::abs(theCurCoord - MidCoord) >= thePeriod / 2 - theTolOnPeriod)
|
|
{
|
|
if (MidPar - FirstT <= Precision::PConfusion() || LastT - MidPar <= Precision::PConfusion())
|
|
{
|
|
Success = false;
|
|
break;
|
|
}
|
|
if (std::abs(MidCoord - thePrevCoord) >= thePeriod / 2 - theTolOnPeriod)
|
|
{
|
|
LastT = (FirstT + LastT) / 2;
|
|
}
|
|
else
|
|
{
|
|
FirstT = (FirstT + LastT) / 2;
|
|
}
|
|
MidPar = (FirstT + LastT) / 2;
|
|
theC3D->D0(MidPar, MidP3d);
|
|
MidP2d = mySurf->ValueOfUV(MidP3d, myPreci);
|
|
MidCoord = MidP2d.Coord(theIndCoord);
|
|
}
|
|
if (Success)
|
|
{
|
|
// Insert additional point - need to resize arrays
|
|
const int aNewLength = thePoints.Length() + 1;
|
|
ArrayOfPnt aNewPoints(1, aNewLength);
|
|
ArrayOfReal aNewParams(1, aNewLength);
|
|
ArrayOfPnt2d aNewPoints2d(1, aNewLength);
|
|
|
|
for (int i = 1; i < theIndex; i++)
|
|
{
|
|
aNewPoints.SetValue(i, thePoints(i));
|
|
aNewParams.SetValue(i, theParams(i));
|
|
aNewPoints2d.SetValue(i, thePoints2d(i));
|
|
}
|
|
aNewPoints.SetValue(theIndex, MidP3d);
|
|
aNewParams.SetValue(theIndex, MidPar);
|
|
aNewPoints2d.SetValue(theIndex, MidP2d);
|
|
for (int i = theIndex; i <= thePoints.Length(); i++)
|
|
{
|
|
aNewPoints.SetValue(i + 1, thePoints(i));
|
|
aNewParams.SetValue(i + 1, theParams(i));
|
|
aNewPoints2d.SetValue(i + 1, thePoints2d(i));
|
|
}
|
|
|
|
thePoints = std::move(aNewPoints);
|
|
theParams = std::move(aNewParams);
|
|
thePoints2d = std::move(aNewPoints2d);
|
|
theIndex++;
|
|
}
|
|
else
|
|
{
|
|
theToAdjust = true;
|
|
}
|
|
}
|
|
}
|
|
if (theToAdjust)
|
|
{
|
|
theCurCoord = CorrectedCurCoord;
|
|
thePoints2d.ChangeValue(theIndex).SetCoord(theIndCoord, theCurCoord);
|
|
}
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
occ::handle<Geom2d_Curve> ShapeConstruct_ProjectCurveOnSurface::interpolatePCurve(
|
|
const int theNbPnt,
|
|
const ArrayOfPnt2d& thePoints2d,
|
|
const ArrayOfReal& theParams) const
|
|
{
|
|
occ::handle<Geom2d_Curve> aC2D;
|
|
const double theTolerance2d = myPreci / (100 * theNbPnt);
|
|
|
|
try
|
|
{
|
|
OCC_CATCH_SIGNALS
|
|
|
|
// Convert to HArrays for Geom2dAPI_Interpolate
|
|
occ::handle<NCollection_HArray1<gp_Pnt2d>> aPnts2d =
|
|
new NCollection_HArray1<gp_Pnt2d>(1, theNbPnt);
|
|
occ::handle<NCollection_HArray1<double>> aParams = new NCollection_HArray1<double>(1, theNbPnt);
|
|
|
|
for (int i = 1; i <= theNbPnt; i++)
|
|
{
|
|
aPnts2d->SetValue(i, thePoints2d(i));
|
|
aParams->SetValue(i, theParams(i));
|
|
}
|
|
|
|
// Check for coincident points
|
|
ArrayOfPnt2d aTmpPnts2d = thePoints2d;
|
|
ArrayOfReal aTmpParams = theParams;
|
|
double aPreci = theTolerance2d;
|
|
checkPoints2d(aTmpPnts2d, aTmpParams, aPreci);
|
|
|
|
if (aTmpPnts2d.Length() < 2)
|
|
{
|
|
return aC2D;
|
|
}
|
|
|
|
// Rebuild HArrays if points were removed
|
|
if (aTmpPnts2d.Length() != theNbPnt)
|
|
{
|
|
aPnts2d = new NCollection_HArray1<gp_Pnt2d>(1, aTmpPnts2d.Length());
|
|
aParams = new NCollection_HArray1<double>(1, aTmpParams.Length());
|
|
for (int i = 1; i <= aTmpPnts2d.Length(); i++)
|
|
{
|
|
aPnts2d->SetValue(i, aTmpPnts2d(i));
|
|
aParams->SetValue(i, aTmpParams(i));
|
|
}
|
|
}
|
|
|
|
Geom2dAPI_Interpolate myInterPol2d(aPnts2d, aParams, false, aPreci);
|
|
myInterPol2d.Perform();
|
|
if (myInterPol2d.IsDone())
|
|
{
|
|
aC2D = myInterPol2d.Curve();
|
|
}
|
|
}
|
|
catch (Standard_Failure const& anException)
|
|
{
|
|
#ifdef OCCT_DEBUG
|
|
std::cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::interpolatePCurve(): Exception: ";
|
|
anException.Print(std::cout);
|
|
std::cout << std::endl;
|
|
#endif
|
|
(void)anException;
|
|
aC2D.Nullify();
|
|
}
|
|
return aC2D;
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
occ::handle<Geom2d_Curve> ShapeConstruct_ProjectCurveOnSurface::approximatePCurve(
|
|
const ArrayOfPnt2d& thePoints2d,
|
|
const ArrayOfReal& theParams) const
|
|
{
|
|
const double theTolerance2d = myPreci;
|
|
occ::handle<Geom2d_Curve> aC2D;
|
|
|
|
try
|
|
{
|
|
OCC_CATCH_SIGNALS
|
|
|
|
ArrayOfPnt2d aTmpPnts2d = thePoints2d;
|
|
ArrayOfReal aTmpParams = theParams;
|
|
double aPreci = theTolerance2d;
|
|
checkPoints2d(aTmpPnts2d, aTmpParams, aPreci);
|
|
|
|
const int numberPnt = aTmpPnts2d.Length();
|
|
if (numberPnt < 2)
|
|
{
|
|
return aC2D;
|
|
}
|
|
|
|
NCollection_Array1<gp_Pnt> points3d(1, numberPnt);
|
|
NCollection_Array1<double> params(1, numberPnt);
|
|
|
|
for (int i = 1; i <= numberPnt; i++)
|
|
{
|
|
const gp_Pnt2d& pnt2d = aTmpPnts2d(i);
|
|
points3d(i).SetCoord(pnt2d.X(), pnt2d.Y(), 0);
|
|
params(i) = aTmpParams(i);
|
|
}
|
|
|
|
GeomAPI_PointsToBSpline appr(points3d, params, 1, 10, GeomAbs_C1, theTolerance2d);
|
|
const occ::handle<Geom_BSplineCurve>& crv3d = appr.Curve();
|
|
|
|
const int NbPoles = crv3d->NbPoles();
|
|
const NCollection_Array1<gp_Pnt>& poles3d = crv3d->Poles();
|
|
NCollection_Array1<gp_Pnt2d> poles2d(1, NbPoles);
|
|
|
|
for (int i = 1; i <= NbPoles; i++)
|
|
{
|
|
poles2d(i).SetCoord(poles3d(i).X(), poles3d(i).Y());
|
|
}
|
|
|
|
const NCollection_Array1<double>& weights = crv3d->WeightsArray();
|
|
const NCollection_Array1<double>& knots = crv3d->Knots();
|
|
const NCollection_Array1<int>& multiplicities = crv3d->Multiplicities();
|
|
|
|
aC2D = new Geom2d_BSplineCurve(poles2d,
|
|
weights,
|
|
knots,
|
|
multiplicities,
|
|
crv3d->Degree(),
|
|
crv3d->IsPeriodic());
|
|
return aC2D;
|
|
}
|
|
catch (Standard_Failure const& anException)
|
|
{
|
|
#ifdef OCCT_DEBUG
|
|
std::cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::approximatePCurve(): Exception: ";
|
|
anException.Print(std::cout);
|
|
std::cout << std::endl;
|
|
#endif
|
|
(void)anException;
|
|
aC2D.Nullify();
|
|
}
|
|
return aC2D;
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
void ShapeConstruct_ProjectCurveOnSurface::checkPoints(ArrayOfPnt& thePoints,
|
|
ArrayOfReal& theParams,
|
|
double& thePreci) const
|
|
{
|
|
const int firstElem = thePoints.Lower();
|
|
const int lastElem = thePoints.Upper();
|
|
|
|
int nbPntDropped = 0;
|
|
int lastValid = firstElem;
|
|
|
|
NCollection_Array1<int> tmpParam(firstElem, lastElem);
|
|
for (int i = firstElem; i <= lastElem; i++)
|
|
{
|
|
tmpParam.SetValue(i, 1);
|
|
}
|
|
|
|
double DistMin2 = RealLast();
|
|
gp_Pnt Prev = thePoints.Value(lastValid);
|
|
|
|
for (int i = firstElem + 1; i <= lastElem; i++)
|
|
{
|
|
const gp_Pnt& Curr = thePoints.Value(i);
|
|
const double CurDist2 = Prev.SquareDistance(Curr);
|
|
if (CurDist2 < gp::Resolution())
|
|
{
|
|
nbPntDropped++;
|
|
if (i == lastElem)
|
|
{
|
|
tmpParam.SetValue(lastValid, 0);
|
|
}
|
|
else
|
|
{
|
|
tmpParam.SetValue(i, 0);
|
|
}
|
|
}
|
|
else
|
|
{
|
|
if (CurDist2 < DistMin2)
|
|
{
|
|
DistMin2 = CurDist2;
|
|
}
|
|
lastValid = i;
|
|
Prev = Curr;
|
|
}
|
|
}
|
|
|
|
if (DistMin2 < RealLast())
|
|
{
|
|
thePreci = 0.9 * std::sqrt(DistMin2);
|
|
}
|
|
|
|
if (nbPntDropped == 0)
|
|
{
|
|
return;
|
|
}
|
|
|
|
const int newLast = lastElem - nbPntDropped;
|
|
if ((newLast - firstElem + 1) < 2)
|
|
{
|
|
return;
|
|
}
|
|
|
|
ArrayOfPnt newPnts(firstElem, newLast);
|
|
ArrayOfReal newParams(firstElem, newLast);
|
|
int newCurr = firstElem;
|
|
|
|
for (int i = firstElem; i <= lastElem; i++)
|
|
{
|
|
if (tmpParam.Value(i) == 1)
|
|
{
|
|
newPnts.SetValue(newCurr, thePoints.Value(i));
|
|
newParams.SetValue(newCurr, theParams.Value(i));
|
|
newCurr++;
|
|
}
|
|
}
|
|
|
|
thePoints = newPnts;
|
|
theParams = newParams;
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
void ShapeConstruct_ProjectCurveOnSurface::checkPoints2d(ArrayOfPnt2d& thePoints2d,
|
|
ArrayOfReal& theParams,
|
|
double& thePreci) const
|
|
{
|
|
const int firstElem = thePoints2d.Lower();
|
|
const int lastElem = thePoints2d.Upper();
|
|
|
|
int nbPntDropped = 0;
|
|
int lastValid = firstElem;
|
|
|
|
NCollection_Array1<int> tmpParam(firstElem, lastElem);
|
|
for (int i = firstElem; i <= lastElem; i++)
|
|
{
|
|
tmpParam.SetValue(i, 1);
|
|
}
|
|
|
|
double DistMin2 = RealLast();
|
|
gp_Pnt2d Prev = thePoints2d.Value(lastValid);
|
|
|
|
for (int i = firstElem + 1; i <= lastElem; i++)
|
|
{
|
|
const gp_Pnt2d& Curr = thePoints2d.Value(i);
|
|
const double CurDist2 = Prev.SquareDistance(Curr);
|
|
if (CurDist2 < gp::Resolution())
|
|
{
|
|
nbPntDropped++;
|
|
if (i == lastElem)
|
|
{
|
|
tmpParam.SetValue(lastValid, 0);
|
|
}
|
|
else
|
|
{
|
|
tmpParam.SetValue(i, 0);
|
|
}
|
|
}
|
|
else
|
|
{
|
|
if (CurDist2 < DistMin2)
|
|
{
|
|
DistMin2 = CurDist2;
|
|
}
|
|
lastValid = i;
|
|
Prev = Curr;
|
|
}
|
|
}
|
|
|
|
if (DistMin2 < RealLast())
|
|
{
|
|
thePreci = 0.9 * std::sqrt(DistMin2);
|
|
}
|
|
|
|
if (nbPntDropped == 0)
|
|
{
|
|
return;
|
|
}
|
|
|
|
int newLast = lastElem - nbPntDropped;
|
|
if ((newLast - firstElem + 1) < 2)
|
|
{
|
|
// Create minimal length pcurve
|
|
tmpParam.SetValue(firstElem, 1);
|
|
tmpParam.SetValue(lastElem, 1);
|
|
gp_XY lastPnt = thePoints2d.Value(lastElem).XY();
|
|
lastPnt.Add(gp_XY(thePreci, thePreci));
|
|
thePoints2d.ChangeValue(lastElem).SetXY(lastPnt);
|
|
newLast = firstElem + 1;
|
|
}
|
|
|
|
ArrayOfPnt2d newPnts(firstElem, newLast);
|
|
ArrayOfReal newParams(firstElem, newLast);
|
|
int newCurr = firstElem;
|
|
|
|
for (int i = firstElem; i <= lastElem; i++)
|
|
{
|
|
if (tmpParam.Value(i) == 1)
|
|
{
|
|
newPnts.SetValue(newCurr, thePoints2d.Value(i));
|
|
newParams.SetValue(newCurr, theParams.Value(i));
|
|
newCurr++;
|
|
}
|
|
}
|
|
|
|
thePoints2d = std::move(newPnts);
|
|
theParams = std::move(newParams);
|
|
}
|
|
|
|
//=================================================================================================
|
|
|
|
bool ShapeConstruct_ProjectCurveOnSurface::isAnIsoparametric(const int theNbPnt,
|
|
const ArrayOfPnt& thePoints,
|
|
const ArrayOfReal& theParams,
|
|
bool& theIsTypeU,
|
|
bool& theP1OnIso,
|
|
gp_Pnt2d& theValueP1,
|
|
bool& theP2OnIso,
|
|
gp_Pnt2d& theValueP2,
|
|
bool& theIsoPar2d3d,
|
|
occ::handle<Geom_Curve>& theCIso,
|
|
double& theT1,
|
|
double& theT2,
|
|
ArrayOfReal& theParamsOut) const
|
|
{
|
|
try
|
|
{
|
|
OCC_CATCH_SIGNALS
|
|
|
|
constexpr double prec = Precision::Confusion();
|
|
|
|
bool isoParam = false;
|
|
theIsoPar2d3d = false;
|
|
|
|
double U1, U2, V1, V2;
|
|
mySurf->Bounds(U1, U2, V1, V2);
|
|
|
|
if (mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface)))
|
|
{
|
|
occ::handle<Geom_RectangularTrimmedSurface> sTrim =
|
|
occ::down_cast<Geom_RectangularTrimmedSurface>(mySurf->Surface());
|
|
sTrim->Bounds(U1, U2, V1, V2);
|
|
}
|
|
|
|
gp_Pnt pt;
|
|
int mpt[2];
|
|
mpt[0] = mpt[1] = 0;
|
|
double t, tpar[2] = {0.0, 0.0}, isoValue = 0.;
|
|
double mindist2;
|
|
double mind2[2];
|
|
mindist2 = mind2[0] = mind2[1] = 4 * prec * prec;
|
|
|
|
theP1OnIso = false;
|
|
theP2OnIso = false;
|
|
const Bnd_Box* aBox = nullptr;
|
|
|
|
for (int j = 1; j <= 4; j++)
|
|
{
|
|
double isoVal = 0.;
|
|
bool isoU = false;
|
|
occ::handle<Geom_Curve> cI;
|
|
double tt1, tt2;
|
|
|
|
if (j == 1)
|
|
{
|
|
if (Precision::IsInfinite(U1))
|
|
{
|
|
continue;
|
|
}
|
|
cI = mySurf->UIso(U1);
|
|
isoU = true;
|
|
isoVal = U1;
|
|
aBox = &mySurf->GetBoxUF();
|
|
}
|
|
else if (j == 2)
|
|
{
|
|
if (Precision::IsInfinite(U2))
|
|
{
|
|
continue;
|
|
}
|
|
cI = mySurf->UIso(U2);
|
|
isoU = true;
|
|
isoVal = U2;
|
|
aBox = &mySurf->GetBoxUL();
|
|
}
|
|
else if (j == 3)
|
|
{
|
|
if (Precision::IsInfinite(V1))
|
|
{
|
|
continue;
|
|
}
|
|
cI = mySurf->VIso(V1);
|
|
isoU = false;
|
|
isoVal = V1;
|
|
aBox = &mySurf->GetBoxVF();
|
|
}
|
|
else if (j == 4)
|
|
{
|
|
if (Precision::IsInfinite(V2))
|
|
{
|
|
continue;
|
|
}
|
|
cI = mySurf->VIso(V2);
|
|
isoU = false;
|
|
isoVal = V2;
|
|
aBox = &mySurf->GetBoxVL();
|
|
}
|
|
if (cI.IsNull())
|
|
{
|
|
continue;
|
|
}
|
|
|
|
if (isoU)
|
|
{
|
|
tt1 = V1;
|
|
tt2 = V2;
|
|
}
|
|
else
|
|
{
|
|
tt1 = U1;
|
|
tt2 = U2;
|
|
}
|
|
|
|
gp_Pnt ext1, ext2;
|
|
cI->D0(tt1, ext1);
|
|
cI->D0(tt2, ext2);
|
|
|
|
gp_Pnt extmi;
|
|
cI->D0((tt1 + tt2) / 2, extmi);
|
|
if (ext1.IsEqual(ext2, prec) && ext1.IsEqual(extmi, prec))
|
|
{
|
|
continue;
|
|
}
|
|
|
|
bool PtEQext1 = false;
|
|
bool PtEQext2 = false;
|
|
|
|
double currd2[2], tp[2] = {0, 0};
|
|
int mp[2];
|
|
|
|
for (int i = 0; i < 2; i++)
|
|
{
|
|
mp[i] = 0;
|
|
const int k = (i == 0 ? 1 : theNbPnt);
|
|
|
|
currd2[i] = thePoints(k).SquareDistance(ext1);
|
|
if (currd2[i] <= prec * prec && !PtEQext1)
|
|
{
|
|
mp[i] = 1;
|
|
tp[i] = tt1;
|
|
PtEQext1 = true;
|
|
continue;
|
|
}
|
|
|
|
currd2[i] = thePoints(k).SquareDistance(ext2);
|
|
if (currd2[i] <= prec * prec && !PtEQext2)
|
|
{
|
|
mp[i] = 2;
|
|
tp[i] = tt2;
|
|
PtEQext2 = true;
|
|
continue;
|
|
}
|
|
|
|
if (mySurf->Surface()->IsKind(STANDARD_TYPE(Geom_SphericalSurface)) && !isoU)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
if (aBox->IsOut(thePoints(k)))
|
|
{
|
|
continue;
|
|
}
|
|
|
|
double Cf = cI->FirstParameter();
|
|
double Cl = cI->LastParameter();
|
|
if (Precision::IsInfinite(Cf))
|
|
{
|
|
Cf = -1000;
|
|
}
|
|
if (Precision::IsInfinite(Cl))
|
|
{
|
|
Cl = +1000;
|
|
}
|
|
|
|
ShapeAnalysis_Curve sac;
|
|
const double dist = sac.Project(cI, thePoints(k), prec, pt, t, Cf, Cl);
|
|
currd2[i] = dist * dist;
|
|
if ((dist <= prec) && (t >= Cf) && (t <= Cl))
|
|
{
|
|
mp[i] = 3;
|
|
tp[i] = t;
|
|
}
|
|
}
|
|
|
|
if (mp[0] > 0 && mp[1] > 0 && std::abs(tp[0] - tp[1]) < Precision::PConfusion())
|
|
{
|
|
continue;
|
|
}
|
|
|
|
if (mp[0] > 0 && (!theP1OnIso || currd2[0] < mind2[0]))
|
|
{
|
|
theP1OnIso = true;
|
|
mind2[0] = currd2[0];
|
|
if (isoU)
|
|
{
|
|
theValueP1.SetCoord(isoVal, tp[0]);
|
|
}
|
|
else
|
|
{
|
|
theValueP1.SetCoord(tp[0], isoVal);
|
|
}
|
|
}
|
|
|
|
if (mp[1] > 0 && (!theP2OnIso || currd2[1] < mind2[1]))
|
|
{
|
|
theP2OnIso = true;
|
|
mind2[1] = currd2[1];
|
|
if (isoU)
|
|
{
|
|
theValueP2.SetCoord(isoVal, tp[1]);
|
|
}
|
|
else
|
|
{
|
|
theValueP2.SetCoord(tp[1], isoVal);
|
|
}
|
|
}
|
|
|
|
if (mp[0] <= 0 || mp[1] <= 0)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
const double md2 = currd2[0] + currd2[1];
|
|
if (mindist2 <= md2)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
mindist2 = md2;
|
|
mpt[0] = mp[0];
|
|
mpt[1] = mp[1];
|
|
tpar[0] = tp[0];
|
|
tpar[1] = tp[1];
|
|
theIsTypeU = isoU;
|
|
isoValue = isoVal;
|
|
theCIso = cI;
|
|
theT1 = tt1;
|
|
theT2 = tt2;
|
|
}
|
|
|
|
if (mpt[0] > 0 && mpt[1] > 0)
|
|
{
|
|
theP1OnIso = theP2OnIso = true;
|
|
if (theIsTypeU)
|
|
{
|
|
theValueP1.SetCoord(isoValue, tpar[0]);
|
|
theValueP2.SetCoord(isoValue, tpar[1]);
|
|
}
|
|
else
|
|
{
|
|
theValueP1.SetCoord(tpar[0], isoValue);
|
|
theValueP2.SetCoord(tpar[1], isoValue);
|
|
}
|
|
|
|
if (mpt[0] != 3 && mpt[1] != 3)
|
|
{
|
|
theIsoPar2d3d = true;
|
|
for (int i = 2; i < theNbPnt && theIsoPar2d3d; i++)
|
|
{
|
|
if (tpar[1] > tpar[0])
|
|
{
|
|
t = theParams(i);
|
|
}
|
|
else
|
|
{
|
|
t = theT1 + theT2 - theParams(i);
|
|
}
|
|
theCIso->D0(t, pt);
|
|
if (!thePoints(i).IsEqual(pt, prec))
|
|
{
|
|
theIsoPar2d3d = false;
|
|
}
|
|
}
|
|
}
|
|
|
|
if (theIsoPar2d3d)
|
|
{
|
|
isoParam = true;
|
|
}
|
|
else
|
|
{
|
|
double prevParam = tpar[0];
|
|
double Cf, Cl;
|
|
bool isoByDistance = true;
|
|
Cf = theCIso->FirstParameter();
|
|
Cl = theCIso->LastParameter();
|
|
if (Precision::IsInfinite(Cf))
|
|
{
|
|
Cf = -1000;
|
|
}
|
|
if (Precision::IsInfinite(Cl))
|
|
{
|
|
Cl = +1000;
|
|
}
|
|
|
|
ShapeAnalysis_Curve sac;
|
|
for (int i = 2; i < theNbPnt && isoByDistance; i++)
|
|
{
|
|
const double dist =
|
|
sac.NextProject(prevParam, theCIso, thePoints(i), prec, pt, t, Cf, Cl, false);
|
|
prevParam = t;
|
|
theParamsOut(i) = t;
|
|
if ((dist > prec) || (t < Cf) || (t > Cl))
|
|
{
|
|
isoByDistance = false;
|
|
}
|
|
}
|
|
if (isoByDistance)
|
|
{
|
|
isoParam = true;
|
|
}
|
|
}
|
|
}
|
|
return isoParam;
|
|
}
|
|
catch (Standard_Failure const& anException)
|
|
{
|
|
#ifdef OCCT_DEBUG
|
|
std::cout << "Warning: ShapeConstruct_ProjectCurveOnSurface::isAnIsoparametric(): Exception: ";
|
|
anException.Print(std::cout);
|
|
std::cout << std::endl;
|
|
#endif
|
|
(void)anException;
|
|
return false;
|
|
}
|
|
}
|