0024510: Remove unused local variables
[occt.git] / src / LocOpe / LocOpe_SplitDrafts.cxx
1 // Created on: 1996-10-02
2 // Created by: Jacques GOUSSARD
3 // Copyright (c) 1996-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
9 // under the terms of the GNU Lesser General Public 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 #include <LocOpe_SplitDrafts.ixx>
18
19 #include <TopExp_Explorer.hxx>
20 #include <TopTools_DataMapIteratorOfDataMapOfShapeListOfShape.hxx>
21
22 #include <BRep_Builder.hxx>
23 #include <BRep_Tool.hxx>
24
25 #include <BRepTools_Substitution.hxx>
26
27 #include <LocOpe_WiresOnShape.hxx>
28 #include <LocOpe_Spliter.hxx>
29 #include <LocOpe_SplitShape.hxx>
30 #include <LocOpe_FindEdges.hxx>
31 #include <LocOpe_BuildShape.hxx>
32
33
34 #include <TopoDS_Vertex.hxx>
35 #include <TopoDS_Edge.hxx>
36
37 #include <TopTools_IndexedDataMapOfShapeListOfShape.hxx>
38 #include <TopTools_DataMapIteratorOfDataMapOfShapeShape.hxx>
39
40 #include <IntAna_QuadQuadGeo.hxx>
41 #include <IntCurveSurface_HInter.hxx>
42 #include <IntCurveSurface_IntersectionPoint.hxx>
43 #include <IntCurveSurface_IntersectionSegment.hxx>
44 #include <GeomInt_IntSS.hxx>
45 #include <Extrema_ExtPC.hxx>
46 #include <GeomAdaptor_HCurve.hxx>
47 #include <GeomAdaptor_HSurface.hxx>
48 #include <GeomAdaptor_Curve.hxx>
49 #include <GeomAdaptor_Surface.hxx>
50
51 #include <Geom_Surface.hxx>
52 #include <Geom_RectangularTrimmedSurface.hxx>
53 #include <Geom_Plane.hxx>
54 #include <Geom_Curve.hxx>
55 #include <Geom_TrimmedCurve.hxx>
56 #include <Geom_Line.hxx>
57
58 #include <Geom2d_Curve.hxx>
59 #include <Geom2d_Line.hxx>
60 #include <gp_Pnt2d.hxx>
61 #include <gp_Vec2d.hxx>
62
63 #include <GeomFill_Pipe.hxx>
64
65 #include <GProp_GProps.hxx>
66
67 #include <Standard_ConstructionError.hxx>
68
69 #include <TopoDS.hxx>
70 #include <TopExp.hxx>
71 #include <Precision.hxx>
72 #include <BRepGProp.hxx>
73 #include <gp.hxx>
74
75 static Standard_Boolean NewPlane(const TopoDS_Face&,
76                                  const gp_Dir&,
77                                  const gp_Pln&,
78                                  const Standard_Real,
79                                  gp_Pln&,
80                                  gp_Ax1&,
81                                  const Standard_Boolean);
82
83 static void MakeFace(TopoDS_Face&,
84                      TopTools_ListOfShape&);
85
86 static TopoDS_Edge  NewEdge(const TopoDS_Edge&,
87                             const TopoDS_Face&,
88                             const Handle(Geom_Surface)&,
89                             const TopoDS_Vertex&,
90                             const TopoDS_Vertex&);
91
92
93 static Standard_Boolean Contains(const TopTools_ListOfShape&,
94                                  const TopoDS_Shape&);
95
96
97 //=======================================================================
98 //function : Init
99 //purpose  : 
100 //=======================================================================
101
102 void LocOpe_SplitDrafts::Init(const TopoDS_Shape& S)
103 {
104   myShape = S;
105   myResult.Nullify();
106   myMap.Clear();
107 }
108
109 //=======================================================================
110 //function : Perform
111 //purpose  : 
112 //=======================================================================
113
114 void LocOpe_SplitDrafts::Perform(const TopoDS_Face& F,
115                                  const TopoDS_Wire& W,
116                                  const gp_Dir& Extr,
117                                  const gp_Pln& NPl,
118                                  const Standard_Real Angle)
119 {
120   Perform(F,W,Extr,NPl,Angle,Extr,NPl,Angle,Standard_True,Standard_False);
121 }
122
123
124 //=======================================================================
125 //function : Perform
126 //purpose  : 
127 //=======================================================================
128
129 void LocOpe_SplitDrafts::Perform(const TopoDS_Face& F,
130                                  const TopoDS_Wire& W,
131                                  const gp_Dir& Extrg,
132                                  const gp_Pln& NPlg,
133                                  const Standard_Real Angleg,
134                                  const gp_Dir& Extrd,
135                                  const gp_Pln& NPld,
136                                  const Standard_Real Angled,
137                                  const Standard_Boolean ModLeft,
138                                  const Standard_Boolean ModRight)
139
140 {
141   Standard_Integer j ;
142
143   myResult.Nullify();
144   myMap.Clear();
145   if (myShape.IsNull() || F.IsNull() || W.IsNull()) {
146     Standard_NullObject::Raise();
147   }    
148
149   if (!ModLeft && !ModRight) {
150     Standard_ConstructionError::Raise();
151   }
152
153   TopAbs_Orientation OriF = TopAbs_FORWARD;
154
155   Standard_Boolean FinS = Standard_False;
156   TopExp_Explorer exp,exp2;
157   for (exp.Init(myShape,TopAbs_FACE); exp.More(); exp.Next()) {
158     const TopoDS_Shape& fac = exp.Current();
159     TopTools_ListOfShape thelist;
160     myMap.Bind(fac, thelist);
161     if (fac.IsSame(F)) {
162       OriF = fac.Orientation();
163       FinS = Standard_True;
164     }
165   }
166
167   if (!FinS) {
168     cout << "LocOpe_SplitDrafts:!Fins Standard_ConstructionError::Raise()" << endl;
169     Standard_ConstructionError::Raise();
170   }    
171
172   gp_Pln NewPlg,NewPld;
173   gp_Ax1 NormalFg,NormalFd;
174   TopoDS_Shape aLocalFace = F.Oriented(OriF);
175
176   if (!NewPlane(TopoDS::Face(aLocalFace),
177                 Extrg,NPlg,Angleg,NewPlg,NormalFg,ModLeft) ||
178       !NewPlane(TopoDS::Face(aLocalFace),
179                 Extrd,NPld,Angled,NewPld,NormalFd,ModRight)) {
180     //  if (!NewPlane(TopoDS::Face(F.Oriented(OriF)),
181     //          Extrg,NPlg,Angleg,NewPlg,NormalFg,ModLeft) ||
182 //      !NewPlane(TopoDS::Face(F.Oriented(OriF)),
183     //          Extrd,NPld,Angled,NewPld,NormalFd,ModRight)) {
184     return;
185   }
186
187
188   TopTools_ListIteratorOfListOfShape itl;
189   BRep_Builder B;
190
191   Handle(Geom_Surface) NewSg = new Geom_Plane(NewPlg);
192   Handle(Geom_Surface) NewSd = new Geom_Plane(NewPld);
193   Handle(Geom_Line) theLinePipe = new Geom_Line(NormalFg); // ou NormalFd
194   GeomInt_IntSS i2s(NewSg,NewSd,Precision::Confusion());
195
196   TopTools_MapOfShape theMap;
197   Handle(GeomAdaptor_HCurve) HAC = new GeomAdaptor_HCurve;
198   GeomAdaptor_Curve AC;
199   Handle(GeomAdaptor_HSurface) HAS = new GeomAdaptor_HSurface;
200   GeomAdaptor_Surface AS;
201   IntCurveSurface_HInter intcs;
202
203   TopoDS_Wire theW = W;
204   if (i2s.IsDone() && i2s.NbLines() > 0) {
205     // on split le wire" << endl;
206
207     GeomFill_Pipe thePipe;
208     thePipe.GenerateParticularCase(Standard_True);
209     thePipe.Init(theLinePipe,i2s.Line(1));
210     thePipe.Perform(Standard_True);
211
212     Handle(Geom_Surface) Spl = thePipe.Surface();
213     AS.Load(Spl);
214     HAS->Set(AS);
215     
216     LocOpe_SplitShape splw(W);
217
218     for (exp.Init(W,TopAbs_EDGE); exp.More(); exp.Next()) {
219       const TopoDS_Edge& edg = TopoDS::Edge(exp.Current());
220       if (theMap.Add(edg)) {
221         TopLoc_Location Loc;
222         Standard_Real f,l;
223         Handle(Geom_Curve) C = BRep_Tool::Curve(edg,f,l);
224         AC.Load(C);
225         HAC->Set(AC);
226         intcs.Perform(HAC,HAS);
227         if (!intcs.IsDone()) {
228           continue; // voir ce qu`on peut faire de mieux
229         }
230
231         if (intcs.NbSegments() >= 2) {
232           continue; // Not yet implemented...and probably never"
233         }
234
235         if (intcs.NbSegments() == 1) {
236           const IntCurveSurface_IntersectionPoint& P1 = 
237             intcs.Segment(1).FirstPoint();
238           const IntCurveSurface_IntersectionPoint& P2 = 
239             intcs.Segment(1).SecondPoint();
240           const gp_Pnt& pf = P1.Pnt();
241           const gp_Pnt& pl = P2.Pnt();
242           TopoDS_Vertex Vf,Vl;
243           TopExp::Vertices(edg,Vf,Vl);
244           gp_Pnt Pf = BRep_Tool::Pnt(Vf);
245           gp_Pnt Pl = BRep_Tool::Pnt(Vl);
246           Standard_Real Tolf = BRep_Tool::Tolerance(Vf);
247           Standard_Real Toll = BRep_Tool::Tolerance(Vl);
248           Tolf *= Tolf;
249           Toll *= Toll;
250
251           Standard_Real dff = pf.SquareDistance(Pf);
252           Standard_Real dfl = pf.SquareDistance(Pl);
253           Standard_Real dlf = pl.SquareDistance(Pf);
254           Standard_Real dll = pl.SquareDistance(Pl);
255
256           if ((dff <= Tolf && dll <= Toll) ||
257               (dlf <= Tolf && dfl <= Toll)) {
258             continue;
259           }
260           else {
261             // on segmente edg en pf et pl
262             TopoDS_Vertex Vnewf,Vnewl;
263             B.MakeVertex(Vnewf,pf,Precision::Confusion());
264             B.MakeVertex(Vnewl,pl,Precision::Confusion());
265             if (P1.W() >= f && P1.W() <= l &&
266                 P2.W() >= f && P2.W() <= l) {
267               splw.Add(Vnewf,P1.W(),edg);
268               splw.Add(Vnewl,P2.W(),edg);
269             }
270             else {
271               continue;
272             }
273           }
274         }
275         else if (intcs.NbPoints() != 0) {
276           TopoDS_Vertex Vf,Vl;
277           TopExp::Vertices(edg,Vf,Vl);
278           gp_Pnt Pf = BRep_Tool::Pnt(Vf);
279           gp_Pnt Pl = BRep_Tool::Pnt(Vl);
280           Standard_Real Tolf = BRep_Tool::Tolerance(Vf);
281           Standard_Real Toll = BRep_Tool::Tolerance(Vl);
282           Tolf *= Tolf;
283           Toll *= Toll;
284
285           for (Standard_Integer i = 1; i <= intcs.NbPoints(); i++) {
286             const IntCurveSurface_IntersectionPoint& Pi = intcs.Point(i);
287             const gp_Pnt& pi = Pi.Pnt();
288             Standard_Real dif = pi.SquareDistance(Pf);
289             Standard_Real dil = pi.SquareDistance(Pl);
290             if (dif <= Tolf) {
291             }
292             else if (dil <= Toll) {
293             }
294             else {
295               if (Pi.W() >= f && Pi.W() <= l) {
296                 TopoDS_Vertex Vnew;
297                 B.MakeVertex(Vnew,pi,Precision::Confusion());
298                 splw.Add(Vnew,Pi.W(),edg);
299               }
300             }
301           }
302         }
303       }
304     }
305
306     const TopTools_ListOfShape& lres = splw.DescendantShapes(W);
307     if (lres.Extent() != 1) {
308       return;
309     }
310
311     if (!W.IsSame(lres.First())) {
312       theW.Nullify();
313       theW = TopoDS::Wire(lres.First());
314     }
315
316     for (exp.ReInit(); exp.More(); exp.Next()) {
317       if (!myMap.IsBound(exp.Current())) {
318         TopTools_ListOfShape thelist1;
319         myMap.Bind(exp.Current(), thelist1);
320         for (itl.Initialize(splw.DescendantShapes(exp.Current())); 
321              itl.More(); itl.Next()) {
322           myMap(exp.Current()).Append(itl.Value());
323         }
324         for (exp2.Init(exp.Current(),TopAbs_VERTEX);exp2.More();exp2.Next()) {
325           if (!myMap.IsBound(exp2.Current())) {
326             TopTools_ListOfShape thelist2;
327             myMap.Bind(exp2.Current(), thelist2);
328             myMap(exp2.Current()).Append(exp2.Current());
329           }
330         }
331       }
332     }
333   }
334   else {
335     for (exp.Init(W,TopAbs_EDGE); exp.More(); exp.Next()) {
336       if (!myMap.IsBound(exp.Current())) {
337         TopTools_ListOfShape thelist3;
338         myMap.Bind(exp.Current(), thelist3);
339         myMap(exp.Current()).Append(exp.Current());
340         for (exp2.Init(exp.Current(),TopAbs_VERTEX);exp2.More();exp2.Next()) {
341           if (!myMap.IsBound(exp2.Current())) {
342             TopTools_ListOfShape thelist4;
343             myMap.Bind(exp2.Current(), thelist4);
344             myMap(exp2.Current()).Append(exp2.Current());
345           }
346         }
347       }
348     }
349   }
350
351   // On split la face par le wire
352   
353   Handle(LocOpe_WiresOnShape) WonS = new LocOpe_WiresOnShape(myShape);
354   LocOpe_Spliter Spls(myShape);
355   WonS->Bind(theW,F);
356
357 // JAG Le code suivant marchera apres integration de thick0
358 //  LocOpe_FindEdges fined(W,F);
359 //  for (fined.InitIterator(); fined.More(); fined.Next()) {
360 //    WonS->Bind(fined.EdgeFrom(),fined.EdgeTo());
361 //  }
362
363   Spls.Perform(WonS);
364   if (!Spls.IsDone()) {
365     return;
366   }
367
368   TopoDS_Shape Res = Spls.ResultingShape();
369   const TopTools_ListOfShape& theLeft = Spls.DirectLeft();
370
371   // Descendants
372   for (exp.Init(myShape,TopAbs_FACE); exp.More(); exp.Next()) {
373     const TopoDS_Shape& fac = exp.Current();
374     for (itl.Initialize(Spls.DescendantShapes(fac)); itl.More(); itl.Next()) {
375       myMap(fac).Append(itl.Value());
376     }
377   }
378
379   TopTools_DataMapOfShapeShape MapW;
380   for (exp.Init(theW,TopAbs_EDGE); exp.More(); exp.Next()) {
381     if (!MapW.IsBound(exp.Current())) {
382       MapW.Bind(exp.Current(),TopoDS_Shape());
383       for (exp2.Init(exp.Current(),TopAbs_VERTEX); exp2.More(); exp2.Next()) {
384         if (!MapW.IsBound(exp2.Current())) {
385           MapW.Bind(exp2.Current(),TopoDS_Shape());
386         }
387
388       }
389     }
390   }
391   
392
393
394   TopTools_IndexedDataMapOfShapeListOfShape theMapEF;
395   TopExp::MapShapesAndAncestors(Res,TopAbs_EDGE,TopAbs_FACE,theMapEF);
396
397   // On stocke les geometries potentiellement generees par les edges
398   TopTools_IndexedDataMapOfShapeShape MapEV; // genere
399   TopTools_DataMapOfShapeListOfShape MapSg,MapSd; // image a gauche et a droite
400
401   Standard_Integer Nbedges,index;
402   for (itl.Initialize(myMap(F)); itl.More(); itl.Next()) {
403     const TopoDS_Shape& fac = TopoDS::Face(itl.Value());
404     for (exp.Init(fac,TopAbs_EDGE); exp.More(); exp.Next()) {
405       const TopoDS_Shape& edg = exp.Current();
406       if (MapEV.FindIndex(edg) != 0) {
407         continue;
408       }
409       if (MapW.IsBound(edg)) { // edge du wire initial
410         TopLoc_Location Loc;
411         Standard_Real f,l;
412         Handle(Geom_Curve) C = BRep_Tool::Curve(TopoDS::Edge(edg),Loc,f,l);
413         if (C->DynamicType() == STANDARD_TYPE(Geom_TrimmedCurve)) {
414           C = Handle(Geom_TrimmedCurve)::DownCast(C)->BasisCurve();
415         }
416         C = Handle(Geom_Curve)::
417           DownCast(C->Transformed(Loc.Transformation()));
418
419         GeomFill_Pipe thePipe;
420         thePipe.GenerateParticularCase(Standard_True);
421         thePipe.Init(theLinePipe,C);
422         thePipe.Perform(Standard_True);
423
424         Handle(Geom_Surface) thePS = thePipe.Surface();
425         if (thePS->DynamicType() == STANDARD_TYPE(Geom_RectangularTrimmedSurface)) {
426           thePS = Handle(Geom_RectangularTrimmedSurface)::DownCast(thePS)
427             ->BasisSurface();
428         }
429
430         TopoDS_Face NewFace;
431         B.MakeFace(NewFace,thePS,Precision::Confusion());
432         MapEV.Add(edg,NewFace);
433       }
434       else { // on recupere la face.
435         index = theMapEF.FindIndex(edg);
436         if (theMapEF(index).Extent() != 2) {
437           return; // NotDone
438         }
439         TopoDS_Face theFace;
440         if (theMapEF(index).First().IsSame(fac)) {
441           MapEV.Add(edg,theMapEF(index).Last());
442         }
443         else {
444           MapEV.Add(edg,theMapEF(index).First());
445         }
446       }
447     }
448   }
449
450
451   TopTools_DataMapOfShapeShape MapSonS;
452
453   Nbedges = MapEV.Extent();
454   for (index = 1; index <= Nbedges; index++) {
455     for (exp.Init(MapEV.FindKey(index),TopAbs_VERTEX); 
456          exp.More(); exp.Next()) {
457       const TopoDS_Vertex& vtx = TopoDS::Vertex(exp.Current());
458       if (MapEV.FindIndex(vtx)!= 0) {
459         continue;
460       }
461
462       // Localisation du vertex :
463       //    - entre 2 edges d`origine : on recupere l`edge qui n`est 
464       //                                pas dans F
465       //    - entre 2 edges du wire   : droite
466       //    - mixte                   : intersection de surfaces
467       for ( j = 1; j<=Nbedges; j++) {
468         if (j == index) {
469           continue;
470         }
471         for (exp2.Init(MapEV.FindKey(j),TopAbs_VERTEX);
472              exp2.More(); exp2.Next()) {
473           const TopoDS_Shape& vtx2 = exp2.Current();
474           if (vtx2.IsSame(vtx)) {
475             break;
476           }
477         }
478         if (exp2.More()) {
479           break;
480         }
481       }
482       Standard_Integer Choice = 0;
483       const TopoDS_Shape& edg1 = MapEV.FindKey(index);
484       TopoDS_Shape edg2;
485       if (j <= Nbedges) {
486         edg2 = MapEV.FindKey(j);
487       }
488       else {
489         edg2 = edg1;
490       }
491       if (MapW.IsBound(edg1)) {
492         if (j>Nbedges) { // doit correspondre a edge ferme
493           Choice = 2; // droite
494         }
495         else if (MapW.IsBound(MapEV.FindKey(j))) {
496           Choice = 2; // droite
497         }
498         else {
499           Choice = 3; // mixte
500         }
501       }
502       else {
503         if (j>Nbedges) { // doit correspondre a edge ferme
504           Choice = 1; // edge a retrouver
505         }
506         else if (!MapW.IsBound(MapEV.FindKey(j))) {
507           Choice = 1; // edge a retrouver
508         }
509         else {
510           Choice = 3; // mixte
511         }
512       }
513       Handle(Geom_Curve) Newc;
514       Handle(Geom2d_Curve) newCs1,newCs2;
515       Standard_Real knownp=0;
516       TopoDS_Edge Ebind;
517       switch (Choice) {
518       case 1:
519         {
520           for (exp2.Init(Res,TopAbs_EDGE); exp2.More(); exp2.Next()) {
521             if (exp2.Current().IsSame(edg1) || exp2.Current().IsSame(edg2)) {
522               continue;
523             }
524 //          for (TopExp_Explorer exp3(exp2.Current().Oriented(TopAbs_FORWARD),
525             TopExp_Explorer exp3(exp2.Current().Oriented(TopAbs_FORWARD),
526                                       TopAbs_VERTEX) ;
527             for ( ; exp3.More(); exp3.Next())  {
528               if (exp3.Current().IsSame(vtx)) {
529                 break;
530               }
531             }
532             if (exp3.More()) {
533               break;
534             }
535           }
536           if (exp2.More()) {
537             Standard_Real f,l;
538             TopLoc_Location Loc;
539             Newc = BRep_Tool::Curve(TopoDS::Edge(exp2.Current()),Loc,f,l);
540             Newc = Handle(Geom_Curve)::DownCast
541               (Newc->Transformed(Loc.Transformation()));
542             Ebind = TopoDS::Edge(exp2.Current());
543             knownp = BRep_Tool::Parameter(vtx,Ebind);
544           }
545           else { // droite ??? il vaudrait mieux sortir
546             return;
547
548 //          gp_Lin theLine(NormalFg);
549 //          theLine.Translate(NormalF.Location(),BRep_Tool::Pnt(vtx));
550 //          Newc = new Geom_Line(theLine);
551 //          knownp = 0.;
552           }
553         }
554         break;
555       case 2:
556         {
557           gp_Lin theLine(NormalFg);
558           theLine.Translate(NormalFg.Location(),BRep_Tool::Pnt(vtx));
559           Newc = new Geom_Line(theLine);
560           knownp = 0.;
561         }
562         break;
563       case 3:
564         {
565           const TopoDS_Face& F1 = TopoDS::Face(MapEV.FindFromKey(edg1));
566           const TopoDS_Face& F2 = TopoDS::Face(MapEV.FindFromKey(edg2));
567           Handle(Geom_Surface) S1 = BRep_Tool::Surface(F1);
568           Handle(Geom_Surface) S2 = BRep_Tool::Surface(F2);
569           Standard_Boolean AppS1 = Standard_False;
570           Standard_Boolean AppS2 = Standard_False;
571           if (S1->DynamicType() != STANDARD_TYPE(Geom_Plane)) {
572             AppS1 = Standard_True;
573           }
574           if (S2->DynamicType() != STANDARD_TYPE(Geom_Plane)) {
575             AppS2 = Standard_True;
576           }
577           i2s.Perform(S1,S2,Precision::Confusion(),Standard_True,AppS1,AppS2);
578           if (!i2s.IsDone() || i2s.NbLines() <= 0) {
579             return;
580           }
581         
582           Standard_Real pmin=0, Dist2, Dist2Min, Glob2Min = RealLast();
583           GeomAdaptor_Curve TheCurve;
584
585           Standard_Integer i,imin,k;
586           gp_Pnt pv = BRep_Tool::Pnt(vtx);
587           imin = 0;
588           for (i=1; i<= i2s.NbLines(); i++) {
589             TheCurve.Load(i2s.Line(i));
590             Extrema_ExtPC myExtPC(pv,TheCurve);
591
592             if (myExtPC.IsDone()) {
593               gp_Pnt p1b,p2b;
594               Standard_Real thepmin = TheCurve.FirstParameter();
595               myExtPC.TrimmedSquareDistances(Dist2Min,Dist2,p1b,p2b);
596               if (Dist2 < Dist2Min) {
597                 thepmin = TheCurve.LastParameter();
598               }
599               for (k=1; k<=myExtPC.NbExt(); k++) {
600                 Dist2 = myExtPC.SquareDistance(k);
601                 if (Dist2 < Dist2Min) {
602                   Dist2Min = Dist2;
603                   thepmin = myExtPC.Point(k).Parameter();
604                 }
605               }
606               
607               if (Dist2Min  < Glob2Min) {
608                 Glob2Min = Dist2Min;
609                 pmin = thepmin;
610                 imin = i;
611               }
612             }
613           }
614           if (imin == 0) {
615             return;
616           }
617           
618           Newc = i2s.Line(imin);
619           knownp = pmin;
620           if (AppS1) {
621             newCs1 = i2s.LineOnS1(imin);
622           }
623           if (AppS2) {
624             newCs2 = i2s.LineOnS2(imin);
625           }
626         }
627         break;
628       }
629
630
631       // Determination des vertex par intersection sur Plg ou/et Pld
632
633       AC.Load(Newc);
634       HAC->Set(AC);
635       Standard_Integer nbfois = 2;
636       TopoDS_Vertex vtx1,vtx2;
637       Standard_Real p1=0,p2=0;
638       Standard_Boolean IsLeft=Standard_False;
639       if (Choice == 1) { 
640         // edge retrouve : on ne fait qu`une seule intersection
641         // il faut utiliser Plg ou Pld
642
643         Standard_Integer indedgf = theMapEF.FindIndex(edg1);
644         for (itl.Initialize(theMapEF(indedgf)); itl.More(); itl.Next()) {
645           if (Contains(myMap(F),itl.Value())) {
646             if (Contains(theLeft,itl.Value())) {
647               AS.Load(NewSg);
648               IsLeft = Standard_True;
649             }
650             else {
651               AS.Load(NewSd);
652               IsLeft = Standard_False;
653             }
654             
655             nbfois = 1;
656             vtx2 = vtx;
657             p2 = knownp;
658             break;
659           }
660         }
661         if (!itl.More()) {
662           cout << "LocOpe_SplitDrafts: betite probleme "<< endl;
663           return;
664         }
665
666       }
667       else {
668         AS.Load(NewSg);
669       }
670
671       for (Standard_Integer it = 1; it<=nbfois; it++) {
672         if (it == 2) {
673           AS.Load(NewSd);
674         }
675         HAS->Set(AS);
676
677         intcs.Perform(HAC,HAS);
678         if (!intcs.IsDone()) {
679           return; // voir ce qu`on peut faire de mieux
680         }
681         Standard_Integer imin = 1;
682         Standard_Real delta = Abs(knownp - intcs.Point(1).W());
683         for (Standard_Integer i = 2;  i<= intcs.NbPoints(); i++) {
684           Standard_Real newdelta =  Abs(knownp - intcs.Point(i).W());
685           if (newdelta < delta) {
686             imin = i;
687             delta = newdelta;
688           }
689         }
690         if (it == 1) {
691           B.MakeVertex(vtx1,intcs.Point(imin).Pnt(),Precision::Confusion());
692           p1 = intcs.Point(imin).W();
693           knownp = p1;
694         }
695         else {
696           B.MakeVertex(vtx2,intcs.Point(imin).Pnt(),Precision::Confusion());
697           p2 = intcs.Point(imin).W();
698         }
699       }
700       if (Abs(p1-p2) > Precision::PConfusion()) {
701         TopoDS_Edge NewEdge;
702         B.MakeEdge(NewEdge,Newc,Precision::Confusion());
703         if (p1 < p2) {
704           B.Add(NewEdge,vtx1.Oriented(TopAbs_FORWARD));
705           B.Add(NewEdge,vtx2.Oriented(TopAbs_REVERSED));
706         }
707         else {
708           B.Add(NewEdge,vtx1.Oriented(TopAbs_REVERSED));
709           B.Add(NewEdge,vtx2.Oriented(TopAbs_FORWARD));
710         }
711         B.UpdateVertex(vtx1,p1,NewEdge,Precision::Confusion());
712         B.UpdateVertex(vtx2,p2,NewEdge,Precision::Confusion());
713         if (!newCs1.IsNull()) {
714           B.UpdateEdge(NewEdge,newCs1,
715                        TopoDS::Face(MapEV.FindFromKey(edg1)),
716                        Precision::Confusion());
717         }
718
719         if (!newCs2.IsNull()) {
720           B.UpdateEdge(NewEdge,newCs2,
721                        TopoDS::Face(MapEV.FindFromKey(edg2)),
722                        Precision::Confusion());
723         }
724
725       
726         MapEV.Add(vtx,NewEdge);
727
728         if (Choice == 1) {
729           TopoDS_Shape aLocalEdge = Ebind.EmptyCopied();
730           TopoDS_Edge NE = TopoDS::Edge(aLocalEdge);
731 //        TopoDS_Edge NE = TopoDS::Edge(Ebind.EmptyCopied());
732           for (exp2.Init(Ebind,TopAbs_VERTEX); exp2.More(); exp2.Next()) {
733             const TopoDS_Vertex& thevtx = TopoDS::Vertex(exp2.Current());
734             if (thevtx.IsSame(vtx)) {
735               B.Add(NE,vtx1.Oriented(thevtx.Orientation()));
736               B.UpdateVertex(vtx1,p1,NE,Precision::Confusion());
737             }
738             else {
739               B.Add(NE,thevtx);
740               Standard_Real theprm = BRep_Tool::Parameter(thevtx,Ebind);
741               B.UpdateVertex(thevtx,theprm,NE,BRep_Tool::Tolerance(thevtx));
742             }
743           }
744           MapSonS.Bind(Ebind,NE.Oriented(TopAbs_FORWARD));
745           if (IsLeft) {
746             TopTools_ListOfShape thelist5;
747             MapSg.Bind(vtx, thelist5);
748             MapSg(vtx).Append(vtx1);
749           }
750           else {
751             TopTools_ListOfShape thelist6;
752             MapSd.Bind(vtx, thelist6);
753             MapSd(vtx).Append(vtx1);
754           }
755         }
756         else {
757           TopTools_ListOfShape thelist7, thelist8;
758           MapSg.Bind(vtx, thelist7);
759           MapSd.Bind(vtx, thelist8);
760           MapSg(vtx).Append(vtx1);
761           MapSd(vtx).Append(vtx2);
762         }
763       }
764       else {
765         MapEV.Add(vtx,vtx2); // on peut avoir vtx2 = vtx si choix == 1
766         if (Choice == 1) {
767           if (IsLeft) {
768             TopTools_ListOfShape thelist9;
769             MapSg.Bind(vtx, thelist9);
770             MapSg(vtx).Append(vtx);
771           }
772           else {
773             TopTools_ListOfShape thelist10;
774             MapSd.Bind(vtx, thelist10);
775             MapSd(vtx).Append(vtx);
776           }
777         }
778         else {
779           TopTools_ListOfShape thelist11, thelist12;
780           MapSg.Bind(vtx, thelist11);
781           MapSd.Bind(vtx, thelist12);
782           MapSg(vtx).Append(vtx2);
783           MapSd(vtx).Append(vtx2);
784         }
785       }
786     }
787   }
788
789
790   theMap.Clear();
791   for (exp.Init(theW,TopAbs_EDGE); exp.More(); exp.Next()) {
792     const TopoDS_Edge& edg = TopoDS::Edge(exp.Current());
793     if (!theMap.Add(edg)) { // precaution sans doute inutile...
794       continue;
795     }
796     Standard_Integer indedg = MapEV.FindIndex(edg);
797     TopoDS_Face& GenF = TopoDS::Face(MapEV(indedg));
798     TopTools_ListOfShape thelist13, thelist14;
799     MapSg.Bind(edg, thelist13);  // genere a gauche
800     MapSd.Bind(edg, thelist14);  // genere a droite
801     TopoDS_Vertex Vf,Vl;
802     TopoDS_Shape aLocalEdge = edg.Oriented(TopAbs_FORWARD);
803     TopExp::Vertices(TopoDS::Edge(aLocalEdge),Vf,Vl);
804 //    TopExp::Vertices(TopoDS::Edge(edg.Oriented(TopAbs_FORWARD)),Vf,Vl);
805     TopoDS_Shape Gvf = MapEV.FindFromKey(Vf); 
806     TopoDS_Shape Gvl = MapEV.FindFromKey(Vl); 
807
808 /* Le code suivant est OK. On essaie de l`ameliorer
809
810     if (Gvf.ShapeType() == TopAbs_VERTEX &&
811         Gvl.ShapeType() == TopAbs_VERTEX) {
812       // en fait on doit pouvoir avoir 1 face a 2 cotes...
813       if (Gvf.IsSame(Vf)) {
814         MapW(edg) = edg;
815         MapSg(edg).Append(edg.Oriented(TopAbs_FORWARD));
816         MapSd(edg).Append(edg.Oriented(TopAbs_FORWARD));
817       }
818       else {
819         TopoDS_Edge NewEdg = NewEdge(edg,
820                                      GenF,NewSg,
821                                      TopoDS::Vertex(Gvf),
822                                      TopoDS::Vertex(Gvl));
823         if (NewEdg.IsNull()) {
824           return;
825         }
826         MapW(edg) = NewEdg;
827         MapSg(edg).Append(NewEdg);
828         MapSd(edg).Append(NewEdg);
829       }
830     }
831     else if (Gvf.ShapeType() == TopAbs_VERTEX  ||
832              Gvl.ShapeType() == TopAbs_VERTEX) {      // face triangulaire
833       TopoDS_Vertex Vfd,Vld,Vfg,Vlg;
834       if (Gvf.ShapeType() == TopAbs_VERTEX) {
835         Vfg = TopoDS::Vertex(Gvf);
836         Vfd = Vfg;
837         Vlg = TopoDS::Vertex(MapSg(Vl).First());
838         Vld = TopoDS::Vertex(MapSd(Vl).First());
839       }
840       else {
841         Vlg = TopoDS::Vertex(Gvl);
842         Vld = Vlg;
843         Vfg = TopoDS::Vertex(MapSg(Vf).First());
844         Vfd = TopoDS::Vertex(MapSd(Vf).First());
845       }
846
847       TopoDS_Edge NewEdgg = NewEdge(edg,GenF,NewSg,Vfg,Vlg);
848       if (NewEdgg.IsNull()) {
849         return;
850       }
851
852       TopoDS_Edge NewEdgd = NewEdge(edg,GenF,NewSd,Vfd,Vld);
853       if (NewEdgg.IsNull()) {
854         return;
855       }
856       MapSg(edg).Append(NewEdgg);
857       MapSd(edg).Append(NewEdgd);
858
859       TopTools_ListOfShape theedges;
860       theedges.Append(NewEdgg);
861       theedges.Append(NewEdgd);
862       if (Gvf.ShapeType() == TopAbs_EDGE) {
863         theedges.Append(Gvf);
864       }
865       else {//if (Gvl.ShapeType() == TopAbs_EDGE) {
866         theedges.Append(Gvl);
867       }
868       MakeFace(GenF,theedges);
869       MapW(edg) = GenF;
870     }
871     else {
872       // une face a 4 cotes
873       TopoDS_Vertex Vfd,Vld,Vfg,Vlg;
874
875       Vfg = TopoDS::Vertex(MapSg(Vf).First());
876       Vfd = TopoDS::Vertex(MapSd(Vf).First());
877       Vlg = TopoDS::Vertex(MapSg(Vl).First());
878       Vld = TopoDS::Vertex(MapSd(Vl).First());
879       
880       TopoDS_Vertex VVf1,VVl1,VVf2,VVl2;
881       TopExp::Vertices(TopoDS::Edge(Gvf.Oriented(TopAbs_FORWARD)),VVf1,VVl1);
882       TopExp::Vertices(TopoDS::Edge(Gvl.Oriented(TopAbs_FORWARD)),VVf2,VVl2);
883
884       TopoDS_Edge NewEdgg = NewEdge(edg,GenF,NewSg,Vfg,Vlg);
885       if (NewEdgg.IsNull()) {
886         return;
887       }
888       
889       TopoDS_Edge NewEdgd = NewEdge(edg,GenF,NewSd,Vfd,Vld);
890       if (NewEdgd.IsNull()) {
891         return;
892       }
893
894       if ((VVf1.IsSame(Vfg) && VVf2.IsSame(Vlg)) ||
895           (VVf1.IsSame(Vfd) && VVf2.IsSame(Vld))) {
896         // 4 cotes
897         MapSg(edg).Append(NewEdgg);
898         MapSd(edg).Append(NewEdgd);
899         
900         TopTools_ListOfShape theedges;
901         theedges.Append(NewEdgg);
902         theedges.Append(NewEdgd);
903         theedges.Append(Gvf);
904         theedges.Append(Gvl);
905         
906         MakeFace(GenF,theedges);
907         MapW(edg) = GenF;
908       }
909       else {
910 #ifdef DEB
911         cout << "Pb d'analyse" << endl;
912 #endif
913         return;
914       }
915     }
916 */
917     // nouveau code
918
919     TopoDS_Vertex Vfd,Vld,Vfg,Vlg;
920     if (Gvf.ShapeType() == TopAbs_VERTEX) {
921       Vfg = TopoDS::Vertex(Gvf);
922       Vfd = Vfg;
923     }
924     else {
925       Vfg = TopoDS::Vertex(MapSg(Vf).First());
926       Vfd = TopoDS::Vertex(MapSd(Vf).First());
927     }
928     if (Gvl.ShapeType() == TopAbs_VERTEX) {
929       Vlg = TopoDS::Vertex(Gvl);
930       Vld = Vlg;
931     }
932     else {
933       Vlg = TopoDS::Vertex(MapSg(Vl).First());
934       Vld = TopoDS::Vertex(MapSd(Vl).First());
935     }
936
937     TopoDS_Edge NewEdgg = NewEdge(edg,GenF,NewSg,Vfg,Vlg);
938     if (NewEdgg.IsNull()) {
939       return;
940     }
941     
942     TopoDS_Edge NewEdgd = NewEdge(edg,GenF,NewSd,Vfd,Vld);
943     if (NewEdgg.IsNull()) {
944       return;
945     }
946
947     Standard_Boolean isedg = Standard_False;
948     if (Gvf.ShapeType() == TopAbs_VERTEX &&
949         Gvl.ShapeType() == TopAbs_VERTEX) {
950       // edg ou face a 2 cotes
951
952       // Comparaison NewEdgg et NewEdgd
953       Standard_Real fg,lg,fd,ld;
954       Handle(Geom_Curve) Cg = BRep_Tool::Curve(NewEdgg,fg,lg);
955       Handle(Geom_Curve) Cd = BRep_Tool::Curve(NewEdgd,fd,ld);
956       Standard_Real prmg = (fg+lg)/2.;
957       Standard_Real prmd = (fd+ld)/2.;
958       gp_Pnt pg = Cg->Value(prmg);
959       gp_Pnt pd = Cd->Value(prmd);
960       Standard_Real Tol = Max(BRep_Tool::Tolerance(NewEdgg),
961                               BRep_Tool::Tolerance(NewEdgg));
962       if (pg.SquareDistance(pd) <= Tol*Tol) {
963         isedg = Standard_True;
964         // raffinement pour essayer de partager l`edge de depart...
965         Standard_Boolean modified = Standard_True;
966         if (Gvf.IsSame(Vf) && Gvl.IsSame(Vl)) {
967           // Comparaison avec l`edge de depart
968           Cd = BRep_Tool::Curve(edg,fd,ld);
969           prmd = (fd+ld)/2.;
970           pd = Cd->Value(prmd);
971           Tol = Max(BRep_Tool::Tolerance(NewEdgg),
972                     BRep_Tool::Tolerance(edg));
973           if (pg.SquareDistance(pd) <= Tol*Tol) {
974             modified = Standard_False;
975           }
976         }
977
978         if (!modified) {
979           MapW(edg) = edg;
980           MapSg(edg).Append(edg);
981           MapSd(edg).Append(edg);
982         }
983         else {
984           MapW(edg) = NewEdgg;
985           MapSg(edg).Append(NewEdgg);
986           MapSd(edg).Append(NewEdgg);
987         }
988       }
989     }
990
991     if (!isedg) {
992       // face a 2 ou 3 ou 4 cotes
993       MapSg(edg).Append(NewEdgg);
994       MapSd(edg).Append(NewEdgd);
995       
996       TopTools_ListOfShape theedges;
997       theedges.Append(NewEdgg);
998       theedges.Append(NewEdgd);
999       if (Gvf.ShapeType() == TopAbs_EDGE) {
1000         theedges.Append(Gvf);
1001       }
1002       if (Gvl.ShapeType() == TopAbs_EDGE) {
1003         theedges.Append(Gvl);
1004       }
1005       MakeFace(GenF,theedges);
1006       MapW(edg) = GenF;
1007     }
1008
1009   }
1010
1011
1012   TopTools_MapOfShape mapedgadded;
1013   TopTools_ListOfShape thefaces;
1014
1015   for (itl.Initialize(myMap(F)); itl.More(); itl.Next()) {
1016     const TopoDS_Face& fac = TopoDS::Face(itl.Value());
1017     theMap.Clear();
1018     TopoDS_Face DrftFace; // elle est FORWARD
1019     Standard_Boolean IsLeft;
1020     if (Contains(theLeft,fac)) {
1021       B.MakeFace(DrftFace,NewSg,BRep_Tool::Tolerance(fac));
1022       IsLeft = Standard_True;
1023     }
1024     else {
1025       B.MakeFace(DrftFace,NewSd,BRep_Tool::Tolerance(fac));
1026       IsLeft = Standard_False;
1027     }
1028     
1029     TopExp_Explorer exp3;
1030     for (exp3.Init(fac.Oriented(TopAbs_FORWARD),TopAbs_WIRE);
1031          exp3.More(); exp3.Next()) {
1032       const TopoDS_Shape& wir = exp3.Current();
1033       TopoDS_Wire NewWireOnF;
1034       B.MakeWire(NewWireOnF);
1035       for (exp.Init(wir.Oriented(TopAbs_FORWARD),TopAbs_EDGE);
1036            exp.More(); exp.Next()) {
1037         const TopoDS_Edge& edg = TopoDS::Edge(exp.Current());
1038         if (!theMap.Add(edg)) { // precaution sans doute inutile...
1039           continue;
1040         }
1041         if (MapW.IsBound(edg)) { // edge du wire d`origine
1042           TopTools_ListIteratorOfListOfShape itld;
1043           TopAbs_Orientation ored = edg.Orientation();
1044           if (IsLeft) {
1045             itld.Initialize(MapSg(edg));
1046           }
1047           else {
1048             itld.Initialize(MapSd(edg));
1049           }
1050           for (; itld.More(); itld.Next()) {
1051             if (itld.Value().Orientation() == TopAbs_REVERSED) {
1052               ored = TopAbs::Reverse(ored);
1053             }
1054             TopoDS_Shape aLocalEdge = itld.Value().Oriented(ored);
1055             B.Add(NewWireOnF,TopoDS::Edge(aLocalEdge));
1056 //          B.Add(NewWireOnF,TopoDS::Edge(itld.Value().Oriented(ored)));
1057           }
1058         }
1059         else {
1060           Handle(Geom_Surface) NewS;
1061           if (IsLeft) {
1062             NewS = NewSg;
1063           }
1064           else {
1065             NewS = NewSd;
1066           }
1067           Standard_Integer indedg = MapEV.FindIndex(edg);
1068           const TopoDS_Face& GenF = TopoDS::Face(MapEV(indedg));
1069           TopoDS_Vertex Vf,Vl;
1070           TopoDS_Shape aLocalEdge = edg.Oriented(TopAbs_FORWARD);
1071           TopExp::Vertices(TopoDS::Edge(aLocalEdge),Vf,Vl);
1072 //        TopExp::Vertices(TopoDS::Edge(edg.Oriented(TopAbs_FORWARD)),Vf,Vl);
1073           TopoDS_Shape Gvf = MapEV.FindFromKey(Vf); 
1074           TopoDS_Shape Gvl = MapEV.FindFromKey(Vl); 
1075           if (Gvf.ShapeType() == TopAbs_VERTEX &&
1076               Gvl.ShapeType() == TopAbs_VERTEX) {
1077             if (!Gvf.IsSame(Vf) || !Gvl.IsSame(Vl)) {
1078               TopoDS_Edge NewEdg = NewEdge(edg,GenF,NewS,
1079                                            TopoDS::Vertex(Gvf),
1080                                            TopoDS::Vertex(Gvl));
1081               if (NewEdg.IsNull()) {
1082                 return;
1083               }
1084
1085               MapSonS.Bind(edg,NewEdg);
1086
1087               if (NewEdg.Orientation() == TopAbs_REVERSED) {
1088                 NewEdg.Orientation(TopAbs::Reverse(edg.Orientation()));
1089               }
1090               else {
1091                 NewEdg.Orientation(edg.Orientation());
1092               }         
1093               B.Add(NewWireOnF,NewEdg);
1094             }
1095             else { // Frozen???
1096               B.Add(NewWireOnF,edg);
1097             }
1098           }
1099           else {
1100             TopoDS_Vertex Vff,Vll;
1101             if (Gvf.ShapeType() == TopAbs_VERTEX) {
1102               Vff = TopoDS::Vertex(Gvf);
1103             }
1104             else {
1105               if (IsLeft) {
1106                 Vff = TopoDS::Vertex(MapSg(Vf).First());
1107               }
1108               else {
1109                 Vff = TopoDS::Vertex(MapSd(Vf).First());
1110               }
1111             }
1112             if (Gvl.ShapeType() == TopAbs_VERTEX) {
1113               Vll = TopoDS::Vertex(Gvl);
1114             }
1115             else {
1116               if (IsLeft) {
1117                 Vll = TopoDS::Vertex(MapSg(Vl).First());
1118               }
1119               else {
1120                 Vll = TopoDS::Vertex(MapSd(Vl).First());
1121               }
1122             }
1123             
1124             TopoDS_Edge NewEdg = NewEdge(edg,GenF,NewS,Vff,Vll);
1125             if (NewEdg.IsNull()) {
1126               return;
1127             }
1128
1129             if (!MapW.IsBound(Vf) && !MapW.IsBound(Vl)) {
1130               MapSonS.Bind(edg,NewEdg);
1131             }
1132 //          else if (MapW.IsBound(Vf) && MapW.IsBound(Vl)) {
1133               
1134
1135 //          }
1136             else {
1137               if (MapW.IsBound(Vf)) {
1138                 if (Gvf.ShapeType() != TopAbs_EDGE || 
1139                     mapedgadded.Contains(Gvf)) {
1140                   MapSonS.Bind(edg,NewEdg);
1141                 }
1142                 else {
1143                   TopoDS_Wire NewWir;
1144                   B.MakeWire(NewWir);
1145                   B.Add(NewWir,NewEdg);
1146
1147                   TopoDS_Vertex Vf2,Vl2;
1148                   TopExp::Vertices(TopoDS::Edge(Gvf),Vf2,Vl2);
1149
1150                   //TopAbs_Orientation ornw = NewEdg.Orientation();
1151
1152                   // ici bug orientation : voir tspdrft6
1153
1154 //                if ((ornw == TopAbs_FORWARD && Vl2.IsSame(Vff)) ||
1155 //                    (ornw == TopAbs_REVERSED && Vl2.IsSame(Vll))) {
1156                   if (Vl2.IsSame(Vff)) {
1157                     B.Add(NewWir,Gvf.Oriented(TopAbs_FORWARD));
1158                   }
1159                   else {
1160                     B.Add(NewWir,Gvf.Oriented(TopAbs_REVERSED));
1161                   }
1162                   mapedgadded.Add(Gvf);
1163                   MapSonS.Bind(edg,NewWir); // NewWire est FORWARD
1164                 }
1165               }
1166               else {
1167                 if (Gvl.ShapeType() != TopAbs_EDGE ||
1168                     mapedgadded.Contains(Gvl)) {
1169                   MapSonS.Bind(edg,NewEdg);
1170                 }
1171                 else {
1172                   TopoDS_Wire NewWir;
1173                   B.MakeWire(NewWir);
1174                   B.Add(NewWir,NewEdg);
1175
1176                   TopoDS_Vertex Vf2,Vl2;
1177                   TopExp::Vertices(TopoDS::Edge(Gvl),Vf2,Vl2);
1178
1179                   //TopAbs_Orientation ornw = NewEdg.Orientation();
1180
1181                   // ici bug orientation : voir tspdrft6
1182
1183 //                if ((ornw == TopAbs_FORWARD && Vl2.IsSame(Vff)) ||
1184 //                    (ornw == TopAbs_REVERSED && Vl2.IsSame(Vll))) {
1185                   if (Vf2.IsSame(Vll)) {
1186                     B.Add(NewWir,Gvl.Oriented(TopAbs_FORWARD));
1187                   }
1188                   else {
1189                     B.Add(NewWir,Gvl.Oriented(TopAbs_REVERSED));
1190                   }
1191                   mapedgadded.Add(Gvl);
1192                   MapSonS.Bind(edg,NewWir); // NewWire est FORWARD
1193                 }
1194               }
1195             }
1196             if (NewEdg.Orientation() == TopAbs_REVERSED) {
1197               NewEdg.Orientation(TopAbs::Reverse(edg.Orientation()));
1198             }
1199             else {
1200               NewEdg.Orientation(edg.Orientation());
1201             }           
1202             B.Add(NewWireOnF,NewEdg);
1203           }
1204         }
1205       }
1206       B.Add(DrftFace,NewWireOnF.Oriented(wir.Orientation()));
1207     }
1208     thefaces.Append(DrftFace);
1209   }
1210
1211   BRepTools_Substitution theSubs;
1212   TopTools_DataMapIteratorOfDataMapOfShapeShape itdmss;
1213   for (itdmss.Initialize(MapSonS);
1214        itdmss.More(); itdmss.Next()) {
1215     TopTools_ListOfShape lsubs;
1216     for (exp.Init(itdmss.Value(),TopAbs_EDGE); exp.More(); exp.Next()) {
1217       lsubs.Append(exp.Current());
1218     }
1219     theSubs.Substitute(itdmss.Key(),lsubs);
1220   }
1221
1222   // on reconstruit les faces
1223   for (exp.Init(Res,TopAbs_FACE); exp.More(); exp.Next()) {
1224     if (Contains(myMap(F),exp.Current())) {
1225       continue;
1226     }
1227     theSubs.Build(exp.Current());
1228   }
1229
1230   // Stockage des descendants
1231 //  for (TopTools_DataMapIteratorOfDataMapOfShapeListOfShape itdmsls(myMap);
1232   TopTools_DataMapIteratorOfDataMapOfShapeListOfShape itdmsls(myMap) ;
1233   for ( ; itdmsls.More(); itdmsls.Next()) {
1234     if (itdmsls.Key().ShapeType() == TopAbs_EDGE) {
1235       TopTools_ListOfShape thedesc;
1236       theMap.Clear();
1237       for (itl.Initialize(itdmsls.Value());itl.More(); itl.Next()) {
1238         if (theMap.Add(MapW(itl.Value()))) {
1239           thedesc.Append(MapW(itl.Value()));
1240         }
1241       }
1242       myMap(itdmsls.Key()) = thedesc;
1243     }
1244     else if (itdmsls.Key().IsSame(F)) {
1245       myMap(F).Clear();
1246       for (itl.Initialize(thefaces); itl.More(); itl.Next()) {
1247         myMap(F).Append(itl.Value());
1248       }
1249     }
1250     else {
1251       TopTools_ListOfShape thedesc;
1252       theMap.Clear();
1253       for (itl.Initialize(itdmsls.Value());itl.More(); itl.Next()) {
1254         if (theSubs.IsCopied(itl.Value())) {
1255           if (theSubs.Copy(itl.Value()).Extent() != 1) {
1256 #ifdef DEB
1257             cout << "Invalid number of descendant" << endl;
1258 #endif
1259             return;
1260           }
1261           else {
1262             if (theMap.Add(theSubs.Copy(itl.Value()).First())) {
1263               thedesc.Append(theSubs.Copy(itl.Value()).First());
1264             }
1265           }
1266         }
1267         else if (theMap.Add(itl.Value())) {
1268           thedesc.Append(itl.Value());
1269         }
1270       }
1271       myMap(itdmsls.Key()) = thedesc;
1272     }
1273   }
1274
1275   theMap.Clear();
1276   thefaces.Clear();
1277   for (itdmsls.Initialize(myMap);itdmsls.More(); itdmsls.Next()) {
1278     for (itl.Initialize(itdmsls.Value());itl.More(); itl.Next()) {
1279       if (itl.Value().ShapeType() == TopAbs_FACE &&
1280           theMap.Add(itl.Value())) {
1281         thefaces.Append(itl.Value());
1282       }
1283     }
1284   }
1285   LocOpe_BuildShape BS(thefaces);
1286   myResult = BS.Shape();
1287 }
1288
1289
1290 //=======================================================================
1291 //function : Shape
1292 //purpose  : 
1293 //=======================================================================
1294
1295 const TopoDS_Shape& LocOpe_SplitDrafts::Shape () const
1296 {
1297   if (myResult.IsNull()) {
1298     StdFail_NotDone::Raise();
1299   }
1300   return myResult;
1301 }
1302
1303 //=======================================================================
1304 //function : ShapesFromShape
1305 //purpose  : 
1306 //=======================================================================
1307
1308 const TopTools_ListOfShape& LocOpe_SplitDrafts::ShapesFromShape 
1309    (const TopoDS_Shape& S) const
1310 {
1311   if (myResult.IsNull()) {
1312     StdFail_NotDone::Raise();
1313   }
1314   return myMap(S);
1315 }
1316
1317
1318 //=======================================================================
1319 //function : NewPlane
1320 //purpose  : 
1321 //=======================================================================
1322
1323 static Standard_Boolean NewPlane(const TopoDS_Face& F,
1324                                  const gp_Dir& Extr,
1325                                  const gp_Pln& Neutr,
1326                                  const Standard_Real Ang,
1327                                  gp_Pln& Newpl,
1328                                  gp_Ax1& NormalF,
1329                                  const Standard_Boolean Modify)
1330 {
1331
1332
1333   // Determination du nouveau plan incline
1334   Handle(Geom_Surface) S = BRep_Tool::Surface(F);
1335   if (S->DynamicType() == STANDARD_TYPE(Geom_RectangularTrimmedSurface)) {
1336     S = Handle(Geom_RectangularTrimmedSurface)::DownCast(S)->BasisSurface();
1337   }
1338   Handle(Geom_Plane) P = Handle(Geom_Plane)::DownCast(S);
1339   if (P.IsNull()) {
1340     return Standard_False;
1341   }
1342
1343   gp_Pln Plorig = P->Pln();
1344   if (!Modify) {
1345     Newpl = Plorig;
1346     NormalF = Newpl.Axis();
1347     if ((Newpl.Direct() && F.Orientation() == TopAbs_REVERSED) ||
1348        (!Newpl.Direct() && F.Orientation() == TopAbs_FORWARD)) {
1349       NormalF.Reverse();
1350     }
1351     return Standard_True;
1352   }
1353
1354   gp_Ax1 Axe;
1355   Standard_Real Theta;
1356
1357   IntAna_QuadQuadGeo i2pl(Plorig,Neutr,
1358                           Precision::Angular(),Precision::Confusion());
1359   
1360   if (i2pl.IsDone() && i2pl.TypeInter() == IntAna_Line) {
1361     gp_Lin LinInters = i2pl.Line(1);
1362     gp_Dir nx = LinInters.Direction();
1363     NormalF = Plorig.Axis();
1364     gp_Dir ny = NormalF.Direction().Crossed(nx);
1365     Standard_Real a = Extr.Dot(nx);
1366     if (Abs(a) <=1-Precision::Angular()) { 
1367       Standard_Real b = Extr.Dot(ny);
1368       Standard_Real c = Extr.Dot(NormalF.Direction());
1369       Standard_Boolean direct(Plorig.Direct());
1370       TopAbs_Orientation Oris = F.Orientation();
1371       if ((direct && Oris == TopAbs_REVERSED) ||
1372           (!direct && Oris == TopAbs_FORWARD)) {
1373         b = -b;
1374         c = -c;
1375         NormalF.Reverse();
1376       }
1377       Standard_Real denom = Sqrt(1-a*a);
1378       Standard_Real Sina = Sin(Ang);
1379       if (denom>Abs(Sina)) {
1380         Standard_Real phi = ATan2(b/denom,c/denom);
1381         Standard_Real theta0 = ACos(Sina/denom); 
1382         Theta = theta0 - phi;
1383         if (Cos(Theta) <0.) {
1384           Theta = -theta0 -phi;
1385         }
1386         Axe = LinInters.Position();
1387         Newpl = Plorig.Rotated(Axe,Theta);
1388         return Standard_True;
1389       }
1390     }
1391   }
1392   cout << "fin newplane return standard_false" << endl;
1393   return Standard_False;
1394 }
1395
1396
1397 //=======================================================================
1398 //function : MakeFace
1399 //purpose  : 
1400 //=======================================================================
1401
1402 static void MakeFace(TopoDS_Face& F,
1403                      TopTools_ListOfShape& ledg)
1404 {
1405
1406   // ledg est une liste d'edge
1407
1408   BRep_Builder B;
1409
1410   // Verification de l`existence des p-curves. Celles qui manquent
1411   // correspondent necessairement a des isos (et meme des iso u).
1412
1413   Standard_Real f,l;
1414 //  for (TopTools_ListIteratorOfListOfShape itl(ledg); 
1415   TopTools_ListIteratorOfListOfShape itl(ledg) ;
1416   for ( ; itl.More(); itl.Next()) {
1417     TopoDS_Edge& edg = TopoDS::Edge(itl.Value());
1418     Handle(Geom2d_Curve) C2d = BRep_Tool::CurveOnSurface(edg,F,f,l);
1419     if (C2d.IsNull()) {
1420       BRep_Tool::Range(edg,f,l);
1421       TopoDS_Vertex V1,V2;
1422       TopExp::Vertices(edg,V1,V2);
1423         TopTools_ListIteratorOfListOfShape itl2;
1424       for (itl2.Initialize(ledg); 
1425            itl2.More(); itl2.Next()) {
1426         const TopoDS_Edge& edg2 = TopoDS::Edge(itl2.Value());
1427         if (edg2.IsSame(edg)) {
1428           continue;
1429         }
1430         TopoDS_Vertex Vp1,Vp2;
1431         TopExp::Vertices(edg2,Vp1,Vp2);
1432         if (Vp1.IsSame(V1) || Vp2.IsSame(V1) ||
1433             Vp1.IsSame(V2) || Vp2.IsSame(V2)) {
1434           Standard_Real f2,l2;
1435           Handle(Geom2d_Curve) C22d = BRep_Tool::CurveOnSurface(edg2,F,f2,l2);
1436           if (!C22d.IsNull()) {
1437             gp_Pnt2d pt2d;
1438             if (Vp1.IsSame(V1)) {
1439               pt2d = C22d->Value(f2);
1440               pt2d.SetY(pt2d.Y()-f);
1441             }
1442             else if (Vp2.IsSame(V1)) {
1443               pt2d = C22d->Value(l2);
1444               pt2d.SetY(pt2d.Y()-f);
1445             }
1446             else if (Vp1.IsSame(V2)) {
1447               pt2d = C22d->Value(f2);
1448               pt2d.SetY(pt2d.Y()-l);
1449             }
1450             else if (Vp2.IsSame(V2)) {
1451               pt2d = C22d->Value(l2);
1452               pt2d.SetY(pt2d.Y()-l);
1453             }
1454             C2d = new Geom2d_Line(pt2d,gp::DY2d());
1455             B.UpdateEdge(edg,C2d,F,BRep_Tool::Tolerance(edg));
1456             break;
1457           }
1458         }
1459       }
1460       if (C2d.IsNull()) {
1461         cout << "Ca merde violemment" << endl;
1462       }
1463     }
1464   }
1465
1466   TopTools_ListOfShape lwires;
1467   Standard_Boolean alldone = ledg.IsEmpty();
1468   while (!alldone) {
1469     TopoDS_Wire Wnew;
1470     B.MakeWire(Wnew);
1471     TopoDS_Shape aLocalShape = ledg.First();
1472     const TopoDS_Edge& edg = TopoDS::Edge(aLocalShape);
1473 //    const TopoDS_Edge& edg = TopoDS::Edge(ledg.First());
1474     TopoDS_Vertex VFirst,VLast;
1475     if (edg.Orientation() == TopAbs_FORWARD) {
1476       TopExp::Vertices(edg,VFirst,VLast);
1477     }
1478     else {
1479       TopExp::Vertices(edg,VLast,VFirst);
1480     }
1481     B.Add(Wnew,edg);
1482     ledg.RemoveFirst();
1483     // on suppose VFirst et VLast non nuls
1484     Standard_Boolean wdone = (ledg.IsEmpty() || VFirst.IsSame(VLast));
1485     while (!wdone) {
1486       TopoDS_Vertex VF,VL;
1487
1488       TopAbs_Orientation oredg = TopAbs_FORWARD;
1489
1490       for (itl.Initialize(ledg); itl.More(); itl.Next()) {
1491         const TopoDS_Edge& edg2 = TopoDS::Edge(itl.Value());
1492         TopoDS_Shape aLocalShape  = edg2.Oriented(TopAbs_FORWARD);
1493         TopExp::Vertices(TopoDS::Edge(aLocalShape),VF,VL);
1494 //      TopExp::Vertices(TopoDS::Edge(edg2.Oriented(TopAbs_FORWARD)),VF,VL);
1495         if (VF.IsSame(VLast)) {
1496           VLast = VL;
1497           oredg = TopAbs_FORWARD;
1498           break;
1499         }
1500         else if (VL.IsSame(VFirst)) {
1501           VFirst = VF;
1502           oredg = TopAbs_FORWARD;
1503           break;
1504         }
1505         else if (VF.IsSame(VFirst)) {
1506           VFirst = VL;
1507           oredg = TopAbs_REVERSED;
1508           break;
1509         }
1510         else if (VL.IsSame(VLast)) {
1511           VLast = VF;
1512           oredg = TopAbs_REVERSED;
1513           break;
1514         }
1515
1516       }
1517       if (!itl.More()) {
1518         wdone = Standard_True;
1519       }
1520       else {
1521         TopoDS_Shape aLocalShape = itl.Value().Oriented(oredg);
1522         B.Add(Wnew,TopoDS::Edge(aLocalShape));
1523 //      B.Add(Wnew,TopoDS::Edge(itl.Value().Oriented(oredg)));
1524         ledg.Remove(itl);
1525         wdone = (ledg.IsEmpty() || VFirst.IsSame(VLast));
1526       }
1527     }
1528     lwires.Append(Wnew);
1529     alldone = ledg.IsEmpty();
1530   }
1531
1532
1533
1534
1535   F.Orientation(TopAbs_FORWARD);
1536   for (itl.Initialize(lwires); itl.More(); itl.Next()) {
1537     TopoDS_Shape aLocalShape = F.EmptyCopied();
1538     TopoDS_Face NewFace = TopoDS::Face(aLocalShape);
1539 //    TopoDS_Face NewFace = TopoDS::Face(F.EmptyCopied());
1540     B.Add(NewFace,itl.Value());
1541     GProp_GProps GP;
1542     BRepGProp::SurfaceProperties(NewFace,GP);
1543     if (GP.Mass() < 0) {
1544       itl.Value().Reverse();
1545     }
1546   }
1547   if (lwires.Extent() == 1) {
1548     B.Add(F,lwires.First());
1549   }
1550   else {
1551     cout << "Not yet implemented : nbwire >= 2" << endl;
1552   }
1553
1554 }
1555
1556
1557 //=======================================================================
1558 //function : Contains
1559 //purpose  : 
1560 //=======================================================================
1561
1562 static Standard_Boolean Contains(const TopTools_ListOfShape& ll,
1563                                  const TopoDS_Shape& s)
1564 {
1565   TopTools_ListIteratorOfListOfShape itl;
1566   for (itl.Initialize(ll); itl.More(); itl.Next()) {
1567     if (itl.Value().IsSame(s)) {
1568       return Standard_True;
1569     }
1570   }
1571   return Standard_False;
1572 }
1573
1574
1575
1576 //=======================================================================
1577 //function : Contains
1578 //purpose  : 
1579 //=======================================================================
1580
1581 static TopoDS_Edge  NewEdge(const TopoDS_Edge& edg,
1582                             const TopoDS_Face& F,
1583                             const Handle(Geom_Surface)& NewS,
1584                             const TopoDS_Vertex& V1,
1585                             const TopoDS_Vertex& V2)
1586 {
1587   TopoDS_Edge NewEdg;
1588   Handle(Geom_Surface) S1 = BRep_Tool::Surface(F);
1589   Standard_Boolean AppS1 = Standard_False;
1590   if (S1->DynamicType() != STANDARD_TYPE(Geom_Plane)) {
1591     AppS1 = Standard_True;
1592   }
1593
1594
1595   GeomInt_IntSS i2s(S1,NewS,Precision::Confusion(),Standard_True,AppS1);
1596   if (!i2s.IsDone() || i2s.NbLines() <= 0) {
1597     return NewEdg;
1598   }
1599
1600   BRep_Builder B;
1601 //  Standard_Real pmin, Dist, DistMin;
1602   Standard_Real Dist2, Dist2Min;
1603   Standard_Real prmf=0,prml=0;
1604   GeomAdaptor_Curve TheCurve;
1605         
1606   Standard_Integer i,k;
1607   gp_Pnt pvf = BRep_Tool::Pnt(V1);
1608   gp_Pnt pvl = BRep_Tool::Pnt(V2);
1609   for (i=1; i<= i2s.NbLines(); i++) {
1610     TheCurve.Load(i2s.Line(i));
1611     Extrema_ExtPC myExtPC(pvf,TheCurve);
1612     
1613     if (myExtPC.IsDone()) {
1614       gp_Pnt p1b,p2b;
1615       Standard_Real thepmin = TheCurve.FirstParameter();
1616       myExtPC.TrimmedSquareDistances(Dist2Min,Dist2,p1b,p2b);
1617       if (Dist2 < Dist2Min && !TheCurve.IsPeriodic()) {
1618         Dist2Min = Dist2;
1619         thepmin = TheCurve.LastParameter();
1620       }
1621       for (k=1; k<=myExtPC.NbExt(); k++) {
1622         Dist2 = myExtPC.SquareDistance(k);
1623         if (Dist2 < Dist2Min) {
1624           Dist2Min = Dist2;
1625           thepmin = myExtPC.Point(k).Parameter();
1626         }
1627       }
1628       
1629       if (Dist2Min  <= Precision::SquareConfusion()) {
1630         prmf = thepmin;
1631         myExtPC.Perform(pvl);
1632         if (myExtPC.IsDone()) {
1633           thepmin = TheCurve.LastParameter();
1634           myExtPC.TrimmedSquareDistances(Dist2,Dist2Min,p1b,p2b);
1635           if (Dist2 < Dist2Min && !TheCurve.IsClosed()) {
1636             Dist2Min = Dist2;
1637             thepmin = TheCurve.FirstParameter();
1638           }
1639           for (k=1; k<=myExtPC.NbExt(); k++) {
1640             Dist2 = myExtPC.SquareDistance(k);
1641             if (Dist2 < Dist2Min) {
1642               Dist2Min = Dist2;
1643               thepmin = myExtPC.Point(k).Parameter();
1644             }
1645           }
1646           
1647       if (Dist2Min  <= Precision::SquareConfusion()) {
1648             prml = thepmin;
1649             break;
1650           }
1651         }
1652       }
1653     }
1654   }
1655
1656   if (i <= i2s.NbLines()) {
1657     Standard_Boolean rev = Standard_False;
1658     TopoDS_Vertex Vf = V1;
1659     TopoDS_Vertex Vl = V2;
1660     Handle(Geom_Curve) Cimg = i2s.Line(i);
1661     Handle(Geom2d_Curve) Cimg2d;
1662     if (AppS1) {
1663       Cimg2d = i2s.LineOnS1(i);
1664     }
1665
1666     if (Cimg->IsPeriodic()) {
1667
1668       Standard_Real period = Cimg->Period();
1669       Standard_Real imf = Cimg->FirstParameter();
1670       Standard_Real iml = Cimg->LastParameter();
1671
1672       Standard_Real f,l;
1673       BRep_Tool::Range(edg,f,l);
1674       Standard_Real delt = l-f;
1675       Standard_Real delt1 = Abs(prml-prmf);
1676       Standard_Real delt2 = Abs(period-delt1);
1677       
1678       if (delt1 == 0 || delt2 == 0) {
1679 //      prmf = 0;
1680 //      prml = period;
1681         prmf = imf;
1682         prml = iml;
1683       }
1684       else {
1685         if (Abs(delt1-delt) > Abs(delt2-delt)) {
1686           // le bon ecart est delt2...
1687           if (prml > prmf) {
1688             if (prml < iml) {
1689               prmf += period;
1690             }
1691             else {
1692               prml -= period;
1693             }
1694           }
1695           else {
1696             if (prmf < iml) {
1697               prml += period;
1698             }
1699             else {
1700               prmf -= period;
1701             }
1702           }
1703         }
1704         else if (Abs(delt1-delt) < Abs(delt2-delt)) {
1705           if (prmf >= iml && prml >= iml) {
1706             prmf -= period;
1707             prml -= period;
1708           }
1709           else if (prmf <= imf && prml <= imf) {
1710             prmf += period;
1711             prml += period;
1712           }
1713         }
1714         else { // egalite; on priveligie l'ordre f,l
1715           if (prmf > prml) {
1716             prmf -= period;
1717           }
1718           if (prmf >= iml && prml >= iml) {
1719             prmf -= period;
1720             prml -= period;
1721           }
1722           else if (prmf <= imf && prml <= imf) {
1723             prmf += period;
1724             prml += period;
1725           }
1726         }
1727       }
1728 #ifdef DEB
1729       Standard_Real ptol = Precision::PConfusion();
1730       if (prmf < imf - ptol || prmf > iml + ptol ||
1731           prml < imf - ptol || prml > iml + ptol) {
1732         cout << "Ca ne va pas aller" << endl;
1733       }
1734 #endif
1735
1736
1737     }
1738
1739     if (S1->IsUPeriodic()) {
1740
1741       Standard_Real speriod = S1->UPeriod();
1742 //      Standard_Real f,l;
1743       gp_Pnt2d pf,pl;
1744       pf = Cimg2d->Value(prmf);
1745       pl = Cimg2d->Value(prml);
1746
1747       Standard_Real Uf = pf.X();
1748       Standard_Real Ul = pl.X();
1749       Standard_Real ptra = 0.0;
1750
1751       Standard_Real Ustart = Min(Uf,Ul);
1752       while (Ustart < -Precision::PConfusion()) {
1753         Ustart += speriod;
1754         ptra += speriod;
1755       }
1756       while (Ustart > speriod - Precision::PConfusion()) {
1757         Ustart -= speriod;
1758         ptra -= speriod;
1759       }
1760       if (ptra != 0) {
1761         Cimg2d->Translate(gp_Vec2d(ptra,0.));
1762       }
1763
1764
1765     }
1766     if (prmf < prml) {
1767       Vf.Orientation(TopAbs_FORWARD);
1768       Vl.Orientation(TopAbs_REVERSED);
1769     }
1770     else {
1771       Vf.Orientation(TopAbs_REVERSED);
1772       Vl.Orientation(TopAbs_FORWARD);
1773       rev = Standard_True;;
1774     }
1775
1776     B.MakeEdge(NewEdg,Cimg,Precision::Confusion());
1777
1778     B.Add(NewEdg,Vf);
1779     B.Add(NewEdg,Vl);
1780     B.UpdateVertex(Vf,prmf,NewEdg,Precision::Confusion());
1781     B.UpdateVertex(Vl,prml,NewEdg,Precision::Confusion());
1782     if (AppS1) {
1783       B.UpdateEdge(NewEdg,Cimg2d,F,Precision::Confusion());
1784     }
1785
1786
1787     if (rev) {
1788       NewEdg.Orientation(TopAbs_REVERSED);
1789     }
1790   }
1791   return NewEdg;
1792 }