// Created on: 1992-08-06
// Created by: Laurent BUCHARD
// Copyright (c) 1992-1999 Matra Datavision
-// Copyright (c) 1999-2012 OPEN CASCADE SAS
+// Copyright (c) 1999-2014 OPEN CASCADE SAS
//
-// The content of this file is subject to the Open CASCADE Technology Public
-// License Version 6.5 (the "License"). You may not use the content of this file
-// except in compliance with the License. Please obtain a copy of the License
-// at http://www.opencascade.org and read it completely before using this file.
+// This file is part of Open CASCADE Technology software library.
//
-// The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
-// main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
+// 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.
//
-// The Original Code and all software distributed under the License is
-// distributed on an "AS IS" basis, without warranty of any kind, and the
-// Initial Developer hereby disclaims all such warranties, including without
-// limitation, any warranties of merchantability, fitness for a particular
-// purpose or non-infringement. Please see the License for the specific terms
-// and conditions governing the rights and limitations under the License.
+// Alternatively, this file may be used under the terms of Open CASCADE
+// commercial license or contractual agreement.
//----------------------------------------------------------------------
//-- Purposse: Geometric Intersection between two Natural Quadric
//-- If the intersection is not a conic,
//-- analytical methods must be called.
//----------------------------------------------------------------------
-#ifndef DEB
+#ifndef OCCT_DEBUG
#define No_Standard_RangeError
#define No_Standard_OutOfRange
#endif
-#include <IntAna_QuadQuadGeo.ixx>
-#include <IntAna_IntConicQuad.hxx>
-#include <StdFail_NotDone.hxx>
-#include <Standard_DomainError.hxx>
-#include <Standard_OutOfRange.hxx>
-#include <math_DirectPolynomialRoots.hxx>
-
-#include <gp.hxx>
-#include <gp_Pln.hxx>
-#include <gp_Vec.hxx>
-#include <ElSLib.hxx>
#include <ElCLib.hxx>
-
+#include <ElSLib.hxx>
+#include <gp.hxx>
+#include <gp_Circ.hxx>
+#include <gp_Cone.hxx>
+#include <gp_Cylinder.hxx>
#include <gp_Dir.hxx>
-#include <gp_XYZ.hxx>
+#include <gp_Dir2d.hxx>
+#include <gp_Elips.hxx>
+#include <gp_Hypr.hxx>
+#include <gp_Lin.hxx>
+#include <gp_Parab.hxx>
+#include <gp_Pln.hxx>
+#include <gp_Pnt.hxx>
#include <gp_Pnt2d.hxx>
+#include <gp_Sphere.hxx>
+#include <gp_Torus.hxx>
+#include <gp_Vec.hxx>
#include <gp_Vec2d.hxx>
-#include <gp_Dir2d.hxx>
-
+#include <gp_XYZ.hxx>
+#include <IntAna_IntConicQuad.hxx>
+#include <IntAna_QuadQuadGeo.hxx>
+#include <math_DirectPolynomialRoots.hxx>
+#include <Standard_DomainError.hxx>
+#include <Standard_OutOfRange.hxx>
+#include <StdFail_NotDone.hxx>
static
gp_Ax2 DirToAx2(const gp_Pnt& P,const gp_Dir& D);
void RefineDir(gp_Dir& aDir);
//=======================================================================
-//class :
+//class : AxeOperator
//purpose : O p e r a t i o n s D i v e r s e s s u r d e s A x 1
//=======================================================================
class AxeOperator {
AxeOperator(const gp_Ax1& A1,const gp_Ax1& A2);
void Distance(Standard_Real& dist,
- Standard_Real& Param1,
- Standard_Real& Param2);
+ Standard_Real& Param1,
+ Standard_Real& Param2);
gp_Pnt PtIntersect() {
return ptintersect;
protected:
Standard_Real Det33(const Standard_Real a11,
- const Standard_Real a12,
- const Standard_Real a13,
- const Standard_Real a21,
- const Standard_Real a22,
- const Standard_Real a23,
- const Standard_Real a31,
- const Standard_Real a32,
- const Standard_Real a33) {
+ const Standard_Real a12,
+ const Standard_Real a13,
+ const Standard_Real a21,
+ const Standard_Real a22,
+ const Standard_Real a23,
+ const Standard_Real a31,
+ const Standard_Real a32,
+ const Standard_Real a33) {
Standard_Real theReturn =
a11*(a22*a33-a32*a23) - a21*(a12*a33-a32*a13) + a31*(a12*a23-a22*a13) ;
return theReturn ;
//function : AxeOperator::AxeOperator
//purpose :
//=======================================================================
- AxeOperator::AxeOperator(const gp_Ax1& A1,const gp_Ax1& A2)
+AxeOperator::AxeOperator(const gp_Ax1& A1,const gp_Ax1& A2)
{
- myEPSILON_DISTANCE=0.00000000000001;
- myEPSILON_AXES_PARA=0.000000000001;
+ myEPSILON_DISTANCE=1.0e-14;
+ myEPSILON_AXES_PARA=Precision::Angular();
Axe1=A1;
Axe2=A2;
//---------------------------------------------------------------------
}
else {
thedistance = Abs(gp_Vec(perp.Normalized()).Dot(gp_Vec(Axe1.Location(),
- Axe2.Location())));
+ Axe2.Location())));
}
//--- check if Axis are Coplanar
Standard_Real D33;
if(thedistance<myEPSILON_DISTANCE) {
D33=Det33(V1.X(),V1.Y(),V1.Z()
- ,V2.X(),V2.Y(),V2.Z()
- ,P1.X()-P2.X(),P1.Y()-P2.Y(),P1.Z()-P2.Z());
+ ,V2.X(),V2.Y(),V2.Z()
+ ,P1.X()-P2.X(),P1.Y()-P2.Y(),P1.Z()-P2.Z());
if(Abs(D33)<=myEPSILON_DISTANCE) {
thecoplanar=Standard_True;
}
A=(smy * V2.X() - smx * V2.Y())/Det1;
}
else if((Det2!=0.0)
- && ((Abs(Det2) >= Abs(Det1))
- &&(Abs(Det2) >= Abs(Det3)))) {
+ && ((Abs(Det2) >= Abs(Det1))
+ &&(Abs(Det2) >= Abs(Det3)))) {
A=(smz * V2.Y() - smy * V2.Z())/Det2;
}
else {
A=(smz * V2.X() - smx * V2.Z())/Det3;
}
ptintersect.SetCoord( P1.X() + A * V1.X()
- ,P1.Y() + A * V1.Y()
- ,P1.Z() + A * V1.Z());
+ ,P1.Y() + A * V1.Y()
+ ,P1.Z() + A * V1.Z());
}
else {
ptintersect.SetCoord(0,0,0); //-- Pour eviter des FPE
//function : Distance
//purpose :
//=======================================================================
- void AxeOperator::Distance(Standard_Real& dist,Standard_Real& Param1,Standard_Real& Param2)
+void AxeOperator::Distance(Standard_Real& dist,
+ Standard_Real& Param1,
+ Standard_Real& Param2)
{
- gp_Vec O1O2(Axe1.Location(),Axe2.Location()); //-----------------------------
+ gp_Vec O1O2(Axe1.Location(),Axe2.Location());
gp_Dir U1 = Axe1.Direction(); //-- juste pour voir.
gp_Dir U2 = Axe2.Direction();
gp_Dir N = U1.Crossed(U2);
Standard_Real D = Det33(U1.X(),U2.X(),N.X(),
- U1.Y(),U2.Y(),N.Y(),
- U1.Z(),U2.Z(),N.Z());
+ U1.Y(),U2.Y(),N.Y(),
+ U1.Z(),U2.Z(),N.Z());
if(D) {
dist = Det33(U1.X(),U2.X(),O1O2.X(),
- U1.Y(),U2.Y(),O1O2.Y(),
- U1.Z(),U2.Z(),O1O2.Z()) / D;
+ U1.Y(),U2.Y(),O1O2.Y(),
+ U1.Z(),U2.Z(),O1O2.Z()) / D;
Param1 = Det33(O1O2.X(),U2.X(),N.X(),
- O1O2.Y(),U2.Y(),N.Y(),
- O1O2.Z(),U2.Z(),N.Z()) / (-D);
+ O1O2.Y(),U2.Y(),N.Y(),
+ O1O2.Z(),U2.Z(),N.Z()) / (-D);
//------------------------------------------------------------
//-- On resout P1 * Dir1 + P2 * Dir2 + d * N = O1O2
//-- soit : Segment perpendiculaire : O1+P1 D1
//-- O2-P2 D2
Param2 = Det33(U1.X(),O1O2.X(),N.X(),
- U1.Y(),O1O2.Y(),N.Y(),
- U1.Z(),O1O2.Z(),N.Z()) / (D);
+ U1.Y(),O1O2.Y(),N.Y(),
+ U1.Z(),O1O2.Z(),N.Z()) / (D);
}
}
//=======================================================================
//function : IntAna_QuadQuadGeo
//purpose : Empty constructor
//=======================================================================
- IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(void)
+IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(void)
: done(Standard_False),
nbint(0),
typeres(IntAna_Empty),
pt1(0,0,0),
pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
param1(0),
param2(0),
+ param3(0),
+ param4(0),
param1bis(0),
param2bis(0),
myCommonGen(Standard_False),
//function : InitTolerances
//purpose :
//=======================================================================
- void IntAna_QuadQuadGeo::InitTolerances()
+void IntAna_QuadQuadGeo::InitTolerances()
{
- myEPSILON_DISTANCE = 0.00000000000001;
- myEPSILON_ANGLE_CONE = 0.000000000001;
- myEPSILON_MINI_CIRCLE_RADIUS = 0.000000001;
- myEPSILON_CYLINDER_DELTA_RADIUS = 0.0000000000001;
- myEPSILON_CYLINDER_DELTA_DISTANCE= 0.0000001;
- myEPSILON_AXES_PARA = 0.000000000001;
+ myEPSILON_DISTANCE = 1.0e-14;
+ myEPSILON_ANGLE_CONE = Precision::Angular();
+ myEPSILON_MINI_CIRCLE_RADIUS = 0.01*Precision::Confusion();
+ myEPSILON_CYLINDER_DELTA_RADIUS = 1.0e-13;
+ myEPSILON_CYLINDER_DELTA_DISTANCE= Precision::Confusion();
+ myEPSILON_AXES_PARA = Precision::Angular();
}
//=======================================================================
//function : IntAna_QuadQuadGeo
//purpose : Pln Pln
//=======================================================================
- IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Pln& P1,
- const gp_Pln& P2,
- const Standard_Real TolAng,
- const Standard_Real Tol)
+IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Pln& P1,
+ const gp_Pln& P2,
+ const Standard_Real TolAng,
+ const Standard_Real Tol)
: done(Standard_False),
nbint(0),
typeres(IntAna_Empty),
pt1(0,0,0),
pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
param1(0),
param2(0),
+ param3(0),
+ param4(0),
param1bis(0),
param2bis(0),
myCommonGen(Standard_False),
//function : Perform
//purpose :
//=======================================================================
- void IntAna_QuadQuadGeo::Perform (const gp_Pln& P1,
- const gp_Pln& P2,
- const Standard_Real TolAng,
- const Standard_Real Tol)
+void IntAna_QuadQuadGeo::Perform (const gp_Pln& P1,
+ const gp_Pln& P2,
+ const Standard_Real TolAng,
+ const Standard_Real Tol)
{
+ Standard_Real A1, B1, C1, D1, A2, B2, C2, D2, dist1, dist2, aMVD;
+ //
done=Standard_False;
+ param2bis=0.;
//
- param2bis=0.0;
-
- Standard_Real A1 = 0., B1 = 0., C1 = 0., D1 = 0., A2 = 0., B2 = 0., C2 = 0., D2 = 0.;
P1.Coefficients(A1,B1,C1,D1);
P2.Coefficients(A2,B2,C2,D2);
-
- gp_Vec vd(gp_Vec(A1,B1,C1).Crossed(gp_Vec(A2,B2,C2)));
- Standard_Real dist1= A2*P1.Location().X() + B2*P1.Location().Y() + C2*P1.Location().Z() + D2;
- Standard_Real dist2= A1*P2.Location().X() + B1*P2.Location().Y() + C1*P2.Location().Z() + D1;
-
- if(vd.Magnitude() <=TolAng) {
+ //
+ gp_Vec aVN1(A1,B1,C1);
+ gp_Vec aVN2(A2,B2,C2);
+ gp_Vec vd(aVN1.Crossed(aVN2));
+ //
+ const gp_Pnt& aLocP1=P1.Location();
+ const gp_Pnt& aLocP2=P2.Location();
+ //
+ dist1=A2*aLocP1.X() + B2*aLocP1.Y() + C2*aLocP1.Z() + D2;
+ dist2=A1*aLocP2.X() + B1*aLocP2.Y() + C1*aLocP2.Z() + D1;
+ //
+ aMVD=vd.Magnitude();
+ if(aMVD <=TolAng) {
// normalles are collinear - planes are same or parallel
- typeres = (Abs(dist1) <= Tol && Abs(dist2) <= Tol) ? IntAna_Same : IntAna_Empty;
+ typeres = (Abs(dist1) <= Tol && Abs(dist2) <= Tol) ? IntAna_Same
+ : IntAna_Empty;
}
else {
- Standard_Real denom=A1*A2 + B1*B2 + C1*C2;
-
- Standard_Real denom2 = denom*denom;
- Standard_Real ddenom = 1. - denom2;
- //denom = ( Abs(ddenom) <= 1.e-9 ) ? 1.e-9 : ddenom;
- denom = ( Abs(ddenom) <= 1.e-16 ) ? 1.e-16 : ddenom;
+ Standard_Real denom, denom2, ddenom, par1, par2;
+ Standard_Real X1, Y1, Z1, X2, Y2, Z2, aEps;
+ //
+ aEps=1.e-16;
+ denom=A1*A2 + B1*B2 + C1*C2;
+ denom2 = denom*denom;
+ ddenom = 1. - denom2;
+
+ denom = ( Abs(ddenom) <= aEps ) ? aEps : ddenom;
- Standard_Real par1 = dist1/denom;
- Standard_Real par2 = -dist2/denom;
+ par1 = dist1/denom;
+ par2 = -dist2/denom;
- gp_Vec inter1(gp_Vec(A1,B1,C1).Crossed(vd));
- gp_Vec inter2(gp_Vec(A2,B2,C2).Crossed(vd));
+ gp_Vec inter1(aVN1.Crossed(vd));
+ gp_Vec inter2(aVN2.Crossed(vd));
- Standard_Real X1=P1.Location().X() + par1*inter1.X();
- Standard_Real Y1=P1.Location().Y() + par1*inter1.Y();
- Standard_Real Z1=P1.Location().Z() + par1*inter1.Z();
- Standard_Real X2=P2.Location().X() + par2*inter2.X();
- Standard_Real Y2=P2.Location().Y() + par2*inter2.Y();
- Standard_Real Z2=P2.Location().Z() + par2*inter2.Z();
+ X1=aLocP1.X() + par1*inter1.X();
+ Y1=aLocP1.Y() + par1*inter1.Y();
+ Z1=aLocP1.Z() + par1*inter1.Z();
+ X2=aLocP2.X() + par2*inter2.X();
+ Y2=aLocP2.Y() + par2*inter2.Y();
+ Z2=aLocP2.Z() + par2*inter2.Z();
pt1=gp_Pnt((X1+X2)*0.5, (Y1+Y2)*0.5, (Z1+Z2)*0.5);
dir1 = gp_Dir(vd);
typeres = IntAna_Line;
nbint = 1;
-
+ //
+ //-------------------------------------------------------
+ // When the value of the angle between the planes is small
+ // the origin of intersection line is computed with error
+ // [ ~0.0001 ] that can not br considered as small one
+ // e.g.
+ // for {A~=2.e-6, dist1=4.2e-5, dist2==1.e-4} =>
+ // {denom=3.4e-12, par1=12550297.6, par2=32605552.9, etc}
+ // So,
+ // the origin should be refined if it is possible
+ //
+ Standard_Real aTreshAng, aTreshDist;
+ //
+ aTreshAng=2.e-6; // 1.e-4 deg
+ aTreshDist=1.e-12;
+ //
+ if (aMVD < aTreshAng) {
+ Standard_Real aDist1, aDist2;
+ //
+ aDist1=A1*pt1.X() + B1*pt1.Y() + C1*pt1.Z() + D1;
+ aDist2=A2*pt1.X() + B2*pt1.Y() + C2*pt1.Z() + D2;
+ //
+ if (fabs(aDist1)>aTreshDist || fabs(aDist2)>aTreshDist) {
+ Standard_Boolean bIsDone, bIsParallel;
+ IntAna_IntConicQuad aICQ;
+ //
+ // 1.
+ gp_Dir aDN1(aVN1);
+ gp_Lin aL1(pt1, aDN1);
+ //
+ aICQ.Perform(aL1, P1, TolAng, Tol);
+ bIsDone=aICQ.IsDone();
+ if (!bIsDone) {
+ return;
+ }
+ //
+ const gp_Pnt& aPnt1=aICQ.Point(1);
+ //----------------------------------
+ // 2.
+ gp_Dir aDL2(dir1.Crossed(aDN1));
+ gp_Lin aL2(aPnt1, aDL2);
+ //
+ aICQ.Perform(aL2, P2, TolAng, Tol);
+ bIsDone=aICQ.IsDone();
+ if (!bIsDone) {
+ return;
+ }
+ //
+ bIsParallel=aICQ.IsParallel();
+ if (bIsParallel) {
+ return;
+ }
+ //
+ const gp_Pnt& aPnt2=aICQ.Point(1);
+ //
+ pt1=aPnt2;
+ }
+ }
}
done=Standard_True;
}
//function : IntAna_QuadQuadGeo
//purpose : Pln Cylinder
//=======================================================================
- IntAna_QuadQuadGeo::IntAna_QuadQuadGeo( const gp_Pln& P
- ,const gp_Cylinder& Cl
- ,const Standard_Real Tolang
- ,const Standard_Real Tol
- ,const Standard_Real H)
- : done(Standard_False),
- nbint(0),
- typeres(IntAna_Empty),
- pt1(0,0,0),
- pt2(0,0,0),
- param1(0),
- param2(0),
- param1bis(0),
- param2bis(0),
- myCommonGen(Standard_False),
- myPChar(0,0,0)
+IntAna_QuadQuadGeo::IntAna_QuadQuadGeo( const gp_Pln& P
+ ,const gp_Cylinder& Cl
+ ,const Standard_Real Tolang
+ ,const Standard_Real Tol
+ ,const Standard_Real H)
+ : done(Standard_False),
+ nbint(0),
+ typeres(IntAna_Empty),
+ pt1(0,0,0),
+ pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
+ param1(0),
+ param2(0),
+ param3(0),
+ param4(0),
+ param1bis(0),
+ param2bis(0),
+ myCommonGen(Standard_False),
+ myPChar(0,0,0)
{
InitTolerances();
Perform(P,Cl,Tolang,Tol,H);
P.Coefficients(A,B,C,D);
axec.Location().Coord(X,Y,Z);
- dist = A*X + B*Y + C*Z + D; // la distance axe/plan est evaluee a l origine de l axe.
+ // la distance axe/plan est evaluee a l origine de l axe.
+ dist = A*X + B*Y + C*Z + D;
Standard_Real tolang = Tolang;
Standard_Boolean newparams = Standard_False;
Standard_Real dang = Abs( ldv.Angle( npv ) ) - M_PI/2.;
Standard_Real dangle = Abs( dang );
if( dangle > Tolang )
- {
- Standard_Real sinda = Abs( Sin( dangle ) );
- Standard_Real dif = Abs( sinda - Tol );
- if( dif < Tol )
- {
- tolang = sinda * 2.;
- newparams = Standard_True;
- }
- }
+ {
+ Standard_Real sinda = Abs( Sin( dangle ) );
+ Standard_Real dif = Abs( sinda - Tol );
+ if( dif < Tol )
+ {
+ tolang = sinda * 2.;
+ newparams = Standard_True;
+ }
+ }
}
nbint = 0;
if (Abs(Abs(dist)-radius) < Tol)
{
- nbint = 1;
- pt1.SetXYZ(omega);
-
- if( newparams )
- {
- gp_XYZ omegaXYZ(X,Y,Z);
- gp_XYZ omegaXYZtrnsl( omegaXYZ + 100.*axec.Direction().XYZ() );
- Standard_Real Xt,Yt,Zt,distt;
- omegaXYZtrnsl.Coord(Xt,Yt,Zt);
- distt = A*Xt + B*Yt + C*Zt + D;
- gp_XYZ omega1( omegaXYZtrnsl.X()-distt*A, omegaXYZtrnsl.Y()-distt*B, omegaXYZtrnsl.Z()-distt*C );
- gp_Pnt ppt1;
- ppt1.SetXYZ( omega1 );
- gp_Vec vv1(pt1,ppt1);
- gp_Dir dd1( vv1 );
- dir1 = dd1;
- }
- else
- dir1 = axec.Direction();
+ nbint = 1;
+ pt1.SetXYZ(omega);
+
+ if( newparams )
+ {
+ gp_XYZ omegaXYZ(X,Y,Z);
+ gp_XYZ omegaXYZtrnsl( omegaXYZ + 100.*axec.Direction().XYZ() );
+ Standard_Real Xt,Yt,Zt,distt;
+ omegaXYZtrnsl.Coord(Xt,Yt,Zt);
+ distt = A*Xt + B*Yt + C*Zt + D;
+ gp_XYZ omega1(omegaXYZtrnsl.X()-distt*A,
+ omegaXYZtrnsl.Y()-distt*B,
+ omegaXYZtrnsl.Z()-distt*C );
+ gp_Pnt ppt1;
+ ppt1.SetXYZ( omega1 );
+ gp_Vec vv1(pt1,ppt1);
+ gp_Dir dd1( vv1 );
+ dir1 = dd1;
+ }
+ else
+ dir1 = axec.Direction();
}
else if (Abs(dist) < radius)
{
- nbint = 2;
- h = Sqrt(radius*radius - dist*dist);
- axey = axec.Direction().XYZ().Crossed(normp); // axey est normalise
-
- pt1.SetXYZ(omega - h*axey);
- pt2.SetXYZ(omega + h*axey);
-
- if( newparams )
- {
- gp_XYZ omegaXYZ(X,Y,Z);
- gp_XYZ omegaXYZtrnsl( omegaXYZ + 100.*axec.Direction().XYZ() );
- Standard_Real Xt,Yt,Zt,distt,ht;
- omegaXYZtrnsl.Coord(Xt,Yt,Zt);
- distt = A*Xt + B*Yt + C*Zt + D;
- // ht = Sqrt(radius*radius - distt*distt);
- Standard_Real anSqrtArg = radius*radius - distt*distt;
- ht = (anSqrtArg > 0.) ? Sqrt(anSqrtArg) : 0.;
-
- gp_XYZ omega1( omegaXYZtrnsl.X()-distt*A, omegaXYZtrnsl.Y()-distt*B, omegaXYZtrnsl.Z()-distt*C );
- gp_Pnt ppt1,ppt2;
- ppt1.SetXYZ( omega1 - ht*axey);
- ppt2.SetXYZ( omega1 + ht*axey);
- gp_Vec vv1(pt1,ppt1);
- gp_Vec vv2(pt2,ppt2);
- gp_Dir dd1( vv1 );
- gp_Dir dd2( vv2 );
- dir1 = dd1;
- dir2 = dd2;
- }
- else
- {
- dir1 = axec.Direction();
- dir2 = axec.Direction();
- }
+ nbint = 2;
+ h = Sqrt(radius*radius - dist*dist);
+ axey = axec.Direction().XYZ().Crossed(normp); // axey est normalise
+
+ pt1.SetXYZ(omega - h*axey);
+ pt2.SetXYZ(omega + h*axey);
+
+ if( newparams )
+ {
+ gp_XYZ omegaXYZ(X,Y,Z);
+ gp_XYZ omegaXYZtrnsl( omegaXYZ + 100.*axec.Direction().XYZ() );
+ Standard_Real Xt,Yt,Zt,distt,ht;
+ omegaXYZtrnsl.Coord(Xt,Yt,Zt);
+ distt = A*Xt + B*Yt + C*Zt + D;
+ // ht = Sqrt(radius*radius - distt*distt);
+ Standard_Real anSqrtArg = radius*radius - distt*distt;
+ ht = (anSqrtArg > 0.) ? Sqrt(anSqrtArg) : 0.;
+
+ gp_XYZ omega1( omegaXYZtrnsl.X()-distt*A,
+ omegaXYZtrnsl.Y()-distt*B,
+ omegaXYZtrnsl.Z()-distt*C );
+ gp_Pnt ppt1,ppt2;
+ ppt1.SetXYZ( omega1 - ht*axey);
+ ppt2.SetXYZ( omega1 + ht*axey);
+ gp_Vec vv1(pt1,ppt1);
+ gp_Vec vv2(pt2,ppt2);
+ gp_Dir dd1( vv1 );
+ gp_Dir dd2( vv2 );
+ dir1 = dd1;
+ dir2 = dd2;
+ }
+ else
+ {
+ dir1 = axec.Direction();
+ dir2 = axec.Direction();
+ }
}
// else nbint = 0
param1bis = radius;
}
}
- if(typeres == IntAna_Ellipse) {
- if( param1>100000.0*param1bis
- || param1bis>100000.0*param1) {
- done = Standard_False;
- return;
- }
- }
+
done = Standard_True;
}
//=======================================================================
//purpose : Pln Cone
//=======================================================================
IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Pln& P,
- const gp_Cone& Co,
- const Standard_Real Tolang,
- const Standard_Real Tol)
+ const gp_Cone& Co,
+ const Standard_Real Tolang,
+ const Standard_Real Tol)
: done(Standard_False),
nbint(0),
typeres(IntAna_Empty),
pt1(0,0,0),
pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
param1(0),
param2(0),
+ param3(0),
+ param4(0),
param1bis(0),
param2bis(0),
myCommonGen(Standard_False),
//purpose :
//=======================================================================
void IntAna_QuadQuadGeo::Perform(const gp_Pln& P,
- const gp_Cone& Co,
- const Standard_Real Tolang,
- const Standard_Real Tol)
+ const gp_Cone& Co,
+ const Standard_Real Tolang,
+ const Standard_Real Tol)
{
done = Standard_False;
distance = apex.Distance(center);
if (inter.ParamOnConic(1) + Co.RefRadius()/Tan(angl) < 0.) {
- axex.Reverse();
- axey.Reverse();
+ axex.Reverse();
+ axey.Reverse();
}
if (Abs(costa) < Tolang) { // plan parallele a une generatrice
- typeres = IntAna_Parabola;
- nbint = 1;
- deltacenter = distance/2./cosa;
- axex.Normalize();
- pt1.SetXYZ(center.XYZ()-deltacenter*axex);
- dir1 = normp;
- dir2.SetXYZ(axex);
- param1 = deltacenter*sina*sina;
+ typeres = IntAna_Parabola;
+ nbint = 1;
+ deltacenter = distance/2./cosa;
+ axex.Normalize();
+ pt1.SetXYZ(center.XYZ()-deltacenter*axex);
+ dir1 = normp;
+ dir2.SetXYZ(axex);
+ param1 = deltacenter*sina*sina;
}
else if (sint < Tolang) { // plan perpendiculaire a l axe
- typeres = IntAna_Circle;
- nbint = 1;
- pt1 = center;
- dir1 = Co.Position().Direction();
- dir2 = Co.Position().XDirection();
- param1 = apex.Distance(center)*Abs(Tan(angl));
+ typeres = IntAna_Circle;
+ nbint = 1;
+ pt1 = center;
+ dir1 = Co.Position().Direction();
+ dir2 = Co.Position().XDirection();
+ param1 = apex.Distance(center)*Abs(Tan(angl));
}
else if (cost < sina ) {
- typeres = IntAna_Hyperbola;
- nbint = 2;
- axex.Normalize();
-
- deltacenter = sint*sina*sina*distance/(sina*sina - cost*cost);
- pt1.SetXYZ(center.XYZ() - deltacenter*axex);
- pt2 = pt1;
- dir1 = normp;
- dir2.SetXYZ(axex);
- param1 = param2 = cost*sina*cosa*distance /(sina*sina-cost*cost);
- param1bis = param2bis = cost*sina*distance / Sqrt(sina*sina-cost*cost);
+ typeres = IntAna_Hyperbola;
+ nbint = 2;
+ axex.Normalize();
+
+ deltacenter = sint*sina*sina*distance/(sina*sina - cost*cost);
+ pt1.SetXYZ(center.XYZ() - deltacenter*axex);
+ pt2 = pt1;
+ dir1 = normp;
+ dir2.SetXYZ(axex);
+ param1 = param2 = cost*sina*cosa*distance /(sina*sina-cost*cost);
+ param1bis = param2bis = cost*sina*distance / Sqrt(sina*sina-cost*cost);
}
else { // on a alors cost > sina
- typeres = IntAna_Ellipse;
- nbint = 1;
- Standard_Real radius = cost*sina*cosa*distance/(cost*cost-sina*sina);
- deltacenter = sint*sina*sina*distance/(cost*cost-sina*sina);
- axex.Normalize();
- pt1.SetXYZ(center.XYZ() + deltacenter*axex);
- dir1 = normp;
- dir2.SetXYZ(axex);
- param1 = radius;
- param1bis = cost*sina*distance/ Sqrt(cost*cost - sina*sina);
+ typeres = IntAna_Ellipse;
+ nbint = 1;
+ Standard_Real radius = cost*sina*cosa*distance/(cost*cost-sina*sina);
+ deltacenter = sint*sina*sina*distance/(cost*cost-sina*sina);
+ axex.Normalize();
+ pt1.SetXYZ(center.XYZ() + deltacenter*axex);
+ dir1 = normp;
+ dir2.SetXYZ(axex);
+ param1 = radius;
+ param1bis = cost*sina*distance/ Sqrt(cost*cost - sina*sina);
}
}
}
//function : IntAna_QuadQuadGeo
//purpose : Pln Sphere
//=======================================================================
- IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Pln& P,
- const gp_Sphere& S)
+IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Pln& P,
+ const gp_Sphere& S)
: done(Standard_False),
nbint(0),
typeres(IntAna_Empty),
pt1(0,0,0),
pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
param1(0),
param2(0),
+ param3(0),
+ param4(0),
param1bis(0),
param2bis(0),
myCommonGen(Standard_False),
//function : Perform
//purpose :
//=======================================================================
- void IntAna_QuadQuadGeo::Perform( const gp_Pln& P
- ,const gp_Sphere& S)
+void IntAna_QuadQuadGeo::Perform( const gp_Pln& P
+ ,const gp_Sphere& S)
{
done = Standard_False;
//function : IntAna_QuadQuadGeo
//purpose : Cylinder - Cylinder
//=======================================================================
- IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Cylinder& Cyl1,
- const gp_Cylinder& Cyl2,
- const Standard_Real Tol)
+IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Cylinder& Cyl1,
+ const gp_Cylinder& Cyl2,
+ const Standard_Real Tol)
: done(Standard_False),
nbint(0),
typeres(IntAna_Empty),
pt1(0,0,0),
pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
param1(0),
param2(0),
+ param3(0),
+ param4(0),
param1bis(0),
param2bis(0),
myCommonGen(Standard_False),
//function : Perform
//purpose :
//=======================================================================
- void IntAna_QuadQuadGeo::Perform(const gp_Cylinder& Cyl1,
- const gp_Cylinder& Cyl2,
- const Standard_Real Tol)
+void IntAna_QuadQuadGeo::Perform(const gp_Cylinder& Cyl1,
+ const gp_Cylinder& Cyl2,
+ const Standard_Real Tol)
{
done=Standard_True;
//---------------------------- Parallel axes -------------------------
Standard_Real RmR, RmR_Relative;
RmR=(R1>R2)? (R1-R2) : (R2-R1);
{
- Standard_Real Rmax, Rmin;
+ Standard_Real Rmax;
Rmax=(R1>R2)? R1 : R2;
- Rmin=(R1>R2)? R2 : R1;
RmR_Relative=RmR/Rmax;
}
Standard_Real DistA1A2=A1A2.Distance();
- if(A1A2.Parallel()) {
- if(DistA1A2<=Tol) {
- if(RmR<=Tol) {
- typeres=IntAna_Same;
+ if(A1A2.Parallel())
+ {
+ if(DistA1A2<=Tol)
+ {
+ if(RmR<=Tol)
+ {
+ typeres=IntAna_Same;
}
- else {
- typeres=IntAna_Empty;
+ else
+ {
+ typeres=IntAna_Empty;
}
}
- else { //-- DistA1A2 > Tol
+ else
+ { //-- DistA1A2 > Tol
gp_Pnt P1=Cyl1.Location();
gp_Pnt P2t=Cyl2.Location();
gp_Pnt P2;
gp_Dir DirCyl = Cyl1.Position().Direction();
Standard_Real ProjP2OnDirCyl1=gp_Vec(DirCyl).Dot(gp_Vec(P1,P2t));
- P2.SetCoord( P2t.X() - ProjP2OnDirCyl1*DirCyl.X()
- ,P2t.Y() - ProjP2OnDirCyl1*DirCyl.Y()
- ,P2t.Z() - ProjP2OnDirCyl1*DirCyl.Z());
+ //P2 is a projection the location of the 2nd cylinder on the base
+ //of the 1st cylinder
+ P2.SetCoord(P2t.X() - ProjP2OnDirCyl1*DirCyl.X(),
+ P2t.Y() - ProjP2OnDirCyl1*DirCyl.Y(),
+ P2t.Z() - ProjP2OnDirCyl1*DirCyl.Z());
//--
Standard_Real R1pR2=R1+R2;
- if(DistA1A2>(R1pR2+Tol)) {
- typeres=IntAna_Empty;
- nbint=0;
+ if(DistA1A2>(R1pR2+Tol))
+ {
+ typeres=IntAna_Empty;
+ nbint=0;
}
- else if(DistA1A2>(R1pR2)) {
- //-- 1 Tangent line -------------------------------------OK
- typeres=IntAna_Line;
-
- nbint=1;
- dir1=DirCyl;
- Standard_Real R1_R1pR2=R1/R1pR2;
- pt1.SetCoord( P1.X() + R1_R1pR2 * (P2.X()-P1.X())
- ,P1.Y() + R1_R1pR2 * (P2.Y()-P1.Y())
- ,P1.Z() + R1_R1pR2 * (P2.Z()-P1.Z()));
-
+ else if((R1pR2 - DistA1A2) <= RealSmall())
+ {
+ //-- 1 Tangent line -------------------------------------OK
+ typeres=IntAna_Line;
+
+ nbint=1;
+ dir1=DirCyl;
+ Standard_Real R1_R1pR2=R1/R1pR2;
+ pt1.SetCoord( P1.X() + R1_R1pR2 * (P2.X()-P1.X()),
+ P1.Y() + R1_R1pR2 * (P2.Y()-P1.Y()),
+ P1.Z() + R1_R1pR2 * (P2.Z()-P1.Z()));
}
- else if(DistA1A2>RmR) {
- //-- 2 lines ---------------------------------------------OK
- typeres=IntAna_Line;
- nbint=2;
- dir1=DirCyl;
- gp_Vec P1P2(P1,P2);
- gp_Dir DirA1A2=gp_Dir(P1P2);
- gp_Dir Ortho_dir1_P1P2 = dir1.Crossed(DirA1A2);
- dir2=dir1;
- Standard_Real Alpha=0.5*(R1*R1-R2*R2+DistA1A2*DistA1A2)/(DistA1A2);
-
-// Standard_Real Beta = Sqrt(R1*R1-Alpha*Alpha);
- Standard_Real anSqrtArg = R1*R1-Alpha*Alpha;
- Standard_Real Beta = (anSqrtArg > 0.) ? Sqrt(anSqrtArg) : 0.;
-
- if((Beta+Beta)<Tol) {
- nbint=1;
- pt1.SetCoord( P1.X() + Alpha*DirA1A2.X()
- ,P1.Y() + Alpha*DirA1A2.Y()
- ,P1.Z() + Alpha*DirA1A2.Z());
- }
- else {
- pt1.SetCoord( P1.X() + Alpha*DirA1A2.X() + Beta*Ortho_dir1_P1P2.X()
- ,P1.Y() + Alpha*DirA1A2.Y() + Beta*Ortho_dir1_P1P2.Y()
- ,P1.Z() + Alpha*DirA1A2.Z() + Beta*Ortho_dir1_P1P2.Z() );
-
- pt2.SetCoord( P1.X() + Alpha*DirA1A2.X() - Beta*Ortho_dir1_P1P2.X()
- ,P1.Y() + Alpha*DirA1A2.Y() - Beta*Ortho_dir1_P1P2.Y()
- ,P1.Z() + Alpha*DirA1A2.Z() - Beta*Ortho_dir1_P1P2.Z());
- }
+ else if(DistA1A2>RmR)
+ {
+ //-- 2 lines ---------------------------------------------OK
+ typeres=IntAna_Line;
+ nbint=2;
+ dir1=DirCyl;
+ dir2=dir1;
+
+ const Standard_Real aR1R1 = R1*R1;
+
+ /*
+ P1
+ o
+ * | *
+ * O1| *
+ A o-----o-----o B
+ * | *
+ * | *
+ o
+ P2
+
+ Two cylinders have axes collinear. Therefore, problem can be reformulated as
+ to find intersection point of two circles (the bases of the cylinders) on
+ the plane: 1st circle has center P1 and radius R1 (the radius of the
+ 1st cylinder) and 2nd circle has center P2 and radius R2 (the radius of the
+ 2nd cylinder). The plane is the base of the 1st cylinder. Points A and B
+ are intersection point of these circles. Distance P1P2 is equal to DistA1A2.
+ O1 is the intersection point of P1P2 and AB segments.
+
+ At that, if distance AB < Tol we consider that the circles are tangent and
+ has only one intersection point.
+
+ AB = 2*R1*sin(angle AP1P2).
+ Accordingly,
+ AB^2 < Tol^2 => 4*R1*R1*sin(angle AP1P2)^2 < Tol*Tol.
+ */
+
+
+ //Cosine and Square of Sine of the A-P1-P2 angle
+ const Standard_Real aCos = 0.5*(aR1R1-R2*R2+DistA1A2*DistA1A2)/(R1*DistA1A2);
+ const Standard_Real aSin2 = 1-aCos*aCos;
+
+ const Standard_Boolean isTangent =((4.0*aR1R1*aSin2) < Tol*Tol);
+
+ //Normalized vector P1P2
+ const gp_Vec DirA1A2((P2.XYZ() - P1.XYZ())/DistA1A2);
+
+ if(isTangent)
+ {
+ //Intercept the segment from P1 point along P1P2 direction
+ //and having |P1O1| length
+ nbint=1;
+ pt1.SetXYZ(P1.XYZ() + DirA1A2.XYZ()*R1*aCos);
+ }
+ else
+ {
+ //Sine of the A-P1-P2 angle (if aSin2 < 0 then isTangent == TRUE =>
+ //go to another branch)
+ const Standard_Real aSin = sqrt(aSin2);
+
+ //1. Rotate P1P2 to the angle A-P1-P2 relative to P1
+ //(clockwise and anticlockwise for getting
+ //two intersection points).
+ //2. Intercept the segment from P1 along direction,
+ //determined in the preview paragraph and having R1 length
+ const gp_Dir &aXDir = Cyl1.Position().XDirection(),
+ &aYDir = Cyl1.Position().YDirection();
+ const gp_Vec aR1Xdir = R1*aXDir.XYZ(),
+ aR1Ydir = R1*aYDir.XYZ();
+
+ //Source 2D-coordinates of the P1P2 vector normalized
+ //in coordinate system, based on the X- and Y-directions
+ //of the 1st cylinder in the plane of the 1st cylinder base
+ //(P1 is the origin of the coordinate system).
+ const Standard_Real aDx = DirA1A2.Dot(aXDir),
+ aDy = DirA1A2.Dot(aYDir);
+
+ //New coordinate (after rotation) of the P1P2 vector normalized.
+ Standard_Real aNewDx = aDx*aCos - aDy*aSin,
+ aNewDy = aDy*aCos + aDx*aSin;
+ pt1.SetXYZ(P1.XYZ() + aNewDx*aR1Xdir.XYZ() + aNewDy*aR1Ydir.XYZ());
+
+ aNewDx = aDx*aCos + aDy*aSin;
+ aNewDy = aDy*aCos - aDx*aSin;
+ pt2.SetXYZ(P1.XYZ() + aNewDx*aR1Xdir.XYZ() + aNewDy*aR1Ydir.XYZ());
+ }
}
- else if(DistA1A2>(RmR-Tol)) {
- //-- 1 Tangent ------------------------------------------OK
- typeres=IntAna_Line;
- nbint=1;
- dir1=DirCyl;
- Standard_Real R1_RmR=R1/RmR;
-
- if(R1 < R2) R1_RmR = -R1_RmR;
-
- pt1.SetCoord( P1.X() + R1_RmR * (P2.X()-P1.X())
- ,P1.Y() + R1_RmR * (P2.Y()-P1.Y())
- ,P1.Z() + R1_RmR * (P2.Z()-P1.Z()));
+ else if(DistA1A2>(RmR-Tol))
+ {
+ //-- 1 Tangent ------------------------------------------OK
+ typeres=IntAna_Line;
+ nbint=1;
+ dir1=DirCyl;
+ Standard_Real R1_RmR=R1/RmR;
+
+ if(R1 < R2)
+ R1_RmR = -R1_RmR;
+
+ pt1.SetCoord( P1.X() + R1_RmR * (P2.X()-P1.X()),
+ P1.Y() + R1_RmR * (P2.Y()-P1.Y()),
+ P1.Z() + R1_RmR * (P2.Z()-P1.Z()));
}
else {
- nbint=0;
- typeres=IntAna_Empty;
+ nbint=0;
+ typeres=IntAna_Empty;
}
}
}
else { //-- No Parallel Axis ---------------------------------OK
if((RmR_Relative<=myEPSILON_CYLINDER_DELTA_RADIUS)
- && (DistA1A2 <= myEPSILON_CYLINDER_DELTA_DISTANCE)) {
+ && (DistA1A2 <= myEPSILON_CYLINDER_DELTA_DISTANCE))
+ {
//-- PI/2 between the two axis and Intersection
//-- and identical radius
typeres=IntAna_Ellipse;
B=Abs(Sin(0.5*(M_PI-A)));
A=Abs(Sin(0.5*A));
- if(A==0.0 || B==0.0) {
- typeres=IntAna_Same;
- return;
+ if(A==0.0 || B==0.0)
+ {
+ typeres=IntAna_Same;
+ return;
}
-
gp_Vec dircyl1(DirCyl1);gp_Vec dircyl2(DirCyl2);
dir1 = gp_Dir(dircyl1.Added(dircyl2));
dir2 = gp_Dir(dircyl1.Subtracted(dircyl2));
-
+
param2 = Cyl1.Radius() / A;
param1 = Cyl1.Radius() / B;
param2bis= param1bis = Cyl1.Radius();
- if(param1 < param1bis) {
- A=param1; param1=param1bis; param1bis=A;
+ if(param1 < param1bis)
+ {
+ A=param1;
+ param1=param1bis;
+ param1bis=A;
}
- if(param2 < param2bis) {
- A=param2; param2=param2bis; param2bis=A;
+
+ if(param2 < param2bis)
+ {
+ A=param2;
+ param2=param2bis;
+ param2bis=A;
}
}
- else {
- if(Abs(DistA1A2-Cyl1.Radius()-Cyl2.Radius())<Tol) {
- typeres = IntAna_Point;
- Standard_Real d,p1,p2;
-
- gp_Dir D1 = Cyl1.Axis().Direction();
- gp_Dir D2 = Cyl2.Axis().Direction();
- A1A2.Distance(d,p1,p2);
- gp_Pnt P = Cyl1.Axis().Location();
- gp_Pnt P1(P.X() - p1*D1.X(),
- P.Y() - p1*D1.Y(),
- P.Z() - p1*D1.Z());
- P = Cyl2.Axis().Location();
- gp_Pnt P2(P.X() - p2*D2.X(),
- P.Y() - p2*D2.Y(),
- P.Z() - p2*D2.Z());
- gp_Vec P1P2(P1,P2);
- D1=gp_Dir(P1P2);
- p1=Cyl1.Radius();
- pt1.SetCoord(P1.X() + p1*D1.X(),
- P1.Y() + p1*D1.Y(),
- P1.Z() + p1*D1.Z());
- nbint = 1;
+ else
+ {
+ if(Abs(DistA1A2-Cyl1.Radius()-Cyl2.Radius())<Tol)
+ {
+ typeres = IntAna_Point;
+ Standard_Real d,p1,p2;
+
+ gp_Dir D1 = Cyl1.Axis().Direction();
+ gp_Dir D2 = Cyl2.Axis().Direction();
+ A1A2.Distance(d,p1,p2);
+ gp_Pnt P = Cyl1.Axis().Location();
+ gp_Pnt P1(P.X() - p1*D1.X(),
+ P.Y() - p1*D1.Y(),
+ P.Z() - p1*D1.Z());
+
+ P = Cyl2.Axis().Location();
+
+ gp_Pnt P2(P.X() - p2*D2.X(),
+ P.Y() - p2*D2.Y(),
+ P.Z() - p2*D2.Z());
+
+ gp_Vec P1P2(P1,P2);
+ D1=gp_Dir(P1P2);
+ p1=Cyl1.Radius();
+
+ pt1.SetCoord(P1.X() + p1*D1.X(),
+ P1.Y() + p1*D1.Y(),
+ P1.Z() + p1*D1.Z());
+ nbint = 1;
}
- else {
- typeres=IntAna_NoGeometricSolution;
+ else
+ {
+ typeres=IntAna_NoGeometricSolution;
}
}
}
//function : IntAna_QuadQuadGeo
//purpose : Cylinder - Cone
//=======================================================================
- IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Cylinder& Cyl,
- const gp_Cone& Con,
- const Standard_Real Tol)
+IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Cylinder& Cyl,
+ const gp_Cone& Con,
+ const Standard_Real Tol)
: done(Standard_False),
nbint(0),
typeres(IntAna_Empty),
pt1(0,0,0),
pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
param1(0),
param2(0),
+ param3(0),
+ param4(0),
param1bis(0),
param2bis(0),
myCommonGen(Standard_False),
//purpose :
//=======================================================================
void IntAna_QuadQuadGeo::Perform(const gp_Cylinder& Cyl,
- const gp_Cone& Con,
- const Standard_Real )
+ const gp_Cone& Con,
+ const Standard_Real )
{
done=Standard_True;
AxeOperator A1A2(Cyl.Axis(),Con.Axis());
Standard_Real dist=Cyl.Radius()/(Tan(Con.SemiAngle()));
gp_Dir dir=Cyl.Position().Direction();
pt1.SetCoord( Pt.X() + dist*dir.X()
- ,Pt.Y() + dist*dir.Y()
- ,Pt.Z() + dist*dir.Z());
+ ,Pt.Y() + dist*dir.Y()
+ ,Pt.Z() + dist*dir.Z());
pt2.SetCoord( Pt.X() - dist*dir.X()
- ,Pt.Y() - dist*dir.Y()
- ,Pt.Z() - dist*dir.Z());
+ ,Pt.Y() - dist*dir.Y()
+ ,Pt.Z() - dist*dir.Z());
dir1=dir2=dir;
param1=param2=Cyl.Radius();
nbint=2;
//purpose : Cylinder - Sphere
//=======================================================================
IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Cylinder& Cyl,
- const gp_Sphere& Sph,
- const Standard_Real Tol)
+ const gp_Sphere& Sph,
+ const Standard_Real Tol)
: done(Standard_False),
nbint(0),
typeres(IntAna_Empty),
pt1(0,0,0),
pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
param1(0),
param2(0),
+ param3(0),
+ param4(0),
param1bis(0),
param2bis(0),
myCommonGen(Standard_False),
//purpose :
//=======================================================================
void IntAna_QuadQuadGeo::Perform( const gp_Cylinder& Cyl
- ,const gp_Sphere& Sph
- ,const Standard_Real)
+ ,const gp_Sphere& Sph
+ ,const Standard_Real)
{
done=Standard_True;
gp_Pnt Pt=Sph.Location();
dir1 = dir2 = dir;
typeres=IntAna_Circle;
pt1.SetCoord( Pt.X() + dist*dir.X()
- ,Pt.Y() + dist*dir.Y()
- ,Pt.Z() + dist*dir.Z());
+ ,Pt.Y() + dist*dir.Y()
+ ,Pt.Z() + dist*dir.Z());
nbint=1;
param1 = Cyl.Radius();
if(dist>RealEpsilon()) {
- pt2.SetCoord( Pt.X() - dist*dir.X()
- ,Pt.Y() - dist*dir.Y()
- ,Pt.Z() - dist*dir.Z());
- param2=Cyl.Radius();
- nbint=2;
+ pt2.SetCoord( Pt.X() - dist*dir.X()
+ ,Pt.Y() - dist*dir.Y()
+ ,Pt.Z() - dist*dir.Z());
+ param2=Cyl.Radius();
+ nbint=2;
}
}
}
//purpose : Cone - Cone
//=======================================================================
IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Cone& Con1,
- const gp_Cone& Con2,
- const Standard_Real Tol)
+ const gp_Cone& Con2,
+ const Standard_Real Tol)
: done(Standard_False),
nbint(0),
typeres(IntAna_Empty),
pt1(0,0,0),
pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
param1(0),
param2(0),
+ param3(0),
+ param4(0),
param1bis(0),
param2bis(0),
myCommonGen(Standard_False),
//purpose :
//=======================================================================
void IntAna_QuadQuadGeo::Perform(const gp_Cone& Con1,
- const gp_Cone& Con2,
- const Standard_Real Tol)
+ const gp_Cone& Con2,
+ const Standard_Real Tol)
{
done=Standard_True;
//
if(Abs(tg1-tg2)>myEPSILON_ANGLE_CONE) {
if (fabs(d) < TOL_APEX_CONF) {
- typeres = IntAna_Point;
- nbint = 1;
- pt1 = P;
- return;
+ typeres = IntAna_Point;
+ nbint = 1;
+ pt1 = P;
+ return;
}
x=(d*tg2)/(tg1+tg2);
pt1.SetCoord( P.X() + x*D.X()
- ,P.Y() + x*D.Y()
- ,P.Z() + x*D.Z());
+ ,P.Y() + x*D.Y()
+ ,P.Z() + x*D.Z());
param1=Abs(x*tg1);
x=(d*tg2)/(tg2-tg1);
pt2.SetCoord( P.X() + x*D.X()
- ,P.Y() + x*D.Y()
- ,P.Z() + x*D.Z());
+ ,P.Y() + x*D.Y()
+ ,P.Z() + x*D.Z());
param2=Abs(x*tg1);
dir1 = dir2 = D;
nbint=2;
}
else {
if (fabs(d) < TOL_APEX_CONF) {
- typeres=IntAna_Same;
+ typeres=IntAna_Same;
}
else {
- typeres=IntAna_Circle;
- nbint=1;
- x=d*0.5;
- pt1.SetCoord( P.X() + x*D.X()
- ,P.Y() + x*D.Y()
- ,P.Z() + x*D.Z());
- param1 = Abs(x * tg1);
- dir1 = D;
+ typeres=IntAna_Circle;
+ nbint=1;
+ x=d*0.5;
+ pt1.SetCoord( P.X() + x*D.X()
+ ,P.Y() + x*D.Y()
+ ,P.Z() + x*D.Z());
+ param1 = Abs(x * tg1);
+ dir1 = D;
}
}
} //-- fin A1A2.Same
Standard_Real DistA1A2=A1A2.Distance();
gp_Dir DA1=Con1.Position().Direction();
gp_Vec O1O2(Con1.Apex(),Con2.Apex());
- Standard_Real O1O2_DA1=gp_Vec(DA1).Dot(O1O2);
-
- gp_Vec O1_Proj_A2(O1O2.X()-O1O2_DA1*DA1.X(),
- O1O2.Y()-O1O2_DA1*DA1.Y(),
- O1O2.Z()-O1O2_DA1*DA1.Z());
+ gp_Dir O1O2n(O1O2); // normalization of the vector before projection
+ Standard_Real O1O2_DA1=gp_Vec(DA1).Dot(gp_Vec(O1O2n));
+
+ gp_Vec O1_Proj_A2(O1O2n.X()-O1O2_DA1*DA1.X(),
+ O1O2n.Y()-O1O2_DA1*DA1.Y(),
+ O1O2n.Z()-O1O2_DA1*DA1.Z());
gp_Dir DB1=gp_Dir(O1_Proj_A2);
-
+
Standard_Real yO1O2=O1O2.Dot(gp_Vec(DA1));
Standard_Real ABSTG1 = Abs(tg1);
Standard_Real X2 = (DistA1A2/ABSTG1 - yO1O2)*0.5;
Standard_Real X1 = X2+yO1O2;
gp_Pnt P1(Con1.Apex().X() + X1*( DA1.X() + ABSTG1*DB1.X()),
- Con1.Apex().Y() + X1*( DA1.Y() + ABSTG1*DB1.Y()),
- Con1.Apex().Z() + X1*( DA1.Z() + ABSTG1*DB1.Z()));
+ Con1.Apex().Y() + X1*( DA1.Y() + ABSTG1*DB1.Y()),
+ Con1.Apex().Z() + X1*( DA1.Z() + ABSTG1*DB1.Z()));
gp_Pnt MO1O2(0.5*(Con1.Apex().X()+Con2.Apex().X()),
- 0.5*(Con1.Apex().Y()+Con2.Apex().Y()),
- 0.5*(Con1.Apex().Z()+Con2.Apex().Z()));
+ 0.5*(Con1.Apex().Y()+Con2.Apex().Y()),
+ 0.5*(Con1.Apex().Z()+Con2.Apex().Z()));
gp_Vec P1MO1O2(P1,MO1O2);
gp_Dir DA1_X_DB1=DA1.Crossed(DB1);
IntAna_QuadQuadGeo INTER_QUAD_PLN(gp_Pln(P1,OrthoPln),Con1,Tol,Tol);
if(INTER_QUAD_PLN.IsDone()) {
switch(INTER_QUAD_PLN.TypeInter()) {
- case IntAna_Ellipse: {
- typeres=IntAna_Ellipse;
- gp_Elips E=INTER_QUAD_PLN.Ellipse(1);
- pt1 = E.Location();
- dir1 = E.Position().Direction();
- dir2 = E.Position().XDirection();
- param1 = E.MajorRadius();
- param1bis = E.MinorRadius();
- nbint = 1;
- break;
+ case IntAna_Ellipse: {
+ typeres=IntAna_Ellipse;
+ gp_Elips E=INTER_QUAD_PLN.Ellipse(1);
+ pt1 = E.Location();
+ dir1 = E.Position().Direction();
+ dir2 = E.Position().XDirection();
+ param1 = E.MajorRadius();
+ param1bis = E.MinorRadius();
+ nbint = 1;
+ break;
}
case IntAna_Circle: {
- typeres=IntAna_Circle;
- gp_Circ C=INTER_QUAD_PLN.Circle(1);
- pt1 = C.Location();
- dir1 = C.Position().XDirection();
- dir2 = C.Position().YDirection();
- param1 = C.Radius();
- nbint = 1;
- break;
+ typeres=IntAna_Circle;
+ gp_Circ C=INTER_QUAD_PLN.Circle(1);
+ pt1 = C.Location();
+ dir1 = C.Position().XDirection();
+ dir2 = C.Position().YDirection();
+ param1 = C.Radius();
+ nbint = 1;
+ break;
}
case IntAna_Hyperbola: {
- typeres=IntAna_Hyperbola;
- gp_Hypr H=INTER_QUAD_PLN.Hyperbola(1);
- pt1 = pt2 = H.Location();
- dir1 = H.Position().Direction();
- dir2 = H.Position().XDirection();
- param1 = param2 = H.MajorRadius();
- param1bis = param2bis = H.MinorRadius();
- nbint = 2;
- break;
+ typeres=IntAna_Hyperbola;
+ gp_Hypr H=INTER_QUAD_PLN.Hyperbola(1);
+ pt1 = pt2 = H.Location();
+ dir1 = H.Position().Direction();
+ dir2 = H.Position().XDirection();
+ param1 = param2 = H.MajorRadius();
+ param1bis = param2bis = H.MinorRadius();
+ nbint = 2;
+ break;
}
case IntAna_Line: {
- typeres=IntAna_Line;
- gp_Lin H=INTER_QUAD_PLN.Line(1);
- pt1 = pt2 = H.Location();
- dir1 = dir2 = H.Position().Direction();
- param1 = param2 = 0.0;
- param1bis = param2bis = 0.0;
- nbint = 2;
- break;
+ typeres=IntAna_Line;
+ gp_Lin H=INTER_QUAD_PLN.Line(1);
+ pt1 = pt2 = H.Location();
+ dir1 = dir2 = H.Position().Direction();
+ param1 = param2 = 0.0;
+ param1bis = param2bis = 0.0;
+ nbint = 2;
+ break;
}
default:
- typeres=IntAna_NoGeometricSolution;
+ typeres=IntAna_NoGeometricSolution;
}
}
}// else if((Abs(tg1-tg2)<EPSILON_ANGLE_CONE) && (A1A2.Parallel()))
aQApex1=Con1.Apex();
aD3Ax1=aAx1.Direction();
aQA1.SetCoord(aQApex1.X()+aD1*aD3Ax1.X(),
- aQApex1.Y()+aD1*aD3Ax1.Y(),
- aQApex1.Z()+aD1*aD3Ax1.Z());
+ aQApex1.Y()+aD1*aD3Ax1.Y(),
+ aQApex1.Z()+aD1*aD3Ax1.Z());
//
aDx=aD3Ax1.Dot(aAx2.Direction());
if (aDx<0.) {
aD2=aD1*sqrt((1.+aTgBeta1*aTgBeta1)/(1.+aTgBeta2*aTgBeta2));
//
aQA2.SetCoord(aQApex1.X()+aD2*aD3Ax2.X(),
- aQApex1.Y()+aD2*aD3Ax2.Y(),
- aQApex1.Z()+aD2*aD3Ax2.Z());
+ aQApex1.Y()+aD2*aD3Ax2.Y(),
+ aQApex1.Z()+aD2*aD3Ax2.Z());
//
gp_Pln aPln1(aQA1, aD3Ax1);
gp_Pln aPln2(aQA2, aD3Ax2);
//
aIntr.Perform(aPln1, aPln2, Tol, Tol);
- if (!aIntr.IsDone()) {
+ if (!aIntr.IsDone() || 0 == aIntr.NbSolutions()) {
iRet=-1; // just in case. it must not be so
typeres=IntAna_NoGeometricSolution;
return;
if(INTER_QUAD_PLN.IsDone()) {
switch(INTER_QUAD_PLN.TypeInter()) {
- case IntAna_Ellipse: {
- typeres=IntAna_Ellipse;
- gp_Elips E=INTER_QUAD_PLN.Ellipse(1);
- pt1 = E.Location();
- dir1 = E.Position().Direction();
- dir2 = E.Position().XDirection();
- param1 = E.MajorRadius();
- param1bis = E.MinorRadius();
- nbint = 1;
- break;
+ case IntAna_Ellipse: {
+ typeres=IntAna_Ellipse;
+ gp_Elips E=INTER_QUAD_PLN.Ellipse(1);
+ pt1 = E.Location();
+ dir1 = E.Position().Direction();
+ dir2 = E.Position().XDirection();
+ param1 = E.MajorRadius();
+ param1bis = E.MinorRadius();
+ nbint = 1;
+ break;
}
case IntAna_Circle: {
- typeres=IntAna_Circle;
- gp_Circ C=INTER_QUAD_PLN.Circle(1);
- pt1 = C.Location();
- dir1 = C.Position().XDirection();
- dir2 = C.Position().YDirection();
- param1 = C.Radius();
- nbint = 1;
- break;
+ typeres=IntAna_Circle;
+ gp_Circ C=INTER_QUAD_PLN.Circle(1);
+ pt1 = C.Location();
+ dir1 = C.Position().XDirection();
+ dir2 = C.Position().YDirection();
+ param1 = C.Radius();
+ nbint = 1;
+ break;
}
case IntAna_Parabola: {
- typeres=IntAna_Parabola;
- gp_Parab Prb=INTER_QUAD_PLN.Parabola(1);
- pt1 = Prb.Location();
- dir1 = Prb.Position().Direction();
- dir2 = Prb.Position().XDirection();
- param1 = Prb.Focal();
- nbint = 1;
- break;
+ typeres=IntAna_Parabola;
+ gp_Parab Prb=INTER_QUAD_PLN.Parabola(1);
+ pt1 = Prb.Location();
+ dir1 = Prb.Position().Direction();
+ dir2 = Prb.Position().XDirection();
+ param1 = Prb.Focal();
+ nbint = 1;
+ break;
}
case IntAna_Hyperbola: {
- typeres=IntAna_Hyperbola;
- gp_Hypr H=INTER_QUAD_PLN.Hyperbola(1);
- pt1 = pt2 = H.Location();
- dir1 = H.Position().Direction();
- dir2 = H.Position().XDirection();
- param1 = param2 = H.MajorRadius();
- param1bis = param2bis = H.MinorRadius();
- nbint = 2;
- break;
+ typeres=IntAna_Hyperbola;
+ gp_Hypr H=INTER_QUAD_PLN.Hyperbola(1);
+ pt1 = pt2 = H.Location();
+ dir1 = H.Position().Direction();
+ dir2 = H.Position().XDirection();
+ param1 = param2 = H.MajorRadius();
+ param1bis = param2bis = H.MinorRadius();
+ nbint = 2;
+ break;
}
default:
- typeres=IntAna_NoGeometricSolution;
+ typeres=IntAna_NoGeometricSolution;
}
}
}
//purpose : Sphere - Cone
//=======================================================================
IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Sphere& Sph,
- const gp_Cone& Con,
- const Standard_Real Tol)
+ const gp_Cone& Con,
+ const Standard_Real Tol)
: done(Standard_False),
nbint(0),
typeres(IntAna_Empty),
pt1(0,0,0),
pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
param1(0),
param2(0),
+ param3(0),
+ param4(0),
param1bis(0),
param2bis(0),
myCommonGen(Standard_False),
//purpose :
//=======================================================================
void IntAna_QuadQuadGeo::Perform(const gp_Sphere& Sph,
- const gp_Cone& Con,
- const Standard_Real)
+ const gp_Cone& Con,
+ const Standard_Real)
{
//
//-- tga = y / (x+dApexSphCenter)
Standard_Real tgatga = tga * tga;
math_DirectPolynomialRoots Eq( 1.0+tgatga
- ,2.0*tgatga*dApexSphCenter
- ,-Rad*Rad + dApexSphCenter*dApexSphCenter*tgatga);
+ ,2.0*tgatga*dApexSphCenter
+ ,-Rad*Rad + dApexSphCenter*dApexSphCenter*tgatga);
if(Eq.IsDone()) {
Standard_Integer nbsol=Eq.NbSolutions();
if(nbsol==0) {
- typeres=IntAna_Empty;
+ typeres=IntAna_Empty;
}
else {
- typeres=IntAna_Circle;
- if(nbsol>=1) {
- Standard_Real x = Eq.Value(1);
- Standard_Real dApexSphCenterpx = dApexSphCenter+x;
- nbint=1;
- pt1.SetCoord( ConApex.X() + (dApexSphCenterpx) * ConDir.X()
- ,ConApex.Y() + (dApexSphCenterpx) * ConDir.Y()
- ,ConApex.Z() + (dApexSphCenterpx) * ConDir.Z());
- param1 = tga * dApexSphCenterpx;
- param1 = Abs(param1);
- dir1 = ConDir;
- if(param1<=myEPSILON_MINI_CIRCLE_RADIUS) {
- typeres=IntAna_PointAndCircle;
- param1=0.0;
- }
- }
- if(nbsol>=2) {
- Standard_Real x=Eq.Value(2);
- Standard_Real dApexSphCenterpx = dApexSphCenter+x;
- nbint=2;
- pt2.SetCoord( ConApex.X() + (dApexSphCenterpx) * ConDir.X()
- ,ConApex.Y() + (dApexSphCenterpx) * ConDir.Y()
- ,ConApex.Z() + (dApexSphCenterpx) * ConDir.Z());
- param2 = tga * dApexSphCenterpx;
- param2 = Abs(param2);
- dir2=ConDir;
- if(param2<=myEPSILON_MINI_CIRCLE_RADIUS) {
- typeres=IntAna_PointAndCircle;
- param2=0.0;
- }
- }
+ typeres=IntAna_Circle;
+ if(nbsol>=1) {
+ Standard_Real x = Eq.Value(1);
+ Standard_Real dApexSphCenterpx = dApexSphCenter+x;
+ nbint=1;
+ pt1.SetCoord( ConApex.X() + (dApexSphCenterpx) * ConDir.X()
+ ,ConApex.Y() + (dApexSphCenterpx) * ConDir.Y()
+ ,ConApex.Z() + (dApexSphCenterpx) * ConDir.Z());
+ param1 = tga * dApexSphCenterpx;
+ param1 = Abs(param1);
+ dir1 = ConDir;
+ if(param1<=myEPSILON_MINI_CIRCLE_RADIUS) {
+ typeres=IntAna_PointAndCircle;
+ param1=0.0;
+ }
+ }
+ if(nbsol>=2) {
+ Standard_Real x=Eq.Value(2);
+ Standard_Real dApexSphCenterpx = dApexSphCenter+x;
+ nbint=2;
+ pt2.SetCoord( ConApex.X() + (dApexSphCenterpx) * ConDir.X()
+ ,ConApex.Y() + (dApexSphCenterpx) * ConDir.Y()
+ ,ConApex.Z() + (dApexSphCenterpx) * ConDir.Z());
+ param2 = tga * dApexSphCenterpx;
+ param2 = Abs(param2);
+ dir2=ConDir;
+ if(param2<=myEPSILON_MINI_CIRCLE_RADIUS) {
+ typeres=IntAna_PointAndCircle;
+ param2=0.0;
+ }
+ }
}
}
else {
//purpose : Sphere - Sphere
//=======================================================================
IntAna_QuadQuadGeo::IntAna_QuadQuadGeo( const gp_Sphere& Sph1
- ,const gp_Sphere& Sph2
- ,const Standard_Real Tol)
+ ,const gp_Sphere& Sph2
+ ,const Standard_Real Tol)
: done(Standard_False),
nbint(0),
typeres(IntAna_Empty),
pt1(0,0,0),
pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
param1(0),
param2(0),
+ param3(0),
+ param4(0),
param1bis(0),
param2bis(0),
myCommonGen(Standard_False),
//purpose :
//=======================================================================
void IntAna_QuadQuadGeo::Perform(const gp_Sphere& Sph1,
- const gp_Sphere& Sph2,
- const Standard_Real Tol)
+ const gp_Sphere& Sph2,
+ const Standard_Real Tol)
{
done=Standard_True;
gp_Pnt O1=Sph1.Location();
Standard_Real t2;
if(R1==Rmax) t2=(R1 + (R2 + dO1O2)) * 0.5;
else t2=(-R1+(dO1O2-R2))*0.5;
-
+
pt1.SetCoord( O1.X() + t2*Dir.X()
- ,O1.Y() + t2*Dir.Y()
- ,O1.Z() + t2*Dir.Z());
+ ,O1.Y() + t2*Dir.Y()
+ ,O1.Z() + t2*Dir.Z());
}
else {
//-----------------------------------------------------------------
//--
//-----------------------------------------------------------------
if((dO1O2 > (R1+R2+Tol)) || (Rmax > (dO1O2+Rmin+Tol))) {
- typeres=IntAna_Empty;
+ typeres=IntAna_Empty;
}
else {
- //---------------------------------------------------------------
- //--
- //--
- //---------------------------------------------------------------
- Standard_Real Alpha=0.5*(R1*R1-R2*R2+dO1O2*dO1O2)/(dO1O2);
- Standard_Real Beta = R1*R1-Alpha*Alpha;
- Beta = (Beta>0.0)? Sqrt(Beta) : 0.0;
-
- if(Beta<= myEPSILON_MINI_CIRCLE_RADIUS) {
- typeres = IntAna_Point;
- Alpha = (R1 + (dO1O2 - R2)) * 0.5;
- }
- else {
- typeres = IntAna_Circle;
- dir1 = Dir;
- param1 = Beta;
- }
- pt1.SetCoord( O1.X() + Alpha*Dir.X()
- ,O1.Y() + Alpha*Dir.Y()
- ,O1.Z() + Alpha*Dir.Z());
-
- nbint=1;
+ //---------------------------------------------------------------
+ //--
+ //--
+ //---------------------------------------------------------------
+ Standard_Real Alpha=0.5*(R1*R1-R2*R2+dO1O2*dO1O2)/(dO1O2);
+ Standard_Real Beta = R1*R1-Alpha*Alpha;
+ Beta = (Beta>0.0)? Sqrt(Beta) : 0.0;
+
+ if(Beta<= myEPSILON_MINI_CIRCLE_RADIUS) {
+ typeres = IntAna_Point;
+ Alpha = (R1 + (dO1O2 - R2)) * 0.5;
+ }
+ else {
+ typeres = IntAna_Circle;
+ dir1 = Dir;
+ param1 = Beta;
+ }
+ pt1.SetCoord( O1.X() + Alpha*Dir.X()
+ ,O1.Y() + Alpha*Dir.Y()
+ ,O1.Z() + Alpha*Dir.Z());
+
+ nbint=1;
}
}
}
}
+
+//=======================================================================
+//function : IntAna_QuadQuadGeo
+//purpose : Plane - Torus
+//=======================================================================
+IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Pln& Pln,
+ const gp_Torus& Tor,
+ const Standard_Real Tol)
+: done(Standard_False),
+ nbint(0),
+ typeres(IntAna_Empty),
+ pt1(0,0,0),
+ pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
+ param1(0),
+ param2(0),
+ param3(0),
+ param4(0),
+ param1bis(0),
+ param2bis(0),
+ myCommonGen(Standard_False),
+ myPChar(0,0,0)
+{
+ InitTolerances();
+ Perform(Pln,Tor,Tol);
+}
+//=======================================================================
+//function : Perform
+//purpose :
+//=======================================================================
+void IntAna_QuadQuadGeo::Perform(const gp_Pln& Pln,
+ const gp_Torus& Tor,
+ const Standard_Real Tol)
+{
+ done = Standard_True;
+ //
+ Standard_Real aRMin, aRMaj;
+ //
+ aRMin = Tor.MinorRadius();
+ aRMaj = Tor.MajorRadius();
+ if (aRMin >= aRMaj) {
+ typeres = IntAna_NoGeometricSolution;
+ return;
+ }
+ //
+ const gp_Ax1 aPlnAx = Pln.Axis();
+ const gp_Ax1 aTorAx = Tor.Axis();
+ //
+ Standard_Boolean bParallel, bNormal;
+ //
+ bParallel = aTorAx.IsParallel(aPlnAx, myEPSILON_AXES_PARA);
+ bNormal = !bParallel ? aTorAx.IsNormal(aPlnAx, myEPSILON_AXES_PARA) : Standard_False;
+ if (!bNormal && !bParallel) {
+ typeres = IntAna_NoGeometricSolution;
+ return;
+ }
+ //
+ Standard_Real aDist;
+ //
+ gp_Pnt aTorLoc = aTorAx.Location();
+ if (bParallel) {
+ Standard_Real aDt, X, Y, Z, A, B, C, D, aDR, aTolNum;
+ //
+ aTolNum=myEPSILON_CYLINDER_DELTA_RADIUS;
+ //
+ Pln.Coefficients(A,B,C,D);
+ aTorLoc.Coord(X,Y,Z);
+ aDist = A*X + B*Y + C*Z + D;
+ //
+ aDR=Abs(aDist) - aRMin;
+ if (aDR > aTolNum) {
+ typeres=IntAna_Empty;
+ return;
+ }
+ //
+ if (Abs(aDR) < aTolNum) {
+ aDist=aRMin;
+ }
+ //
+ typeres = IntAna_Circle;
+ //
+ pt1.SetCoord(X - aDist*A, Y - aDist*B, Z - aDist*C);
+ aDt = Sqrt(Abs(aRMin*aRMin - aDist*aDist));
+ param1 = aRMaj + aDt;
+ dir1 = aTorAx.Direction();
+ nbint = 1;
+ if ((aDR < -aTolNum) && (aDt > Tol)) {
+ pt2 = pt1;
+ param2 = aRMaj - aDt;
+ dir2 = dir1;
+ nbint = 2;
+ }
+ }
+ //
+ else {
+ aDist = Pln.Distance(aTorLoc);
+ if (aDist > myEPSILON_DISTANCE) {
+ typeres = IntAna_NoGeometricSolution;
+ return;
+ }
+ //
+ typeres = IntAna_Circle;
+ param2 = param1 = aRMin;
+ dir2 = dir1 = aPlnAx.Direction();
+ nbint = 2;
+ //
+ gp_Dir aDir = aTorAx.Direction()^dir1;
+ pt1.SetXYZ(aTorLoc.XYZ() + aRMaj*aDir.XYZ());
+ pt2.SetXYZ(aTorLoc.XYZ() - aRMaj*aDir.XYZ());
+ }
+}
+
+//=======================================================================
+//function : IntAna_QuadQuadGeo
+//purpose : Cylinder - Torus
+//=======================================================================
+IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Cylinder& Cyl,
+ const gp_Torus& Tor,
+ const Standard_Real Tol)
+: done(Standard_False),
+ nbint(0),
+ typeres(IntAna_Empty),
+ pt1(0,0,0),
+ pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
+ param1(0),
+ param2(0),
+ param3(0),
+ param4(0),
+ param1bis(0),
+ param2bis(0),
+ myCommonGen(Standard_False),
+ myPChar(0,0,0)
+{
+ InitTolerances();
+ Perform(Cyl,Tor,Tol);
+}
+//=======================================================================
+//function : Perform
+//purpose :
+//=======================================================================
+void IntAna_QuadQuadGeo::Perform(const gp_Cylinder& Cyl,
+ const gp_Torus& Tor,
+ const Standard_Real Tol)
+{
+ done = Standard_True;
+ //
+ Standard_Real aRMin, aRMaj;
+ //
+ aRMin = Tor.MinorRadius();
+ aRMaj = Tor.MajorRadius();
+ if (aRMin >= aRMaj) {
+ typeres = IntAna_NoGeometricSolution;
+ return;
+ }
+ //
+ const gp_Ax1 aCylAx = Cyl.Axis();
+ const gp_Ax1 aTorAx = Tor.Axis();
+ //
+ const gp_Lin aLin(aTorAx);
+ const gp_Pnt aLocCyl = Cyl.Location();
+ //
+ if (!aTorAx.IsParallel(aCylAx, myEPSILON_AXES_PARA) ||
+ (aLin.Distance(aLocCyl) > myEPSILON_DISTANCE)) {
+ typeres = IntAna_NoGeometricSolution;
+ return;
+ }
+ //
+ Standard_Real aRCyl;
+ //
+ aRCyl = Cyl.Radius();
+ if (((aRCyl + Tol) < (aRMaj - aRMin)) || ((aRCyl - Tol) > (aRMaj + aRMin))) {
+ typeres = IntAna_Empty;
+ return;
+ }
+ //
+ typeres = IntAna_Circle;
+ //
+ Standard_Real aDist = Sqrt(Abs(aRMin*aRMin - (aRCyl-aRMaj)*(aRCyl-aRMaj)));
+ gp_XYZ aTorLoc = aTorAx.Location().XYZ();
+ //
+ dir1 = aTorAx.Direction();
+ pt1.SetXYZ(aTorLoc + aDist*dir1.XYZ());
+ param1 = aRCyl;
+ nbint = 1;
+ if ((aDist > Tol) && (aRCyl > (aRMaj - aRMin)) &&
+ (aRCyl < (aRMaj + aRMin))) {
+ dir2 = dir1;
+ pt2.SetXYZ(aTorLoc - aDist*dir2.XYZ());
+ param2 = param1;
+ nbint = 2;
+ }
+}
+
+//=======================================================================
+//function : IntAna_QuadQuadGeo
+//purpose : Cone - Torus
+//=======================================================================
+IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Cone& Con,
+ const gp_Torus& Tor,
+ const Standard_Real Tol)
+: done(Standard_False),
+ nbint(0),
+ typeres(IntAna_Empty),
+ pt1(0,0,0),
+ pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
+ param1(0),
+ param2(0),
+ param3(0),
+ param4(0),
+ param1bis(0),
+ param2bis(0),
+ myCommonGen(Standard_False),
+ myPChar(0,0,0)
+{
+ InitTolerances();
+ Perform(Con,Tor,Tol);
+}
+//=======================================================================
+//function : Perform
+//purpose :
+//=======================================================================
+void IntAna_QuadQuadGeo::Perform(const gp_Cone& Con,
+ const gp_Torus& Tor,
+ const Standard_Real Tol)
+{
+ done = Standard_True;
+ //
+ Standard_Real aRMin, aRMaj;
+ //
+ aRMin = Tor.MinorRadius();
+ aRMaj = Tor.MajorRadius();
+ if (aRMin >= aRMaj) {
+ typeres = IntAna_NoGeometricSolution;
+ return;
+ }
+ //
+ const gp_Ax1 aConAx = Con.Axis();
+ const gp_Ax1 aTorAx = Tor.Axis();
+ //
+ const gp_Lin aLin(aTorAx);
+ const gp_Pnt aConApex = Con.Apex();
+ //
+ if (!aTorAx.IsParallel(aConAx, myEPSILON_AXES_PARA) ||
+ (aLin.Distance(aConApex) > myEPSILON_DISTANCE)) {
+ typeres = IntAna_NoGeometricSolution;
+ return;
+ }
+ //
+ Standard_Real anAngle, aDist, aParam[4], aDt;
+ Standard_Integer i;
+ gp_Pnt aTorLoc, aPCT, aPN, aPt[4];
+ gp_Dir aDir[4];
+ //
+ anAngle = Con.SemiAngle();
+ aTorLoc = aTorAx.Location();
+ //
+ aPN.SetXYZ(aTorLoc.XYZ() + aRMaj*Tor.YAxis().Direction().XYZ());
+ gp_Dir aDN (gp_Vec(aTorLoc, aPN));
+ gp_Ax1 anAxCLRot(aConApex, aDN);
+ gp_Lin aConL = aLin.Rotated(anAxCLRot, anAngle);
+ gp_Dir aDL = aConL.Position().Direction();
+ gp_Dir aXDir = Tor.XAxis().Direction();
+ //
+ typeres = IntAna_Empty;
+ //
+ for (i = 0; i < 2; ++i) {
+ if (i) {
+ aXDir.Reverse();
+ }
+ aPCT.SetXYZ(aTorLoc.XYZ() + aRMaj*aXDir.XYZ());
+ //
+ aDist = aConL.Distance(aPCT);
+ if (aDist > aRMin+Tol) {
+ continue;
+ }
+ //
+ typeres = IntAna_Circle;
+ //
+ gp_XYZ aPh = aPCT.XYZ() - aDist*aConL.Normal(aPCT).Direction().XYZ();
+ aDt = Sqrt(Abs(aRMin*aRMin - aDist*aDist));
+ //
+ gp_Pnt aP;
+ gp_XYZ aDVal = aDt*aDL.XYZ();
+ aP.SetXYZ(aPh + aDVal);
+ aParam[nbint] = aLin.Distance(aP);
+ aPt[nbint].SetXYZ(aP.XYZ() - aParam[nbint]*aXDir.XYZ());
+ aDir[nbint] = aTorAx.Direction();
+ ++nbint;
+ if ((aDist < aRMin) && (aDt > Tol)) {
+ aP.SetXYZ(aPh - aDVal);
+ aParam[nbint] = aLin.Distance(aP);
+ aPt[nbint].SetXYZ(aP.XYZ() - aParam[nbint]*aXDir.XYZ());
+ aDir[nbint] = aDir[nbint-1];
+ ++nbint;
+ }
+ }
+ //
+ for (i = 0; i < nbint; ++i) {
+ switch (i) {
+ case 0:{
+ pt1 = aPt[i];
+ param1 = aParam[i];
+ dir1 = aDir[i];
+ break;
+ }
+ case 1:{
+ pt2 = aPt[i];
+ param2 = aParam[i];
+ dir2 = aDir[i];
+ break;
+ }
+ case 2:{
+ pt3 = aPt[i];
+ param3 = aParam[i];
+ dir3 = aDir[i];
+ break;
+ }
+ case 3:{
+ pt4 = aPt[i];
+ param4 = aParam[i];
+ dir4 = aDir[i];
+ break;
+ }
+ default:
+ break;
+ }
+ }
+}
+
+//=======================================================================
+//function : IntAna_QuadQuadGeo
+//purpose : Sphere - Torus
+//=======================================================================
+IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Sphere& Sph,
+ const gp_Torus& Tor,
+ const Standard_Real Tol)
+: done(Standard_False),
+ nbint(0),
+ typeres(IntAna_Empty),
+ pt1(0,0,0),
+ pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
+ param1(0),
+ param2(0),
+ param3(0),
+ param4(0),
+ param1bis(0),
+ param2bis(0),
+ myCommonGen(Standard_False),
+ myPChar(0,0,0)
+{
+ InitTolerances();
+ Perform(Sph,Tor,Tol);
+}
+//=======================================================================
+//function : Perform
+//purpose :
+//=======================================================================
+void IntAna_QuadQuadGeo::Perform(const gp_Sphere& Sph,
+ const gp_Torus& Tor,
+ const Standard_Real Tol)
+{
+ done = Standard_True;
+ //
+ Standard_Real aRMin, aRMaj;
+ //
+ aRMin = Tor.MinorRadius();
+ aRMaj = Tor.MajorRadius();
+ if (aRMin >= aRMaj) {
+ typeres = IntAna_NoGeometricSolution;
+ return;
+ }
+ //
+ const gp_Ax1 aTorAx = Tor.Axis();
+ const gp_Lin aLin(aTorAx);
+ const gp_Pnt aSphLoc = Sph.Location();
+ //
+ if (aLin.Distance(aSphLoc) > myEPSILON_DISTANCE) {
+ typeres = IntAna_NoGeometricSolution;
+ return;
+ }
+ //
+ Standard_Real aRSph, aDist;
+ gp_Pnt aTorLoc;
+ //
+ gp_Dir aXDir = Tor.XAxis().Direction();
+ aTorLoc.SetXYZ(aTorAx.Location().XYZ() + aRMaj*aXDir.XYZ());
+ aRSph = Sph.Radius();
+ //
+ gp_Vec aVec12(aTorLoc, aSphLoc);
+ aDist = aVec12.Magnitude();
+ if (((aDist - Tol) > (aRMin + aRSph)) ||
+ ((aDist + Tol) < Abs(aRMin - aRSph))) {
+ typeres = IntAna_Empty;
+ return;
+ }
+ //
+ typeres = IntAna_Circle;
+ //
+ Standard_Real anAlpha, aBeta;
+ //
+ anAlpha = 0.5*(aRMin*aRMin - aRSph*aRSph + aDist*aDist ) / aDist;
+ aBeta = Sqrt(Abs(aRMin*aRMin - anAlpha*anAlpha));
+ //
+ gp_Dir aDir12(aVec12);
+ gp_XYZ aPh = aTorLoc.XYZ() + anAlpha*aDir12.XYZ();
+ gp_Dir aDC = Tor.YAxis().Direction()^aDir12;
+ //
+ gp_Pnt aP;
+ gp_XYZ aDVal = aBeta*aDC.XYZ();
+ aP.SetXYZ(aPh + aDVal);
+ param1 = aLin.Distance(aP);
+ pt1.SetXYZ(aP.XYZ() - param1*aXDir.XYZ());
+ dir1 = aTorAx.Direction();
+ nbint = 1;
+ if ((aDist < (aRSph + aRMin)) && (aDist > Abs(aRSph - aRMin)) &&
+ (aDVal.Modulus() > Tol)) {
+ aP.SetXYZ(aPh - aDVal);
+ param2 = aLin.Distance(aP);
+ pt2.SetXYZ(aP.XYZ() - param2*aXDir.XYZ());
+ dir2 = dir1;
+ nbint = 2;
+ }
+}
+
+//=======================================================================
+//function : IntAna_QuadQuadGeo
+//purpose : Torus - Torus
+//=======================================================================
+IntAna_QuadQuadGeo::IntAna_QuadQuadGeo(const gp_Torus& Tor1,
+ const gp_Torus& Tor2,
+ const Standard_Real Tol)
+: done(Standard_False),
+ nbint(0),
+ typeres(IntAna_Empty),
+ pt1(0,0,0),
+ pt2(0,0,0),
+ pt3(0,0,0),
+ pt4(0,0,0),
+ param1(0),
+ param2(0),
+ param3(0),
+ param4(0),
+ param1bis(0),
+ param2bis(0),
+ myCommonGen(Standard_False),
+ myPChar(0,0,0)
+{
+ InitTolerances();
+ Perform(Tor1,Tor2,Tol);
+}
+//=======================================================================
+//function : Perform
+//purpose :
+//=======================================================================
+void IntAna_QuadQuadGeo::Perform(const gp_Torus& Tor1,
+ const gp_Torus& Tor2,
+ const Standard_Real Tol)
+{
+ done = Standard_True;
+ //
+ Standard_Real aRMin1, aRMin2, aRMaj1, aRMaj2;
+ //
+ aRMin1 = Tor1.MinorRadius();
+ aRMaj1 = Tor1.MajorRadius();
+ aRMin2 = Tor2.MinorRadius();
+ aRMaj2 = Tor2.MajorRadius();
+ //
+ const gp_Ax1& anAx1 = Tor1.Axis();
+ const gp_Ax1& anAx2 = Tor2.Axis();
+ //
+ const gp_Pnt& aLoc1 = anAx1.Location();
+ const gp_Pnt& aLoc2 = anAx2.Location();
+ //
+ gp_Lin aL1(anAx1);
+ if (!anAx1.IsParallel(anAx2, myEPSILON_AXES_PARA) ||
+ (aL1.Distance(aLoc2) > myEPSILON_DISTANCE)) {
+ typeres = IntAna_NoGeometricSolution;
+ return;
+ }
+ //
+ if (aLoc1.IsEqual(aLoc2, Tol) &&
+ (Abs(aRMin1 - aRMin2) <= Tol) &&
+ (Abs(aRMaj1 - aRMaj2) <= Tol)) {
+ typeres = IntAna_Same;
+ return;
+ }
+ //
+ if (aRMin1 >= aRMaj1 || aRMin2 >= aRMaj2) {
+ typeres = IntAna_NoGeometricSolution;
+ return;
+ }
+ //
+ Standard_Real aDist;
+ gp_Pnt aP1, aP2;
+ //
+ gp_Dir aXDir1 = Tor1.XAxis().Direction();
+ aP1.SetXYZ(aLoc1.XYZ() + aRMaj1*aXDir1.XYZ());
+ aP2.SetXYZ(aLoc2.XYZ() + aRMaj2*aXDir1.XYZ());
+ //
+ gp_Vec aV12(aP1, aP2);
+ aDist = aV12.Magnitude();
+ if (((aDist - Tol) > (aRMin1 + aRMin2)) ||
+ ((aDist + Tol) < Abs(aRMin1 - aRMin2))) {
+ typeres = IntAna_Empty;
+ return;
+ }
+ //
+ typeres = IntAna_Circle;
+ //
+ Standard_Real anAlpha, aBeta;
+ //
+ anAlpha = 0.5*(aRMin1*aRMin1 - aRMin2*aRMin2 + aDist*aDist ) / aDist;
+ aBeta = Sqrt(Abs(aRMin1*aRMin1 - anAlpha*anAlpha));
+ //
+ gp_Dir aDir12(aV12);
+ gp_XYZ aPh = aP1.XYZ() + anAlpha*aDir12.XYZ();
+ gp_Dir aDC = Tor1.YAxis().Direction()^aDir12;
+ //
+ gp_Pnt aP;
+ gp_XYZ aDVal = aBeta*aDC.XYZ();
+ aP.SetXYZ(aPh + aDVal);
+ param1 = aL1.Distance(aP);
+ pt1.SetXYZ(aP.XYZ() - param1*aXDir1.XYZ());
+ dir1 = anAx1.Direction();
+ nbint = 1;
+ if ((aDist < (aRMin1 + aRMin2)) && (aDist > Abs(aRMin1 - aRMin2)) &&
+ aDVal.Modulus() > Tol) {
+ aP.SetXYZ(aPh - aDVal);
+ param2 = aL1.Distance(aP);
+ pt2.SetXYZ(aP.XYZ() - param2*aXDir1.XYZ());
+ dir2 = dir1;
+ nbint = 2;
+ }
+}
+
//=======================================================================
//function : Point
//purpose : Returns a Point
//=======================================================================
gp_Pnt IntAna_QuadQuadGeo::Point(const Standard_Integer n) const
{
- if(!done) { StdFail_NotDone::Raise(); }
- if(n>nbint || n<1) { Standard_DomainError::Raise(); }
+ if(!done) { throw StdFail_NotDone(); }
+ if(n>nbint || n<1) { throw Standard_DomainError(); }
if(typeres==IntAna_PointAndCircle) {
- if(n!=1) { Standard_DomainError::Raise(); }
+ if(n!=1) { throw Standard_DomainError(); }
if(param1==0.0) return(pt1);
return(pt2);
}
//=======================================================================
gp_Lin IntAna_QuadQuadGeo::Line(const Standard_Integer n) const
{
- if(!done) { StdFail_NotDone::Raise(); }
+ if(!done) { throw StdFail_NotDone(); }
if((n>nbint) || (n<1) || (typeres!=IntAna_Line)) {
- Standard_DomainError::Raise();
+ throw Standard_DomainError();
}
if(n==1) { return(gp_Lin(pt1,dir1)); }
else { return(gp_Lin(pt2,dir2)); }
//=======================================================================
gp_Circ IntAna_QuadQuadGeo::Circle(const Standard_Integer n) const
{
- if(!done) { StdFail_NotDone::Raise(); }
+ if(!done) { throw StdFail_NotDone(); }
if(typeres==IntAna_PointAndCircle) {
- if(n!=1) { Standard_DomainError::Raise(); }
+ if(n!=1) { throw Standard_DomainError(); }
if(param2==0.0) return(gp_Circ(DirToAx2(pt1,dir1),param1));
return(gp_Circ(DirToAx2(pt2,dir2),param2));
}
else if((n>nbint) || (n<1) || (typeres!=IntAna_Circle)) {
- Standard_DomainError::Raise();
+ throw Standard_DomainError();
}
- if(n==1) { return(gp_Circ(DirToAx2(pt1,dir1),param1)); }
- else { return(gp_Circ(DirToAx2(pt2,dir2),param2)); }
+ if (n==1) { return(gp_Circ(DirToAx2(pt1,dir1),param1));}
+ else if (n==2) { return(gp_Circ(DirToAx2(pt2,dir2),param2));}
+ else if (n==3) { return(gp_Circ(DirToAx2(pt3,dir3),param3));}
+ else { return(gp_Circ(DirToAx2(pt4,dir4),param4));}
}
//=======================================================================
//=======================================================================
gp_Elips IntAna_QuadQuadGeo::Ellipse(const Standard_Integer n) const
{
- if(!done) { StdFail_NotDone::Raise(); }
+ if(!done) { throw StdFail_NotDone(); }
if((n>nbint) || (n<1) || (typeres!=IntAna_Ellipse)) {
- Standard_DomainError::Raise();
+ throw Standard_DomainError();
}
if(n==1) {
gp_Parab IntAna_QuadQuadGeo::Parabola(const Standard_Integer n) const
{
if(!done) {
- StdFail_NotDone::Raise();
+ throw StdFail_NotDone();
}
if (typeres!=IntAna_Parabola) {
- Standard_DomainError::Raise();
+ throw Standard_DomainError();
}
if((n>nbint) || (n!=1)) {
- Standard_OutOfRange::Raise();
+ throw Standard_OutOfRange();
}
return(gp_Parab(gp_Ax2( pt1
- ,dir1
- ,dir2)
- ,param1));
+ ,dir1
+ ,dir2)
+ ,param1));
}
//=======================================================================
//function : Hyperbola
gp_Hypr IntAna_QuadQuadGeo::Hyperbola(const Standard_Integer n) const
{
if(!done) {
- StdFail_NotDone::Raise();
+ throw StdFail_NotDone();
}
if((n>nbint) || (n<1) || (typeres!=IntAna_Hyperbola)) {
- Standard_DomainError::Raise();
+ throw Standard_DomainError();
}
if(n==1) {
return(gp_Hypr(gp_Ax2( pt1
- ,dir1
- ,dir2)
- ,param1,param1bis));
+ ,dir1
+ ,dir2)
+ ,param1,param1bis));
}
else {
return(gp_Hypr(gp_Ax2( pt2
- ,dir1
- ,dir2.Reversed())
- ,param2,param2bis));
+ ,dir1
+ ,dir2.Reversed())
+ ,param2,param2bis));
}
}
//=======================================================================
m=(aC[k]>0.);
aNum=(m)? aC[k] : -aC[k];
if (aNum>aR1 && aNum<aR2) {
- if (m) {
- aC[k]=1.;
- }
- else {
- aC[k]=-1.;
- }
- //
- aC[(k+1)%3]=0.;
- aC[(k+2)%3]=0.;
- break;
+ if (m) {
+ aC[k]=1.;
+ }
+ else {
+ aC[k]=-1.;
+ }
+ //
+ aC[(k+1)%3]=0.;
+ aC[(k+2)%3]=0.;
+ break;
}
}
aDir.SetCoord(aC[0], aC[1], aC[2]);
projects / occt.git / blobdiff