summaryrefslogtreecommitdiff
path: root/src/IntAna2d/IntAna2d_AnaIntersection_2.cxx
blob: 5fc149a26b59446cd018afdf7b4c9b28a5e66955 (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
#include <IntAna2d_AnaIntersection.jxx>

#include <gp_Vec2d.hxx>

void IntAna2d_AnaIntersection::Perform (const gp_Circ2d& C1,
					const gp_Circ2d& C2) {

  done=Standard_False;
  Standard_Real d=C1.Location().Distance(C2.Location());
  Standard_Real R1=C1.Radius();
  Standard_Real R2=C2.Radius();
  Standard_Real sum=R1+R2;
  Standard_Real dif=Abs(R1-R2);
 

  if (d<=RealEpsilon()) {                 // Cercle concentriques
    para=Standard_True;
    nbp=0;
    if (dif<=RealEpsilon()) {             // Cercles confondus
      empt=Standard_False;
      iden=Standard_True;
    }
    else {                               // Cercles paralleles
      empt=Standard_True;
      iden=Standard_False;
    }
  }
  else if ((d-sum)>Epsilon(sum)) {        // Cercles exterieurs l un a l autre
                                         // et No solution
    empt=Standard_True;
    para=Standard_False;
    iden=Standard_False;
    nbp=0;
  }
  else if (Abs(d-sum)<=Epsilon(sum)) {    // Cercles exterieurs et tangents
    empt=Standard_False;
    para=Standard_False;
    iden=Standard_False;
    nbp=1;
    gp_Vec2d ax(C1.Location(),C2.Location());
    gp_Vec2d Ox1(C1.XAxis().Direction());
    gp_Vec2d Ox2(C2.XAxis().Direction());

    Standard_Real XS = ( C1.Location().X()*R2 + C2.Location().X()*R1 ) / sum;
    Standard_Real YS = ( C1.Location().Y()*R2 + C2.Location().Y()*R1 ) / sum;
    Standard_Real ang1=Ox1.Angle(ax);                     // Resultat entre -PI et +PI
    Standard_Real ang2=Ox2.Angle(ax) + PI;
    if (ang1<0) {ang1=2*PI+ang1;}                // On revient entre 0 et 2PI
    lpnt[0].SetValue(XS,YS,ang1,ang2);
  }
  else if (((sum-d)>Epsilon(d)) && ((d-dif)>Epsilon(d))) {
    empt=Standard_False;
    para=Standard_False;
    iden=Standard_False;
    nbp=2;
    gp_Vec2d ax(C1.Location(),C2.Location());
    gp_Vec2d Ox1(C1.XAxis().Direction());
    gp_Vec2d Ox2(C2.XAxis().Direction());
    Standard_Real ref1=Ox1.Angle(ax);                       // Resultat entre -PI et +PI
    Standard_Real ref2=Ox2.Angle(ax);                       // Resultat entre -PI et +PI

    Standard_Real l1=(d*d + R1*R1 -R2*R2)/(2.0*d);
    Standard_Real h= Sqrt(R1*R1-l1*l1);

    Standard_Real XS1= C1.Location().X() + l1*ax.X()/d - h*ax.Y()/d;
    Standard_Real YS1= C1.Location().Y() + l1*ax.Y()/d + h*ax.X()/d;

    Standard_Real XS2= C1.Location().X() + l1*ax.X()/d + h*ax.Y()/d;
    Standard_Real YS2= C1.Location().Y() + l1*ax.Y()/d - h*ax.X()/d;

    Standard_Real sint1=h /R1;
    Standard_Real cost1=l1/R1;

    Standard_Real sint2=h     /R2;
    Standard_Real cost2=(l1-d)/R2;

    Standard_Real ang1,ang2;

    // ang1 et ang2 correspondent aux solutions avec sinus positif
    // si l'axe de reference est l'axe des centres C1C2
    // On prend l'arccos entre pi/2 et 3pi/2, l'arcsin sinon.

    if (Abs(cost1)<=0.707) {
      ang1=ACos(cost1); 
    }
    else {
      ang1=ASin(sint1);
      if (cost1<0.0) {ang1=PI-ang1;}
    }
    if (Abs(cost2)<=0.707) {
      ang2=ACos(cost2);
    }
    else {
      ang2=ASin(sint2);
      if (cost2<0.0) {ang2=PI-ang2;}
    }
    Standard_Real ang11=ref1+ang1;
    Standard_Real ang21=ref2+ang2;
    Standard_Real ang12=ref1-ang1;
    Standard_Real ang22=ref2-ang2;
    if (ang11<0.) {
      ang11=2*PI+ang11;
    }
    else if (ang11>=2*PI) {
      ang11=ang11-2*PI;
    }
    if (ang21<0.) {
      ang21=2*PI+ang21;
    }
    else if (ang21>=2*PI) {
      ang21=ang21-2*PI;
    }
    if (ang12<0.) {
      ang12=2*PI+ang12;
    }
    else if (ang12>=2*PI) {
      ang12=ang12-2*PI;
    }
    if (ang22<0.) {
      ang22=2*PI+ang22;
    }
    else if (ang22>=2*PI) {
      ang22=ang22-2*PI;
    }
    lpnt[0].SetValue(XS1,YS1,ang11,ang21);
    lpnt[1].SetValue(XS2,YS2,ang12,ang22);
  }
  else if (Abs(d-dif)<=Epsilon(d)) {    // Cercles tangents interieurs
    empt=Standard_False;
    para=Standard_False;
    iden=Standard_False;
    nbp=1;
    gp_Vec2d ax(C1.Location(),C2.Location());
    gp_Vec2d Ox1(C1.XAxis().Direction());
    gp_Vec2d Ox2(C2.XAxis().Direction());
    Standard_Real ang1=Ox1.Angle(ax);                       // Resultat entre -PI et +PI
    Standard_Real ang2=Ox2.Angle(ax);
    if (ang1<0) {ang1=2*PI+ang1;}                  // On revient entre 0 et 2PI
    if (ang2<0) {ang2=2*PI+ang2;}                  // On revient entre 0 et 2PI
    Standard_Real XS = ( C1.Location().X()*R2 - C2.Location().X()*R1 ) / (R2 - R1);
    Standard_Real YS = ( C1.Location().Y()*R2 - C2.Location().Y()*R1 ) / (R2 - R1);
    lpnt[0].SetValue(XS,YS,ang1,ang2);
  }
  else {                           // On doit avoir d<dif-Resol et d<>0 donc
                                   // 1 cercle dans l autre et no solution
    empt=Standard_True;
    para=Standard_False;
    iden=Standard_False;
    nbp=0;
  }
  done=Standard_True;
}