0022792: Globally defined symbol PI conflicts with VTK definition (Intel compiler)
[occt.git] / src / IGESGeom / IGESGeom_ConicArc.cxx
1 //--------------------------------------------------------------------
2 //
3 //  File Name : IGESGeom_ConicArc.cxx
4 //  Date      :
5 //  Author    : CKY / Contract Toubro-Larsen
6 //  Copyright : MATRA-DATAVISION 1993
7 //
8 //--------------------------------------------------------------------
9 //#59 rln 29.12.98 PRO17015
10
11 #include <IGESGeom_ConicArc.ixx>
12 #include <gp_Dir2d.hxx>
13 #include <gp_GTrsf.hxx>
14
15
16     IGESGeom_ConicArc::IGESGeom_ConicArc ()    {  }
17
18
19     void  IGESGeom_ConicArc::Init
20   (const Standard_Real A, const Standard_Real B,
21    const Standard_Real C, const Standard_Real D,  const Standard_Real E,
22    const Standard_Real F, const Standard_Real ZT, const gp_XY& aStart,
23    const gp_XY& anEnd)
24 {
25   theA = A;
26   theB = B;
27   theC = C;
28   theD = D;
29   theE = E;
30   theF = F;
31   theZT = ZT;
32   theStart = aStart;
33   theEnd   = anEnd;
34
35   Standard_Integer fn = FormNumber();
36   if (fn == 0) fn = ComputedFormNumber();
37   InitTypeAndForm(104,fn);
38 }
39
40     Standard_Boolean  IGESGeom_ConicArc::OwnCorrect ()
41 {
42   Standard_Integer cfn = ComputedFormNumber();
43   if (FormNumber() == cfn) return Standard_False;
44   InitTypeAndForm(104,cfn);
45   return Standard_True;
46 }
47
48     void  IGESGeom_ConicArc::Equation
49   (Standard_Real& A, Standard_Real& B, Standard_Real& C, Standard_Real& D,
50    Standard_Real& E, Standard_Real& F) const
51 {
52   A = theA;
53   B = theB;
54   C = theC;
55   D = theD;
56   E = theE;
57   F = theF;
58 }
59
60     Standard_Real  IGESGeom_ConicArc::ZPlane () const
61 {
62   return theZT;
63 }
64
65     gp_Pnt2d  IGESGeom_ConicArc::StartPoint () const
66 {
67   gp_Pnt2d start(theStart.X(), theStart.Y());
68   return start;
69 }
70
71     gp_Pnt  IGESGeom_ConicArc::TransformedStartPoint () const
72 {
73   gp_XYZ start(theStart.X(), theStart.Y(), theZT);
74   if (HasTransf()) Location().Transforms(start);
75   gp_Pnt transStart(start);
76   return transStart;
77 }
78
79     gp_Pnt2d  IGESGeom_ConicArc::EndPoint () const
80 {
81   gp_Pnt2d end(theEnd.X(), theEnd.Y());
82   return end;
83 }
84
85     gp_Pnt  IGESGeom_ConicArc::TransformedEndPoint () const
86 {
87   gp_XYZ end(theEnd.X(), theEnd.Y(), theZT);
88   if (HasTransf()) Location().Transforms(end);
89   gp_Pnt transEnd(end);
90   return transEnd;
91 }
92
93     Standard_Integer  IGESGeom_ConicArc::ComputedFormNumber () const
94 {
95   Standard_Real eps,eps2,eps4;
96   eps = 1.E-08;  eps2 = eps*eps; eps4 = eps2*eps2;//#59 rln
97   Standard_Real Q1 = theA   * (theC*theF   - theE*theE/4.)
98     + theB/2. * (theE*theD/4. - theB*theF/2.)
99       + theD/2. * (theB*theE/4. - theC*theD/2.);
100   Standard_Real Q2 = theA*theC - theB*theB/4;
101   Standard_Real Q3 = theA + theC;
102
103 //  Resultats
104   //#59 rln 29.12.98 PRO17015 face#67, ellipse
105   //each Qi has its own dimension:
106   //[Q1] = L^-4, [Q2]=L^-4, [Q3]=L^-2
107   if (Q2       >  eps4 && Q1*Q3    < 0   ) return 1;    // Ellipse
108   if (Q2       < -eps4 && Abs (Q1) > eps4) return 2;    // Hyperbola
109   if (Abs (Q2) <= eps4 && Abs (Q1) > eps4) return 3;    // Parabola
110   return 0;
111 }
112
113     Standard_Boolean  IGESGeom_ConicArc::IsFromParabola () const
114 {
115   Standard_Integer fn = FormNumber();
116   if (fn == 0)     fn = ComputedFormNumber();
117   return (fn == 3);
118 }
119
120     Standard_Boolean  IGESGeom_ConicArc::IsFromEllipse () const
121 {
122   Standard_Integer fn = FormNumber();
123   if (fn == 0)     fn = ComputedFormNumber();
124   return (fn == 1);
125 }
126
127     Standard_Boolean  IGESGeom_ConicArc::IsFromHyperbola () const
128 {
129   Standard_Integer fn = FormNumber();
130   if (fn == 0)     fn = ComputedFormNumber();
131   return (fn == 2);
132 }
133
134     Standard_Boolean  IGESGeom_ConicArc::IsClosed () const
135 {
136   return ((theStart.X() == theEnd.X()) && (theStart.Y() == theEnd.Y()));
137 }
138
139     gp_Dir IGESGeom_ConicArc::Axis () const
140 {
141   gp_Dir axis(0.0 , 0.0, 1.0);
142   return axis;
143 }
144
145 //    Valeurs calculees
146
147     gp_Dir IGESGeom_ConicArc::TransformedAxis () const
148 {
149   gp_XYZ axis(0.0 , 0.0, 1.0);
150   if (!HasTransf()) return gp_Dir(axis);
151   gp_GTrsf loc = Location();
152   loc.SetTranslationPart (gp_XYZ(0.,0.,0.));
153   loc.Transforms(axis);
154   return gp_Dir(axis);
155 }
156
157
158     void  IGESGeom_ConicArc::Definition
159   (gp_Pnt& Center, gp_Dir& MainAxis,
160    Standard_Real& Rmin, Standard_Real& Rmax) const
161 {
162   Standard_Real Xcen,Ycen, Xax,Yax;
163   ComputedDefinition (Xcen,Ycen, Xax,Yax, Rmin,Rmax);
164   Center.SetCoord (Xcen,Ycen,theZT);
165   MainAxis.SetCoord (Xax,Yax,0.);
166 }
167
168     void  IGESGeom_ConicArc::TransformedDefinition
169   (gp_Pnt& Center, gp_Dir& MainAxis,
170    Standard_Real& Rmin, Standard_Real& Rmax) const
171 {
172   if (!HasTransf()) {
173     Definition (Center,MainAxis,Rmin,Rmax);
174     return;
175   }
176   Standard_Real Xcen,Ycen, Xax,Yax;
177   ComputedDefinition (Xcen,Ycen, Xax,Yax, Rmin,Rmax);
178   gp_GTrsf loc = Location();
179   gp_XYZ cen  (Xcen,Ycen,theZT);
180   gp_XYZ axis (Xax, Yax, 0.);
181   loc.Transforms (cen);
182   loc.SetTranslationPart (gp_XYZ(0.,0.,0.));
183   loc.Transforms (axis);
184   Center.SetCoord   (cen.X(), cen.Y(), cen.Z() );
185   MainAxis.SetCoord (axis.X(),axis.Y(),axis.Z());
186 }
187
188
189     void  IGESGeom_ConicArc::ComputedDefinition
190   (Standard_Real& Xcen, Standard_Real& Ycen,
191    Standard_Real& Xax,  Standard_Real& Yax,
192    Standard_Real& Rmin, Standard_Real& Rmax) const
193 {
194   Standard_Real a,b,c,d,e,f;
195   //  conic : a*x2 + 2*b*x*y + c*y2 + 2*d*x + 2*e*y + f = 0.
196   Equation (a,b,c,d,e,f);
197   b = b/2.;  d = d/2.;  e = e/2.;    // chgt de variable
198
199   Standard_Real eps = 1.E-08;  // ?? comme ComputedForm
200
201   if (IsFromParabola()) {
202     Rmin = Rmax = -1.;  // rayons : yena pas
203     if ( (Abs(a) <= eps) && (Abs(b) <= eps)) {
204       Xcen = (f*c - e*e) /c /d /2.;
205       Ycen = e/c;
206       Standard_Real focal = -d/c;
207       Xax  = (focal >= 0 ? 1. : -1.);
208       Yax  = 0.;
209       Rmin = Rmax = Abs(focal);
210     } 
211     else {
212       Standard_Real ss = a+c;
213       Standard_Real cc =   - (a*d+b*e) / ss;
214       Standard_Real dd = d + (c*d - b*e) / ss;
215       Standard_Real fc =     (a*e - b*d) / ss;
216       Standard_Real ee = e + fc;
217
218       Standard_Real dn = a*ee - dd*b;
219       Xcen = ( cc*ee + f*b) / dn;
220       Ycen = (-cc*dd - f*a) / dn;
221
222       Standard_Real teta = M_PI/2.;
223       if (Abs(b) > eps) teta = ATan (-a/b);
224       if (fc < 0) teta += M_PI;
225       Xax  = Cos(teta);
226       Yax  = Sin(teta);
227
228       Rmin = Rmax = Abs(fc)/sqrt(a*a+b*b)/2.;
229     }
230   } 
231
232   else {
233     //   -> Conique a centre, cas general
234     //  On utilise les Determinants des matrices :
235     //               | a b d |
236     //  gdet (3x3) = | b c e |  et pdet (2X2) = | a b |
237     //               | d e f |                  | b c |
238     
239     Standard_Real gdet = a*c*f + 2*b*d*e - c*d*d - a*e*e - b*b*f;
240     Standard_Real pdet = a*c - b*b;
241     
242     Xcen = (b*e - c*d) / pdet;
243     Ycen = (b*d - a*e) / pdet;
244
245     Standard_Real term1 = a-c;
246     Standard_Real term2 = 2*b;
247     Standard_Real cos2t;
248     Standard_Real auxil;
249
250     if (Abs(term2)<= eps && Abs(term1)<= eps) {
251       cos2t = 1.;
252       auxil = 0.;
253     }
254     else {
255       Standard_Real t2d = term2/term1; //skl 28.12.2001
256       cos2t = 1./sqrt(1+t2d*t2d);
257       auxil = sqrt (term1*term1 + term2*term2);
258     }
259     
260     Standard_Real cost = sqrt ( (1+cos2t)/2. );
261     Standard_Real sint = sqrt ( (1-cos2t)/2. );
262
263     Standard_Real aprim = (a+c+auxil)/2.;
264     Standard_Real cprim = (a+c-auxil)/2.;
265     
266     term1 = -gdet/(aprim*pdet);
267     term2 = -gdet/(cprim*pdet);
268       
269     if (IsFromEllipse()) {
270       Xax = cost;
271       Yax = sint;
272       Rmin = sqrt ( term1);
273       Rmax = sqrt ( term2);
274       if(Rmax<Rmin){    //skl 28.12.2001
275         Rmax = sqrt ( term1);
276         Rmin = sqrt ( term2);
277       }
278     }
279     else if (term1 <= eps){
280       Xax  = -sint;
281       Yax  =  cost;
282       Rmin = sqrt (-term1);
283       Rmax = sqrt (term2);
284     } 
285     else {
286       Xax  =  cost;
287       Yax  =  sint;
288       Rmin = sqrt (-term2);
289       Rmax = sqrt (term1);
290     }
291   }
292 }