0031661: Modeling Data - Exception when projecting parabola or hyperbola to plane

ProjLib/ProjLib_ProjectOnPlane.cxx - formatting

0031661: Modeling Data - Algorithm crashes when projecting parabola or hyperbola to plane

ProjLib/ProjLib_ProjectOnPlane.cxx - building of analytical parabola and hyperbola is added
bugs/moddata_3/bug31661_* - new test cases are added

src/ProjLib/ProjLib_ProjectOnPlane.cxx

@@ -47,8 +47,14 @@
 #include <GeomLib_Tool.hxx>
 #include <math_Jacobi.hxx>
 #include <math_Matrix.hxx>
+#include <gce_MakeParab.hxx>
+#include <gce_MakeDir.hxx>
+#include <LProp3d_CLProps.hxx>
+#include <math_Function.hxx>
+#include <math_BrentMinimum.hxx>
 
-
+const Standard_Real aParabolaLimit = 20000.;
+const Standard_Real aHyperbolaLimit = 10.;
 
 //=======================================================================
 //function : OnPlane_Value
@@ -269,10 +275,14 @@ public :
   }
 
   Standard_Real FirstParameter() const
-    {return myCurve->FirstParameter();}
+  {
+    return myCurve->FirstParameter();
+  }
 
   Standard_Real LastParameter() const
-    {return myCurve->LastParameter();}
+  {
+    return myCurve->LastParameter();
+  }
 
   Standard_Boolean Value(const Standard_Real   theT,
     NCollection_Array1<gp_Pnt2d>& /*thePnt2d*/,
@@ -292,6 +302,35 @@ public :
 };
 
 
+//=======================================================================
+//  class  : ProjLib_MaxCurvature
+//purpose  : Use to search apex of parabola or hyperbola, which is its projection  
+//           on a plane. Apex is point with maximal curvature
+//=======================================================================
+
+class ProjLib_MaxCurvature : public math_Function
+
+{
+
+public:
+
+  ProjLib_MaxCurvature(LProp3d_CLProps& theProps):
+    myProps(&theProps)
+  {
+  }
+
+  virtual Standard_Boolean Value(const Standard_Real X, Standard_Real& F)
+  {
+    myProps->SetParameter(X);
+    F = -myProps->Curvature();
+    return Standard_True;
+  }
+  
+private:
+
+  LProp3d_CLProps* myProps;
+ 
+};
 
 
 //=====================================================================//
