From 21e336f8174e88960005ec35ff5dd154e6bce78d Mon Sep 17 00:00:00 2001 From: nbv Date: Mon, 12 Nov 2018 17:00:06 +0300 Subject: [PATCH] 0030365: Create tool to compute deviation between any 2D-curve and some its segment The tool, DRAW-command for its check and corresponding test cases have been created. See documentation about new tool in GeomLib_Tool.hxx file. --- src/GeomLib/GeomLib_Tool.cxx | 354 +++++++++++++++--- src/GeomLib/GeomLib_Tool.hxx | 103 +++-- .../GeomliteTest_API2dCommands.cxx | 99 +++++ tests/lowalgos/2ddeviation/A1 | 24 ++ tests/lowalgos/2ddeviation/A2 | 24 ++ tests/lowalgos/2ddeviation/A3 | 22 ++ tests/lowalgos/grids.list | 3 +- 7 files changed, 539 insertions(+), 90 deletions(-) create mode 100644 tests/lowalgos/2ddeviation/A1 create mode 100644 tests/lowalgos/2ddeviation/A2 create mode 100644 tests/lowalgos/2ddeviation/A3 diff --git a/src/GeomLib/GeomLib_Tool.cxx b/src/GeomLib/GeomLib_Tool.cxx index 4996027d3e..93db74fc23 100644 --- a/src/GeomLib/GeomLib_Tool.cxx +++ b/src/GeomLib/GeomLib_Tool.cxx @@ -13,66 +13,17 @@ // Alternatively, this file may be used under the terms of Open CASCADE // commercial license or contractual agreement. +#include #include -#include #include #include #include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include #include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include #include #include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include +#include +#include // The functions Parameter(s) are used to compute parameter(s) of point // on curves and surfaces. The main rule is that tested point must lied @@ -92,7 +43,7 @@ static const Standard_Real PARTOLERANCE = 1.e-9; Standard_Boolean GeomLib_Tool::Parameter(const Handle(Geom_Curve)& Curve, const gp_Pnt& Point, const Standard_Real MaxDist, - Standard_Real& U) + Standard_Real& U) { if( Curve.IsNull() ) return Standard_False; // @@ -140,8 +91,8 @@ Standard_Boolean GeomLib_Tool::Parameter(const Handle(Geom_Curve)& Curve, Standard_Boolean GeomLib_Tool::Parameters(const Handle(Geom_Surface)& Surface, const gp_Pnt& Point, const Standard_Real MaxDist, - Standard_Real& U, - Standard_Real& V) + Standard_Real& U, + Standard_Real& V) { if( Surface.IsNull() ) return Standard_False; // @@ -192,7 +143,7 @@ Standard_Boolean GeomLib_Tool::Parameters(const Handle(Geom_Surface)& Surface, Standard_Boolean GeomLib_Tool::Parameter(const Handle(Geom2d_Curve)& Curve, const gp_Pnt2d& Point, const Standard_Real MaxDist, - Standard_Real& U) + Standard_Real& U) { if( Curve.IsNull() ) return Standard_False; // @@ -225,3 +176,294 @@ Standard_Boolean GeomLib_Tool::Parameter(const Handle(Geom2d_Curve)& Curve, return Standard_True; } + +//! Target function to compute deviation of the source +//! 2D-curve. It is one-variate function. Its parameter +//! is a parameter on the curve. Deviation is a maximal +//! distance between any point in the curve and the given line. +class FuncSolveDeviation : public math_MultipleVarFunction +{ +public: + + //! Constructor. Initializes the curve and the line + //! going through two given points. + FuncSolveDeviation(const Geom2dAdaptor_Curve& theCurve, + const gp_XY& thePf, + const gp_XY& thePl) : myCurve(theCurve), + myPRef(thePf) + { + myDirRef = thePl - thePf; + mySqMod = myDirRef.SquareModulus(); + myIsValid = (mySqMod > Precision::SquarePConfusion()); + } + + //! Compute additional parameters depending on the argument + //! of *this + void UpdateFields(const Standard_Real theParam) + { + myCurve.D0(theParam, myPointOnCurve); + const gp_XY aVt = myPointOnCurve.XY() - myPRef; + myVecCurvLine = aVt.Dot(myDirRef)*myDirRef / mySqMod - aVt; + } + + //! Returns value of *this (square deviation) and its 1st and 2nd derivative. + void ValueAndDerives(const Standard_Real theParam, Standard_Real& theVal, + Standard_Real& theD1, Standard_Real& theD2) + { + gp_Vec2d aD1; + gp_Vec2d aD2; + myCurve.D2(theParam, myPointOnCurve, aD1, aD2); + + const gp_XY aVt = myPointOnCurve.XY() - myPRef; + theVal = aVt.Crossed(myDirRef); + theD1 = aD1.Crossed(myDirRef); + theD2 = 2.0*(theD1*theD1 + theVal*aD2.Crossed(myDirRef)); + theD1 *= 2.0*theVal; + theVal *= theVal / mySqMod; + } + + //! Returns TRUE if the function has been initializes correctly. + inline Standard_Boolean IsValid() const + { + return myIsValid; + } + + //! Returns number of variables + virtual Standard_Integer NbVariables() const Standard_OVERRIDE + { + return 1; + } + + //! Returns last computed Point in the given curve. + //! Its value will be recomputed after calling UpdateFields(...) method, + //! which sets this point correspond to the input parameter. + inline const gp_Pnt2d& PointOnCurve() const + { + return myPointOnCurve; + } + + //! Returns last computed vector directed from some point on the curve + //! to the given line. This vector is correspond to the found deviation. + //! Its value will be recomputed after calling UpdateFields(...) method, + //! which set this vector correspond to the input parameter. + inline const gp_Vec2d& VecCurveLine() const + { + return myVecCurvLine; + } + + //! Returns the given line + inline void GetLine(gp_Lin2d* const theLine) const + { + if (theLine == nullptr) + return; + + theLine->SetDirection(myDirRef); + theLine->SetLocation(myPRef); + } + + //! Returns value of *this (square deviation) + virtual Standard_Boolean Value(const math_Vector& thePrm, + Standard_Real& theVal) Standard_OVERRIDE + { + Standard_Real aD1; + Standard_Real aD2; + ValueAndDerives(thePrm.Value(thePrm.Lower()), theVal, aD1, aD2); + theVal = -theVal; + return Standard_True; + } + + //! Always returns 0. It is used for compatibility with the parent class. + virtual Standard_Integer GetStateNumber() Standard_OVERRIDE + { + return 0; + } + +private: + + //! The curve + Geom2dAdaptor_Curve myCurve; + + //! Square modulus of myDirRef (it is constant) + Standard_Real mySqMod; + + //! TRUE if *this is initialized correctly + Standard_Boolean myIsValid; + + //! Sets the given line + gp_XY myPRef, myDirRef; + + //! Last computed point in the curve + gp_Pnt2d myPointOnCurve; + + //! Always directed from myPointOnCurve to the line + gp_Vec2d myVecCurvLine; +}; + +//======================================================================= +//function : ComputeDevation +//purpose : Computes parameter on curve (*thePrmOnCurve) where maximal deviation +// (maximal value of correspond function FuncSolveDeviation) is obtained. +// Returns the (positive) value of deflection. Returns negative value if +// the deflection cannot be computed. +// The returned parameter (in case of successful) will always be in +// the range [theFPar, theLPar]. +// +// ALGORITHM! +// The point is looked for where 1st derivative of the function +// FuncSolveDeviation is equal to 0. It is made by iterative formula: +// +// U(n+1)=U(n) - D1/D2, +// +// where D1 and D2 are 1st and 2nd derivative of the function, computed in +// the point U(n). U(0) = theStartParameter. +// +// ATTENTION!!! +// 1. The computed value is heavily dependent on the theStartParameter. +// E.g. the method can reach a minimal value of the function instead of its +// maximal value. Or it can reach some local (not global) maximum. +// 2. Recommend value of theStartParameter can be found with +// 2nd method GeomLib_Tool::ComputeDevation (see below). +//======================================================================= +Standard_Real GeomLib_Tool::ComputeDevation(const Geom2dAdaptor_Curve& theCurve, + const Standard_Real theFPar, + const Standard_Real theLPar, + const Standard_Real theStartParameter, + const Standard_Integer theNbIters, + Standard_Real* const thePrmOnCurve, + gp_Pnt2d* const thePtOnCurve, + gp_Vec2d* const theVecCurvLine, + gp_Lin2d* const theLine) +{ + // Computed maximal deflection + Standard_Real aSqDefl = -1.0; + + if ((theStartParameter < theFPar) || (theStartParameter > theLPar)) + { + return aSqDefl; + } + + const gp_Pnt2d aPf(theCurve.Value(theFPar)); + const gp_Pnt2d aPl(theCurve.Value(theLPar)); + + FuncSolveDeviation aFunc(theCurve, aPf.XY(), aPl.XY()); + + if (!aFunc.IsValid()) + { + return aSqDefl; + } + + aFunc.GetLine(theLine); + + const Standard_Real aTolDefl = Precision::PConfusion(); + + Standard_Real aD1 = 0.0; + Standard_Real aD2 = 0.0; + Standard_Real aU0 = theStartParameter; + Standard_Real aUmax = theStartParameter; + aFunc.ValueAndDerives(aU0, aSqDefl, aD1, aD2); + for (Standard_Integer anItr = 1; anItr <= theNbIters; anItr++) + { + if (Abs(aD2) < Precision::PConfusion()) + { + break; + } + + const Standard_Real aDelta = aD1 / aD2; + const Standard_Real aU1 = aU0 - aDelta; + + if ((aU1 < theFPar) || (aU1 > theLPar)) + break; + + Standard_Real aSqD = aSqDefl; + aFunc.ValueAndDerives(aU1, aSqD, aD1, aD2); + if (aSqD > aSqDefl) + { + aUmax = aU1; + const Standard_Real aDD = aSqDefl > 0.0 ? + Abs(Sqrt(aSqD) - Sqrt(aSqDefl)) : Sqrt(aSqD); + aSqDefl = aSqD; + if (aDD < aTolDefl) + break; + } + + if (Abs(aU0 - aU1) < Precision::PConfusion()) + break; + + aU0 = aU1; + } + + if (aSqDefl < 0.0) + return aSqDefl; + + if (thePrmOnCurve) + *thePrmOnCurve = aUmax; + + if ((thePtOnCurve != nullptr) || (theVecCurvLine != nullptr)) + { + aFunc.UpdateFields(aUmax); + + if (thePtOnCurve != nullptr) + thePtOnCurve->SetXY(aFunc.PointOnCurve().XY()); + + if (theVecCurvLine != nullptr) + theVecCurvLine->SetXY(aFunc.VecCurveLine().XY()); + } + + return Sqrt(aSqDefl); +} + +//======================================================================= +//function : ComputeDevation +//purpose : Computes parameter on curve (*thePrmOnCurve) where maximal deviation +// (maximal value of correspond function FuncSolveDeviation) is obtained. +// Returns the (positive) value of deflection. Returns negative value if +// the deflection cannot be computed. +// The returned parameter (in case of successful) will always be in +// the range [theFPar, theLPar]. +// +// math_PSO Algorithm is used. +// +// ATTENTION!!! +// The main idea of this method is to obtain rough value of returned +// parameter (fast but not precisely). The found value can be precised +// by 2nd method GeomLib_Tool::ComputeDevation(...) (see above). +//======================================================================= +Standard_Real GeomLib_Tool::ComputeDevation(const Geom2dAdaptor_Curve& theCurve, + const Standard_Real theFPar, + const Standard_Real theLPar, + const Standard_Integer theNbSubIntervals, + const Standard_Integer theNbIters, + Standard_Real* const thePrmOnCurve) +{ + // Computed maximal deflection + Standard_Real aSqDefl = -1.0; + + const gp_Pnt2d aPf(theCurve.Value(theFPar)); + const gp_Pnt2d aPl(theCurve.Value(theLPar)); + + FuncSolveDeviation aFunc(theCurve, aPf.XY(), aPl.XY()); + + if (!aFunc.IsValid()) + { + return aSqDefl; + } + + const math_Vector aFPar(1, 1, theFPar); + const math_Vector aLPar(1, 1, theLPar); + const math_Vector aStep(1, 1, (theLPar - theFPar) / (10.0*theNbSubIntervals)); + math_Vector anOutputPnt(1, 1, theFPar); + math_PSO aMPSO(&aFunc, aFPar, aLPar, aStep, theNbSubIntervals, theNbIters); + + aSqDefl = RealLast(); + aMPSO.Perform(aStep, aSqDefl, anOutputPnt, theNbIters); + + if (aSqDefl == RealLast()) + { + return -1.0; + } + + if (thePrmOnCurve) + *thePrmOnCurve = anOutputPnt(1); + + return Sqrt(Abs(aSqDefl)); +} \ No newline at end of file diff --git a/src/GeomLib/GeomLib_Tool.hxx b/src/GeomLib/GeomLib_Tool.hxx index 5746a6d155..fbbd55e143 100644 --- a/src/GeomLib/GeomLib_Tool.hxx +++ b/src/GeomLib/GeomLib_Tool.hxx @@ -22,63 +22,100 @@ #include #include + class Geom_Curve; -class gp_Pnt; class Geom_Surface; class Geom2d_Curve; +class Geom2dAdaptor_Curve; +class gp_Lin2d; +class gp_Pnt; class gp_Pnt2d; - +class gp_Vec2d; //! Provides various methods with Geom2d and Geom curves and surfaces. //! The methods of this class compute the parameter(s) of a given point on a -//! curve or a surface. To get the valid result the point must be located rather close +//! curve or a surface. To get the valid result the point must be located rather close //! to the curve (surface) or at least to allow getting unambiguous result //! (do not put point at center of circle...), -//! but choice of "trust" distance between curve/surface and point is -//! responcibility of user (parameter MaxDist). +//! but choice of "trust" distance between curve/surface and point is +//! responsibility of user (parameter MaxDist). //! Return FALSE if the point is beyond the MaxDist //! limit or if computation fails. -class GeomLib_Tool +class GeomLib_Tool { public: DEFINE_STANDARD_ALLOC - + //! Extracts the parameter of a 3D point lying on a 3D curve + //! or at a distance less than the MaxDist value. + Standard_EXPORT static Standard_Boolean Parameter(const Handle(Geom_Curve)& Curve, const gp_Pnt& Point, const Standard_Real MaxDist, Standard_Real& U); - //! Extracts the parameter of a 3D point lying on a 3D curve - //! or at a distance less than the MaxDist value. - Standard_EXPORT static Standard_Boolean Parameter (const Handle(Geom_Curve)& Curve, const gp_Pnt& Point, const Standard_Real MaxDist, Standard_Real& U); - //! Extracts the parameter of a 3D point lying on a surface //! or at a distance less than the MaxDist value. - Standard_EXPORT static Standard_Boolean Parameters (const Handle(Geom_Surface)& Surface, const gp_Pnt& Point, const Standard_Real MaxDist, Standard_Real& U, Standard_Real& V); - + Standard_EXPORT static Standard_Boolean Parameters(const Handle(Geom_Surface)& Surface, const gp_Pnt& Point, const Standard_Real MaxDist, Standard_Real& U, Standard_Real& V); + //! Extracts the parameter of a 2D point lying on a 2D curve //! or at a distance less than the MaxDist value. - Standard_EXPORT static Standard_Boolean Parameter (const Handle(Geom2d_Curve)& Curve, const gp_Pnt2d& Point, const Standard_Real MaxDist, Standard_Real& U); - - - + Standard_EXPORT static Standard_Boolean Parameter(const Handle(Geom2d_Curve)& Curve, const gp_Pnt2d& Point, const Standard_Real MaxDist, Standard_Real& U); + + //! Computes parameter in theCurve (*thePrmOnCurve) where maximal deviation + //! between theCurve and the linear segment joining its points with + //! the parameters theFPar and theLPar is obtained. + //! Returns the (positive) value of deviation. Returns negative value if + //! the deviation cannot be computed. + //! The returned parameter (in case of successful) will always be in + //! the range [theFPar, theLPar]. + //! Iterative method is used for computation. So, theStartParameter is + //! needed to be set. Recommend value of theStartParameter can be found with + //! the overloaded method. + //! Additionally, following values can be returned (optionally): + //! - thePtOnCurve - the point on curve where maximal deviation is achieved; + //! - thePrmOnCurve - the parameter of thePtOnCurve; + //! - theVecCurvLine - the vector along which is computed (this vector is always + //! perpendicular theLine); + //! - theLine - the linear segment joining the point of theCurve having parameters + //! theFPar and theLPar. + Standard_EXPORT static + Standard_Real ComputeDevation(const Geom2dAdaptor_Curve& theCurve, + const Standard_Real theFPar, + const Standard_Real theLPar, + const Standard_Real theStartParameter, + const Standard_Integer theNbIters = 100, + Standard_Real* const thePrmOnCurve = nullptr, + gp_Pnt2d* const thePtOnCurve = nullptr, + gp_Vec2d* const theVecCurvLine = nullptr, + gp_Lin2d* const theLine = nullptr); + + //! Computes parameter in theCurve (*thePrmOnCurve) where maximal deviation + //! between theCurve and the linear segment joining its points with + //! the parameters theFPar and theLPar is obtained. + //! Returns the (positive) value of deviation. Returns negative value if + //! the deviation cannot be computed. + //! The returned parameter (in case of successful) will always be in + //! the range [theFPar, theLPar]. + //! theNbSubIntervals defines discretization of the given interval [theFPar, theLPar] + //! to provide better search condition. This value should be chosen taking into + //! account complexity of the curve in considered interval. E.g. if there are many + //! oscillations of the curve in the interval then theNbSubIntervals mus be + //! great number. However, the greater value of theNbSubIntervals the slower the + //! algorithm will compute. + //! theNbIters sets number of iterations. + //! ATTENTION!!! + //! This algorithm cannot compute deviation precisely (so, there is no point in + //! setting big value of theNbIters). But it can give some start point for + //! the overloaded method. + Standard_EXPORT static + Standard_Real ComputeDevation(const Geom2dAdaptor_Curve& theCurve, + const Standard_Real theFPar, + const Standard_Real theLPar, + const Standard_Integer theNbSubIntervals, + const Standard_Integer theNbIters = 10, + Standard_Real * const thePrmOnCurve = nullptr); protected: - - - - private: - - - - - }; - - - - - - -#endif // _GeomLib_Tool_HeaderFile +#endif // _GeomLib_Tool_HeaderFile \ No newline at end of file diff --git a/src/GeomliteTest/GeomliteTest_API2dCommands.cxx b/src/GeomliteTest/GeomliteTest_API2dCommands.cxx index 89e5f30bea..300e59abd3 100644 --- a/src/GeomliteTest/GeomliteTest_API2dCommands.cxx +++ b/src/GeomliteTest/GeomliteTest_API2dCommands.cxx @@ -30,6 +30,7 @@ #include #include #include +#include #include #include #include @@ -548,7 +549,99 @@ static Standard_Integer intconcon(Draw_Interpretor& di, Standard_Integer n, cons return 0; } +//======================================================================= +//function : deviation +//purpose : +//======================================================================= +static Standard_Integer deviation(Draw_Interpretor& theDI, Standard_Integer theNArg, const char** theArgv) +{ + if (theNArg < 3) + { + theDI << "Use: deviation2d result curve [U0]\n"; + return 1; + } + + const Handle(Geom2d_Curve) aC = DrawTrSurf::GetCurve2d(theArgv[2]); + + if (aC.IsNull()) + { + theDI << "Error: " << theArgv[2] << " is not a 2D-curve.\n"; + return 1; + } + Geom2dAdaptor_Curve anAC(aC); + + Standard_Integer aNbInterv = 2; + Standard_Real aU0 = RealLast(); + + for (Standard_Integer aCurrArg = 3; aCurrArg < theNArg; aCurrArg++) + { + if (!strcmp(theArgv[aCurrArg], "-i")) + { + aU0 = Draw::Atof(theArgv[++aCurrArg]); + } + else if (!strcmp(theArgv[aCurrArg], "-d")) + { + aNbInterv = Draw::Atoi(theArgv[++aCurrArg]); + } + else + { + theDI << "Error: Wrong option " << theArgv[aCurrArg] << "\n"; + return 1; + } + } + + const Standard_Real aU1 = anAC.FirstParameter(); + const Standard_Real aU2 = anAC.LastParameter(); + + Standard_Real aRetCurvParam = aU0; + gp_Pnt2d aPtOnCurv; + gp_Vec2d aRetVec; + gp_Lin2d aLinSegm; + + Standard_Real aDefl = RealLast(); + + if (aU0 == RealLast()) + { + aDefl = GeomLib_Tool::ComputeDevation(anAC, aU1, aU2, + aNbInterv, 10, &aU0); + + if (aDefl < 0.0) + { + theDI << "Error: Cannot compute deviation on interval.\n"; + return 0; + } + } + + aDefl = GeomLib_Tool::ComputeDevation(anAC, aU1, aU2, aU0, 100, + &aRetCurvParam, &aPtOnCurv, + &aRetVec, &aLinSegm); + + if (aDefl < 0.0) + { + theDI << "Error: Cannot compute a deviation!\n"; + return 0; + } + + char aBuff[1000]; + theDI << "Computed value is: " << aDefl << "\n"; + + Sprintf(aBuff, "%s_pnt", theArgv[1]); + DrawTrSurf::Set(aBuff, aPtOnCurv); + theDI << "From point " << aBuff << " (with parameter " << aRetCurvParam << ") to "; + Sprintf(aBuff, "%s_lin", theArgv[1]); + + Handle(Geom2d_Curve) aLine = new Geom2d_Line(aLinSegm); + DrawTrSurf::Set(aBuff, aLine); + theDI << "the line " << aBuff << ".\n"; + Sprintf(aBuff, "%s_norm", theArgv[1]); + aLine = new Geom2d_Line(aPtOnCurv, aRetVec); + aLine = new Geom2d_TrimmedCurve(aLine, 0.0, aDefl); + DrawTrSurf::Set(aBuff, aLine); + theDI << "The deflection is measured along the line " << aBuff << ".\n"; + + return 0; +} void GeomliteTest::API2dCommands(Draw_Interpretor& theCommands) { @@ -588,4 +681,10 @@ void GeomliteTest::API2dCommands(Draw_Interpretor& theCommands) intersect_ana,g); theCommands.Add("intconcon", "intersect conic curve1 and conic curve2 using IntAna", __FILE__, intconcon, g); + + theCommands.Add("2ddeviation", "deviation2d result curve [-i U0] [-d N]\n" + "\"-i\" sets an initial parameter for computation by iterative method;\n" + "\"-d\" sets number of sub-intervals for searching. Default val is 2.", + __FILE__, deviation, g); + } diff --git a/tests/lowalgos/2ddeviation/A1 b/tests/lowalgos/2ddeviation/A1 new file mode 100644 index 0000000000..b0c66cc8e5 --- /dev/null +++ b/tests/lowalgos/2ddeviation/A1 @@ -0,0 +1,24 @@ +set ExpectDeviation 0.36597801294402493 + +restore [locate_data_file OCC538.brep] r +smallview -2D- +pcurve r + +trim cc r_3 1.5704826035188950 3.1409652070377900 +don cc +2dfit + +set log_result1 [2ddeviation res cc -d 8] +regexp {Computed value is: +([-0-9.+eE]+)} $log_result1 full aDev1 +checkreal FoundDeviation $aDev1 $ExpectDeviation 1.0e-9 0.0 +checkview -screenshot -2d -path ${imagedir}/${test_image}_1.png + +reverse cc +don cc +2dfit + +set log_result2 [2ddeviation res cc -d 8] +regexp {Computed value is: +([-0-9.+eE]+)} $log_result2 full aDev2 +checkview -screenshot -2d -path ${imagedir}/${test_image}_2.png +checkreal FoundDeviation $aDev2 $ExpectDeviation 1.0e-9 0.0 + diff --git a/tests/lowalgos/2ddeviation/A2 b/tests/lowalgos/2ddeviation/A2 new file mode 100644 index 0000000000..48b9a1b55e --- /dev/null +++ b/tests/lowalgos/2ddeviation/A2 @@ -0,0 +1,24 @@ +set ExpectDeviation 0.019156616993488952 + +restore [locate_data_file bug28211_cface7.brep] r +smallview -2D- +pcurve r + +trim cc r_4 0 0.970848497621606 +don cc +2dfit + +set log_result1 [2ddeviation res cc -d 8] +regexp {Computed value is: +([-0-9.+eE]+)} $log_result1 full aDev1 +checkreal FoundDeviation $aDev1 $ExpectDeviation 1.0e-9 0.0 +checkview -screenshot -2d -path ${imagedir}/${test_image}_1.png + +reverse cc +don cc +2dfit + +set log_result2 [2ddeviation res cc -d 8] +regexp {Computed value is: +([-0-9.+eE]+)} $log_result2 full aDev2 +checkview -screenshot -2d -path ${imagedir}/${test_image}_2.png +checkreal FoundDeviation $aDev2 $ExpectDeviation 1.0e-9 0.0 + diff --git a/tests/lowalgos/2ddeviation/A3 b/tests/lowalgos/2ddeviation/A3 new file mode 100644 index 0000000000..e678f899a3 --- /dev/null +++ b/tests/lowalgos/2ddeviation/A3 @@ -0,0 +1,22 @@ +set ExpectDeviation 199.9999995 + +circle cc 0 0 100 +trim cc cc 0.0001 2*pi-0.0001 + +smallview -2D- +don cc +2dfit + +set log_result1 [2ddeviation res cc] +regexp {Computed value is: +([-0-9.+eE]+)} $log_result1 full aDev1 +checkreal FoundDeviation $aDev1 $ExpectDeviation 1.0e-9 0.0 +checkview -screenshot -2d -path ${imagedir}/${test_image}_1.png + +reverse cc +don cc +2dfit + +set log_result2 [2ddeviation res cc] +regexp {Computed value is: +([-0-9.+eE]+)} $log_result2 full aDev2 +checkreal FoundDeviation $aDev2 $ExpectDeviation 1.0e-9 0.0 +checkview -screenshot -2d -path ${imagedir}/${test_image}_2.png diff --git a/tests/lowalgos/grids.list b/tests/lowalgos/grids.list index 398ff7c2b2..68f1679350 100644 --- a/tests/lowalgos/grids.list +++ b/tests/lowalgos/grids.list @@ -4,4 +4,5 @@ 004 extcc 005 2dgcc 006 intss -007 classifier \ No newline at end of file +007 classifier +008 2ddeviation \ No newline at end of file -- 2.39.5