0023300: ShapeFix_EdgeProjAux.cxx: line 283 - identical sub-expression in 'if-clause'
[occt.git] / src / ShapeFix / ShapeFix_EdgeProjAux.cxx
1 // Copyright (c) 1999-2012 OPEN CASCADE SAS
2 //
3 // The content of this file is subject to the Open CASCADE Technology Public
4 // License Version 6.5 (the "License"). You may not use the content of this file
5 // except in compliance with the License. Please obtain a copy of the License
6 // at http://www.opencascade.org and read it completely before using this file.
7 //
8 // The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
9 // main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
10 //
11 // The Original Code and all software distributed under the License is
12 // distributed on an "AS IS" basis, without warranty of any kind, and the
13 // Initial Developer hereby disclaims all such warranties, including without
14 // limitation, any warranties of merchantability, fitness for a particular
15 // purpose or non-infringement. Please see the License for the specific terms
16 // and conditions governing the rights and limitations under the License.
17
18 //:r5 abv 06.04.99: ec_turbine-A.stp, #4313: protect against null curve
19 //    abv 09.04.99  S4136: add parameter preci (to eliminate BRepAPI::Precision)
20 #include <ShapeFix_EdgeProjAux.ixx>
21
22 #include <Standard_ErrorHandler.hxx>
23 #include <Standard_Failure.hxx>
24
25 #include <ElCLib.hxx>
26
27 #include <Adaptor3d_CurveOnSurface.hxx>
28 #include <Geom2dAdaptor_Curve.hxx>
29 #include <Geom2dAdaptor_HCurve.hxx>
30 #include <GeomAdaptor_Surface.hxx>
31 #include <GeomAdaptor_HSurface.hxx>
32
33 #include <Geom2d_Curve.hxx>
34 #include <Geom2d_BSplineCurve.hxx>
35 #include <Geom2d_Line.hxx>
36 #include <Geom_Curve.hxx>
37 #include <Geom_Surface.hxx> 
38 #include <Geom_SphericalSurface.hxx>
39
40 #include <BRep_Tool.hxx>
41 #include <Precision.hxx>
42 #include <TopoDS_Vertex.hxx>
43 #include <TopExp.hxx>
44
45 #include <ShapeAnalysis.hxx>
46 #include <ShapeAnalysis_Curve.hxx>
47 #include <ShapeAnalysis_Edge.hxx>
48 #include <ShapeAnalysis_Surface.hxx>
49
50 #include <Extrema_ExtPC.hxx>
51 #include <gp_Pnt.hxx>
52
53 //=======================================================================
54 //function : ShapeFix_EdgeProjAux
55 //purpose  : 
56 //=======================================================================
57
58 ShapeFix_EdgeProjAux::ShapeFix_EdgeProjAux ()
59 {
60   myFirstDone = myLastDone = Standard_False;
61 }
62
63 //=======================================================================
64 //function : ShapeFix_EdgeProjAux
65 //purpose  : 
66 //=======================================================================
67
68 ShapeFix_EdgeProjAux::ShapeFix_EdgeProjAux (const TopoDS_Face& F,
69                                             const TopoDS_Edge& E)
70 {
71   Init (F, E);
72 }
73
74 //=======================================================================
75 //function : Init
76 //purpose  : 
77 //=======================================================================
78
79 void ShapeFix_EdgeProjAux::Init (const TopoDS_Face& F,
80                                  const TopoDS_Edge& E)
81 {
82   myFace = F;
83   myEdge = E;
84   myFirstDone = myLastDone = Standard_False;
85 }
86
87 //=======================================================================
88 //function : Compute
89 //purpose  : 
90 //=======================================================================
91
92 void ShapeFix_EdgeProjAux::Compute (const Standard_Real preci) 
93 {
94   myFirstDone = myLastDone = Standard_False;
95     
96   // Project Point3d on Surface
97   // TEMPORARY Call ShapeFix_EdgeProjAux
98
99   Init2d(preci);
100   if (IsFirstDone() && IsLastDone()) {
101     Standard_Real U1 = FirstParam();
102     Standard_Real U2 = LastParam();
103     if (U1>=U2) {
104 #ifdef DEBUG
105       cout << "Parametres inverses ... " << endl;
106 #endif
107       Standard_Real tmp = U1;
108       U1 = U2; U2 = tmp;
109     }
110     myFirstParam = U1;
111     myFirstDone  = Standard_True;
112     myLastParam  = U2;
113     myLastDone   = Standard_True;
114   }
115 }
116
117 //=======================================================================
118 //function : IsFirstDone
119 //purpose  : 
120 //=======================================================================
121
122 Standard_Boolean ShapeFix_EdgeProjAux::IsFirstDone() const
123 {
124   return myFirstDone;
125 }
126
127 //=======================================================================
128 //function : IsLastDone
129 //purpose  : 
130 //=======================================================================
131
132 Standard_Boolean ShapeFix_EdgeProjAux::IsLastDone() const
133 {
134   return myLastDone;
135 }
136
137 //=======================================================================
138 //function : FirstParam
139 //purpose  : 
140 //=======================================================================
141
142 Standard_Real ShapeFix_EdgeProjAux::FirstParam() const
143 {
144   return myFirstParam;
145 }
146
147 //=======================================================================
148 //function : LastParam
149 //purpose  : 
150 //=======================================================================
151
152 Standard_Real ShapeFix_EdgeProjAux::LastParam() const
153 {
154   return myLastParam;
155 }
156
157 //=======================================================================
158 //function : IsIso
159 //purpose  : 
160 //=======================================================================
161
162 Standard_Boolean ShapeFix_EdgeProjAux::IsIso (const Handle(Geom2d_Curve)& /*theCurve2d*/)
163 {
164   // Until an ISO is recognized by Adaptor3d_Curve
165 /*  
166   if (theCurve2d->IsKind(STANDARD_TYPE(Geom2d_Line))) {
167     Handle(Geom2d_Line) theLine2d = Handle(Geom2d_Line)::DownCast(theCurve2d);
168     gp_Dir2d theDir2d = theLine2d->Direction();
169     
170     gp_Dir2d dir1(0.,1.);
171     gp_Dir2d dir2(0.,-1.);
172     
173     return (theDir2d.IsEqual(dir1,Precision::Angular()) ||
174             theDir2d.IsEqual(dir2,Precision::Angular()) ||
175             theDir2d.IsNormal(dir1,Precision::Angular()) ||
176             theDir2d.IsNormal(dir2,Precision::Angular()) );
177   }
178 */
179   return Standard_False;
180 }
181
182 // ----------------------------------------------------------------------------
183 // static  : FindParameterWithExt
184 // Purpose : Computes the trimming parameter of Pt1 on COnS
185 // ----------------------------------------------------------------------------
186
187 static Standard_Boolean FindParameterWithExt (const gp_Pnt& Pt1, 
188                                               const Adaptor3d_CurveOnSurface& COnS,
189                                               const Standard_Real Uinf,
190                                               const Standard_Real Usup, 
191                                               const Standard_Real preci,
192                                               Standard_Real& w1)
193 {
194   try {  // et allez donc !
195     OCC_CATCH_SIGNALS
196     Extrema_ExtPC myExtPC (Pt1, COnS, Uinf, Usup, preci);
197
198     //ShapeFix_ExtPCOnS myExtPCOnS1 = 
199     //ShapeFix_ExtPCOnS(Pt1, COnS, Uinf, Usup, preci);
200   
201     if (myExtPC.IsDone()) {
202       Standard_Integer NbExt1 = myExtPC.NbExt();
203       for (Standard_Integer i=1; i<=NbExt1; i++) {
204         if (myExtPC.IsMin(i)) {
205           //Standard_Real dist = myExtPC.Value(i); //szv#4:S4163:12Mar99 debug mode only
206           w1   = myExtPC.Point(i).Parameter();
207         }
208       }
209       return Standard_True;
210     }
211     else return Standard_False;
212   }  // end try
213   catch(Standard_Failure) {
214 #ifdef DEB //:s5
215     cout << "Warning: ShapeFix_EdgeProjAux, FindParameterWithExt(): Exception: ";
216     Standard_Failure::Caught()->Print(cout); cout << endl;
217 #endif
218     return Standard_False;
219   }
220   return Standard_False;  // normalement, on n y passe jamais
221 }
222
223 //=======================================================================
224 //function : Init2d
225 //purpose  : 
226 //=======================================================================
227
228 void ShapeFix_EdgeProjAux::Init2d (const Standard_Real preci) 
229 {
230   Standard_Real cl, cf;
231
232   // Extract Geometries
233   Handle(Geom_Surface) theSurface = BRep_Tool::Surface(myFace);
234   Handle(Geom2d_Curve) theCurve2d = BRep_Tool::CurveOnSurface(myEdge, myFace, cf, cl);
235   if ( theCurve2d.IsNull() ) return; //:r5 abv 6 Apr 99:  ec_turbine-A.stp, #4313
236
237   TopoDS_Vertex V1,V2;
238   TopExp::Vertices(myEdge, V1, V2);
239   gp_Pnt Pt1,Pt2;
240   // pdn 28.12.98: r_39-db.stp #605: use ends of 3d curve instead of vertices 
241   ShapeAnalysis_Edge sae;
242   Standard_Real a,b;
243   Handle(Geom_Curve) C3d;
244   if(sae.Curve3d(myEdge,C3d,a,b,Standard_False)) {
245     Pt1 = C3d->Value(a);
246     Pt2 = C3d->Value(b);
247   } 
248   else {
249     Pt1 = BRep_Tool::Pnt(V1);
250     Pt2 = BRep_Tool::Pnt(V2);
251   }
252 //:S4136  Standard_Real preci = BRepAPI::Precision();
253   //pdn to manage degenerated case
254   if (V1.IsSame(V2)) {
255     Handle(ShapeAnalysis_Surface) stsu = new ShapeAnalysis_Surface (theSurface);
256     gp_Pnt2d aPt1,aPt2;
257     Standard_Real firstpar,lastpar;
258     if (stsu->DegeneratedValues(Pt1,preci,aPt1,aPt2,firstpar,lastpar)){
259       
260       if(theCurve2d->IsKind(STANDARD_TYPE(Geom2d_Line))) {
261         if (aPt1.IsEqual(theCurve2d->Value(firstpar),preci) && 
262             aPt2.IsEqual(theCurve2d->Value(lastpar),preci)){
263           myFirstParam = firstpar;
264           myLastParam  = lastpar;
265           myFirstDone = myLastDone = Standard_True;
266           return;
267         }
268       } 
269 #ifdef DEBUG
270       else cout <<"Other type of deg curve"<<endl;
271 #endif
272       
273     }
274   }
275   
276   Standard_Boolean parU = Standard_False, parV = Standard_False;
277   GeomAdaptor_Surface          SA     = GeomAdaptor_Surface(theSurface);
278   Handle(GeomAdaptor_HSurface) myHSur = new GeomAdaptor_HSurface(SA);
279   
280   cf = theCurve2d->FirstParameter();
281   cl = theCurve2d->LastParameter();
282   //pdn cutting pcurve by suface bounds
283   if (Precision::IsInfinite(cf)||Precision::IsInfinite(cl)) {
284     if(theCurve2d->IsKind(STANDARD_TYPE(Geom2d_Line))) {
285       Standard_Real uf,ul,vf,vl;
286       theSurface->Bounds(uf,ul,vf,vl);
287       if(!Precision::IsInfinite(uf)&&!Precision::IsInfinite(ul)&&
288          !Precision::IsInfinite(vf)&&!Precision::IsInfinite(vl)) {
289         Standard_Real cfi,cli;
290         Handle(Geom2d_Line) lin = Handle(Geom2d_Line)::DownCast(theCurve2d);
291         gp_Pnt2d pnt = lin->Location();
292         gp_Dir2d dir = lin->Direction();
293         if (dir.Y()==0) {
294           parU = Standard_True;
295           cfi = (uf-pnt.X())/dir.X();
296           cli = (ul-pnt.X())/dir.X();
297         } 
298         else if (dir.X()==0) {
299           parV = Standard_True;
300           cfi = (vf-pnt.Y())/dir.Y();
301           cli = (vl-pnt.Y())/dir.Y();
302         } 
303         else {//common case
304           Standard_Real xfi, xli, yfi, yli;
305           xfi = (uf-pnt.X())/dir.X();
306           xli = (ul-pnt.X())/dir.X();
307           yfi = (vf-pnt.Y())/dir.Y();
308           yli = (vl-pnt.Y())/dir.Y();
309           if (dir.X()*dir.Y() > 0) {
310             cfi = (Abs(xli-xfi) < Abs(xli-yfi)? xfi : yfi);
311             cli = (Abs(xfi-xli) < Abs(xfi-yli)? xli : yli);
312           } else {
313             cfi = (Abs(xli-xfi) < Abs(xli-yli)? xfi : yli);
314             cli = (Abs(yli-xli) < Abs(yli-yfi)? xli : yfi);
315           }
316         }
317         if (cfi < cli) { cf = cfi; cl = cli; }
318         else { cf = cli; cl = cfi; }
319       } 
320       else if(!Precision::IsInfinite(uf)&&!Precision::IsInfinite(ul)){
321         Handle(Geom2d_Line) lin = Handle(Geom2d_Line)::DownCast(theCurve2d);
322         gp_Dir2d dir = lin->Direction();
323         if (dir.X()!=0) {
324           if (dir.Y()==0) parU = Standard_True;
325           gp_Pnt2d pnt = lin->Location(); //szv#4:S4163:12Mar99 moved
326           Standard_Real cfi = (uf-pnt.X())/dir.X();
327           Standard_Real cli = (ul-pnt.X())/dir.X();
328           if (cfi < cli) { cf = cfi; cl = cli; }
329           else { cf = cli; cl = cfi; }
330         }
331         else {
332           cf=-10000;
333           cl= 10000;
334         }
335       }
336       else {
337         cf=-10000;
338         cl= 10000;
339         //pdn not cutted by bounds
340 #ifdef DEBUG
341         cout<<"Infinite Surface"<<endl;
342 #endif  
343       }
344     }
345     else {
346       //pdn not linear case not managed
347       cf=-10000;
348       cl= 10000;
349 #ifdef DEBUG     
350       cout<<"Some infinite curve"<<endl;
351 #endif 
352     }
353   }
354   
355   Geom2dAdaptor_Curve          CA     = Geom2dAdaptor_Curve(theCurve2d,cf,cl);
356   Handle(Geom2dAdaptor_HCurve) myHCur = new Geom2dAdaptor_HCurve(CA);
357   
358   Adaptor3d_CurveOnSurface COnS = Adaptor3d_CurveOnSurface(myHCur, myHSur);
359
360   // ----------------------------------------------
361   // --- topological limit == geometric limit ? ---
362   // ----------------------------------------------
363   Standard_Real Uinf = COnS.FirstParameter();
364   Standard_Real Usup = COnS.LastParameter();
365   
366   Standard_Real w1,w2;
367   ShapeAnalysis_Curve sac;
368   gp_Pnt pnt;
369   Standard_Real dist = sac.Project(COnS,Pt1,preci,pnt,w1,Standard_False);
370   if (dist > preci) return;
371   dist = sac.Project(COnS,Pt2,preci,pnt,w2,Standard_False);
372   if (dist > preci) return;
373   
374   myFirstParam = w1;
375   myLastParam  = w2;
376   myFirstDone = myLastDone = Standard_True;
377   if ( myFirstParam == Uinf && myLastParam == Usup ) return;
378   if ( myFirstParam == Usup && myLastParam == Uinf ) {
379     myFirstParam = theCurve2d->ReversedParameter(Usup);
380     myLastParam  = theCurve2d->ReversedParameter(Uinf);
381     theCurve2d->Reverse();
382 #ifdef DEB
383     cout << "Warning: ShapeFix_EdgeProjAux: pcurve reversed" << endl;
384 #endif
385     return;
386   }
387   //:abv 29.08.01: SAT: fix for closed case
388   if ( COnS.Value(Uinf).Distance ( COnS.Value(Usup) ) < Precision::Confusion() ) {
389     // 18.11.2002 SKL OCC630 compare values with tolerance Precision::PConfusion() instead of "=="
390     if ( Abs(myFirstParam-Uinf) < ::Precision::PConfusion() &&
391          Abs(myLastParam-Uinf) < ::Precision::PConfusion() )
392         myLastParam = w2 = Usup;
393     // 18.11.2002 SKL OCC630 compare values with tolerance Precision::PConfusion() instead of "=="
394     else if ( Abs(myFirstParam-Usup) < ::Precision::PConfusion() && 
395               Abs(myLastParam-Usup) < ::Precision::PConfusion() )
396         myFirstParam = w1 = Uinf;
397   }
398    
399   //pdn adjust parameters in periodic case
400   if(parU || parV) {
401     Standard_Real uf,ul,vf,vl;
402     theSurface->Bounds(uf,ul,vf,vl);
403     Standard_Real period = (parU ? ul-uf : vl-vf);
404     w1+=ShapeAnalysis::AdjustToPeriod(w1,0,period);
405     myFirstParam = w1;
406     w2+=ShapeAnalysis::AdjustToPeriod(w2,0,period);
407     myLastParam = w2;
408     Handle(Geom_Curve) C3d1;
409     if(!sae.Curve3d (myEdge, C3d1, cf, cl, Standard_False )) {
410       UpdateParam2d(theCurve2d);
411       return;
412     }
413     gp_Pnt mid = C3d1->Value((cf+cl)/2);
414     Standard_Real wmid;
415     sac.Project(COnS,mid,preci,pnt,wmid,Standard_False);
416     wmid+=ShapeAnalysis::AdjustToPeriod(wmid,0,period);
417     if(w1>w2) {
418       if(w2 > wmid) myFirstParam -= period;
419       else if (w1 > wmid)
420         UpdateParam2d(theCurve2d);
421       else {
422         myLastParam+=period;
423 #ifdef DEBUG
424         cout <<" Added"<<endl;
425 #endif  
426       }
427     }
428     else {
429       if(w1 > wmid) {
430         myLastParam -=period;
431         UpdateParam2d(theCurve2d);
432 #ifdef DEBUG
433         cout <<" Added & Inverted"<<endl;
434 #endif  
435       } else if (w2 < wmid) {
436         myFirstParam += period;
437         UpdateParam2d(theCurve2d);
438       }
439     }
440   }
441   UpdateParam2d(theCurve2d);
442   return;
443 }
444
445 //=======================================================================
446 //function : Init3d
447 //purpose  : 
448 //=======================================================================
449
450 void ShapeFix_EdgeProjAux::Init3d (const Standard_Real preci) 
451 {
452   Standard_Real cl, cf;
453
454   // Extract Geometries
455   Handle(Geom_Surface) theSurface = BRep_Tool::Surface(myFace);
456   Handle(Geom2d_Curve) theCurve2d = BRep_Tool::CurveOnSurface(myEdge, myFace, cf, cl);
457   if ( theCurve2d.IsNull() ) return; //:r5 abv 6 Apr 99:  ec_turbine-A.stp, #4313
458   TopoDS_Vertex V1,V2;
459   
460   V1 = TopExp::FirstVertex(myEdge);
461   V2 = TopExp::LastVertex(myEdge);
462   gp_Pnt Pt1 = BRep_Tool::Pnt(V1);
463   gp_Pnt Pt2 = BRep_Tool::Pnt(V2);
464
465
466   GeomAdaptor_Surface          SA     = GeomAdaptor_Surface(theSurface);  
467   Handle(GeomAdaptor_HSurface) myHSur = new GeomAdaptor_HSurface(SA);
468
469   Geom2dAdaptor_Curve          CA     = Geom2dAdaptor_Curve(theCurve2d);
470   Handle(Geom2dAdaptor_HCurve) myHCur = new Geom2dAdaptor_HCurve(CA);
471
472   Adaptor3d_CurveOnSurface COnS = Adaptor3d_CurveOnSurface(myHCur, myHSur);
473   
474 //:S4136  Standard_Real preci = BRepAPI::Precision();
475   Standard_Real Uinf = theCurve2d->FirstParameter();
476   Standard_Real Usup = theCurve2d->LastParameter();
477
478   // ----------------------------------------------
479   // --- topological limit == geometric limit ? ---
480   // ----------------------------------------------
481   
482   if (theCurve2d->IsKind(STANDARD_TYPE(Geom2d_BoundedCurve))) {
483     
484     gp_Pnt Pdeb = COnS.Value(Uinf);
485     gp_Pnt Pfin = COnS.Value(Usup);
486     
487     //szv#4:S4163:12Mar99 optimized
488     if ( Pdeb.IsEqual(Pt1, preci) && Pfin.IsEqual(Pt2, preci) ) {
489       myFirstParam = Uinf;
490       myLastParam  = Usup;
491       myFirstDone = myLastDone = Standard_True;
492       return;
493     }
494   }
495
496   // ------------------------------------------
497   // --- The CurveOnSurface is not infinite ---
498   // ---          Try with Extrema          ---
499   // ------------------------------------------
500
501   Standard_Real w1 = COnS.FirstParameter();
502   Standard_Real w2 = COnS.LastParameter();
503
504   if ((!Precision::IsInfinite(w1) &&
505        !Precision::IsInfinite(w2) &&
506        theCurve2d->Continuity() != GeomAbs_C0) ||
507       IsIso(theCurve2d))  {
508
509     //szv#4:S4163:12Mar99 optimized
510     if ( FindParameterWithExt(Pt1, COnS, Uinf, Usup, preci, w1) &&
511          FindParameterWithExt(Pt2, COnS, Uinf, Usup, preci, w2) ) {
512       myFirstParam = w1;
513       myLastParam  = w2;
514       UpdateParam2d(theCurve2d);
515       myFirstDone = myLastDone = Standard_True;
516       return;
517     }
518   }
519   myFirstDone = myLastDone = Standard_True;
520 }
521
522 //=======================================================================
523 //function : UpdateParam2d
524 //purpose  : 
525 //=======================================================================
526
527 void ShapeFix_EdgeProjAux::UpdateParam2d (const Handle(Geom2d_Curve)& theCurve2d)
528 {
529   if (myFirstParam < myLastParam)  return;
530
531   Standard_Real cf = theCurve2d->FirstParameter();
532   Standard_Real cl = theCurve2d->LastParameter();
533 //:S4136  Standard_Real preci = BRepAPI::Precision();
534   Standard_Real preci2d = Precision::PConfusion(); //:S4136: Parametric(preci, 0.01);
535
536   // 15.11.2002 PTV OCC966
537   if (ShapeAnalysis_Curve::IsPeriodic(theCurve2d)) {
538     ElCLib::AdjustPeriodic(cf,cl,preci2d,myFirstParam,myLastParam);
539   }
540   else if (theCurve2d->IsClosed()) {
541     //szv#4:S4163:12Mar99 optimized
542     if      ( Abs ( myFirstParam - cl ) <= preci2d ) myFirstParam = cf;
543     else if ( Abs ( myLastParam  - cf ) <= preci2d ) myLastParam  = cl;
544     else {
545 #ifdef DEBUG
546       cout << "Error : curve 2d range crossing non periodic curve origin";
547       cout <<  endl;
548 #endif
549       // add fail result;
550       return;
551     }
552   }
553   // the curve is closed within the 'file' 2D tolerance 
554   else if (theCurve2d->IsKind(STANDARD_TYPE(Geom2d_BSplineCurve))) {
555     Handle(Geom2d_BSplineCurve) aBSpline2d =
556       Handle(Geom2d_BSplineCurve)::DownCast(theCurve2d);
557     if (aBSpline2d->StartPoint().Distance(aBSpline2d->EndPoint()) <= preci2d) {
558       if      ( Abs ( myFirstParam - cl ) <= preci2d ) myFirstParam = cf;
559       else if ( Abs ( myLastParam  - cf ) <= preci2d ) myLastParam  = cl;
560     }
561   }
562   else {
563 #ifdef DEBUG
564     cout << "Warning : non increasing parameters for 2d curve." << endl;
565     cout << "          update parameter 2d uncertain." << endl;
566 #endif
567     Standard_Real tmp1 = myFirstParam, tmp2 = myLastParam;
568     myFirstParam = theCurve2d->ReversedParameter(tmp1);
569     myLastParam = theCurve2d->ReversedParameter(tmp2);
570     theCurve2d->Reverse();
571     //cout<<"Reversed case 2"<<endl;
572   }
573 }