@@ -319,10 +358,20 @@ static void  PerformApprox (const Handle(Adaptor3d_Curve)& C,
 
   Standard_Integer Deg1, Deg2;
   Deg1 = 8; Deg2 = 8;
-  
+  if (C->GetType() == GeomAbs_Parabola)
+  {
+    Deg1 = 2; Deg2 = 2;
+  }
+  Standard_Integer aNbSegm = 100;
+  if (C->GetType() == GeomAbs_Hyperbola)
+  {
+    Deg1 = 14;
+    Deg2 = 14;
+    aNbSegm = 1000;
+  }
   Approx_FitAndDivide Fit(Deg1, Deg2, Precision::Approximation(),
     Precision::PApproximation(), Standard_True);
-  Fit.SetMaxSegments(100);
+  Fit.SetMaxSegments(aNbSegm); 
   Fit.Perform(F);
   if (!Fit.IsAllApproximated())
   {
@@ -348,9 +397,11 @@ static void  PerformApprox (const Handle(Adaptor3d_Curve)& C,
   TColStd_Array1OfReal Knots(1, NbCurves + 1);  //Noeuds de la BSpline
 
   Standard_Integer Compt = 1;
+  Standard_Real anErrMax = 0., anErr3d, anErr2d;
   for (i = 1; i <= Fit.NbMultiCurves(); i++) {
     Fit.Parameters(i, Knots(i), Knots(i + 1));
-    
+    Fit.Error(i, anErr3d, anErr2d);
+    anErrMax = Max(anErrMax, anErr3d);
     AppParCurves_MultiCurve MC = Fit.Value(i);   //Charge la Ieme Curve
     TColgp_Array1OfPnt LocalPoles(1, MC.Degree() + 1);//Recupere les poles
     MC.Curve(1, LocalPoles);
@@ -389,6 +440,17 @@ static void  PerformApprox (const Handle(Adaptor3d_Curve)& C,
   Mults.SetValue(NbKnots, MaxDeg + 1);
   BSplineCurvePtr =
     new Geom_BSplineCurve(Poles, Knots, Mults, MaxDeg, Standard_False);
+
+  //Try to smooth
+  Standard_Integer m1 = MaxDeg - 1;
+  for (i = 2; i < NbKnots; ++i)
+  {
+    if (BSplineCurvePtr->Multiplicity(i) == MaxDeg)
+    {
+      BSplineCurvePtr->RemoveKnot(i, m1, anErrMax);
+    }
+  }
+
 }
 
 
@@ -537,6 +599,7 @@ void ProjLib_ProjectOnPlane::Load(const Handle(Adaptor3d_Curve)&    C,
   Handle(Geom_Circle)    GeomCirclePtr;
   Handle(Geom_Ellipse)   GeomEllipsePtr;
   Handle(Geom_Hyperbola) GeomHyperbolaPtr;
+  Handle(Geom_Parabola) GeomParabolaPtr;
 
   gp_Lin   aLine;
   gp_Elips Elips;
@@ -806,7 +869,6 @@ void ProjLib_ProjectOnPlane::Load(const Handle(Adaptor3d_Curve)&    C,
   break;
   case GeomAbs_Parabola:
   {
-      myKeepParam = Standard_True;
     //     P(u) = O + (u*u)/(4*f) * Xc + u * Yc
     // ==> Q(u) = f(P(u)) 
     //          = f(O) + (u*u)/(4*f) * f(Xc) + u * f(Yc)
@@ -815,50 +877,51 @@ void ProjLib_ProjectOnPlane::Load(const Handle(Adaptor3d_Curve)&    C,
     gp_Ax2   AxeRef = Parab.Position();
     gp_Vec Xc = ProjectVec(myPlane, myDirection, gp_Vec(AxeRef.XDirection()));
     gp_Vec Yc = ProjectVec(myPlane, myDirection, gp_Vec(AxeRef.YDirection()));
-
-      // fix for case when no one action is done. 28.03.2002
-      Standard_Boolean alocalIsDone = Standard_False;
-
-      if ( Abs( Yc.Magnitude() - 1.) < Precision::Confusion()) {
     gp_Pnt P = ProjectPnt(myPlane, myDirection, AxeRef.Location());
-	if ( Xc.Magnitude() < Precision::Confusion()) {
+
+    myIsApprox = Standard_False;
+
+    if ((Abs(Yc.Magnitude() - 1.) < Precision::Confusion()) &&
+      (Xc.Magnitude() < Precision::Confusion()))
+    {
       myType = GeomAbs_Line;
       aLine = gp_Lin(P, gp_Dir(Yc));
-          
       GeomLinePtr = new Geom_Line(aLine);
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:29 2002 Begin
-	  GeomAdaptor_Curve aGACurve(GeomLinePtr);
-	  myResult = new GeomAdaptor_Curve(aGACurve);
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:30 2002 End
-          alocalIsDone = Standard_True;
     }
     else if (Xc.IsNormal(Yc, Precision::Angular())) {
       myType = GeomAbs_Parabola;
       Standard_Real F = Parab.Focal() / Xc.Magnitude();
-	  gp_Parab aParab = gp_Parab( gp_Ax2(P,Xc^Yc,Xc), F);
+      gp_Parab aProjParab = gp_Parab(gp_Ax2(P, Xc^Yc, Xc), F);
-          Handle(Geom_Parabola) GeomParabolaPtr =
+      GeomParabolaPtr =
-	    new Geom_Parabola(aParab) ;
+        new Geom_Parabola(aProjParab);
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:29 2002 Begin
-	  GeomAdaptor_Curve aGACurve(GeomParabolaPtr);
-	  myResult = new GeomAdaptor_Curve(aGACurve);
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:30 2002 End
-          alocalIsDone = Standard_True;
     }
-      }
+    else if (Yc.Magnitude() < Precision::Confusion() || 
-      if (!alocalIsDone)/*else*/ {
+      Yc.IsParallel(Xc, Precision::Angular()))
+    {
       myIsApprox = Standard_True;
-	myType     = GeomAbs_BSplineCurve;
+    }
-        PerformApprox(myCurve,myPlane,myDirection,ApproxCurve);
+    else if(!myKeepParam)
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:29 2002 Begin
+    {
-	GeomAdaptor_Curve aGACurve(ApproxCurve);
+      // Try building parabola with help of apex position
-	myResult = new GeomAdaptor_Curve(aGACurve);
+      myIsApprox = !BuildParabolaByApex(GeomParabolaPtr);
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:30 2002 End
+    }
+    else
+    {
+      myIsApprox = Standard_True;
+    }
+
+    if (!myIsApprox)
+    {
+      GetTrimmedResult(GeomParabolaPtr);
+    }
+    else
+    {
+      BuildByApprox(aParabolaLimit);
     }
   }
   break;
   case GeomAbs_Hyperbola:
   {
-      myKeepParam = Standard_True;
     //     P(u) = O + R1 * Cosh(u) * Xc + R2 * Sinh(u) * Yc
     // ==> Q(u) = f(P(u)) 
     //          = f(O) + R1 * Cosh(u) * f(Xc) + R2 * Sinh(u) * f(Yc)
@@ -871,6 +934,7 @@ void ProjLib_ProjectOnPlane::Load(const Handle(Adaptor3d_Curve)&    C,
     Standard_Real aR1 = Hypr.MajorRadius();
     Standard_Real aR2 = Hypr.MinorRadius();
     gp_Dir Z = myPlane.Direction();
+    myIsApprox = Standard_False;
 
     if (Xc.Magnitude() < Precision::Confusion()) {
       myType = GeomAbs_Hyperbola;
@@ -878,10 +942,6 @@ void ProjLib_ProjectOnPlane::Load(const Handle(Adaptor3d_Curve)&    C,
       Hypr = gp_Hypr(gp_Ax2(P, Z, X), 0., aR2 * Yc.Magnitude());
       GeomHyperbolaPtr =
         new Geom_Hyperbola(Hypr);
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:29 2002 Begin
-	GeomAdaptor_Curve aGACurve(GeomHyperbolaPtr);
-	myResult = new GeomAdaptor_Curve(aGACurve);
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:30 2002 End
     }
     else if (Yc.Magnitude() < Precision::Confusion()) {
       myType = GeomAbs_Hyperbola;
@@ -889,10 +949,6 @@ void ProjLib_ProjectOnPlane::Load(const Handle(Adaptor3d_Curve)&    C,
         gp_Hypr(gp_Ax2(P, Z, gp_Dir(Xc)), aR1 * Xc.Magnitude(), 0.);
       GeomHyperbolaPtr =
         new Geom_Hyperbola(Hypr);
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:29 2002 Begin
-	GeomAdaptor_Curve aGACurve(GeomHyperbolaPtr);
-	myResult = new GeomAdaptor_Curve(aGACurve);
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:30 2002 End
     }
     else if (Xc.IsNormal(Yc, Precision::Angular())) {
       myType = GeomAbs_Hyperbola;
@@ -900,19 +956,27 @@ void ProjLib_ProjectOnPlane::Load(const Handle(Adaptor3d_Curve)&    C,
         aR1 * Xc.Magnitude(), aR2 * Yc.Magnitude());
       GeomHyperbolaPtr =
         new Geom_Hyperbola(Hypr);
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:29 2002 Begin
-	GeomAdaptor_Curve aGACurve(GeomHyperbolaPtr);
-	myResult = new GeomAdaptor_Curve(aGACurve);
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:30 2002 End
     }
-      else {
+    else if (Yc.Magnitude() < Precision::Confusion() ||
+      Yc.IsParallel(Xc, Precision::Angular()))
+    {
       myIsApprox = Standard_True;
-	myType     = GeomAbs_BSplineCurve;
+    }
-        PerformApprox(myCurve,myPlane,myDirection,ApproxCurve);
+    else if(!myKeepParam)
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:29 2002 Begin
+    {
-	GeomAdaptor_Curve aGACurve(ApproxCurve);
+      myIsApprox = !BuildHyperbolaByApex(GeomHyperbolaPtr);
-	myResult = new GeomAdaptor_Curve(aGACurve);
+    }
-//  Modified by Sergey KHROMOV - Tue Jan 29 16:57:30 2002 End
+    else
+    {
+      myIsApprox = Standard_True;
+    }
+    if ( !myIsApprox )
+    {
+      GetTrimmedResult(GeomHyperbolaPtr);
+    }
+    else
+    {
+      BuildByApprox(aHyperbolaLimit);
     }
   }
   break;
@@ -1458,7 +1522,7 @@ Handle(Geom_BezierCurve)  ProjLib_ProjectOnPlane::Bezier() const
 }
 
 //=======================================================================
-//function : Bezier
+//function : BSpline
 //purpose  : 
 //=======================================================================
 
@@ -1470,3 +1534,222 @@ Handle(Geom_BSplineCurve)  ProjLib_ProjectOnPlane::BSpline() const
   return myResult->BSpline();
 }
 
+//=======================================================================
+//function : GetTrimmedResult
+//purpose  : 
+//=======================================================================
+
+void  ProjLib_ProjectOnPlane::GetTrimmedResult(const Handle(Geom_Curve)& theProjCurve)
+{
+  gp_Lin aLin;
+  gp_Parab aParab;
+  gp_Hypr aHypr;
+  if (myType == GeomAbs_Line)
+  {
+    aLin = Handle(Geom_Line)::DownCast(theProjCurve)->Lin();
+  }
+  else if (myType == GeomAbs_Parabola)
+  {
+    aParab = Handle(Geom_Parabola)::DownCast(theProjCurve)->Parab();
+  }
+  else if (myType == GeomAbs_Hyperbola)
+  {
+    aHypr = Handle(Geom_Hyperbola)::DownCast(theProjCurve)->Hypr();
+  }
+
+  myFirstPar = theProjCurve->FirstParameter();
+  myLastPar = theProjCurve->LastParameter();
+  if (!Precision::IsInfinite(myCurve->FirstParameter()))
+  {
+    gp_Pnt aP = myCurve->Value(myCurve->FirstParameter());
+    aP = ProjectPnt(myPlane, myDirection, aP);
+    if (myType == GeomAbs_Line)
+    {
+      myFirstPar = ElCLib::Parameter(aLin, aP);
+    }
+    else if (myType == GeomAbs_Parabola)
+    {
+      myFirstPar = ElCLib::Parameter(aParab, aP);
+    }
+    else if (myType == GeomAbs_Hyperbola)
+    {
+      myFirstPar = ElCLib::Parameter(aHypr, aP);
+    }
+    else
+    {
+      GeomLib_Tool::Parameter(theProjCurve, aP, Precision::Confusion(), myFirstPar);
+    }
+  }
+  if (!Precision::IsInfinite(myCurve->LastParameter()))
+  {
+    gp_Pnt aP = myCurve->Value(myCurve->LastParameter());
+    aP = ProjectPnt(myPlane, myDirection, aP);
+    if (myType == GeomAbs_Line)
+    {
+      myLastPar = ElCLib::Parameter(aLin, aP);
+    }
+    else if (myType == GeomAbs_Parabola)
+    {
+      myLastPar = ElCLib::Parameter(aParab, aP);
+    }
+    else if (myType == GeomAbs_Hyperbola)
+    {
+      myLastPar = ElCLib::Parameter(aHypr, aP);
+    }
+    else
+    {
+      GeomLib_Tool::Parameter(theProjCurve, aP, Precision::Confusion(), myLastPar);
+    }
+  }
+  myResult = new GeomAdaptor_Curve(theProjCurve, myFirstPar, myLastPar);
+
+}
+
+//=======================================================================
+//function : BuildParabolaByApex
+//purpose  : 
+//=======================================================================
+
+Standard_Boolean ProjLib_ProjectOnPlane::BuildParabolaByApex(Handle(Geom_Curve)& theGeomParabolaPtr)
+{
+  //
+  //Searching parabola apex as point with maximal curvature
+  Standard_Real aF = myCurve->Parabola().Focal();
+  GeomAbs_CurveType aCurType = myType;
+  myType = GeomAbs_OtherCurve; //To provide correct calculation of derivativesb by projection for
+                               //copy of instance;
+  Handle(Adaptor3d_Curve) aProjCrv = ShallowCopy();
+  myType = aCurType;
+  LProp3d_CLProps aProps(aProjCrv, 2, Precision::Confusion());
+  ProjLib_MaxCurvature aMaxCur(aProps);
+  math_BrentMinimum aSolver(Precision::PConfusion());
+  aSolver.Perform(aMaxCur, -10.*aF, 0., 10.*aF);
+
+  if (!aSolver.IsDone())
+  {
+    return Standard_False;
+  }
+  
+  Standard_Real aT;
+  aT = aSolver.Location();
+  aProps.SetParameter(aT);
+  gp_Pnt aP0 = aProps.Value();
+  gp_Vec aDY = aProps.D1();
+  gp_Dir anYDir(aDY);
+  gp_Dir anXDir;
+  Standard_Real aCurv = aProps.Curvature();
+  if (Precision::IsInfinite(aCurv) || aCurv < Precision::Confusion())
+  {
+    return Standard_False;
+  }
+  aProps.Normal(anXDir);
+  //
+  gp_Lin anXLine(aP0, anXDir);
+  gp_Pnt aP1 = Value(aT + 10.*aF);
+  //
+  Standard_Real anX = ElCLib::LineParameter(anXLine.Position(), aP1);
+  Standard_Real anY = anXLine.Distance(aP1);
+  Standard_Real aNewF = anY * anY / 4. / anX;
+  gp_Dir anN = anXDir^anYDir;
+  gp_Ax2 anA2(aP0, anN, anXDir);
+  gce_MakeParab aMkParab(anA2, aNewF);
+  if (!aMkParab.IsDone())
+  {
+    return Standard_False;
+  }
+
+  gp_Parab aProjParab = aMkParab.Value();
+  
+  myType = GeomAbs_Parabola;
+  theGeomParabolaPtr = new Geom_Parabola(aProjParab);
+  //GetTrimmedResult(theGeomParabolaPtr);
+
+  return Standard_True;
+}
+
+//=======================================================================
+//function : BuildHyperbolaByApex
+//purpose  : 
+//=======================================================================
+
+Standard_Boolean ProjLib_ProjectOnPlane::BuildHyperbolaByApex(Handle(Geom_Curve)& theGeomHyperbolaPtr)
+{
+  //Try to build hyperbola with help of apex position
+  GeomAbs_CurveType aCurType = myType;
+  myType = GeomAbs_OtherCurve; //To provide correct calculation of derivativesb by projection for
+                               //copy of instance;
+  Handle(Adaptor3d_Curve) aProjCrv = ShallowCopy();
+  myType = aCurType;
+  //Searching hyperbola apex as point with maximal curvature
+  LProp3d_CLProps aProps(aProjCrv, 2, Precision::Confusion());
+  ProjLib_MaxCurvature aMaxCur(aProps);
+  math_BrentMinimum aSolver(Precision::PConfusion());
+  aSolver.Perform(aMaxCur, -5., 0., 5.);
+
+  if (aSolver.IsDone())
+  {
+    Standard_Real aT;
+    aT = aSolver.Location();
+    aProps.SetParameter(aT);
+    Standard_Real aCurv = aProps.Curvature();
+    if (Precision::IsInfinite(aCurv) || aCurv < Precision::Confusion())
+    {
+      return Standard_False;
+    }
+    else
+    {
+      gp_Hypr Hypr = myCurve->Hyperbola();
+      gp_Ax2 AxeRef = Hypr.Position();
+      gp_Pnt P = ProjectPnt(myPlane, myDirection, AxeRef.Location());
+      gp_Dir Z = myPlane.Direction();
+      gp_Pnt aP0 = aProps.Value();
+      gp_Dir anXDir = gce_MakeDir(P, aP0);
+      gp_Dir anYDir = gce_MakeDir(aProps.D1());
+      //
+      Standard_Real aMajRad = P.Distance(aP0);
+      gp_Pnt aP1 = Value(aT + 1.);
+      gp_Vec aV(P, aP1);
+      Standard_Real anX = aV * anXDir;
+      Standard_Real anY = aV * anYDir;
+      Standard_Real aMinRad = anY / Sqrt(anX * anX / aMajRad / aMajRad - 1.);
+      gp_Ax2 anA2(P, Z, anXDir);
+      gp_Hypr anHypr(anA2, aMajRad, aMinRad);
+      theGeomHyperbolaPtr =
+        new Geom_Hyperbola(anHypr);
+      myType = GeomAbs_Hyperbola;
+    }
+  }
+  else
+  {
+    return Standard_False;
+  }
+  return Standard_True;
+}
+
+//=======================================================================
+//function : BuilByApprox
+//purpose  : 
+//=======================================================================
+
+void ProjLib_ProjectOnPlane::BuildByApprox(const Standard_Real theLimitParameter)
+{
+  myType = GeomAbs_BSplineCurve;
+  Handle(Geom_BSplineCurve)  anApproxCurve;
+  if (Precision::IsInfinite(myCurve->FirstParameter()) ||
+    Precision::IsInfinite(myCurve->LastParameter()))
+  {
+    //To avoid exception in approximation
+    Standard_Real f = Max(-theLimitParameter, myCurve->FirstParameter());
+    Standard_Real l = Min(theLimitParameter, myCurve->LastParameter());
+    Handle(Adaptor3d_Curve) aTrimCurve = myCurve->Trim(f, l, Precision::Confusion());
+    PerformApprox(aTrimCurve, myPlane, myDirection, anApproxCurve);
+  }
+  else
+  {
+    PerformApprox(myCurve, myPlane, myDirection, anApproxCurve);
+  }
+  myFirstPar = anApproxCurve->FirstParameter();
+  myLastPar = anApproxCurve->LastParameter();
+  GeomAdaptor_Curve aGACurve(anApproxCurve);
+  myResult = new GeomAdaptor_Curve(aGACurve);
+}

src/ProjLib/ProjLib_ProjectOnPlane.hxx

@@ -199,9 +199,12 @@ public:
 
 protected:
 
+  void GetTrimmedResult(const Handle(Geom_Curve)& theProjCurve);
 
+  Standard_Boolean BuildParabolaByApex(Handle(Geom_Curve)& theGeomParabolaPtr);
+  Standard_Boolean BuildHyperbolaByApex(Handle(Geom_Curve)& theGeomParabolaPtr);
 
-
+  void BuildByApprox(const Standard_Real theLimitParameter);
 
 private:
 

tests/bugs/moddata_3/bug31661_1

@@ -0,0 +1,26 @@
+puts ""
+puts "=========================================================================="
+puts "OCC31661: Modeling Data - Algorithm crashes when projecting parabola or hyperbola to plane"
+puts "=========================================================================="
+puts ""
+
+parabola p 0 0 0  1 1 1  2 0 -2  10
+plane pln 0 0 0  0 0 1
+projonplane r p pln 0
+
+if {![regexp {Parabola} [dump r]]} {
+  puts "ERROR: Projected curve is not a parabola"
+}
+
+trim p p -100 100
+projonplane rp p pln 0
+
+if {![regexp {Parabola} [dump rp]]} {
+  puts "ERROR: Projected curve is not a parabola"
+}
+
+checklength rp -l 408.40363195229503
+bounds rp t1 t2
+checkreal t1 [dval t1] -91.077748943768597 1.e-7 1.e-7
+checkreal t2 [dval t2] 72.221567418462357 1.e-7 1.e-7
+

tests/bugs/moddata_3/bug31661_2

@@ -0,0 +1,26 @@
+puts ""
+puts "=========================================================================="
+puts "OCC31661: Modeling Data - Algorithm crashes when projecting parabola or hyperbola to plane"
+puts "=========================================================================="
+puts ""
+
+hyperbola h 0 0 0  1 1 1  2 0 -2  10 10
+plane pln 0 0 0  0 0 1
+projonplane rh h pln 0
+
+if {![regexp {Hyperbola} [dump rh]]} {
+  puts "ERROR: Projected curve is not a hyperbola"
+}
+
+trim h h -5 5
+projonplane rh h pln 0
+
+if {![regexp {Hyperbola} [dump rh]]} {
+  puts "ERROR: Projected curve is not a hyperbola"
+}
+
+checklength rh -l 1664.3732976598988
+bounds rh t1 t2
+checkreal t1 [dval t1] -5.23179933356147 1.e-7 1.e-7
+checkreal t2 [dval t2] 4.76820064934972 1.e-7 1.e-7
+