2f49afbb34f157778f729ad626683a133bd9a238
[occt.git] / src / Draft / Draft_Modification_1.cxx
1 // Created on: 1994-12-02
2 // Created by: Jacques GOUSSARD
3 // Copyright (c) 1994-1999 Matra Datavision
4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
5 //
6 // This file is part of Open CASCADE Technology software library.
7 //
8 // This library is free software; you can redistribute it and/or modify it under
9 // the terms of the GNU Lesser General Public License version 2.1 as published
10 // by the Free Software Foundation, with special exception defined in the file
11 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12 // distribution for complete text of the license and disclaimer of any warranty.
13 //
14 // Alternatively, this file may be used under the terms of Open CASCADE
15 // commercial license or contractual agreement.
16
17
18 #include <Adaptor3d_CurveOnSurface.hxx>
19 #include <Adaptor3d_HCurveOnSurface.hxx>
20 #include <Adaptor3d_SurfaceOfLinearExtrusion.hxx>
21 #include <Approx_CurveOnSurface.hxx>
22 #include <BRep_Builder.hxx>
23 #include <BRep_Tool.hxx>
24 #include <BRepExtrema_ExtCF.hxx>
25 #include <BRepExtrema_ExtPC.hxx>
26 #include <BRepLib_MakeFace.hxx>
27 #include <BRepTools.hxx>
28 #include <Draft_DataMapIteratorOfDataMapOfEdgeEdgeInfo.hxx>
29 #include <Draft_DataMapIteratorOfDataMapOfFaceFaceInfo.hxx>
30 #include <Draft_DataMapIteratorOfDataMapOfVertexVertexInfo.hxx>
31 #include <Draft_EdgeInfo.hxx>
32 #include <Draft_FaceInfo.hxx>
33 #include <Draft_Modification.hxx>
34 #include <Draft_VertexInfo.hxx>
35 #include <ElCLib.hxx>
36 #include <ElSLib.hxx>
37 #include <Extrema_ExtPC.hxx>
38 #include <Geom2d_BezierCurve.hxx>
39 #include <Geom2d_BSplineCurve.hxx>
40 #include <Geom2d_Curve.hxx>
41 #include <Geom2d_Line.hxx>
42 #include <Geom2d_TrimmedCurve.hxx>
43 #include <Geom2dAdaptor_HCurve.hxx>
44 #include <Geom2dAPI_ProjectPointOnCurve.hxx>
45 #include <Geom2dConvert.hxx>
46 #include <Geom2dConvert_CompCurveToBSplineCurve.hxx>
47 #include <Geom_BSplineCurve.hxx>
48 #include <Geom_Circle.hxx>
49 #include <Geom_Conic.hxx>
50 #include <Geom_ConicalSurface.hxx>
51 #include <Geom_Curve.hxx>
52 #include <Geom_CylindricalSurface.hxx>
53 #include <Geom_ElementarySurface.hxx>
54 #include <Geom_Ellipse.hxx>
55 #include <Geom_Hyperbola.hxx>
56 #include <Geom_Line.hxx>
57 #include <Geom_Parabola.hxx>
58 #include <Geom_Plane.hxx>
59 #include <Geom_RectangularTrimmedSurface.hxx>
60 #include <Geom_Surface.hxx>
61 #include <Geom_SurfaceOfLinearExtrusion.hxx>
62 #include <Geom_TrimmedCurve.hxx>
63 #include <GeomAdaptor_Curve.hxx>
64 #include <GeomAdaptor_HCurve.hxx>
65 #include <GeomAdaptor_HSurface.hxx>
66 #include <GeomAdaptor_Surface.hxx>
67 #include <GeomAPI_ProjectPointOnCurve.hxx>
68 #include <GeomAPI_ProjectPointOnSurf.hxx>
69 #include <GeomConvert_CompCurveToBSplineCurve.hxx>
70 #include <GeomInt_IntSS.hxx>
71 #include <GeomProjLib.hxx>
72 #include <gp.hxx>
73 #include <gp_Circ.hxx>
74 #include <gp_Dir.hxx>
75 #include <gp_Elips.hxx>
76 #include <gp_Hypr.hxx>
77 #include <gp_Lin.hxx>
78 #include <gp_Parab.hxx>
79 #include <gp_Pln.hxx>
80 #include <gp_Pnt.hxx>
81 #include <gp_Vec.hxx>
82 #include <IntAna_IntConicQuad.hxx>
83 #include <IntAna_QuadQuadGeo.hxx>
84 #include <IntCurveSurface_HInter.hxx>
85 #include <IntCurveSurface_IntersectionPoint.hxx>
86 #include <Precision.hxx>
87 #include <ProjLib_CompProjectedCurve.hxx>
88 #include <ProjLib_HCompProjectedCurve.hxx>
89 #include <Standard_ConstructionError.hxx>
90 #include <Standard_DomainError.hxx>
91 #include <Standard_Failure.hxx>
92 #include <Standard_NoSuchObject.hxx>
93 #include <Standard_NotImplemented.hxx>
94 #include <StdFail_NotDone.hxx>
95 #include <TColgp_Array1OfPnt2d.hxx>
96 #include <TopExp.hxx>
97 #include <TopExp_Explorer.hxx>
98 #include <TopLoc_Location.hxx>
99 #include <TopoDS.hxx>
100 #include <TopoDS_Edge.hxx>
101 #include <TopoDS_Face.hxx>
102 #include <TopoDS_Shape.hxx>
103 #include <TopoDS_Vertex.hxx>
104 #include <TopTools_ListIteratorOfListOfShape.hxx>
105 #include <TopTools_MapIteratorOfMapOfShape.hxx>
106 #include <TopTools_MapOfShape.hxx>
107
108 static Standard_Boolean Choose(const Draft_DataMapOfFaceFaceInfo&,
109                                Draft_DataMapOfEdgeEdgeInfo&,
110                                const TopoDS_Vertex&,
111                                Draft_VertexInfo&,
112                                GeomAdaptor_Curve&,
113                                GeomAdaptor_Surface&);
114
115 static Standard_Real Parameter(const Handle(Geom_Curve)&,
116                                const gp_Pnt&,
117                                Standard_Integer&);
118
119 static Standard_Real SmartParameter(Draft_EdgeInfo&,
120                                     const Standard_Real EdgeTol,
121                                     const gp_Pnt&,
122                                     const Standard_Integer,
123                                     const Handle(Geom_Surface)&,
124                                     const Handle(Geom_Surface)&);
125
126 static TopAbs_Orientation Orientation(const TopoDS_Shape&,
127                                       const TopoDS_Face&);
128
129 static Standard_Boolean FindRotation(const gp_Pln&,
130                                      const TopAbs_Orientation,
131                                      const gp_Dir&,
132                                      const Standard_Real,
133                                      const gp_Pln&,
134                                      gp_Ax1&,
135                                      Standard_Real&);
136
137
138 //=======================================================================
139 //function : InternalAdd
140 //purpose  : 
141 //=======================================================================
142
143 Standard_Boolean Draft_Modification::InternalAdd(const TopoDS_Face& F,
144                                                  const gp_Dir& Direction,
145                                                  const Standard_Real Angle,
146                                                  const gp_Pln& NeutralPlane,
147                                                  const Standard_Boolean Flag)
148 {
149
150   if (myFMap.IsBound(F)) {
151     return (badShape.IsNull());
152   }
153
154   TopAbs_Orientation oris = Orientation(myShape,F);
155   TopLoc_Location Lo;
156   //gp_Dir NewDirection = Direction;
157   //Standard_Real NewAngle = Angle;
158   Handle(Geom_Surface) S = BRep_Tool::Surface(F,Lo);
159   S = Handle(Geom_Surface)::DownCast(S->Transformed(Lo.Transformation()));
160   if (S->DynamicType() == STANDARD_TYPE(Geom_RectangularTrimmedSurface)) {
161     S = Handle(Geom_RectangularTrimmedSurface)::DownCast(S)->BasisSurface();
162   }
163   Handle(Geom_Surface) NewS;
164   Handle(Geom_Circle) theCircle;
165
166   Standard_Boolean postponed = (Flag == Standard_False);
167   if (postponed) {
168     Handle(Standard_Type) typS = S->DynamicType();
169     if (typS == STANDARD_TYPE(Geom_CylindricalSurface) ||
170         typS == STANDARD_TYPE(Geom_SurfaceOfLinearExtrusion)) {
171       gp_Circ Cir;
172       if (typS == STANDARD_TYPE(Geom_CylindricalSurface)) {
173         gp_Cylinder cyl = 
174           Handle(Geom_CylindricalSurface)::DownCast(S)->Cylinder();
175         gp_Ax1 axcyl = cyl.Axis();
176         Cir = ElSLib::CylinderVIso( cyl.Position(), cyl.Radius(), 0.);
177         gp_Vec VV(cyl.Location(),NeutralPlane.Location());
178         Cir.Translate(VV.Dot(axcyl.Direction())*axcyl.Direction());
179       }
180       else {
181         Handle(Geom_Curve) Cbas = 
182           Handle(Geom_SurfaceOfLinearExtrusion)::DownCast(S)->BasisCurve();
183         gp_Dir theDirextr = 
184           Handle(Geom_SurfaceOfLinearExtrusion)::DownCast(S)->Direction();
185
186         if (Cbas->IsKind(STANDARD_TYPE(Geom_TrimmedCurve))) {
187           Cbas = Handle(Geom_TrimmedCurve)::DownCast(Cbas)->BasisCurve();
188         }
189         if (Cbas->IsKind(STANDARD_TYPE(Geom_Circle))) {
190           Cir = Handle(Geom_Circle)::DownCast(Cbas)->Circ();
191           gp_Dir dircir = Cir.Axis().Direction();
192           if (!Direction.IsParallel(dircir,Precision::Angular())) {
193             badShape = F;
194             errStat = Draft_FaceRecomputation;
195             return Standard_False;
196           }
197         }
198         else {
199           badShape = F;
200           errStat = Draft_FaceRecomputation;
201           return Standard_False;
202         }
203
204         gp_Ax3 Axis = NeutralPlane.Position();
205         Standard_Real L =
206           gp_Vec(Cir.Location(),Axis.Location()).
207             Dot(Axis.Direction());
208         Standard_Real Cos = theDirextr.Dot(Axis.Direction());
209         gp_Vec VV = ( L / Cos) * theDirextr;
210         Cir.Translate(VV);
211       }
212
213       theCircle = new Geom_Circle(Cir);
214
215     }
216     else {
217       postponed = Standard_False;
218     }
219   }
220
221
222   if (!postponed) {
223     NewS = NewSurface(S,oris,Direction,Angle,NeutralPlane); 
224     if (NewS.IsNull()) {
225       badShape = F;
226       errStat = Draft_FaceRecomputation;
227       return Standard_False;
228     }
229     // To avoid some problems with infinite restrictions
230     const Handle(Standard_Type)& typs = NewS->DynamicType();
231     if (typs == STANDARD_TYPE(Geom_CylindricalSurface) ||
232         typs == STANDARD_TYPE(Geom_ConicalSurface)) {
233       Standard_Real umin,umax,vmin,vmax;
234       BRepTools::UVBounds(F,umin,umax,vmin,vmax);
235       if (!Precision::IsNegativeInfinite(vmin) &&
236           !Precision::IsPositiveInfinite(vmax)) {
237         Standard_Real deltav = 10.*(vmax-vmin);
238         if(typs == STANDARD_TYPE(Geom_CylindricalSurface)) {
239           vmin = vmin - deltav;
240           vmax = vmax + deltav;
241         }
242         else {
243           gp_Cone Co = Handle(Geom_ConicalSurface)::DownCast(NewS)->Cone();
244           Standard_Real Vapex = - Co.RefRadius()/Sin(Co.SemiAngle());
245           if (vmin < Vapex) { // vmax should not exceed Vapex
246             if (vmax + deltav > Vapex) {
247               vmax = Vapex;
248               vmin = vmin - 10.*(vmax - vmin);
249               // JAG debug to avoid apex
250               vmax = vmax-Precision::Confusion();
251             }
252             else {
253               vmin = vmin - deltav;
254               vmax = vmax + deltav;
255             }
256           }
257           else { // Vapex <= vmin < vmax
258             if (vmin - deltav < Vapex) {
259               vmin = Vapex;
260               vmax = vmax + 10.*(vmax - vmin);
261               // JAG debug to avoid apex
262               vmin = vmin+Precision::Confusion();
263             }
264             else {
265               vmin = vmin - deltav;
266               vmax = vmax + deltav;
267             }
268           }
269         }
270         NewS = new Geom_RectangularTrimmedSurface(NewS,0.,2.*M_PI,vmin,vmax);
271       }
272     }
273   }
274
275   if (postponed || S != NewS) {
276     Draft_FaceInfo FI(NewS,Standard_True);
277     FI.RootFace(curFace);
278     myFMap.Bind(F,FI);
279     if (postponed) {
280       myFMap(F).ChangeCurve() = theCircle;
281     }
282   }    
283
284   TopExp_Explorer expl(F,TopAbs_EDGE);
285   TopTools_MapOfShape MapOfE;
286   while (expl.More() && badShape.IsNull()) {
287     const TopoDS_Edge& edg = TopoDS::Edge(expl.Current());
288     if (!myEMap.IsBound(edg)) {
289       Standard_Boolean addedg  = Standard_False;
290       Standard_Boolean addface = Standard_False;
291       TopoDS_Face OtherF;
292       //        if (BRep_Tool::IsClosed(edg,F)) {
293       if (BRepTools::IsReallyClosed(edg,F)) {
294         addedg = Standard_True;
295         addface = Standard_False;
296       }
297       else {
298         // Find the other face containing the edge.
299         TopTools_ListIteratorOfListOfShape it;
300         it.Initialize(myEFMap.FindFromKey(edg));
301         Standard_Integer nbother = 0;
302         while (it.More()) {
303           if (!it.Value().IsSame(F)) {
304             if (OtherF.IsNull()) {
305               OtherF = TopoDS::Face(it.Value());
306             }
307             nbother++;
308           }
309           it.Next();
310         }         
311         if (nbother >=2) {
312           badShape = edg;
313           errStat = Draft_EdgeRecomputation;
314         }
315         else if (! OtherF.IsNull() && 
316                  BRep_Tool::Continuity(edg,F,OtherF) >= GeomAbs_G1) {
317           addface= Standard_True;
318           addedg = Standard_True;
319         }
320         else if (nbother == 0) {
321           //        badShape = F;
322         }
323       }
324       if (addedg) {
325         if (postponed) {
326           myFMap(F).Add(OtherF);
327         }
328         Standard_Real f,l;
329         TopLoc_Location L;
330         Handle(Geom_Curve) C = BRep_Tool::Curve(edg,L,f,l);
331         C = Handle(Geom_Curve)::DownCast(C->Transformed(L));
332         if (C->DynamicType() == STANDARD_TYPE(Geom_TrimmedCurve)) {
333           C = Handle(Geom_TrimmedCurve)::DownCast(C)->BasisCurve();
334         }
335         Handle(Geom_Curve) NewC;
336         Draft_EdgeInfo EInf(Standard_True);     
337         if(postponed) {
338           EInf.Add(F);
339           EInf.Add(OtherF);
340
341           // find fixed point 
342           Handle(Geom_Line) aLocalGeom = Handle(Geom_Line)::DownCast(C);
343           if (aLocalGeom.IsNull()) {
344             badShape = edg;
345             errStat = Draft_EdgeRecomputation;
346           }
347           else {
348             gp_Lin lin = aLocalGeom->Lin();
349             IntAna_IntConicQuad ilipl(lin,NeutralPlane,Precision::Angular());
350             if (ilipl.IsDone() && ilipl.NbPoints() != 0){
351               EInf.Tangent(ilipl.Point(1));
352             }
353             else {
354               badShape = edg;
355               errStat = Draft_EdgeRecomputation;
356             }
357           }
358         }
359         else {
360           NewC = NewCurve(C,S,oris,Direction,Angle,NeutralPlane, Flag);
361           if (NewC.IsNull()) {
362             badShape = edg;
363             errStat = Draft_EdgeRecomputation;
364           }
365         }
366
367         Handle(Geom_TrimmedCurve) T = Handle(Geom_TrimmedCurve)::DownCast(NewC);
368         if (!T.IsNull()) NewC = T->BasisCurve();
369         EInf.ChangeGeometry() = NewC;
370
371         EInf.RootFace(curFace);
372         myEMap.Bind(edg,EInf);
373         MapOfE.Add(edg);
374         if (addface) {
375           Standard_Boolean Fl = Flag;
376           Handle(Geom_Surface) alocalSurface = BRep_Tool::Surface(OtherF,Lo);
377           if (alocalSurface->DynamicType() == 
378               STANDARD_TYPE(Geom_RectangularTrimmedSurface)) {
379             alocalSurface = Handle(Geom_RectangularTrimmedSurface)::
380               DownCast(alocalSurface)->BasisSurface();
381           }
382           Handle(Standard_Type) typS = alocalSurface->DynamicType();
383           if (typS == STANDARD_TYPE(Geom_CylindricalSurface) || 
384               typS == STANDARD_TYPE(Geom_SurfaceOfLinearExtrusion)) {
385             if ( myFMap.IsBound(F)) {
386               if ( Flag == Standard_False && !postponed) {
387                 myFMap.UnBind(F);
388                 TopTools_MapIteratorOfMapOfShape itm(MapOfE);
389                 for ( ; itm.More(); itm.Next())
390                   myEMap.UnBind(TopoDS::Edge(itm.Key()));
391               }
392             }
393           } 
394           InternalAdd(OtherF,Direction,Angle,NeutralPlane, Fl);
395         }
396       }
397     }
398     expl.Next();
399   }
400   return (badShape.IsNull());
401 }
402     
403
404 //=======================================================================
405 //function : Propagate
406 //purpose  : 
407 //=======================================================================
408
409 Standard_Boolean Draft_Modification::Propagate ()
410 {
411
412   if (!badShape.IsNull()) return Standard_False;
413
414   Draft_DataMapIteratorOfDataMapOfFaceFaceInfo itf(myFMap);
415
416   // Set all edges and vertices of modified faces
417   TopoDS_Face F;
418   TopoDS_Edge E;
419   TopoDS_Vertex V;
420   TopExp_Explorer editer;
421   TopExp_Explorer vtiter;
422
423   while (itf.More()) {
424     const TopoDS_Face& Fc = itf.Key();
425
426     // Exploration of the edges of the face
427     editer.Init(Fc,TopAbs_EDGE);
428     while (editer.More()) {
429       E = TopoDS::Edge(editer.Current());
430
431       if (!myEMap.IsBound(E)) {
432         Draft_EdgeInfo EInf(Standard_True);     
433         myEMap.Bind(E,EInf);
434       }
435       myEMap(E).Add(Fc);
436
437       // Exploration of the vertices of the edge
438       vtiter.Init(E,TopAbs_VERTEX);
439       while (vtiter.More()) {
440         V = TopoDS::Vertex(vtiter.Current());
441         if (!myVMap.IsBound(V)) {
442           Draft_VertexInfo VInf;
443           myVMap.Bind(V,VInf);
444         }
445
446         myVMap(V).Add(E);
447         myVMap(V).ChangeParameter(E) = BRep_Tool::Parameter(V, E);
448         vtiter.Next();
449       }
450       editer.Next();
451     }
452     itf.Next();
453   }
454
455
456   TopExp_Explorer anc;
457   Standard_Boolean found;
458
459   // Set edges containing modified vertices.
460
461   Draft_DataMapIteratorOfDataMapOfVertexVertexInfo itv(myVMap);
462
463   while (itv.More()) {
464     const TopoDS_Vertex& Vt = itv.Key();
465
466     // Exploration of the ancestors of the vertex
467     anc.Init(myShape,TopAbs_EDGE);
468
469     while (anc.More()) {
470       E = TopoDS::Edge(anc.Current());
471       vtiter.Init(E,TopAbs_VERTEX);
472       found = Standard_False;
473       while (vtiter.More()) {
474         if (Vt.IsSame(TopoDS::Vertex(vtiter.Current()))) {
475           found = Standard_True;
476           break;
477         }
478         vtiter.Next();
479       }
480       if (found) {
481         if (!myEMap.IsBound(E)) {
482           Draft_EdgeInfo EInf(Standard_False);
483           myEMap.Bind(E,EInf);
484         }
485         myVMap(Vt).Add(E);
486         myVMap(Vt).ChangeParameter(E) = BRep_Tool::Parameter(Vt, E);
487       }
488       anc.Next();
489     }
490     itv.Next();
491   }
492
493
494   // Set faces containing modified edges
495
496   Draft_DataMapIteratorOfDataMapOfEdgeEdgeInfo ite(myEMap);
497
498   while (ite.More()) {
499     const TopoDS_Edge& Ed = ite.Key();
500     TopTools_ListIteratorOfListOfShape it;
501     for (it.Initialize(myEFMap.FindFromKey(Ed)); it.More(); it.Next()) {
502       F = TopoDS::Face(it.Value());
503       if (!myFMap.IsBound(F)) {
504         TopLoc_Location L;
505         Handle(Geom_Surface) S = BRep_Tool::Surface(F,L);
506         Handle(Geom_Surface) NewS = 
507           Handle(Geom_Surface)::DownCast(S->Transformed(L.Transformation()));
508         
509         const Handle(Standard_Type)& typs = S->DynamicType();
510         if (typs == STANDARD_TYPE(Geom_CylindricalSurface) ||
511             typs == STANDARD_TYPE(Geom_ConicalSurface)) {
512           Standard_Real umin,umax,vmin,vmax;
513           BRepTools::UVBounds(F,umin,umax,vmin,vmax);
514           if (!Precision::IsNegativeInfinite(vmin) &&
515                 !Precision::IsPositiveInfinite(vmax)) {
516             Standard_Real deltav = 10.*(vmax-vmin);
517             vmin = vmin - deltav;
518             vmax = vmax + deltav;
519             NewS = 
520               new Geom_RectangularTrimmedSurface(NewS,0.,2.*M_PI,vmin,vmax);
521           }
522         }
523         
524         Draft_FaceInfo FInf(NewS,Standard_False);
525         myFMap.Bind(F,FInf);
526       }
527       myEMap(Ed).Add(F);
528     }
529     ite.Next();
530   }
531
532   //  Try to add faces for free borders...
533   // JAG 09.11.95
534   ite.Initialize(myEMap);
535   for (; ite.More(); ite.Next()) {
536     Draft_EdgeInfo& Einf = myEMap(ite.Key());
537     if (Einf.NewGeometry() && Einf.Geometry().IsNull() && 
538         Einf.SecondFace().IsNull()) {
539       
540       TopLoc_Location Loc;
541       Handle(Geom_Surface) S1 = BRep_Tool::Surface(Einf.FirstFace(),Loc);
542       S1 = Handle(Geom_Surface)::
543         DownCast(S1->Transformed(Loc.Transformation()));
544       Handle(Geom_Surface) S2;
545       
546       Standard_Real f,l;
547       Handle(Geom_Curve) C = BRep_Tool::Curve(ite.Key(),Loc,f,l);
548       C = Handle(Geom_Curve)::DownCast(C->Transformed(Loc.Transformation()));
549       if (C->DynamicType() == STANDARD_TYPE(Geom_TrimmedCurve)) {
550         C = Handle(Geom_TrimmedCurve)::DownCast(C)->BasisCurve();
551       }
552       if (!S1->IsKind(STANDARD_TYPE(Geom_Plane))) {
553         if (C->IsKind(STANDARD_TYPE(Geom_Conic))) {
554           gp_Ax3 thePl(Handle(Geom_Conic)::DownCast(C)->Position());
555           S2 = new Geom_Plane(thePl);
556         }
557         else if (C->DynamicType() == STANDARD_TYPE(Geom_Line)) {
558           gp_Ax1 axis;
559           if (S1->DynamicType()== STANDARD_TYPE(Geom_RectangularTrimmedSurface)) {
560             axis = Handle(Geom_ElementarySurface)::DownCast
561               (Handle(Geom_RectangularTrimmedSurface)::DownCast(S1)->
562                BasisSurface())->Axis();
563           }
564           else {
565             axis = Handle(Geom_ElementarySurface)::DownCast(S1)->Axis();
566           }
567           gp_Vec they(axis.Location(), C->Value(0.));
568           gp_Dir axz(axis.Direction().Crossed(they));
569           S2=new Geom_Plane(gp_Ax3(axis.Location(),axz,axis.Direction()));
570
571         }
572         else {
573           badShape = TopoDS::Edge(ite.Key());
574           errStat = Draft_EdgeRecomputation;
575           break; // leave from for
576         }
577       }
578       else { // on the plane
579         Draft_DataMapIteratorOfDataMapOfVertexVertexInfo anewitv(myVMap);
580         while (anewitv.More()) {
581           Draft_VertexInfo& Vinf = myVMap(anewitv.Key());
582           for (Vinf.InitEdgeIterator();Vinf.MoreEdge();Vinf.NextEdge()) {
583             if (Vinf.Edge().IsSame(ite.Key())) {
584               break;
585             }
586           }
587           if (Vinf.MoreEdge()) {
588             for (Vinf.InitEdgeIterator();Vinf.MoreEdge();Vinf.NextEdge()) {
589               const TopoDS_Edge& edg = Vinf.Edge();
590               if (!edg.IsSame(ite.Key())) {
591                 if (!myEMap(edg).FirstFace().IsSame(Einf.FirstFace()) &&
592                     (myEMap(edg).SecondFace().IsNull() || 
593                      !myEMap(edg).SecondFace().IsSame(Einf.FirstFace()))) {
594                   break;
595                 }
596               }
597             }
598             if (Vinf.MoreEdge()) {
599               Handle(Geom_Curve) C2 = BRep_Tool::Curve(Vinf.Edge(), Loc,f,l);
600               Handle(GeomAdaptor_HCurve) HCur;
601               gp_Vec Direc;
602               C2 = Handle(Geom_Curve)::DownCast
603                 (C2->Transformed(Loc.Transformation()));
604               if (C2->DynamicType() == STANDARD_TYPE(Geom_TrimmedCurve)) {
605                 C2 = Handle(Geom_TrimmedCurve)::DownCast(C2)->BasisCurve();
606               }
607               if (C->DynamicType() == STANDARD_TYPE(Geom_Line)) {
608                 Direc = Handle(Geom_Line)::DownCast(C)->Lin().Direction();
609                 HCur = new GeomAdaptor_HCurve(C2);
610
611               }
612               else if (C2->DynamicType() == STANDARD_TYPE(Geom_Line)) {
613                 Direc = Handle(Geom_Line)::DownCast(C2)->Lin().Direction();
614                 HCur = new GeomAdaptor_HCurve(C);
615               }
616               else {
617                 badShape = TopoDS::Edge(ite.Key());
618                 errStat = Draft_EdgeRecomputation;
619                 break; // leave from while
620               }
621               Adaptor3d_SurfaceOfLinearExtrusion SLE(HCur,Direc);
622               switch(SLE.GetType()){
623
624               case GeomAbs_Plane :
625                 {
626                   S2 = new Geom_Plane(SLE.Plane());
627                 }
628                 break;
629               case GeomAbs_Cylinder :
630                 {
631                   S2 =   new Geom_CylindricalSurface(SLE.Cylinder());
632                 }
633                 break;
634               default :
635                 {
636                   S2 = new Geom_SurfaceOfLinearExtrusion(HCur->ChangeCurve().
637                                                          Curve(),
638                                                          Direc);
639                 }
640                 break;
641               }
642               
643             }
644             else {
645               badShape = TopoDS::Edge(ite.Key());
646               errStat = Draft_EdgeRecomputation;
647               break; // leave from while
648             }
649             break;
650           }
651           anewitv.Next();
652         }
653       }
654
655       if (badShape.IsNull()) {
656         BRep_Builder B;
657         TopoDS_Face TheNewFace;
658         B.MakeFace(TheNewFace,S2,Precision::Confusion());
659         Einf.Add(TheNewFace);
660         Draft_FaceInfo FI(S2,Standard_False);
661         myFMap.Bind(TheNewFace,FI);
662       }
663       else {
664         break; // leave from for
665       }
666       // Fin JAG 09.11.95
667     }
668   }
669   return (badShape.IsNull());
670 }
671
672
673
674
675 //=======================================================================
676 //function : Perform
677 //purpose  : 
678 //=======================================================================
679
680 void Draft_Modification::Perform ()
681 {
682   if (!badShape.IsNull())  Standard_ConstructionError::Raise();
683
684   if (!myComp) {
685     myComp = Standard_True;
686     if (!Propagate()) {
687       return;
688     }
689
690     // Calculate eventual faces
691
692     Draft_DataMapIteratorOfDataMapOfFaceFaceInfo itf(myFMap);
693     while (itf.More()) {
694       Draft_FaceInfo& Finf = myFMap(itf.Key());
695       if (Finf.NewGeometry() && Finf.Geometry().IsNull()) {
696         const TopoDS_Face& F1 = Finf.FirstFace();
697         const TopoDS_Face& F2 = Finf.SecondFace();
698
699         if (F1.IsNull() || F2.IsNull()) {
700           errStat = Draft_FaceRecomputation;
701           badShape = TopoDS::Face(itf.Key());
702           return;
703         }
704         Handle(Geom_Surface) S1 = myFMap(F1).Geometry();
705         Handle(Geom_Surface) S2 = myFMap(F2).Geometry();
706         if (S1.IsNull() || S2.IsNull()) {
707           errStat = Draft_FaceRecomputation;
708           badShape = TopoDS::Face(itf.Key());
709           return;
710         }
711         if (S1->DynamicType() == STANDARD_TYPE(Geom_RectangularTrimmedSurface)) {
712           S1 = Handle(Geom_RectangularTrimmedSurface)::
713             DownCast(S1)->BasisSurface();
714         }
715         if (S2->DynamicType() == STANDARD_TYPE(Geom_RectangularTrimmedSurface)) {
716           S2 = Handle(Geom_RectangularTrimmedSurface)::
717             DownCast(S2)->BasisSurface();
718         }
719         Handle(Geom_Plane) P1 = Handle(Geom_Plane)::DownCast(S1);
720         Handle(Geom_Plane) P2 = Handle(Geom_Plane)::DownCast(S2);
721         if (P1.IsNull() || P2.IsNull()) {
722           errStat = Draft_FaceRecomputation;
723           badShape = TopoDS::Face(itf.Key());
724           return;
725         }
726         gp_Pln pp1 = P1->Pln();
727         gp_Pln pp2 = P2->Pln();
728         IntAna_QuadQuadGeo i2p(pp1,pp2,
729                                Precision::Angular(),Precision::Confusion());
730         if (!i2p.IsDone() || i2p.TypeInter() != IntAna_Line) {
731           errStat = Draft_FaceRecomputation;
732           badShape = TopoDS::Face(itf.Key());
733           return;
734         }
735
736         gp_Dir extrdir = i2p.Line(1).Direction();
737
738         // Preserve the same direction as the base face
739         Handle(Geom_Surface) RefSurf = 
740           BRep_Tool::Surface(TopoDS::Face(itf.Key()));
741         if (RefSurf->DynamicType() == STANDARD_TYPE(Geom_RectangularTrimmedSurface)) {
742           RefSurf = 
743             Handle(Geom_RectangularTrimmedSurface)::DownCast(RefSurf)
744               ->BasisSurface();
745         }
746         gp_Dir DirRef;
747
748         if ( RefSurf->DynamicType() == STANDARD_TYPE(Geom_CylindricalSurface)) {
749           gp_Ax3 AxeRef = 
750             Handle(Geom_CylindricalSurface)::DownCast(RefSurf)
751               ->Cylinder().Position();
752           DirRef = AxeRef.Direction();
753         }
754         else if (RefSurf->DynamicType() == 
755                  STANDARD_TYPE(Geom_SurfaceOfLinearExtrusion)) {
756           DirRef = 
757             Handle(Geom_SurfaceOfLinearExtrusion)::DownCast(RefSurf)
758               ->Direction();
759         }
760         
761         if (extrdir.Dot(DirRef) < 0.) extrdir.Reverse();
762
763         // it is possible to accelerate speed by storing the info during
764         // InternalAdd --> modification of FaceInfo to preserve the circle
765
766         Handle(Geom_Circle) CCir = 
767           Handle(Geom_Circle)::DownCast(Finf.Curve());
768         Handle(Geom_Surface) NewS = 
769           new Geom_SurfaceOfLinearExtrusion(CCir, extrdir);    
770
771         Standard_Real umin, umax, vmin, vmax;
772         BRepTools::UVBounds(TopoDS::Face(itf.Key()),umin,umax,vmin,vmax);
773         if (!Precision::IsNegativeInfinite(vmin) &&
774             !Precision::IsPositiveInfinite(vmax)) {
775           Standard_Real deltav = 2.*(vmax-vmin);
776           vmin = vmin - deltav;
777           vmax = vmax + deltav;
778         }
779
780         // very temporary
781         else {
782           vmax = 300;
783           vmin = -300;
784         }
785
786         NewS = new Geom_RectangularTrimmedSurface(NewS,0.,1.9*M_PI,vmin,vmax);
787         Finf.ChangeGeometry() = NewS;
788       }
789       itf.Next();
790     }
791     
792     // Calculate new edges.
793
794     Handle(Geom_Surface) S1,S2;
795     Handle(Geom_Curve) C, newC;
796     Standard_Real f,l;
797     TopLoc_Location L;
798     Draft_DataMapIteratorOfDataMapOfEdgeEdgeInfo ite(myEMap);
799
800     while (ite.More()) {
801       Draft_EdgeInfo& Einf = myEMap(ite.Key());
802
803       const TopoDS_Edge& theEdge = TopoDS::Edge(ite.Key());
804
805       C = BRep_Tool::Curve(theEdge,L,f,l);
806       C = Handle(Geom_Curve)::DownCast(C->Transformed(L.Transformation()));
807
808       if (Einf.NewGeometry() && Einf.Geometry().IsNull()) {
809         gp_Pnt ptfixe;
810         if (!Einf.IsTangent(ptfixe)) {
811           const TopoDS_Face& FirstFace = Einf.FirstFace();
812           const TopoDS_Face& SecondFace = Einf.SecondFace();
813
814           S1 = myFMap(FirstFace).Geometry();
815           S2 = myFMap(SecondFace).Geometry();
816
817           Standard_Integer detrompeur = 0;
818
819           // Return FirstVertex and the tangent at this point.
820           TopoDS_Vertex FV = TopExp::FirstVertex(theEdge);
821           TopoDS_Vertex LV = TopExp::LastVertex(theEdge);
822           Standard_Real pmin = 0.;
823           Standard_Real prmfv = BRep_Tool::Parameter(FV,ite.Key());
824           Standard_Real prmlv = BRep_Tool::Parameter(LV,ite.Key());
825           gp_Pnt pfv, plv;
826           gp_Vec d1fv,d1lv, newd1;
827           C->D1(prmfv,pfv,d1fv);
828           C->D1(prmlv,plv,d1lv);
829
830           Standard_Real TolF1 = BRep_Tool::Tolerance (FirstFace);
831           Standard_Real TolF2 = BRep_Tool::Tolerance (SecondFace);
832
833           //Pass the tolerance of the face to project
834           GeomAPI_ProjectPointOnSurf proj1 (pfv, S1, TolF1);
835           GeomAPI_ProjectPointOnSurf proj2 (plv, S1, TolF1);
836           GeomAPI_ProjectPointOnSurf proj3 (pfv, S2, TolF2);
837           GeomAPI_ProjectPointOnSurf proj4 (plv, S2, TolF2);
838
839           if (proj1.IsDone () && proj2.IsDone ()) {
840             if(proj1.LowerDistance()<= Precision::Confusion() &&
841                proj2.LowerDistance()<= Precision::Confusion())  {
842               detrompeur = 1;
843             }
844           }
845
846           if (proj3.IsDone () && proj4.IsDone ()) {
847             if(proj3.LowerDistance() <= Precision::Confusion() &&
848                proj4.LowerDistance() <= Precision::Confusion())  {
849               detrompeur = 2;
850             }
851           }
852           
853           gp_Dir TheDirExtr;
854           gp_Ax3 Axis;
855           Handle(Geom_Curve) TheNewCurve;
856           Standard_Boolean KPart = Standard_False;
857
858           if ( S1->DynamicType() == STANDARD_TYPE(Geom_RectangularTrimmedSurface)) {
859             S1 = Handle(Geom_RectangularTrimmedSurface)::
860               DownCast(S1)->BasisSurface();
861           }
862           if ( S2->DynamicType() == STANDARD_TYPE(Geom_RectangularTrimmedSurface)) {
863             S2 = Handle(Geom_RectangularTrimmedSurface)::
864               DownCast(S2)->BasisSurface();
865           }
866
867           Standard_Boolean PC1 = Standard_True; // KPart on S1
868           if (S1->DynamicType() == STANDARD_TYPE(Geom_SurfaceOfLinearExtrusion) &&
869               S2->DynamicType() == STANDARD_TYPE(Geom_Plane) ) {
870             KPart = Standard_True;
871             Axis = Handle(Geom_Plane)::DownCast(S2)->Position();
872             TheNewCurve = Handle(Geom_SurfaceOfLinearExtrusion)::
873               DownCast(S1)->BasisCurve();
874             TheDirExtr = Handle(Geom_SurfaceOfLinearExtrusion)::
875               DownCast(S1)->Direction();
876           }
877           else if (S2->DynamicType() == STANDARD_TYPE(Geom_SurfaceOfLinearExtrusion) &&
878                    S1->DynamicType() == STANDARD_TYPE(Geom_Plane) ) {
879             KPart = Standard_True;
880             PC1 = Standard_False;
881             Axis = Handle(Geom_Plane)::DownCast(S1)->Position();
882             TheNewCurve = Handle(Geom_SurfaceOfLinearExtrusion)::
883               DownCast(S2)->BasisCurve();
884             TheDirExtr = Handle(Geom_SurfaceOfLinearExtrusion)::
885               DownCast(S2)->Direction();
886           }
887           Handle(Geom_Circle) aCirc ;
888           if ( KPart) {  // very temporary on circles !!!
889             aCirc = Handle(Geom_Circle)::DownCast(TheNewCurve);
890             if (aCirc.IsNull())
891               KPart = Standard_False;
892             else
893               {
894                 gp_Dir AxofCirc = aCirc->Position().Direction();
895                 if (AxofCirc.IsParallel(Axis.Direction(),Precision::Angular()))
896                   KPart = Standard_True;
897                 else
898                   KPart = Standard_False;
899               }
900           }
901
902           Standard_Integer imin;
903           GeomInt_IntSS i2s;
904           if ( KPart) {
905             // direct calculation of NewC
906             Standard_Real aLocalReal =
907               gp_Vec(aCirc->Circ().Location(),Axis.Location()).
908                 Dot(Axis.Direction());
909             Standard_Real Cos = TheDirExtr.Dot(Axis.Direction());
910             gp_Vec VV = ( aLocalReal / Cos) * TheDirExtr;
911             newC = Handle(Geom_Curve)::DownCast(TheNewCurve->Translated(VV));
912             // it is possible to calculate PCurve
913             Handle(Geom2d_Line) L2d 
914               = new Geom2d_Line(gp_Pnt2d(0.,aLocalReal/Cos),
915                                 gp::DX2d());
916
917             if ( PC1) 
918               Einf.ChangeFirstPC() = L2d;
919             else
920               Einf.ChangeSecondPC() = L2d;
921           }
922           else {
923             S1 = myFMap(Einf.FirstFace()).Geometry();
924             S2 = myFMap(Einf.SecondFace()).Geometry();
925
926
927             // PCurves are not calculated immediately for 2 reasons:
928             // 1 - If ProjLib should make an Approx, it is stupid to approximate the 
929             //     entire intersection curve.
930             // 2 - Additionally, if YaRev, there is a risk to not be SameRange.
931             i2s.Perform(S1,S2,Precision::Confusion(),
932                         Standard_True,Standard_False,Standard_False);
933             
934             if (!i2s.IsDone() || i2s.NbLines() <= 0) {
935               errStat = Draft_EdgeRecomputation;
936               badShape = TopoDS::Edge(ite.Key());
937               return;
938             }
939             
940             Standard_Real Dist2, Dist2Min = 0., Glob2Min = RealLast();
941             GeomAdaptor_Curve TheCurve;
942             
943             Standard_Integer i,j; //,jmin;
944             
945             if (i2s.Line(1)->DynamicType() != STANDARD_TYPE(Geom_BSplineCurve))
946               {
947                 imin = 0;
948                 for (i=1; i<= i2s.NbLines(); i++) {
949                   TheCurve.Load(i2s.Line(i));
950                   Extrema_ExtPC myExtPC(pfv,TheCurve);
951                   
952                   Standard_Real locpmin = 0.;
953                   if (myExtPC.IsDone()) {
954                     if(myExtPC.NbExt() >= 1) {
955                       Dist2Min = myExtPC.SquareDistance(1);
956                       locpmin = myExtPC.Point(1).Parameter();
957                     }
958             if(myExtPC.NbExt() == 2 && Dist2Min > Precision::SquareConfusion()) {
959                       //to avoid incorrectly choosing the image 
960                       //of the first vertex of the initial edge
961                       Standard_Real d1_2 = myExtPC.SquareDistance(1);
962                       Standard_Real d2_2 = myExtPC.SquareDistance(2);
963                       if(d1_2 > 1.21*d2_2) {
964                         Dist2Min = myExtPC.SquareDistance(2);
965                         locpmin = myExtPC.Point(2).Parameter();
966                       }
967                       else if(d2_2 > 1.21*d1_2) {
968                         Dist2Min = myExtPC.SquareDistance(1);
969                         locpmin = myExtPC.Point(1).Parameter();
970                       }
971                       else {
972                         Standard_Real pfvpar = myExtPC.Point(1).Parameter();
973                         Standard_Real plvpar = myExtPC.Point(2).Parameter();
974                         newC = i2s.Line(i);
975                         
976                         gp_Pnt pfvprim, plvprim;
977                         
978                         newC->D0(pfvpar,pfvprim);
979                         newC->D0(plvpar,plvprim);
980                         
981                         Handle(Geom_Surface) theSurf;
982                         if(detrompeur == 1) {
983                           if(S1->DynamicType() == 
984                              STANDARD_TYPE(Geom_RectangularTrimmedSurface))  
985                             S1 = Handle(Geom_RectangularTrimmedSurface)::
986                               DownCast(S1)->BasisSurface();
987                           theSurf = S1;
988                           
989                         }
990                         else if(detrompeur == 2) {
991                           if(S2->DynamicType() == 
992                              STANDARD_TYPE(Geom_RectangularTrimmedSurface)) 
993                             S2 = Handle(Geom_RectangularTrimmedSurface)::
994                               DownCast(S2)->BasisSurface();
995                           theSurf = S2;             
996                         }
997                         if(detrompeur != 0 && detrompeur != 4) {
998                           Standard_Real ul = 0., vl = 0., uf = 0., vf = 0.;
999                           Standard_Real ufprim = 0., ulprim = 0., vfprim = 0., vlprim = 0.;
1000                           
1001                           if(theSurf->DynamicType() == STANDARD_TYPE(Geom_Plane)) {     
1002                             gp_Pln pl = Handle(Geom_Plane)::DownCast(S2)->Pln();
1003                             ElSLib::Parameters(pl, plv, ul, vl);
1004                             ElSLib::Parameters(pl, pfv, uf, vf);
1005                             ElSLib::Parameters(pl, plvprim, ulprim, vlprim);
1006                             ElSLib::Parameters(pl, pfvprim, ufprim, vfprim);
1007                           }
1008                           else if(theSurf->DynamicType() == 
1009                                   STANDARD_TYPE(Geom_CylindricalSurface)) {
1010                             gp_Cylinder cy = Handle(Geom_CylindricalSurface)
1011                               ::DownCast(S2)->Cylinder();
1012                             ElSLib::Parameters(cy, plv, ul, vl);
1013                             ElSLib::Parameters(cy, pfv, uf, vf);
1014                             ElSLib::Parameters(cy, plvprim, ulprim, vlprim);
1015                             ElSLib::Parameters(cy, pfvprim, ufprim, vfprim);
1016                           }
1017                           else detrompeur = 4;
1018                           
1019                           if(detrompeur == 1 || detrompeur == 2) {
1020                             gp_Vec2d v1((ul-ufprim), (vl-vfprim));
1021                             gp_Vec2d norm((vf-vfprim), (ufprim-uf)); 
1022                             gp_Vec2d v2((ulprim-ufprim), (vlprim-vfprim));
1023                             if( (v1.Dot(norm))*(v2.Dot(norm)) < 0) {
1024                               Dist2Min = myExtPC.SquareDistance(2);
1025                               locpmin = myExtPC.Point(2).Parameter();
1026                             }
1027                           }
1028                         }
1029                       }
1030                     }
1031                     if (myExtPC.NbExt() == 1 || myExtPC.NbExt() > 2 || detrompeur ==4) {
1032                       Dist2Min = myExtPC.SquareDistance(1);
1033                       locpmin = myExtPC.Point(1).Parameter();
1034                       for (j=2; j<=myExtPC.NbExt(); j++) {
1035                         Dist2 = myExtPC.SquareDistance(j);
1036                         if (Dist2 < Dist2Min) {
1037                           Dist2Min = Dist2;
1038                           locpmin = myExtPC.Point(j).Parameter();
1039                         }
1040                       }
1041                     }
1042                     else if(myExtPC.NbExt() < 1){
1043                       Standard_Real dist1_2,dist2_2;
1044                       gp_Pnt p1b,p2b;
1045                       myExtPC.TrimmedSquareDistances(dist1_2,dist2_2,p1b,p2b);
1046                       if (dist1_2 < dist2_2) {
1047                         Dist2Min = dist1_2;
1048                         locpmin = TheCurve.FirstParameter();
1049                       }
1050                       else {
1051                         Dist2Min = dist2_2;
1052                         locpmin = TheCurve.LastParameter();
1053                       }
1054                     }
1055                     
1056                     if (Dist2Min  < Glob2Min) {
1057                       Glob2Min = Dist2Min;
1058                       imin = i;
1059                       pmin = locpmin;
1060                     }
1061                   }
1062                 }
1063                 if (imin == 0) {
1064                   errStat = Draft_EdgeRecomputation;
1065                   badShape = TopoDS::Edge(ite.Key());
1066                   return;
1067                 }
1068                 
1069                 newC = i2s.Line(imin);
1070                 
1071                 newC->D1(pmin,pfv,newd1);
1072                 Standard_Boolean YaRev = d1fv.Dot(newd1) <0.; 
1073                 
1074                 if (YaRev)
1075                   newC->Reverse();
1076                 
1077                 if (i2s.HasLineOnS1(imin)) {
1078                   Einf.ChangeFirstPC() = i2s.LineOnS1(imin);
1079                   if ( YaRev) 
1080                     Einf.ChangeFirstPC()->Reverse();
1081                 }
1082                 
1083                 if (i2s.HasLineOnS2(imin)) {
1084                   Einf.ChangeSecondPC() = i2s.LineOnS2(imin);
1085                   if ( YaRev) 
1086                     Einf.ChangeSecondPC()->Reverse();
1087                 }
1088               } // if (i2s.Line(1)->DynamicType() != STANDARD_TYPE(Geom_BSplineCurve))
1089             else // i2s.Line(1) is BSplineCurve
1090               {
1091                 //Find the first curve to glue
1092                 TColGeom_SequenceOfCurve Candidates;
1093                 if (S1->DynamicType() == STANDARD_TYPE(Geom_CylindricalSurface) ||
1094                     S1->DynamicType() == STANDARD_TYPE(Geom_ConicalSurface)) 
1095                   {
1096                     for (i = 1; i <= i2s.NbLines(); i++)
1097                       {
1098                         Handle( Geom_Curve ) aCurve = i2s.Line(i);
1099                         gp_Pnt Pnt = aCurve->Value( aCurve->FirstParameter() );
1100                         GeomAPI_ProjectPointOnSurf projector( Pnt, S1, Precision::Confusion() );
1101                         Standard_Real U, V;
1102                         projector.LowerDistanceParameters( U, V );
1103                         if (Abs(U) <= Precision::Confusion() || Abs(U-2.*M_PI) <= Precision::Confusion())
1104                           Candidates.Append( aCurve );
1105                         else
1106                           {
1107                             Pnt = aCurve->Value( aCurve->LastParameter() );
1108                             projector.Init( Pnt, S1, Precision::Confusion() );
1109                             projector.LowerDistanceParameters( U, V );
1110                             if (Abs(U) <= Precision::Confusion() || Abs(U-2.*M_PI) <= Precision::Confusion())
1111                               {
1112                                 aCurve->Reverse();
1113                                 Candidates.Append( aCurve );
1114                               }
1115                           }
1116                       }
1117
1118                     if(Candidates.Length() == 0) 
1119                       {
1120 //                      errStat = Draft_EdgeRecomputation;
1121 //                      badShape = TopoDS::Edge(ite.Key());
1122 //                      return;
1123                         for (i = 1; i <= i2s.NbLines(); i++)
1124                           Candidates.Append( i2s.Line(i) );
1125                       }
1126                   }
1127                 else 
1128                   {
1129                     for (i = 1; i <= i2s.NbLines(); i++)
1130                       Candidates.Append( i2s.Line(i) );
1131                   }
1132                 
1133                 Handle( Geom_Curve ) FirstCurve;
1134                 if (Candidates.Length() > 1)
1135                   {
1136                     Dist2Min = RealLast();
1137                     for (i = 1; i <= Candidates.Length(); i++)
1138                       {
1139                         Handle( Geom_Curve ) aCurve = Candidates(i);
1140                         gp_Pnt Pnt = aCurve->Value( aCurve->FirstParameter() );
1141                         Dist2 = Pnt.SquareDistance( pfv );
1142                         if (Dist2 < Dist2Min)
1143                           {
1144                             Dist2Min = Dist2;
1145                             FirstCurve = aCurve;
1146                           }
1147                       }
1148                   }
1149                 else
1150                   FirstCurve = Candidates(1);
1151
1152                 //Glueing
1153                 TColGeom_SequenceOfCurve Curves;
1154                 for (i = 1; i <= i2s.NbLines(); i++)
1155                   if (FirstCurve != i2s.Line(i))
1156                     Curves.Append( i2s.Line(i) );
1157
1158                 TColGeom_SequenceOfCurve ToGlue;
1159                 gp_Pnt EndPoint = FirstCurve->Value( FirstCurve->LastParameter() );
1160                 Standard_Boolean added = Standard_True;
1161                 while (added)
1162                   {
1163                     added = Standard_False;
1164                     for (i = 1; i <= Curves.Length(); i++)
1165                       {
1166                         Handle( Geom_Curve ) aCurve = Curves(i);
1167                         gp_Pnt pfirst, plast;
1168                         pfirst = aCurve->Value( aCurve->FirstParameter() );
1169                         plast = aCurve->Value( aCurve->LastParameter() );
1170                         if (pfirst.Distance( EndPoint ) <= Precision::Confusion())
1171                           {
1172                             ToGlue.Append( aCurve );
1173                             EndPoint = plast;
1174                             Curves.Remove(i);
1175                             added = Standard_True;
1176                             break;
1177                           }
1178                         if (plast.Distance( EndPoint ) <= Precision::Confusion())
1179                           {
1180                             aCurve->Reverse();
1181                             ToGlue.Append( aCurve );
1182                             EndPoint = pfirst;
1183                             Curves.Remove(i);
1184                             added = Standard_True;
1185                             break;
1186                           }
1187                       }
1188                   }
1189
1190                 if (FirstCurve.IsNull()) {
1191                   errStat = Draft_EdgeRecomputation;
1192                   badShape = TopoDS::Edge(ite.Key());
1193                   return;
1194                 }
1195                 
1196                 GeomConvert_CompCurveToBSplineCurve Concat( Handle(Geom_BSplineCurve)::DownCast(FirstCurve) );
1197                 for (i = 1; i <= ToGlue.Length(); i++)
1198                   Concat.Add( Handle(Geom_BSplineCurve)::DownCast(ToGlue(i)), Precision::Confusion(), Standard_True );
1199
1200                 newC = Concat.BSplineCurve();
1201                 
1202                 TheCurve.Load( newC );
1203                 Extrema_ExtPC myExtPC( pfv, TheCurve );
1204                 Dist2Min = RealLast();
1205                 for (i = 1; i <= myExtPC.NbExt(); i++)
1206                   {
1207                     if (myExtPC.IsMin(i))
1208                     {
1209                       Dist2 = myExtPC.SquareDistance(i);
1210                       if (Dist2 < Dist2Min)
1211                       {
1212                         Dist2Min = Dist2;
1213                         pmin = myExtPC.Point(i).Parameter();
1214                       }
1215                     }
1216                   }
1217                 newC->D1(pmin,pfv,newd1);
1218                 Standard_Boolean YaRev = d1fv.Dot(newd1) < 0.; 
1219                 
1220                 if (YaRev)
1221                   newC->Reverse();
1222                 /*
1223                 if (i2s.HasLineOnS1(imin)) {
1224                   Einf.ChangeFirstPC() = i2s.LineOnS1(imin);
1225                   if ( YaRev) 
1226                     Einf.ChangeFirstPC()->Reverse();
1227                 }
1228                 
1229                 if (i2s.HasLineOnS2(imin)) {
1230                   Einf.ChangeSecondPC() = i2s.LineOnS2(imin);
1231                   if ( YaRev) 
1232                     Einf.ChangeSecondPC()->Reverse();
1233                 }
1234                 */
1235               } // else: i2s.NbLines() > 2 && S1 is Cylinder or Cone
1236             
1237             Einf.Tolerance(Max(Einf.Tolerance(), i2s.TolReached3d()));
1238           }  // End step KPart
1239         }
1240         else { // case of tangency
1241           const TopoDS_Face& F1 = Einf.FirstFace();
1242           const TopoDS_Face& F2 = Einf.SecondFace();
1243
1244           Handle(Geom_Surface) aLocalS1 = myFMap(F1).Geometry();
1245           Handle(Geom_Surface) aLocalS2 = myFMap(F2).Geometry();
1246           if (aLocalS1.IsNull() || aLocalS2.IsNull()) {
1247             errStat = Draft_EdgeRecomputation;
1248             badShape = TopoDS::Edge(ite.Key());
1249             return;
1250           }
1251           if (aLocalS1->DynamicType() == STANDARD_TYPE(Geom_RectangularTrimmedSurface)) {
1252             aLocalS1 = Handle(Geom_RectangularTrimmedSurface)::
1253             DownCast(aLocalS1)->BasisSurface();
1254           }
1255           if (aLocalS2->DynamicType() == STANDARD_TYPE(Geom_RectangularTrimmedSurface)) {
1256             aLocalS2 = Handle(Geom_RectangularTrimmedSurface)::
1257             DownCast(aLocalS2)->BasisSurface();
1258           }
1259
1260           gp_Dir dirextr;
1261           //Standard_Boolean dirfound = Standard_False;
1262           if (aLocalS1->DynamicType() == STANDARD_TYPE(Geom_CylindricalSurface)) {
1263             gp_Cylinder cyl = 
1264               Handle(Geom_CylindricalSurface)::DownCast(aLocalS1)->Cylinder();
1265             dirextr = cyl.Axis().Direction();
1266             //dirfound = Standard_True;
1267             // see direction...
1268
1269           }
1270           else if (aLocalS1->DynamicType() == STANDARD_TYPE(Geom_SurfaceOfLinearExtrusion)) {
1271             dirextr = Handle(Geom_SurfaceOfLinearExtrusion)::
1272               DownCast(aLocalS1)->Direction();
1273             //dirfound = Standard_True;
1274             // see direction...
1275
1276             // Here it is possible to calculate PCurve.
1277             Handle(Geom_SurfaceOfLinearExtrusion) SEL = 
1278               Handle(Geom_SurfaceOfLinearExtrusion)::DownCast(aLocalS1);
1279             Handle(Geom_Circle) GCir = 
1280               Handle(Geom_Circle)::DownCast(SEL->BasisCurve());
1281             if ( !GCir.IsNull()) {
1282               Standard_Real U = ElCLib::Parameter(GCir->Circ(),ptfixe);
1283               Handle(Geom2d_Line) PC1 = 
1284                 new Geom2d_Line(gp_Pnt2d(U,0.),gp::DY2d());
1285               Einf.ChangeFirstPC() = PC1;
1286             }
1287           }
1288
1289           else if (aLocalS2->DynamicType() == STANDARD_TYPE(Geom_CylindricalSurface)) {
1290             gp_Cylinder cyl = 
1291               Handle(Geom_CylindricalSurface)::DownCast(aLocalS2)->Cylinder();
1292             dirextr = cyl.Axis().Direction();
1293             // dirfound = Standard_True;
1294             // see direction...
1295
1296           }
1297           else if (aLocalS2->DynamicType() == STANDARD_TYPE(Geom_SurfaceOfLinearExtrusion)) {
1298             dirextr = Handle(Geom_SurfaceOfLinearExtrusion)::
1299               DownCast(aLocalS2)->Direction();
1300             // dirfound = Standard_True;
1301             // see direction...
1302
1303             // Here it is possible to calculate PCurve.
1304             Handle(Geom_SurfaceOfLinearExtrusion) SEL = 
1305               Handle(Geom_SurfaceOfLinearExtrusion)::DownCast(aLocalS2);
1306             Handle(Geom_Circle) GCir = 
1307               Handle(Geom_Circle)::DownCast(SEL->BasisCurve());
1308             if ( !GCir.IsNull()) {
1309               Standard_Real U = ElCLib::Parameter(GCir->Circ(),ptfixe);
1310               Handle(Geom2d_Line) PC2 = 
1311                 new Geom2d_Line(gp_Pnt2d(U,0.),gp::DY2d());
1312               Einf.ChangeSecondPC() = PC2;
1313             }
1314           }
1315           newC = new Geom_Line(ptfixe,dirextr);
1316
1317           gp_Pnt pfv;
1318           gp_Vec d1fv,newd1;
1319           C->D1(0.,pfv,d1fv);
1320           newC->D1(0.,pfv,newd1);
1321           Standard_Boolean YaRev = d1fv.Dot(newd1) <0.; 
1322           if (YaRev) {
1323             newC->Reverse();
1324             if(!Einf.FirstPC().IsNull()) {
1325               Einf.ChangeFirstPC()->Reverse();
1326             }
1327             if(!Einf.SecondPC().IsNull()) {
1328               Einf.ChangeSecondPC()->Reverse();
1329             }
1330           }
1331         }
1332
1333         Handle(Geom_TrimmedCurve) T = Handle(Geom_TrimmedCurve)::DownCast(newC);
1334         if (!T.IsNull()) newC = T->BasisCurve();
1335         Einf.ChangeGeometry() = newC;
1336       }
1337       else if (!Einf.NewGeometry()){
1338         // set existing curve 3D
1339         Handle(Geom_TrimmedCurve) T = Handle(Geom_TrimmedCurve)::DownCast(C);
1340         if (!T.IsNull()) C = T->BasisCurve();
1341         Einf.ChangeGeometry() = C;
1342       }
1343       ite.Next();
1344     }
1345
1346     // Calculate new vertices.
1347
1348     Draft_DataMapIteratorOfDataMapOfVertexVertexInfo itv(myVMap);
1349
1350     Handle(GeomAdaptor_HCurve)   HAC = new GeomAdaptor_HCurve;
1351     Handle(GeomAdaptor_HSurface) HAS = new GeomAdaptor_HSurface;
1352     
1353     while (itv.More()) {
1354       GeomAdaptor_Curve    AC;
1355       GeomAdaptor_Surface  AS;
1356
1357       Draft_VertexInfo& Vinf = myVMap(itv.Key());
1358       if (!Choose(myFMap,myEMap,itv.Key(),Vinf,AC,AS)) {
1359
1360 // no concerted edge => alignment of two consecutive edges.
1361         gp_Pnt pvt;
1362         Vinf.ChangeGeometry() = pvt;
1363         Vinf.InitEdgeIterator();
1364         if (Vinf.MoreEdge()) {
1365           const TopoDS_Edge& Edg1 = Vinf.Edge();
1366           //const Draft_EdgeInfo& Einf1 = myEMap(Edg1);
1367           Draft_EdgeInfo& Einf1 = myEMap(Edg1);
1368           gp_Pnt vtori = BRep_Tool::Pnt(itv.Key());
1369           //Einf1.Geometry()->D0(Vinf.Parameter(Edg1), pvt);
1370           GeomAPI_ProjectPointOnCurve Projector( vtori, Einf1.Geometry() ); //patch
1371           pvt = Projector.NearestPoint();
1372
1373 #ifdef OCCT_DEBUG
1374           static Standard_Integer VertexRecomp = 1;
1375           if (VertexRecomp!=0) {
1376             cout << "pori :" << vtori.X() << " " << vtori.Y() << " " << vtori.Z() << endl;
1377             cout << "  Edg 1 :" << Vinf.Parameter(Edg1) << endl;
1378             cout << "pvt :" << pvt.X() << " " << pvt.Y() << " " << pvt.Z() << endl;
1379           }
1380 #endif
1381
1382           Standard_Real dion=pvt.SquareDistance(vtori);
1383           Vinf.NextEdge();
1384           if (Vinf.MoreEdge()) {
1385             const TopoDS_Edge& Edg2 = Vinf.Edge();
1386             //const Draft_EdgeInfo& Einf2 = myEMap(Edg2);
1387             Draft_EdgeInfo& Einf2 = myEMap(Edg2);
1388 //          Standard_Real f;
1389             gp_Pnt opvt;
1390             Einf2.Geometry()->D0(Vinf.Parameter(Edg2), opvt);
1391
1392 #ifdef OCCT_DEBUG
1393           if (VertexRecomp!=0) {
1394             cout << "  Edg 2 :" << Vinf.Parameter(Vinf.Edge()) << endl;
1395             cout << "opvt " << opvt.X() << " " << opvt.Y() << " " << opvt.Z() << endl;
1396           }
1397 #endif
1398
1399             if (opvt.SquareDistance(vtori) < dion) {
1400               pvt = opvt;
1401             }
1402             //Vinf.ChangeParameter(Edg2) = Parameter(Einf2.Geometry(), pvt);
1403             Standard_Integer done;
1404             Standard_Real param = Parameter(Einf2.Geometry(), pvt, done);
1405             if (done != 0)
1406               {
1407                 S1 = myFMap(Einf2.FirstFace()).Geometry();
1408                 S2 = myFMap(Einf2.SecondFace()).Geometry();
1409                 Vinf.ChangeParameter(Edg2) = SmartParameter( Einf2, BRep_Tool::Tolerance(Edg2), pvt, done, S1, S2 );
1410               }
1411             else
1412               Vinf.ChangeParameter(Edg2) = param;
1413           }
1414
1415           Vinf.ChangeGeometry() = pvt;
1416           //Vinf.ChangeParameter(Edg1) = Parameter(Einf1.Geometry(), pvt);
1417           Standard_Integer done;
1418           Standard_Real param = Parameter(Einf1.Geometry(), pvt, done);
1419           if (done != 0)
1420             {
1421               S1 = myFMap(Einf1.FirstFace()).Geometry();
1422               S2 = myFMap(Einf1.SecondFace()).Geometry();
1423               Vinf.ChangeParameter(Edg1) = SmartParameter( Einf1, BRep_Tool::Tolerance(Edg1), pvt, done, S1, S2 );
1424             }
1425           else
1426             Vinf.ChangeParameter(Edg1) = param;
1427           itv.Next();
1428           continue;
1429         }
1430
1431
1432         errStat = Draft_VertexRecomputation;
1433         badShape = TopoDS::Vertex(itv.Key());
1434         return;
1435       }
1436
1437       IntCurveSurface_HInter myintcs;
1438       HAC->Set(AC);
1439       HAS->Set(AS);
1440       myintcs.Perform(HAC,HAS);
1441
1442       if (!myintcs.IsDone()) {
1443         errStat = Draft_VertexRecomputation;
1444         badShape = TopoDS::Vertex(itv.Key());
1445         return;
1446       }
1447
1448       gp_Pnt vtori = BRep_Tool::Pnt(itv.Key());
1449       gp_Pnt pvt;
1450
1451       Standard_Integer nbsol = myintcs.NbPoints();
1452       if (nbsol <= 0)
1453         {
1454           Extrema_ExtCS extr( AC, AS, Precision::PConfusion(), Precision::PConfusion() );
1455
1456           if(!extr.IsDone() || extr.NbExt() == 0) {
1457             errStat = Draft_VertexRecomputation;
1458             badShape = TopoDS::Vertex(itv.Key());
1459             return;
1460           }
1461             
1462
1463           Standard_Real disref = RealLast();
1464           Standard_Integer iref = 0;
1465           Extrema_POnCurv Pc;
1466           Extrema_POnSurf Ps;
1467           for (Standard_Integer i = 1; i <= extr.NbExt(); i++)
1468             {
1469               extr.Points( i, Pc, Ps );
1470               Standard_Real distemp = Pc.Value().SquareDistance(vtori);
1471               if ( distemp < disref)
1472                 {
1473                   disref = distemp;
1474                   iref = i;
1475                 }
1476             } 
1477           extr.Points( iref, Pc, Ps );
1478           pvt = Pc.Value();
1479         }
1480       else
1481         {
1482           Standard_Real disref = RealLast();
1483           Standard_Integer iref = 0;
1484           for (Standard_Integer i = 1; i <= nbsol; i++)
1485             {
1486               Standard_Real distemp = myintcs.Point(i).Pnt().SquareDistance(vtori);
1487               if ( distemp < disref)
1488                 {
1489                   disref = distemp;
1490                   iref = i;
1491                 }
1492             } 
1493           pvt = myintcs.Point(iref).Pnt();
1494         }
1495           
1496       Vinf.ChangeGeometry() = pvt;
1497
1498       for (Vinf.InitEdgeIterator();Vinf.MoreEdge(); Vinf.NextEdge()) {
1499         const TopoDS_Edge& Edg = Vinf.Edge();
1500         //const Draft_EdgeInfo& Einf = myEMap(Edg);
1501         Draft_EdgeInfo& Einf = myEMap(Edg);
1502         //Vinf.ChangeParameter(Edg) = Parameter(Einf.Geometry(),pvt);
1503         Standard_Integer done;
1504         Standard_Real param = Parameter(Einf.Geometry(), pvt, done);
1505         if (done != 0)
1506           {
1507             S1 = myFMap(Einf.FirstFace()).Geometry();
1508             S2 = myFMap(Einf.SecondFace()).Geometry();
1509             Vinf.ChangeParameter(Edg) = SmartParameter( Einf, BRep_Tool::Tolerance(Edg), pvt, done, S1, S2 );
1510           }
1511         else
1512           Vinf.ChangeParameter(Edg) = param;
1513       }
1514       itv.Next();
1515     }
1516   }
1517
1518   // small loop of validation/protection
1519
1520   for (Draft_DataMapIteratorOfDataMapOfEdgeEdgeInfo ite(myEMap);
1521        ite.More(); ite.Next()) {
1522     const TopoDS_Edge& edg = TopoDS::Edge(ite.Key());
1523
1524     TopoDS_Vertex Vf,Vl;
1525     TopExp::Vertices(edg,Vf,Vl);
1526     if (edg.Orientation() == TopAbs_REVERSED) {
1527       Vf.Reverse();
1528       Vl.Reverse();
1529     }
1530     Standard_Real pf,pl,tolerance;
1531     if (!NewParameter(Vf,edg,pf,tolerance)) {
1532       pf = BRep_Tool::Parameter(Vf,edg);
1533     }
1534     if (!NewParameter(Vl,edg,pl,tolerance)) {
1535       pl = BRep_Tool::Parameter(Vl,edg);
1536     }
1537     if (pl <= pf) {
1538 //      const Handle(Geom_Curve) gc=ite.Value().Geometry();
1539 //      if (!gc.IsNull()) {
1540 //      pl = gc->LastParameter();
1541 //      pf = gc->FirstParameter();
1542 //      }
1543       Handle( Geom_Curve ) theCurve = myEMap(edg).Geometry();
1544       if (theCurve->IsClosed())
1545         {
1546           // pf >= pl
1547           Standard_Real FirstPar = theCurve->FirstParameter(), LastPar = theCurve->LastParameter();
1548           Standard_Real pconf = Precision::PConfusion();
1549           if (Abs( pf - LastPar ) <= pconf)
1550             pf = FirstPar;
1551           else if (Abs( pl - FirstPar ) <= pconf)
1552             pl = LastPar;
1553
1554           if(pl <= pf) {
1555             pl += (LastPar-FirstPar);
1556           }
1557
1558         }
1559       if (pl <= pf) {
1560         errStat = Draft_EdgeRecomputation;
1561         badShape = TopoDS::Edge(ite.Key());
1562         return;
1563       }
1564     }
1565     if (myVMap.IsBound( Vf ))
1566       myVMap(Vf).ChangeParameter(edg) = pf;
1567     if (myVMap.IsBound( Vl ))
1568       myVMap(Vl).ChangeParameter(edg) = pl;
1569   }
1570 }
1571
1572
1573
1574 //=======================================================================
1575 //function : NewSurface
1576 //purpose  : 
1577 //=======================================================================
1578
1579 Handle(Geom_Surface) Draft_Modification::NewSurface
1580   (const Handle(Geom_Surface)& S,
1581    const TopAbs_Orientation Oris,
1582    const gp_Dir& Direction,
1583    const Standard_Real Angle,
1584    const gp_Pln& NeutralPlane) 
1585 {
1586   Handle(Geom_Surface) NewS;
1587
1588   Handle(Standard_Type) TypeS = S->DynamicType();
1589
1590   if (TypeS == STANDARD_TYPE(Geom_Plane)) {
1591     gp_Pln Pl = Handle(Geom_Plane)::DownCast(S)->Pln();
1592     gp_Ax1 Axe;
1593     Standard_Real Theta;
1594     if (FindRotation(Pl,Oris,Direction,Angle,NeutralPlane,Axe,Theta)) {
1595       if ( Abs(Theta) > Precision::Angular()) {
1596         NewS = Handle(Geom_Surface)::DownCast(S->Rotated(Axe,Theta));
1597       }
1598       else {
1599         NewS = S;
1600       }
1601     }
1602   }
1603   else if (TypeS == STANDARD_TYPE(Geom_CylindricalSurface)) {
1604     Standard_Real testdir = Direction.Dot(NeutralPlane.Axis().Direction());
1605     if (Abs(testdir) <= 1.-Precision::Angular()) {      
1606 #ifdef OCCT_DEBUG
1607     cout << "NewSurfaceCyl:Draft_Direction_and_Neutral_Perpendicular" << endl;
1608 #endif
1609       return NewS;      
1610     }     
1611     gp_Cylinder Cy = Handle(Geom_CylindricalSurface)::DownCast(S)->Cylinder();     
1612     testdir = Direction.Dot(Cy.Axis().Direction());
1613     if (Abs(testdir) <= 1.-Precision::Angular()) {
1614 #ifdef OCCT_DEBUG
1615     cout << "NewSurfaceCyl:Draft_Direction_and_Cylinder_Perpendicular" << endl;
1616 #endif
1617       return NewS;
1618     }
1619     if (Abs(Angle) > Precision::Angular())
1620     {
1621       IntAna_QuadQuadGeo i2s;
1622       i2s.Perform(NeutralPlane,Cy,Precision::Angular(),Precision::Confusion());
1623       Standard_Boolean isIntDone = i2s.IsDone();
1624
1625       if(i2s.TypeInter() == IntAna_Ellipse)
1626       {
1627         const gp_Elips anEl = i2s.Ellipse(1);
1628         const Standard_Real aMajorR = anEl.MajorRadius();
1629         const Standard_Real aMinorR = anEl.MinorRadius();
1630         isIntDone = (aMajorR < 100000.0 * aMinorR);
1631       }
1632
1633       if (!isIntDone || i2s.TypeInter() != IntAna_Circle) {
1634 #ifdef OCCT_DEBUG
1635     cout << "NewSurfaceCyl:Draft_Intersection_Neutral_Cylinder_NotDone" << endl;
1636 #endif
1637         return NewS;
1638       } 
1639       gp_Ax3 axcone = Cy.Position();
1640       // Pb : Where is the material???
1641       Standard_Real alpha = Angle;
1642       Standard_Boolean direct(axcone.Direct());
1643       if ((direct && Oris == TopAbs_REVERSED) ||
1644           (!direct && Oris == TopAbs_FORWARD)) {
1645         alpha = -alpha;
1646       }
1647       
1648       gp_Pnt Center = i2s.Circle(1).Location();
1649       if (testdir <0.) {
1650         alpha = -alpha;
1651       }
1652       Standard_Real Z = ElCLib::LineParameter(Cy.Axis(),Center);
1653       Standard_Real Rad = Cy.Radius()+Z*Tan(alpha);
1654       if (Rad < 0.) {
1655         Rad = -Rad;
1656       }
1657       else {
1658         alpha = -alpha;
1659       }
1660       gp_Cone co(axcone,alpha,Rad);
1661       NewS = new Geom_ConicalSurface(co);
1662     }
1663     else {
1664       NewS = S;
1665     }
1666   }
1667   else if (TypeS == STANDARD_TYPE(Geom_ConicalSurface)) {
1668     
1669     Standard_Real testdir = Direction.Dot(NeutralPlane.Axis().Direction());
1670     if (Abs(testdir) <= 1.-Precision::Angular()) {      
1671 #ifdef OCCT_DEBUG
1672     cout << "NewSurfaceCone:Draft_Direction_and_Neutral_Perpendicular" << endl;
1673 #endif
1674       return NewS;      
1675     }   
1676     
1677     gp_Cone Co1 = Handle(Geom_ConicalSurface)::DownCast(S)->Cone();
1678     
1679     testdir = Direction.Dot(Co1.Axis().Direction());
1680     if (Abs(testdir) <= 1.-Precision::Angular()) {
1681 #ifdef OCCT_DEBUG
1682     cout << "NewSurfaceCone:Draft_Direction_and_Cone_Perpendicular" << endl;
1683 #endif
1684       return NewS;
1685     }
1686
1687
1688     IntAna_QuadQuadGeo i2s;
1689     i2s.Perform(NeutralPlane,Co1,Precision::Angular(),Precision::Confusion());
1690     if (!i2s.IsDone() || i2s.TypeInter() != IntAna_Circle) {
1691 #ifdef OCCT_DEBUG
1692     cout << "NewSurfaceCone:Draft_Intersection_Neutral_Conical_NotDone" << endl;
1693 #endif
1694       return NewS;
1695     }
1696     gp_Ax3 axcone = Co1.Position();
1697     // Pb : Where is the material???
1698     Standard_Real alpha = Angle;
1699     Standard_Boolean direct(axcone.Direct());
1700     if ((direct && Oris == TopAbs_REVERSED) ||
1701         (!direct && Oris == TopAbs_FORWARD)) {
1702       alpha = -alpha;
1703     }
1704
1705     gp_Pnt Center = i2s.Circle(1).Location();
1706     if (Abs(Angle) > Precision::Angular()) {
1707       if (testdir <0.) {
1708         alpha = -alpha;
1709       }
1710       Standard_Real Z = ElCLib::LineParameter(Co1.Axis(),Center);
1711       Standard_Real Rad = i2s.Circle(1).Radius()+Z*Tan(alpha);
1712       if (Rad < 0.) {
1713         Rad = -Rad;
1714       }
1715       else {
1716         alpha = -alpha;
1717       }
1718       if (Abs(alpha-Co1.SemiAngle()) < Precision::Angular()) {
1719         NewS = S;
1720       }
1721       else {
1722         gp_Cone co(axcone,alpha,Rad);
1723         NewS = new Geom_ConicalSurface(co);
1724       }
1725     }
1726     else {
1727       NewS = new
1728         Geom_CylindricalSurface(gp_Cylinder(axcone,i2s.Circle(1).Radius()));
1729     }
1730   }
1731   else {
1732 #ifdef OCCT_DEBUG
1733     cout << "NewSurface:Draft_SurfNotYetImplemented" << endl;
1734 #endif
1735   }
1736   return NewS;
1737
1738
1739
1740 //=======================================================================
1741 //function : NewCurve
1742 //purpose  : 
1743 //=======================================================================
1744
1745 Handle(Geom_Curve) Draft_Modification::NewCurve
1746   (const Handle(Geom_Curve)& C,
1747    const Handle(Geom_Surface)& S,
1748    const TopAbs_Orientation Oris,
1749    const gp_Dir& Direction,
1750    const Standard_Real Angle,
1751    const gp_Pln& NeutralPlane,
1752    const Standard_Boolean )
1753
1754 {
1755   Handle(Geom_Curve) NewC;
1756
1757   Handle(Standard_Type) TypeS = S->DynamicType();
1758
1759   if (TypeS == STANDARD_TYPE(Geom_Plane)) {
1760     gp_Pln Pl = Handle(Geom_Plane)::DownCast(S)->Pln();
1761     gp_Ax1 Axe;
1762     Standard_Real Theta;
1763     if (FindRotation(Pl,Oris,Direction,Angle,NeutralPlane,Axe,Theta)) {
1764       if ( Abs(Theta) > Precision::Angular()) {
1765         NewC = Handle(Geom_Curve)::DownCast(C->Rotated(Axe,Theta));
1766       }
1767       else {
1768         NewC = C;
1769       }
1770     }
1771     return NewC;
1772   }
1773
1774
1775   if (C->DynamicType() != STANDARD_TYPE(Geom_Line)) {
1776     return NewC;
1777   }
1778
1779
1780   gp_Lin lin = Handle(Geom_Line)::DownCast(C)->Lin();
1781 //  Standard_Real testdir = Direction.Dot(lin.Direction());
1782 //  if (Abs(testdir) <= 1.-Precision::Angular()) {
1783 //    return NewC;
1784 //  }
1785   gp_Dir Norm;
1786   if (TypeS == STANDARD_TYPE(Geom_CylindricalSurface)) {
1787     Standard_Real U,V;
1788     gp_Vec d1u,d1v;
1789     gp_Pnt pbid;
1790     gp_Cylinder Cy = Handle(Geom_CylindricalSurface)::DownCast(S)->Cylinder();
1791     ElSLib::Parameters(Cy,lin.Location(),U,V);
1792     ElSLib::D1(U,V,Cy,pbid,d1u,d1v);
1793     Norm = d1u.Crossed(d1v);
1794   }
1795   else if (TypeS == STANDARD_TYPE(Geom_ConicalSurface)) {
1796     Standard_Real U,V;
1797     gp_Vec d1u,d1v;
1798     gp_Pnt pbid;
1799     gp_Cone Co = Handle(Geom_ConicalSurface)::DownCast(S)->Cone();
1800     ElSLib::Parameters(Co,lin.Location(),U,V);
1801     ElSLib::D1(U,V,Co,pbid,d1u,d1v);
1802     Norm = d1u.Crossed(d1v);
1803   }
1804
1805   IntAna_IntConicQuad ilipl(lin,NeutralPlane,Precision::Angular());
1806   if (ilipl.IsDone() && ilipl.NbPoints() != 0){
1807     if (Oris == TopAbs_REVERSED) {
1808       Norm.Reverse();
1809     }
1810     gp_Ax1 axrot(ilipl.Point(1), Norm.Crossed(Direction));
1811     gp_Lin lires = gp_Lin(gp_Ax1(ilipl.Point(1),Direction)).
1812       Rotated(axrot,Angle);
1813     if (lires.Direction().Dot(lin.Direction()) < 0.) {
1814       lires.Reverse();
1815     }
1816     NewC = new Geom_Line(lires);
1817   }
1818   return NewC;
1819 }
1820
1821
1822 //=======================================================================
1823 //function : Choose
1824 //purpose  : 
1825 //=======================================================================
1826
1827 static Standard_Boolean Choose(const Draft_DataMapOfFaceFaceInfo& theFMap,
1828                                Draft_DataMapOfEdgeEdgeInfo& theEMap,
1829                                const TopoDS_Vertex& Vtx,
1830                                Draft_VertexInfo& Vinf,
1831                                GeomAdaptor_Curve& AC,
1832                                GeomAdaptor_Surface& AS)
1833 {
1834   gp_Vec tgref; 
1835   Vinf.InitEdgeIterator();
1836
1837   // Find a regular edge with null SecondFace
1838   while (Vinf.MoreEdge()) {
1839     const TopoDS_Edge& E1 = Vinf.Edge();
1840     const Draft_EdgeInfo& Einf1 = theEMap(E1);
1841     if (Einf1.SecondFace().IsNull()) {
1842       break;
1843     }
1844     else {
1845       GeomAbs_Shape te = BRep_Tool::Continuity(E1,Einf1.FirstFace(),
1846                                                   Einf1.SecondFace());
1847       if (te >= GeomAbs_G1) {
1848         break;
1849       }
1850     }
1851     Vinf.NextEdge();
1852   }
1853   if (!Vinf.MoreEdge()) { // take the first edge
1854     Vinf.InitEdgeIterator();
1855   }
1856
1857   const TopoDS_Edge& Eref = Vinf.Edge();
1858   //const Draft_EdgeInfo& Einf = theEMap(Eref);
1859   Draft_EdgeInfo& Einf = theEMap(Eref);
1860
1861   AC.Load(Einf.Geometry());
1862
1863   Standard_Real f,l,prm;
1864   TopLoc_Location Loc;
1865   Handle(Geom_Curve) C = BRep_Tool::Curve(Eref,Loc,f,l);
1866   C = Handle(Geom_Curve)::DownCast(C->Transformed(Loc.Transformation()));
1867   gp_Pnt ptbid;
1868   //prm = Parameter(C,BRep_Tool::Pnt(Vtx));
1869   Standard_Integer done;
1870   Standard_Real param = Parameter( C, BRep_Tool::Pnt(Vtx), done );
1871   if (done != 0)
1872     {
1873       Handle( Geom_Surface ) S1 = theFMap(Einf.FirstFace()).Geometry();
1874       Handle( Geom_Surface ) S2 = theFMap(Einf.SecondFace()).Geometry();
1875       prm = SmartParameter( Einf, BRep_Tool::Tolerance(Eref), BRep_Tool::Pnt(Vtx), done, S1, S2 );
1876     }
1877   else
1878     prm = param;
1879   C->D1(prm,ptbid,tgref);
1880
1881
1882   Vinf.InitEdgeIterator();
1883   while (Vinf.MoreEdge()) {
1884     // Find a non tangent edge
1885     const TopoDS_Edge& Edg = Vinf.Edge();
1886     if (!Edg.IsSame(Eref)) {
1887       //const Draft_EdgeInfo& Einfo = theEMap(Edg);
1888       Draft_EdgeInfo& Einfo = theEMap(Edg);
1889       if (!Einfo.SecondFace().IsNull() &&
1890           BRep_Tool::Continuity(Edg,Einfo.FirstFace(),Einfo.SecondFace()) 
1891           <= GeomAbs_C0) {
1892         C = BRep_Tool::Curve(Edg,Loc,f,l);
1893         C = Handle(Geom_Curve)::DownCast(C->Transformed(Loc.Transformation()));
1894         //prm = Parameter(C,BRep_Tool::Pnt(Vtx));
1895         Standard_Integer anewdone;
1896         Standard_Real anewparam = Parameter( C, BRep_Tool::Pnt(Vtx), anewdone );
1897         if (anewdone != 0)
1898           {
1899             Handle( Geom_Surface ) S1 = theFMap(Einfo.FirstFace()).Geometry();
1900             Handle( Geom_Surface ) S2 = theFMap(Einfo.SecondFace()).Geometry();
1901             prm = SmartParameter( Einfo, BRep_Tool::Tolerance(Edg), BRep_Tool::Pnt(Vtx), anewdone, S1, S2 );
1902           }
1903         else
1904           prm = anewparam;
1905         gp_Vec tg;
1906         C->D1(prm,ptbid,tg);
1907         if (tg.CrossMagnitude(tgref) > Precision::Confusion()) {
1908           break;
1909         }
1910       }
1911     }
1912     Vinf.NextEdge();
1913   }
1914   if (!Vinf.MoreEdge()) {
1915     return Standard_False;
1916   }
1917
1918   const Draft_EdgeInfo& Einf2 = theEMap(Vinf.Edge());
1919   if (!Einf.SecondFace().IsNull()) {
1920
1921     if (Einf2.FirstFace().IsSame(Einf.FirstFace()) ||
1922         Einf2.FirstFace().IsSame(Einf.SecondFace())) {
1923       AS.Load(theFMap(Einf2.SecondFace()).Geometry());
1924     }
1925     else {
1926       AS.Load(theFMap(Einf2.FirstFace()).Geometry());
1927     }
1928   }
1929   else {
1930     if (Einf2.FirstFace().IsSame(Einf.FirstFace())) {
1931       AS.Load(theFMap(Einf2.SecondFace()).Geometry());
1932     }
1933     else {
1934       AS.Load(theFMap(Einf2.FirstFace()).Geometry());
1935     }
1936   }
1937   return Standard_True;
1938 }
1939
1940
1941 //=======================================================================
1942 //function : Parameter
1943 //purpose  : 
1944 //=======================================================================
1945
1946 static Standard_Real Parameter(const Handle(Geom_Curve)& C,
1947                                const gp_Pnt& P,
1948                                Standard_Integer& done)
1949 {
1950   done = 0;
1951   Handle(Geom_Curve) cbase = C;
1952   Handle(Standard_Type) ctyp = C->DynamicType();
1953   if (ctyp == STANDARD_TYPE(Geom_TrimmedCurve)) {
1954     cbase = Handle(Geom_TrimmedCurve)::DownCast(C)->BasisCurve();
1955     ctyp = cbase->DynamicType();
1956   }
1957   Standard_Real param;
1958   if (ctyp == STANDARD_TYPE(Geom_Line)) {
1959     param = ElCLib::Parameter(Handle(Geom_Line)::DownCast(cbase)->Lin(),P);
1960   }
1961   else if (ctyp == STANDARD_TYPE(Geom_Circle)) {
1962     param = ElCLib::Parameter(Handle(Geom_Circle)::DownCast(cbase)->Circ(),P);
1963     if (Abs(2.*M_PI-param) <=Epsilon(2.*M_PI)) {
1964       param = 0.;
1965     }
1966   }
1967   else if (ctyp == STANDARD_TYPE(Geom_Ellipse)) {
1968     param = ElCLib::Parameter(Handle(Geom_Ellipse)::DownCast(cbase)->Elips(),P);
1969     if (Abs(2.*M_PI-param) <=Epsilon(2.*M_PI)) {
1970       param = 0.;
1971     }
1972   }
1973   else if (ctyp == STANDARD_TYPE(Geom_Parabola)) {
1974     param = ElCLib::Parameter(Handle(Geom_Parabola)::DownCast(cbase)->Parab(),P);
1975   }
1976   else if (ctyp == STANDARD_TYPE(Geom_Hyperbola)) {
1977     param = ElCLib::Parameter(Handle(Geom_Hyperbola)::DownCast(cbase)->Hypr(),P);
1978   }
1979   else {
1980     GeomAdaptor_Curve TheCurve(C);
1981     Extrema_ExtPC myExtPC(P,TheCurve);
1982     if (!myExtPC.IsDone()) {
1983       Standard_Failure::Raise("Draft_Modification_1::Parameter: ExtremaPC not done.");
1984     }
1985     if (myExtPC.NbExt() >= 1) {
1986       Standard_Real Dist2, Dist2Min = myExtPC.SquareDistance(1);
1987       Standard_Integer j, jmin = 1;
1988       for (j = 2; j <= myExtPC.NbExt(); j++) {
1989         Dist2 = myExtPC.SquareDistance(j);
1990         if (Dist2 < Dist2Min) {
1991           Dist2Min = Dist2;
1992           jmin = j;
1993         }
1994       }
1995       param = myExtPC.Point(jmin).Parameter();
1996     }
1997     else {
1998       Standard_Real dist1_2,dist2_2;
1999       gp_Pnt p1b,p2b;
2000       myExtPC.TrimmedSquareDistances(dist1_2,dist2_2,p1b,p2b);
2001       if (dist1_2 < dist2_2) {
2002         done = -1;
2003         param = TheCurve.FirstParameter();
2004       }
2005       else {
2006         done = 1;
2007         param = TheCurve.LastParameter();
2008       }
2009     }
2010
2011     if (cbase->IsPeriodic()) {
2012       Standard_Real Per  = cbase->Period();
2013       Standard_Real Tolp = Precision::Parametric(Precision::Confusion());  
2014       if (Abs(Per-param) <= Tolp) {
2015         param = 0.;
2016       }
2017     }
2018   }
2019   return param;
2020 }
2021
2022 //=======================================================================
2023 //function : SmartParameter
2024 //purpose  : 
2025 //=======================================================================
2026
2027 static Standard_Real SmartParameter(Draft_EdgeInfo& Einf,
2028                                     const Standard_Real EdgeTol,
2029                                     const gp_Pnt& Pnt,
2030                                     const Standard_Integer sign,
2031                                     const Handle(Geom_Surface)& S1,
2032                                     const Handle(Geom_Surface)& S2)
2033 {
2034   Handle( Geom2d_Curve ) NewC2d;
2035   Standard_Real Tol = Precision::Confusion();
2036   Standard_Real Etol = EdgeTol;
2037
2038   Handle( Geom2d_Curve ) pcu1 = Einf.FirstPC();
2039   Handle( Geom2d_Curve ) pcu2 = Einf.SecondPC();
2040
2041   if (pcu1.IsNull())
2042     {
2043       Handle( Geom_Curve ) theCurve = Einf.Geometry();
2044       pcu1 = GeomProjLib::Curve2d( theCurve, theCurve->FirstParameter(), theCurve->LastParameter(), S1, Etol );
2045       Einf.ChangeFirstPC() = pcu1;
2046     }
2047   if (pcu2.IsNull())
2048     {
2049       Handle( Geom_Curve ) theCurve = Einf.Geometry();
2050       pcu2 = GeomProjLib::Curve2d( theCurve, theCurve->FirstParameter(), theCurve->LastParameter(), S2, Etol );
2051       Einf.ChangeSecondPC() = pcu2;
2052     }
2053
2054   GeomAPI_ProjectPointOnSurf Projector( Pnt, S1 );
2055   Standard_Real U, V;
2056   Projector.LowerDistanceParameters( U, V );
2057   
2058   NewC2d = Einf.FirstPC();
2059   if (NewC2d->DynamicType() == STANDARD_TYPE(Geom2d_TrimmedCurve))
2060     NewC2d = (Handle(Geom2d_TrimmedCurve)::DownCast(NewC2d))->BasisCurve();
2061   
2062   gp_Pnt2d P2d( U, V );
2063   Geom2dAPI_ProjectPointOnCurve Projector2d( P2d, NewC2d );
2064   if (Projector2d.NbPoints() == 0 || Projector2d.LowerDistance() > Tol)
2065     {
2066       Handle( Geom2d_BSplineCurve ) BCurve;
2067       if (NewC2d->DynamicType() != STANDARD_TYPE(Geom2d_BSplineCurve))
2068         BCurve = Geom2dConvert::CurveToBSplineCurve( NewC2d );
2069       else
2070         BCurve = Handle( Geom2d_BSplineCurve )::DownCast( NewC2d );
2071       if (sign == -1)
2072         {
2073           TColgp_Array1OfPnt2d PntArray( 1, 2 );
2074           PntArray(1) = P2d;
2075           PntArray(2) = BCurve->Pole(1);
2076           Handle( Geom2d_BezierCurve ) Patch = new Geom2d_BezierCurve( PntArray );
2077           Geom2dConvert_CompCurveToBSplineCurve Concat( BCurve, Convert_QuasiAngular );
2078           Concat.Add( Patch, Tol, Standard_False );
2079           BCurve = Concat.BSplineCurve();
2080         }
2081       else
2082         {
2083           TColgp_Array1OfPnt2d PntArray( 1, 2 );
2084           PntArray(1) = BCurve->Pole( BCurve->NbPoles() );
2085           PntArray(2) = P2d;
2086           Handle( Geom2d_BezierCurve ) Patch = new Geom2d_BezierCurve( PntArray );
2087           Geom2dConvert_CompCurveToBSplineCurve Concat( BCurve, Convert_QuasiAngular );
2088           Concat.Add( Patch, Tol, Standard_True );
2089           BCurve = Concat.BSplineCurve();
2090         }
2091       NewC2d = BCurve;
2092     }
2093   Einf.ChangeFirstPC() = NewC2d;
2094   Handle( Geom2dAdaptor_HCurve ) hcur = new Geom2dAdaptor_HCurve( NewC2d );
2095   Handle( GeomAdaptor_HSurface ) hsur = new GeomAdaptor_HSurface( S1 );
2096   Adaptor3d_CurveOnSurface cons( hcur, hsur );
2097   Handle( Adaptor3d_HCurveOnSurface ) hcons = new Adaptor3d_HCurveOnSurface( cons );
2098   Handle( GeomAdaptor_HSurface ) hsur2 = new GeomAdaptor_HSurface( S2 );
2099   ProjLib_CompProjectedCurve ProjCurve( hsur2, hcons, Tol, Tol );
2100   Handle(ProjLib_HCompProjectedCurve) HProjector = new ProjLib_HCompProjectedCurve();
2101   HProjector->Set( ProjCurve );
2102   Standard_Real Udeb, Ufin;
2103   ProjCurve.Bounds(1, Udeb, Ufin);
2104   Standard_Integer MaxSeg = 20 + HProjector->NbIntervals(GeomAbs_C3);
2105   Approx_CurveOnSurface appr( HProjector, hsur2, Udeb, Ufin, Tol,
2106                               GeomAbs_C1, 10, MaxSeg, 
2107                               Standard_False, Standard_False );
2108   Einf.ChangeSecondPC() = appr.Curve2d();
2109   Einf.ChangeGeometry() = appr.Curve3d();
2110   Einf.SetNewGeometry( Standard_True );
2111
2112   if (sign == -1)
2113     return Einf.Geometry()->FirstParameter();
2114   else
2115     return Einf.Geometry()->LastParameter();
2116
2117 }
2118
2119 //=======================================================================
2120 //function : Orientation
2121 //purpose  : 
2122 //=======================================================================
2123
2124 static TopAbs_Orientation Orientation(const TopoDS_Shape& S,
2125                                       const TopoDS_Face& F)
2126 {
2127 //
2128 // change porting NT
2129 //
2130   TopExp_Explorer expl ;
2131   expl.Init(S,
2132             TopAbs_FACE) ;
2133   while (expl.More()) {
2134     if (TopoDS::Face(expl.Current()).IsSame(F)) {
2135       return expl.Current().Orientation();
2136     }
2137     expl.Next();
2138   }
2139   return TopAbs_FORWARD;
2140 }
2141
2142
2143 //=======================================================================
2144 //function : FindRotation
2145 //purpose  : 
2146 //=======================================================================
2147
2148 static Standard_Boolean FindRotation(const gp_Pln& Pl,
2149                                      const TopAbs_Orientation Oris,
2150                                      const gp_Dir& Direction,
2151                                      const Standard_Real Angle,
2152                                      const gp_Pln& NeutralPlane,
2153                                      gp_Ax1& Axe,
2154                                      Standard_Real& theta)
2155 {                                    
2156   IntAna_QuadQuadGeo i2pl(Pl,NeutralPlane,
2157                           Precision::Angular(),Precision::Confusion());
2158   
2159   if (i2pl.IsDone() && i2pl.TypeInter() == IntAna_Line) {
2160     gp_Lin li = i2pl.Line(1);
2161     // Try to turn around this line
2162     gp_Dir nx = li.Direction();
2163     gp_Dir ny = Pl.Axis().Direction().Crossed(nx);
2164     Standard_Real a = Direction.Dot(nx);
2165     if (Abs(a) <=1-Precision::Angular()) { 
2166       Standard_Real b = Direction.Dot(ny);
2167       Standard_Real c = Direction.Dot(Pl.Axis().Direction());
2168       Standard_Boolean direct(Pl.Position().Direct());
2169       if ((direct && Oris == TopAbs_REVERSED) ||
2170           (!direct && Oris == TopAbs_FORWARD)) {
2171         b = -b;
2172         c = -c;
2173       }
2174       Standard_Real denom = Sqrt(1-a*a);
2175       Standard_Real Sina = Sin(Angle);
2176       if (denom>Abs(Sina)) {
2177         Standard_Real phi = ATan2(b/denom,c/denom);
2178         Standard_Real theta0 = ACos(Sina/denom); 
2179         theta = theta0 - phi;
2180         if (Cos(theta) <0.) {
2181           theta = -theta0 -phi;
2182         }
2183         //  modified by NIZHNY-EAP Tue Nov 16 15:51:38 1999 ___BEGIN___
2184         while (Abs(theta)>M_PI) {
2185           theta = theta + M_PI*(theta<0 ? 1 : -1);
2186         }
2187         //  modified by NIZHNY-EAP Tue Nov 16 15:53:32 1999 ___END___
2188         Axe = li.Position();
2189         return Standard_True;
2190       }
2191     }
2192   }
2193   return Standard_False;
2194 }
2195