gp_Dir DirCyl = Cyl1.Position().Direction();
Standard_Real ProjP2OnDirCyl1=gp_Vec(DirCyl).Dot(gp_Vec(P1,P2t));
+ //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());
typeres=IntAna_Empty;
nbint=0;
}
- else if(DistA1A2>(R1pR2))
+ else if((R1pR2 - DistA1A2) <= RealSmall())
{
//-- 1 Tangent line -------------------------------------OK
typeres=IntAna_Line;
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.;
+ 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.
+ */
+
- if((Beta+Beta)<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.SetCoord( P1.X() + Alpha*DirA1A2.X()
- ,P1.Y() + Alpha*DirA1A2.Y()
- ,P1.Z() + Alpha*DirA1A2.Z());
+ pt1.SetXYZ(P1.XYZ() + DirA1A2.XYZ()*R1*aCos);
}
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());
+ {
+ //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))
--- /dev/null
+puts "========"
+puts "OCC25782"
+puts "========"
+puts ""
+######################################################
+# The result of intersection between two cylinders is incorrect
+# Algorithm must find one curves only
+######################################################
+
+set GoodNbCurv 1
+
+restore [locate_data_file bug25782_fz19.brep] b1
+restore [locate_data_file bug25782_fz53.brep] b2
+
+mksurface s1 b1
+mksurface s2 b2
+
+intersect res s1 s2
+
+set che [whatis res]
+set ind [string first "3d curve" $che]
+if {${ind} >= 0} {
+ #Only variable "res" exists
+
+ copy res res_1
+}
+
+set ic 1
+set AllowRepeate 1
+while { $AllowRepeate != 0 } {
+ set che [whatis res_$ic]
+ set ind [string first "3d curve" $che]
+ if {${ind} < 0} {
+ set AllowRepeate 0
+ } else {
+ dlog reset
+ dlog on
+ xdistcs res_$ic s1 0 100 10
+ set Log1 [dlog get]
+ set List1 [split ${Log1} {TD= \t\n}]
+ set Tolerance 1.0e-7
+ set Limit_Tol 1.0e-7
+ set D_good 0.
+ checkList ${List1} ${Tolerance} ${D_good} ${Limit_Tol}
+
+ dlog reset
+ dlog on
+ xdistcs res_$ic s2 0 100 10
+ set Log1 [dlog get]
+ set List1 [split ${Log1} {TD= \t\n}]
+ set Tolerance 1.0e-7
+ set Limit_Tol 1.0e-7
+ set D_good 0.
+ checkList ${List1} ${Tolerance} ${D_good} ${Limit_Tol}
+
+ incr ic
+ }
+}
+
+if {[expr {$ic - 1}] == $GoodNbCurv} {
+ puts "OK: Curve Number is good!"
+} else {
+ puts "Error: Curves Number is bad!"
+}
+
+set log [bopcurves b1 b2]
+
+regexp {Tolerance Reached=+([-0-9.+eE]+)\n+([-0-9.+eE]+)} ${log} full Toler NbCurv
+set MaxTol 1.e-7
+if {${Toler} > ${MaxTol}} {
+ puts "Error: Tolerance is too big!"
+}
+
+if {$NbCurv != $GoodNbCurv} {
+ puts "Error: Curves Number is bad!"
+}
+
+for {set i 1} {$i <= ${NbCurv}} {incr i} {
+ bounds c_$i U1 U2
+
+ if {[dval U2-U1] < 1.0e-9} {
+ puts "Error: Wrong curve's range!"
+ }
+
+ dlog reset
+ dlog on
+ xdistcs c_$i s1 U1 U2 10
+ set Log2 [dlog get]
+ set List2 [split ${Log2} {TD= \t\n}]
+ set Tolerance 1.0e-7
+ set Limit_Tol 1.0e-7
+ set D_good 0.
+ checkList ${List2} ${Tolerance} ${D_good} ${Limit_Tol}
+
+ dlog reset
+ dlog on
+ xdistcs c_$i s2 U1 U2 10
+ set Log2 [dlog get]
+ set List2 [split ${Log2} {TD= \t\n}]
+ set Tolerance 1.0e-7
+ set Limit_Tol 1.0e-7
+ set D_good 0.
+ checkList ${List2} ${Tolerance} ${D_good} ${Limit_Tol}
+}
--- /dev/null
+puts "========"
+puts "OCC25782"
+puts "========"
+puts ""
+######################################################
+# The result of intersection between two cylinders is incorrect
+######################################################
+
+set GoodNbCurv 2
+
+cylinder s1 0 0 0 12 35 47 5
+cylinder s2 3 2 8 12 35 47 4
+
+set bug_info [intersect res s1 s2]
+
+set che [whatis res]
+set ind [string first "3d curve" $che]
+if {${ind} >= 0} {
+ #Only variable "res" exists
+
+ copy res res_1
+}
+
+if {[llength ${bug_info}] != $GoodNbCurv} {
+ puts "Error: The result of intersection between two cylinders is incorrect"
+}
+
+set Tolerance 1.e-7
+set D_good 0.
+set Limit_Tol 1.0e-7
+
+set ic 1
+set AllowRepeate 1
+while { $AllowRepeate != 0 } {
+ set che [whatis res_$ic]
+ set ind [string first "3d curve" $che]
+ if {${ind} < 0} {
+ set AllowRepeate 0
+ } else {
+ if { [regexp {\*\*\nLine} [dump res_$ic]] } {
+ #puts "OK : Correct intersection"
+ } else {
+ puts "Error : Bad intersection"
+ }
+
+ dlog reset
+ dlog on
+ xdistcs res_$ic s1 0 100 10
+ set Log1 [dlog get]
+ set List1 [split ${Log1} {TD= \t\n}]
+ set Tolerance 1.0e-7
+ set Limit_Tol 1.0e-7
+ set D_good 0.
+ checkList ${List1} ${Tolerance} ${D_good} ${Limit_Tol}
+
+ dlog reset
+ dlog on
+ xdistcs res_$ic s2 0 100 10
+ set Log1 [dlog get]
+ set List1 [split ${Log1} {TD= \t\n}]
+ set Tolerance 1.0e-7
+ set Limit_Tol 1.0e-7
+ set D_good 0.
+ checkList ${List1} ${Tolerance} ${D_good} ${Limit_Tol}
+
+ incr ic
+ }
+}