src/ShapeFix/ShapeFix_Wire_1.cxx

   1 // Copyright (c) 1999-2014 OPEN CASCADE SAS
   2 //
   3 // This file is part of Open CASCADE Technology software library.
   4 //
   5 // This library is free software; you can redistribute it and/or modify it under
   6 // the terms of the GNU Lesser General Public License version 2.1 as published
   7 // by the Free Software Foundation, with special exception defined in the file
   8 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
   9 // distribution for complete text of the license and disclaimer of any warranty.
  10 //
  11 // Alternatively, this file may be used under the terms of Open CASCADE
  12 // commercial license or contractual agreement.
  13
  14 // szv 19.08.99: new methods for fixing gaps between edges (3d curves and pcurves)
  15 #include <ShapeFix_Wire.hxx>
  16 #include <Standard_ErrorHandler.hxx>
  17 #include <Standard_Failure.hxx>
  18
  19 #include <Precision.hxx>
  20 #include <Bnd_Box2d.hxx>
  21 #include <Geom_Curve.hxx>
  22 #include <Geom2d_Curve.hxx>
  23 #include <Geom2d_Line.hxx>
  24
  25 #include <IntRes2d_SequenceOfIntersectionPoint.hxx>
  26 #include <IntRes2d_IntersectionPoint.hxx>
  27 #include <TColgp_SequenceOfPnt.hxx>
  28
  29 #include <TopoDS.hxx>
  30 #include <TopoDS_Vertex.hxx>
  31 #include <TopoDS_Edge.hxx>
  32 #include <TopTools_HSequenceOfShape.hxx>
  33
  34 #include <BRep_Tool.hxx>
  35 #include <BRep_Builder.hxx>
  36
  37 #include <ShapeExtend.hxx>
  38 #include <ShapeBuild_Edge.hxx>
  39 #include <ShapeBuild_Vertex.hxx>
  40 #include <ShapeAnalysis_Curve.hxx>
  41 #include <ShapeAnalysis_Edge.hxx>
  42 #include <ShapeAnalysis_Surface.hxx>
  43 #include <ShapeAnalysis.hxx>
  44 #include <GeomConvert_CompCurveToBSplineCurve.hxx>
  45 #include <Geom_TrimmedCurve.hxx>
  46 #include <gp_Pln.hxx>
  47 #include <GeomAPI.hxx>
  48 #include <TColgp_Array1OfPnt.hxx>
  49 #include <TColStd_Array1OfReal.hxx>
  50 #include <TColStd_Array1OfInteger.hxx>
  51 #include <Geom_BSplineCurve.hxx>
  52 #include <Geom_SphericalSurface.hxx> //S4135
  53 #include <Geom2d_BSplineCurve.hxx>
  54 #include <GeomAdaptor_Curve.hxx>
  55 #include <GeomAdaptor_HSurface.hxx>
  56 #include <GeomAdaptor_Surface.hxx>
  57 #include <TopTools_Array1OfShape.hxx>
  58 #include <BRepTools.hxx>
  59 #include <Bnd_Array1OfBox2d.hxx>
  60 #include <BndLib_Add2dCurve.hxx>
  61 #include <Geom2dAdaptor_Curve.hxx>
  62 #include <Geom2dConvert.hxx>
  63 #include <Geom2d_TrimmedCurve.hxx>
  64 #include <ShapeBuild_ReShape.hxx>
  65
  66 //szv
  67 #include <TColgp_Array1OfPnt2d.hxx>
  68 #include <Geom2d_Circle.hxx>
  69 #include <Geom2d_Ellipse.hxx>
  70 #include <Geom2d_Parabola.hxx>
  71 #include <Geom2d_Hyperbola.hxx>
  72 #include <Geom2d_OffsetCurve.hxx>
  73 #include <Geom2dInt_GInter.hxx>
  74 #include <IntRes2d_Domain.hxx>
  75 #include <IntRes2d_IntersectionSegment.hxx>
  76 #include <Geom2dAPI_ExtremaCurveCurve.hxx>
  77 #include <Geom2dAPI_ProjectPointOnCurve.hxx>
  78 #include <Geom2dAdaptor_HCurve.hxx>
  79 #include <Approx_Curve2d.hxx>
  80 #include <Geom2dConvert.hxx>
  81
  82 #include <Geom_Line.hxx>
  83 #include <Geom_Circle.hxx>
  84 #include <Geom_Ellipse.hxx>
  85 #include <Geom_Parabola.hxx>
  86 #include <Geom_Hyperbola.hxx>
  87 #include <Geom_OffsetCurve.hxx>
  88 #include <GeomAPI_ExtremaCurveCurve.hxx>
  89 #include <GeomAPI_ProjectPointOnCurve.hxx>
  90 #include <GeomAdaptor_HCurve.hxx>
  91 #include <Approx_Curve3d.hxx>
  92 #include <GeomConvert.hxx>
  93 #include <TopoDS_Iterator.hxx>
  94 #include <ShapeFix_ShapeTolerance.hxx>
  95 #include <ShapeAnalysis_TransferParametersProj.hxx>
  96 //=======================================================================
  97 //function : FixGaps3d
  98 //purpose  :
  99 //=======================================================================
 100
 101 Standard_Boolean ShapeFix_Wire::FixGaps3d ()
 102 {
 103   myStatusGaps3d = ShapeExtend::EncodeStatus ( ShapeExtend_OK );
 104 // if ( !IsReady() ) return Standard_False;
 105
 106   Standard_Integer i, start = ( myClosedMode ? 1 : 2 );
 107   if (myFixGapsByRanges)
 108   {
 109     for ( i = start; i <= NbEdges(); i++ )
 110     {
 111       FixGap3d ( i );
 112       myStatusGaps3d |= myLastFixStatus;
 113     }
 114   }
 115   for ( i = start; i <= NbEdges(); i++ )
 116   {
 117     FixGap3d ( i, Standard_True );
 118     myStatusGaps3d |= myLastFixStatus;
 119   }
 120
 121   return StatusGaps3d ( ShapeExtend_DONE );
 122 }
 123
 124 //=======================================================================
 125 //function : FixGaps2d
 126 //purpose  :
 127 //=======================================================================
 128
 129  Standard_Boolean ShapeFix_Wire::FixGaps2d ()
 130 {
 131   myStatusGaps2d = ShapeExtend::EncodeStatus ( ShapeExtend_OK );
 132 //  if ( !IsReady() ) return Standard_False;
 133
 134   Standard_Integer i, start = ( myClosedMode ? 1 : 2 );
 135   if (myFixGapsByRanges)
 136   {
 137     for ( i = start; i <= NbEdges(); i++ )
 138     {
 139       FixGap2d ( i );
 140       myStatusGaps2d |= myLastFixStatus;
 141     }
 142   }
 143   for ( i = start; i <= NbEdges(); i++ )
 144   {
 145     FixGap2d ( i, Standard_True );
 146     myStatusGaps2d |= myLastFixStatus;
 147   }
 148
 149   return StatusGaps2d ( ShapeExtend_DONE );
 150 }
 151
 152 //=======================================================================
 153 //function : FixGap3d
 154 //purpose  :
 155 //=======================================================================
 156
 157 static Standard_Real AdjustOnPeriodic3d (const Handle(Geom_Curve)& c,
 158                                          const Standard_Boolean takefirst,
 159                                          const Standard_Real first,
 160                                          const Standard_Real last,
 161                                          const Standard_Real param)
 162 {
 163   // 15.11.2002 PTV OCC966
 164   if (ShapeAnalysis_Curve::IsPeriodic(c))
 165   {
 166     Standard_Real T = c->Period();
 167     Standard_Real shift = -IntegerPart(first/T)*T; if (first<0.) shift += T;
 168     Standard_Real sfirst = first+shift, slast = last+shift;
 169     if ( takefirst && (param>slast) && (param>sfirst)) return param-T-shift;
 170     if (!takefirst && (param<slast) && (param<sfirst)) return param+T-shift;
 171   }
 172   return param;
 173 }
 174
 175  Standard_Boolean ShapeFix_Wire::FixGap3d (const Standard_Integer num,
 176                                            const Standard_Boolean convert)
 177 {
 178   myLastFixStatus = ShapeExtend::EncodeStatus( ShapeExtend_OK );
 179 //  if ( !IsReady() ) return Standard_False;
 180
 181   //=============
 182   // First phase: analysis whether the problem (gap) exists
 183   //=============
 184
 185   if (Context().IsNull()) SetContext(new ShapeBuild_ReShape);
 186
 187   Standard_Real preci = Precision();
 188
 189   Handle(ShapeExtend_WireData) sbwd = WireData();
 190   Standard_Integer n2 = ( num >0 ? num  : sbwd->NbEdges() );
 191   Standard_Integer n1 = ( n2  >1 ? n2-1 : sbwd->NbEdges() );
 192 //smh#8
 193   TopoDS_Shape tmp1 = Context()->Apply(sbwd->Edge(n1)),
 194                tmp2 = Context()->Apply(sbwd->Edge(n2));
 195   TopoDS_Edge E1 = TopoDS::Edge(tmp1),
 196               E2 = TopoDS::Edge(tmp2);
 197 //  TopoDS_Face face = myAnalyzer->Face(); // comment by enk
 198
 199   // Retrieve curves on edges
 200   Standard_Real cfirst1, clast1, cfirst2, clast2;
 201   Handle(Geom_Curve) C1, C2;
 202   ShapeAnalysis_Edge SAE;
 203   if (!SAE.Curve3d(E1,C1,cfirst1,clast1) ||
 204       !SAE.Curve3d(E2,C2,cfirst2,clast2))
 205   {
 206     myLastFixStatus |= ShapeExtend::EncodeStatus( ShapeExtend_FAIL1 );
 207     return Standard_False;
 208   }
 209
 210   // Check gap in 3d space
 211   gp_Pnt cpnt1 = C1->Value(clast1), cpnt2 = C2->Value(cfirst2);
 212   Standard_Real gap = cpnt1.Distance(cpnt2);
 213   if (!convert && gap<=preci) return Standard_False;
 214
 215   //=============
 216   // Second phase: collecting data necessary for further analysis
 217   //=============
 218
 219   Standard_Boolean reversed1 = (E1.Orientation()==TopAbs_REVERSED),
 220                    reversed2 = (E2.Orientation()==TopAbs_REVERSED);
 221
 222   TopoDS_Vertex V1 = SAE.LastVertex(E1), V2 = SAE.FirstVertex(E2);
 223   gp_Pnt vpnt = (V1.IsSame(V2))? BRep_Tool::Pnt(V1) :
 224     gp_Pnt((BRep_Tool::Pnt(V1).XYZ()+BRep_Tool::Pnt(V2).XYZ())*0.5);
 225
 226   Standard_Real first1, last1, first2, last2;
 227   if (reversed1)
 228   {
 229     first1 = clast1;  last1 = cfirst1;
 230   }
 231   else
 232   {
 233     first1 = cfirst1; last1 = clast1;
 234   }
 235   if (reversed2)
 236   {
 237     first2 = clast2;  last2 = cfirst2;
 238   }
 239   else
 240   {
 241     first2 = cfirst2; last2 = clast2;
 242   }
 243
 244   Handle(Geom_Curve) c1 = C1, c2 = C2;
 245
 246   // Extract basic curves from trimmed and offset
 247   Standard_Boolean basic = Standard_False;
 248   Standard_Boolean trimmed1 = Standard_False, offset1 = Standard_False;
 249   gp_XYZ offval1(0.,0.,0.);
 250   while (!basic)
 251   {
 252     if (c1->IsKind(STANDARD_TYPE(Geom_TrimmedCurve)))
 253     {
 254       c1 = Handle(Geom_TrimmedCurve)::DownCast(c1)->BasisCurve();
 255       trimmed1 = Standard_True;
 256     }
 257     else if (c1->IsKind(STANDARD_TYPE(Geom_OffsetCurve)))
 258     {
 259       Handle(Geom_OffsetCurve) oc = Handle(Geom_OffsetCurve)::DownCast(c1);
 260       c1 = oc->BasisCurve();
 261       offval1 += oc->Offset()*oc->Direction().XYZ();
 262       offset1 = Standard_True;
 263     }
 264     else basic = Standard_True;
 265   }
 266   basic = Standard_False;
 267   Standard_Boolean trimmed2 = Standard_False, offset2 = Standard_False;
 268   gp_XYZ offval2(0.,0.,0.);
 269   while (!basic)
 270   {
 271     if (c2->IsKind(STANDARD_TYPE(Geom_TrimmedCurve)))
 272     {
 273       c2 = Handle(Geom_TrimmedCurve)::DownCast(c2)->BasisCurve();
 274       trimmed2 = Standard_True;
 275     }
 276     else if (c2->IsKind(STANDARD_TYPE(Geom_OffsetCurve)))
 277     {
 278       Handle(Geom_OffsetCurve) oc = Handle(Geom_OffsetCurve)::DownCast(c2);
 279       c2 = oc->BasisCurve();
 280       offval2 += oc->Offset()*oc->Direction().XYZ();
 281       offset2 = Standard_True;
 282     }
 283     else basic = Standard_True;
 284   }
 285   // Restore offset curves
 286   if (offset1) c1 = new Geom_OffsetCurve(c1,offval1.Modulus(),gp_Dir(offval1));
 287   if (offset2) c2 = new Geom_OffsetCurve(c2,offval2.Modulus(),gp_Dir(offval2));
 288
 289   Standard_Boolean done1 = Standard_False, done2 = Standard_False;
 290
 291   if (convert)
 292   {
 293
 294     Handle(Geom_BSplineCurve) bsp1, bsp2;
 295     Handle(Geom_Curve) c;
 296     Standard_Real first, last;
 297
 298     // iterate on curves
 299     Standard_Integer nbcurv = (n1==n2? 1 : 2);
 300     for (Standard_Integer j=1; j<=nbcurv; j++)
 301     {
 302       //Standard_Boolean trim = Standard_False;  // skl
 303       if (j==1)
 304       {
 305         if (cpnt1.Distance(vpnt)<preci)
 306         {
 307           if (n1==n2)
 308           {
 309             if (cpnt2.Distance(vpnt)<preci) continue;
 310           }
 311           else continue;
 312         }
 313         c = c1; first = first1; last = last1; /*trim = trimmed1;*/ // skl
 314       }
 315       else
 316       {
 317         if (cpnt2.Distance(vpnt)<preci) continue;
 318         c = c2; first = first2; last = last2; /*trim = trimmed2;*/ // skl
 319       }
 320
 321       Handle(Geom_BSplineCurve) bsp;
 322
 323       // Convert curve to bspline
 324       if (c->IsKind(STANDARD_TYPE(Geom_BSplineCurve)))
 325       {
 326         bsp = Handle(Geom_BSplineCurve)::DownCast(c->Copy());
 327         // take segment if trim and range differ
 328         Standard_Real fbsp = bsp->FirstParameter(), lbsp = bsp->LastParameter();
 329         Standard_Boolean segment = Standard_False;
 330         if (first>fbsp)
 331         {
 332           fbsp = first; segment = Standard_True;
 333         }
 334         if (last<lbsp)
 335         {
 336           lbsp = last; segment = Standard_True;
 337         }
 338         if (segment)
 339           bsp = GeomConvert::SplitBSplineCurve(bsp,fbsp,lbsp,
 340                                                ::Precision::Confusion());
 341       }
 342       else if (c->IsKind(STANDARD_TYPE(Geom_Conic)))
 343       {
 344         Approx_Curve3d Conv(new GeomAdaptor_HCurve(c,first,last),
 345                             myAnalyzer->Precision(),GeomAbs_C1,9,1000);
 346         if (Conv.IsDone() || Conv.HasResult()) bsp = Conv.Curve();
 347       }
 348       else
 349       {
 350         // Restore trim for pcurve
 351         Handle(Geom_Curve) tc ;
 352         try
 353         {
 354           OCC_CATCH_SIGNALS
 355           // 15.11.2002 PTV OCC966
 356           if(!ShapeAnalysis_Curve::IsPeriodic(c))
 357             tc = new Geom_TrimmedCurve(c,Max(first,c->FirstParameter()),Min(last,c->LastParameter()));
 358           else tc = new Geom_TrimmedCurve(c,first,last);
 359           bsp = GeomConvert::CurveToBSplineCurve(tc);
 360         }
 361         catch (Standard_Failure)
 362         {
 363 #ifdef SHAPEFIX_DEB
 364           cout << "Warning: ShapeFix_Wire_1::FixGap3d: Exception in TrimmedCurve" <<first<<" " <<last<<endl;
 365           Standard_Failure::Caught()->Print(cout); cout << endl;
 366 #endif
 367         }
 368       }
 369
 370       if (j==1) bsp1 = bsp; else bsp2 = bsp;
 371     }
 372
 373     // Take curves ends if could not convert
 374     if (bsp1.IsNull()) vpnt = cpnt1;
 375     else if (bsp2.IsNull()) vpnt = cpnt2;
 376
 377     if (!bsp1.IsNull())
 378     {
 379       if(bsp1->Degree() == 1) bsp1->IncreaseDegree(2); //gka
 380       if (n1==n2)
 381       {
 382         bsp1->SetPole(1,vpnt); bsp1->SetPole(bsp1->NbPoles(),vpnt);
 383       }
 384       else
 385       {
 386         if (reversed1) bsp1->SetPole(1,vpnt);
 387         else bsp1->SetPole(bsp1->NbPoles(),vpnt);
 388       }
 389       first1 = bsp1->FirstParameter(); last1 = bsp1->LastParameter();
 390       c1 = bsp1;
 391       done1 = Standard_True;
 392     }
 393     if (!bsp2.IsNull())
 394     {
 395       if(bsp2->Degree() == 1) bsp2->IncreaseDegree(2); //gka
 396       if (reversed2) bsp2->SetPole(bsp2->NbPoles(),vpnt);
 397       else bsp2->SetPole(1,vpnt);
 398       first2 = bsp2->FirstParameter(); last2 = bsp2->LastParameter();
 399       c2 = bsp2;
 400       done2 = Standard_True;
 401     }
 402   }
 403   else
 404   {
 405
 406     if (n1==n2)
 407     {
 408       if (c1->IsKind(STANDARD_TYPE(Geom_Circle)) ||
 409           c1->IsKind(STANDARD_TYPE(Geom_Ellipse)))
 410       {
 411         Standard_Real diff = M_PI - Abs(clast1-cfirst2)*0.5;
 412         first1 -= diff; last1 += diff;
 413         done1 = Standard_True;
 414       }
 415     }
 416     else
 417     {
 418
 419       // Determine domains for extremal points locating
 420       Standard_Real domfirst1 = first1, domlast1 = last1;
 421       if (c1->IsKind(STANDARD_TYPE(Geom_BSplineCurve)) ||
 422           c1->IsKind(STANDARD_TYPE(Geom_BezierCurve)))
 423       {
 424         domfirst1 = c1->FirstParameter();
 425         domlast1  = c1->LastParameter();
 426       }
 427       else if (c1->IsKind(STANDARD_TYPE(Geom_Line)) ||
 428                c1->IsKind(STANDARD_TYPE(Geom_Parabola)) ||
 429                c1->IsKind(STANDARD_TYPE(Geom_Hyperbola)))
 430       {
 431         Standard_Real diff = domlast1 - domfirst1;
 432         if (reversed1) domfirst1 -= 10.*diff;
 433         else           domlast1 += 10.*diff;
 434       }
 435       else if (c1->IsKind(STANDARD_TYPE(Geom_Circle)) ||
 436                c1->IsKind(STANDARD_TYPE(Geom_Ellipse)))
 437       {
 438         domfirst1 = 0.; domlast1 = 2*M_PI;
 439       }
 440       Standard_Real domfirst2 = first2, domlast2 = last2;
 441       if (c2->IsKind(STANDARD_TYPE(Geom_BSplineCurve)) ||
 442           c2->IsKind(STANDARD_TYPE(Geom_BezierCurve)))
 443       {
 444         domfirst2 = c2->FirstParameter();
 445         domlast2  = c2->LastParameter();
 446       }
 447       else if (c2->IsKind(STANDARD_TYPE(Geom_Line)) ||
 448                c2->IsKind(STANDARD_TYPE(Geom_Parabola)) ||
 449                c2->IsKind(STANDARD_TYPE(Geom_Hyperbola)))
 450       {
 451         Standard_Real diff = domlast2 - domfirst2;
 452         if (reversed2) domlast2 += 10.*diff;
 453         else           domfirst2 -= 10.*diff;
 454       }
 455       else if (c2->IsKind(STANDARD_TYPE(Geom_Circle)) ||
 456                c2->IsKind(STANDARD_TYPE(Geom_Ellipse)))
 457       {
 458         domfirst2 = 0.; domlast2 = 2*M_PI;
 459       }
 460
 461       Standard_Real ipar1 = clast1, ipar2 = cfirst2;
 462
 463       // Try to find projections of vertex point
 464       GeomAPI_ProjectPointOnCurve Proj;
 465       Standard_Real u1 = ipar1, u2 = ipar2;
 466       Proj.Init(vpnt,c1,domfirst1,domlast1);
 467       if (Proj.NbPoints())
 468       {
 469         Standard_Integer index = 1;
 470         Standard_Real dist, mindist=-1.;
 471         for (Standard_Integer i=1; i<=Proj.NbPoints(); i++)
 472         {
 473           dist = vpnt.Distance(Proj.Point(i));
 474           if (mindist>dist || mindist<0.)
 475           {
 476             index = i; mindist = dist;
 477           }
 478           u1 = Proj.Parameter(index);
 479         }
 480       }
 481       Proj.Init(vpnt,c2,domfirst2,domlast2);
 482       if (Proj.NbPoints())
 483       {
 484         Standard_Integer index = 1;
 485         Standard_Real dist, mindist=-1.;
 486         for (Standard_Integer i=1; i<=Proj.NbPoints(); i++)
 487         {
 488           dist = vpnt.Distance(Proj.Point(i));
 489           if (mindist>dist || mindist<0.)
 490           {
 491             index = i; mindist = dist;
 492           }
 493           u2 = Proj.Parameter(index);
 494         }
 495       }
 496       // Ajust parameters on periodic curves
 497       u1 = AdjustOnPeriodic3d(c1,reversed1,first1,last1,u1);
 498       u2 = AdjustOnPeriodic3d(c2,!reversed2,first2,last2,u2);
 499       // Check points to satisfy distance criterium
 500       gp_Pnt p1 = c1->Value(u1), p2 = c2->Value(u2);
 501       if (p1.Distance(p2)<=gap &&
 502           Abs(cfirst1-u1) > ::Precision::PConfusion() &&
 503           Abs(clast2-u2) > ::Precision::PConfusion() &&
 504           (((u1>first1) && (u1<last1)) || ((u2>first2) && (u2<last2)) ||
 505            (cpnt1.Distance(p1)<=gap) || (cpnt2.Distance(p2)<=gap)))
 506       {
 507         ipar1 = u1; ipar2 = u2;
 508         done1 = done2 = Standard_True;
 509       }
 510
 511       // Try to find closest points if nothing yet found
 512       if (!done1)
 513       {
 514
 515         // Recompute domains
 516         if (reversed1)
 517         {
 518           domfirst1 = ipar1; domlast1 = last1;
 519         }
 520         else
 521         {
 522           domfirst1 = first1; domlast1 = ipar1;
 523         }
 524         if (reversed2)
 525         {
 526           domfirst2 = first2; domlast2 = ipar2;
 527         }
 528         else
 529         {
 530           domfirst2 = ipar2; domlast2 = last2;
 531         }
 532
 533         GeomAPI_ExtremaCurveCurve Extr(c1,c2,domfirst1,domlast1,domfirst2,domlast2);
 534         if (Extr.NbExtrema())
 535         {
 536           try
 537           {
 538             OCC_CATCH_SIGNALS
 539             // First find all intersections
 540             gp_Pnt pp1, pp2;
 541             Standard_Integer index1=0, index2=0;
 542             Standard_Real uu1, uu2, pardist, pardist1=-1., pardist2=-1.;
 543             for (Standard_Integer i=1; i<=Extr.NbExtrema(); i++)
 544             {
 545               Extr.Parameters(i,uu1,uu2);
 546               // Ajust parameters on periodic curves
 547               uu1 = AdjustOnPeriodic3d(c1,reversed1,first1,last1,uu1);
 548               uu2 = AdjustOnPeriodic3d(c2,!reversed2,first2,last2,uu2);
 549               pp1 = c1->Value(uu1); pp2 = c2->Value(uu2);
 550               if (pp1.Distance(pp2) < ::Precision::Confusion())
 551               {
 552                 // assume intersection
 553                 pardist = Abs(cfirst1-uu1);
 554                 if (pardist1>pardist || pardist1<0.)
 555                 {
 556                   index1 = i; pardist1 = pardist;
 557                 }
 558                 pardist = Abs(clast2-uu2);
 559                 if (pardist2>pardist || pardist2<0.)
 560                 {
 561                   index2 = i; pardist2 = pardist;
 562                 }
 563               }
 564             }
 565             if (index1!=0 && index2!=0)
 566             {
 567               if (index1!=index2)
 568               {
 569                 // take intersection closer to vertex point
 570                 Extr.Parameters(index1,uu1,uu2);
 571                 pp1 = gp_Pnt((c1->Value(uu1).XYZ()+c2->Value(uu2).XYZ())*0.5);
 572                 Extr.Parameters(index2,uu1,uu2);
 573                 pp2 = gp_Pnt((c1->Value(uu1).XYZ()+c2->Value(uu2).XYZ())*0.5);
 574                 if (pp2.Distance(vpnt) < pp1.Distance(vpnt)) index1 = index2;
 575               }
 576               Extr.Parameters(index1,uu1,uu2);
 577             }
 578             else Extr.LowerDistanceParameters(uu1,uu2);
 579             // Ajust parameters on periodic curves
 580             uu1 = AdjustOnPeriodic3d(c1,reversed1,first1,last1,uu1);
 581             uu2 = AdjustOnPeriodic3d(c2,!reversed2,first2,last2,uu2);
 582             // Check points to satisfy distance criterium
 583             pp1 = c1->Value(uu1), pp2 = c2->Value(uu2);
 584             if (pp1.Distance(pp2)<=gap &&
 585                 Abs(cfirst1-uu1) > ::Precision::PConfusion() &&
 586                 Abs(clast2-uu2) > ::Precision::PConfusion() &&
 587                 (((uu1>first1) && (uu1<last1)) || ((uu2>first2) && (uu2<last2)) ||
 588                  (cpnt1.Distance(pp1)<=gap) || (cpnt2.Distance(pp2)<=gap)))
 589             {
 590               ipar1 = uu1; ipar2 = uu2;
 591               done1 = done2 = Standard_True;
 592             }
 593           }
 594           catch ( Standard_Failure )
 595           {
 596           }
 597         }
 598       }
 599
 600       try
 601       {
 602         OCC_CATCH_SIGNALS
 603         if (done1)
 604         {
 605           if (ipar1==clast1) done1 = Standard_False;
 606           else
 607           {
 608             // Set up new bounds for curve
 609             if (reversed1) first1 = ipar1; else last1 = ipar1;
 610             // Set new trim for old curve
 611             if (trimmed1)
 612             {
 613               // Standard_Real ff1 = c1->FirstParameter();
 614               // Standard_Real ll1 = c1->LastParameter();
 615               c1 = new Geom_TrimmedCurve(c1,first1,last1);
 616             }
 617           }
 618         }
 619         if (done2)
 620         {
 621           if (ipar2==cfirst2) done2 = Standard_False;
 622           else
 623           {
 624             // Set up new bounds for curve
 625             if (reversed2) last2 = ipar2; else first2 = ipar2;
 626             // Set new trim for old curve
 627             if (trimmed2)
 628             {
 629               // Standard_Real ff2 = c2->FirstParameter();
 630               // Standard_Real ll2 = c2->LastParameter();
 631               c2 = new Geom_TrimmedCurve(c2,first2,last2);
 632             }
 633           }
 634         }
 635       }
 636       catch (Standard_Failure)
 637       {
 638 #ifdef SHAPEFIX_DEB
 639         cout << "Warning: ShapeFix_Wire_1::FixGap3d: Exception in TrimmedCurve      :"<<endl;
 640         Standard_Failure::Caught()->Print(cout); cout << endl;
 641 #endif
 642       }
 643     }
 644   }
 645
 646   if (done1 || done2)
 647   {
 648
 649     BRep_Builder B;
 650     ShapeBuild_Edge SBE;
 651     ShapeFix_ShapeTolerance SFST;
 652
 653     // Update vertices
 654     TopoDS_Vertex nullV, newV1;
 655 //smh#8
 656     TopoDS_Shape emptyCopiedV2 = V2.EmptyCopied();
 657     TopoDS_Vertex newV2 = TopoDS::Vertex(emptyCopiedV2);
 658     SFST.SetTolerance(newV2,::Precision::Confusion());
 659     Context()->Replace(V2,newV2);
 660     if (V1.IsSame(V2))
 661 //smh#8
 662       {
 663         TopoDS_Shape tmpV2 = newV2.Oriented(TopAbs_REVERSED);
 664         newV1 = TopoDS::Vertex(tmpV2);
 665       }
 666     else
 667     {
 668 //smh#8
 669       TopoDS_Shape emptyCopied = V1.EmptyCopied();
 670       newV1 = TopoDS::Vertex(emptyCopied);
 671       SFST.SetTolerance(newV1,::Precision::Confusion());
 672       Context()->Replace(V1,newV1);
 673     }
 674
 675     if (done1)
 676     {
 677       // Update first edge
 678       TopoDS_Edge newE1 = SBE.CopyReplaceVertices(E1,nullV,newV1);
 679 //smh#8
 680       TopoDS_Shape tmpE1 = newE1.Oriented(TopAbs_FORWARD);
 681       B.UpdateEdge(TopoDS::Edge(tmpE1),c1,0.);
 682       SBE.SetRange3d(TopoDS::Edge(tmpE1),first1,last1);
 683       SFST.SetTolerance(newE1,::Precision::Confusion(),TopAbs_EDGE);
 684       B.SameRange(newE1,Standard_False);
 685 //      B.SameParameter(newE1,Standard_False);
 686
 687       //To keep NM vertices belonging initial edges
 688       TopoDS_Iterator aItv(E1,Standard_False);
 689       for( ; aItv.More(); aItv.Next()) {
 690         if(aItv.Value().Orientation() == TopAbs_INTERNAL ||
 691            aItv.Value().Orientation() == TopAbs_EXTERNAL) {
 692           TopoDS_Vertex aOldV = TopoDS::Vertex(aItv.Value());
 693           TopoDS_Vertex anewV =  ShapeAnalysis_TransferParametersProj::CopyNMVertex(aOldV,newE1,E1);
 694           B.Add(newE1,anewV);
 695           Context()->Replace(aOldV,anewV);
 696         }
 697       }
 698
 699       Context()->Replace(E1,newE1);
 700       sbwd->Set(newE1,n1);
 701     }
 702
 703     if (done2)
 704     {
 705       // Update second edge
 706       TopoDS_Edge newE2 = SBE.CopyReplaceVertices(E2,newV2,nullV);
 707 //smh#8
 708       TopoDS_Shape tmpE2 = newE2.Oriented(TopAbs_FORWARD);
 709       B.UpdateEdge(TopoDS::Edge(tmpE2),c2,0.);
 710       SBE.SetRange3d(TopoDS::Edge(tmpE2),first2,last2);
 711       SFST.SetTolerance(newE2,::Precision::Confusion(),TopAbs_EDGE);
 712       B.SameRange(newE2,Standard_False);
 713 //      B.SameParameter(newE2,Standard_False);
 714
 715       //To keep NM vertices belonging initial edges
 716       TopoDS_Iterator aItv(E2,Standard_False);
 717       for( ; aItv.More(); aItv.Next()) {
 718         if(aItv.Value().Orientation() == TopAbs_INTERNAL ||
 719            aItv.Value().Orientation() == TopAbs_EXTERNAL) {
 720           TopoDS_Vertex aOldV = TopoDS::Vertex(aItv.Value());
 721           TopoDS_Vertex anewV =  ShapeAnalysis_TransferParametersProj::CopyNMVertex(aOldV,newE2,E2);
 722           B.Add(newE2,anewV);
 723           Context()->Replace(aOldV,anewV);
 724         }
 725       }
 726       Context()->Replace(E2,newE2);
 727       sbwd->Set(newE2,n2);
 728     }
 729
 730     myLastFixStatus |= ShapeExtend::EncodeStatus( ShapeExtend_DONE1 );
 731   }
 732   else
 733     if (convert)
 734       myLastFixStatus |= ShapeExtend::EncodeStatus( ShapeExtend_FAIL2 );
 735
 736   return (done1 || done2);
 737 }
 738
 739 //=======================================================================
 740 //function : FixGap2d
 741 //purpose  :
 742 //=======================================================================
 743
 744 static Standard_Real AdjustOnPeriodic2d (const Handle(Geom2d_Curve)& pc,
 745                                          const Standard_Boolean takefirst,
 746                                          const Standard_Real first,
 747                                          const Standard_Real last,
 748                                          const Standard_Real param)
 749 {
 750   // 15.11.2002 PTV OCC966
 751   if (ShapeAnalysis_Curve::IsPeriodic(pc))
 752   {
 753     Standard_Real T = pc->Period();
 754     Standard_Real shift = -IntegerPart(first/T)*T; if (first<0.) shift += T;
 755     Standard_Real sfirst = first+shift, slast = last+shift;
 756     if ( takefirst && (param>slast) && (param>sfirst)) return param-T-shift;
 757     if (!takefirst && (param<slast) && (param<sfirst)) return param+T-shift;
 758   }
 759   return param;
 760 }
 761
 762  Standard_Boolean ShapeFix_Wire::FixGap2d (const Standard_Integer num,
 763                                            const Standard_Boolean convert)
 764 {
 765   myLastFixStatus = ShapeExtend::EncodeStatus( ShapeExtend_OK );
 766   if ( !IsReady() ) return Standard_False;
 767
 768   //=============
 769   // First phase: analysis whether the problem (gap) exists
 770   //=============
 771
 772   if (Context().IsNull()) SetContext(new ShapeBuild_ReShape);
 773
 774   Standard_Real preci = ::Precision::PConfusion();
 775   //Standard_Real preci = Precision();
 776   //GeomAdaptor_Surface& SA = Analyzer().Surface()->Adaptor()->ChangeSurface();
 777   //preci = Max(SA.UResolution(preci), SA.VResolution(preci));
 778
 779   Handle(ShapeExtend_WireData) sbwd = WireData();
 780   Standard_Integer n2 = ( num >0 ? num  : sbwd->NbEdges() );
 781   Standard_Integer n1 = ( n2  >1 ? n2-1 : sbwd->NbEdges() );
 782 //smh#8
 783   TopoDS_Shape tmp1 = Context()->Apply(sbwd->Edge(n1)),
 784                tmp2 = Context()->Apply(sbwd->Edge(n2));
 785   TopoDS_Edge E1 = TopoDS::Edge(tmp1),
 786               E2 = TopoDS::Edge(tmp2);
 787   TopoDS_Face face = myAnalyzer->Face();
 788
 789   // Retrieve pcurves on edges
 790   Standard_Real cfirst1, clast1, cfirst2, clast2;
 791   Handle(Geom2d_Curve) PC1, PC2;
 792   ShapeAnalysis_Edge SAE;
 793   if (!SAE.PCurve(E1,face,PC1,cfirst1,clast1) ||
 794       !SAE.PCurve(E2,face,PC2,cfirst2,clast2) ||
 795       sbwd->IsSeam(n1) || sbwd->IsSeam(n2))
 796   {
 797     myLastFixStatus |= ShapeExtend::EncodeStatus( ShapeExtend_FAIL1 );
 798     return Standard_False;
 799   }
 800
 801   // Check gap in 2d space
 802   gp_Pnt2d cpnt1 = PC1->Value(clast1), cpnt2 = PC2->Value(cfirst2);
 803   Standard_Real gap = cpnt1.Distance(cpnt2);
 804   if (gap<=preci) return Standard_False;
 805
 806   //=============
 807   // Second phase: collecting data necessary for further analysis
 808   //=============
 809
 810   Standard_Boolean reversed1 = (E1.Orientation()==TopAbs_REVERSED),
 811                    reversed2 = (E2.Orientation()==TopAbs_REVERSED);
 812
 813   Standard_Real first1, last1, first2, last2;
 814   if (reversed1)
 815   {
 816     first1 = clast1;  last1 = cfirst1;
 817   }
 818   else
 819   {
 820     first1 = cfirst1; last1 = clast1;
 821   }
 822   if (reversed2)
 823   {
 824     first2 = clast2;  last2 = cfirst2;
 825   }
 826   else
 827   {
 828     first2 = cfirst2; last2 = clast2;
 829   }
 830
 831   Handle(Geom2d_Curve) pc1 = PC1, pc2 = PC2;
 832
 833   // Extract basic curves from trimmed and offset
 834   Standard_Boolean basic = Standard_False;
 835   Standard_Boolean trimmed1 = Standard_False, offset1 = Standard_False;
 836   Standard_Real offval1 = 0.;
 837   while (!basic)
 838   {
 839     if (pc1->IsKind(STANDARD_TYPE(Geom2d_TrimmedCurve)))
 840     {
 841       pc1 = Handle(Geom2d_TrimmedCurve)::DownCast(pc1)->BasisCurve();
 842       trimmed1 = Standard_True;
 843     }
 844     else if (pc1->IsKind(STANDARD_TYPE(Geom2d_OffsetCurve)))
 845     {
 846       Handle(Geom2d_OffsetCurve) oc = Handle(Geom2d_OffsetCurve)::DownCast(pc1);
 847       pc1 = oc->BasisCurve();
 848       offval1 += oc->Offset();
 849       offset1 = Standard_True;
 850     }
 851     else basic = Standard_True;
 852   }
 853   basic = Standard_False;
 854   Standard_Boolean trimmed2 = Standard_False, offset2 = Standard_False;
 855   Standard_Real offval2 = 0.;
 856   while (!basic)
 857   {
 858     if (pc2->IsKind(STANDARD_TYPE(Geom2d_TrimmedCurve)))
 859     {
 860       pc2 = Handle(Geom2d_TrimmedCurve)::DownCast(pc2)->BasisCurve();
 861       trimmed2 = Standard_True;
 862     }
 863     else if (pc2->IsKind(STANDARD_TYPE(Geom2d_OffsetCurve)))
 864     {
 865       Handle(Geom2d_OffsetCurve) oc = Handle(Geom2d_OffsetCurve)::DownCast(pc2);
 866       pc2 = oc->BasisCurve();
 867       offval2 += oc->Offset();
 868       offset2 = Standard_True;
 869     }
 870     else basic = Standard_True;
 871   }
 872   // Restore offset curves
 873   if (offset1) pc1 = new Geom2d_OffsetCurve(pc1,offval1);
 874   if (offset2) pc2 = new Geom2d_OffsetCurve(pc2,offval2);
 875
 876   Standard_Boolean done1 = Standard_False, done2 = Standard_False;
 877
 878   // Determine same edge case
 879   if (convert)
 880   {
 881
 882     Handle(Geom2d_BSplineCurve) bsp1, bsp2;
 883     Handle(Geom2d_Curve) pc;
 884     Standard_Real first, last;
 885
 886     // iterate on pcurves
 887     Standard_Integer nbcurv = (n1==n2? 1 : 2);
 888     for (Standard_Integer j=1; j<=nbcurv; j++)
 889     {
 890       //Standard_Integer trim = Standard_False; // skl
 891       Handle(Geom2d_BSplineCurve) bsp;
 892
 893       if (j==1)
 894       {
 895         pc = pc1; first = first1; last = last1; /*trim = trimmed1;*/
 896       } // skl
 897       else
 898       {
 899         pc = pc2; first = first2; last = last2; /*trim = trimmed2;*/
 900       } // skl
 901
 902       // Convert pcurve to bspline
 903       if (pc->IsKind(STANDARD_TYPE(Geom2d_BSplineCurve)))
 904       {
 905         bsp = Handle(Geom2d_BSplineCurve)::DownCast(pc->Copy());
 906         // take segment if trim and range differ
 907         Standard_Real fbsp = bsp->FirstParameter(), lbsp = bsp->LastParameter();
 908         Standard_Boolean segment = Standard_False;
 909         if (first>fbsp)
 910         {
 911           fbsp = first; segment = Standard_True;
 912         }
 913         if (last<lbsp)
 914         {
 915           lbsp = last; segment = Standard_True;
 916         }
 917         if (segment)
 918           bsp = Geom2dConvert::SplitBSplineCurve(bsp,fbsp,lbsp,
 919                                                  ::Precision::PConfusion());
 920       }
 921       else if (pc->IsKind(STANDARD_TYPE(Geom2d_Conic)))
 922       {
 923         GeomAdaptor_Surface& AS = myAnalyzer->Surface()->Adaptor3d()->ChangeSurface();
 924         Standard_Real tolu = AS.UResolution(myAnalyzer->Precision()),
 925         tolv = AS.VResolution(myAnalyzer->Precision());
 926         Approx_Curve2d Conv(new Geom2dAdaptor_HCurve(pc,first,last),
 927                             first,last,tolu,tolv,GeomAbs_C1,9,1000);
 928         if (Conv.IsDone() || Conv.HasResult()) bsp = Conv.Curve();
 929       }
 930       else
 931       {
 932         // Restore trim for pcurve
 933         try
 934         {
 935           OCC_CATCH_SIGNALS
 936           Handle(Geom2d_Curve) c;
 937           // 15.11.2002 PTV OCC966
 938           if(!ShapeAnalysis_Curve::IsPeriodic(pc))
 939             c = new Geom2d_TrimmedCurve(pc,Max(first,pc->FirstParameter()),Min(last,pc->LastParameter()));
 940           else
 941             c = new Geom2d_TrimmedCurve(pc,first,last);
 942           bsp = Geom2dConvert::CurveToBSplineCurve(c);
 943         }
 944         catch (Standard_Failure)
 945         {
 946 #ifdef SHAPEFIX_DEB
 947           cout << "Warning: ShapeFix_Wire_1::FixGap2d: Exception in TrimmedCurve2d" <<first<<" " <<last<<endl;
 948           Standard_Failure::Caught()->Print(cout); cout << endl;
 949 #endif
 950         }
 951       }
 952
 953       if (j==1) bsp1 = bsp; else bsp2 = bsp;
 954     }
 955
 956     // Take curves ends if could not convert
 957     gp_Pnt2d mpnt((cpnt1.XY()+cpnt2.XY())*0.5);
 958     if (bsp1.IsNull()) mpnt = cpnt1;
 959     else if (bsp2.IsNull()) mpnt = cpnt2;
 960
 961     if (!bsp1.IsNull())
 962     {
 963       if(bsp1->Degree() == 1) bsp1->IncreaseDegree(2);
 964       if (n1==n2)
 965       {
 966         bsp1->SetPole(1,mpnt); bsp1->SetPole(bsp1->NbPoles(),mpnt);
 967       }
 968       else
 969       {
 970         if (reversed1) bsp1->SetPole(1,mpnt);
 971         else bsp1->SetPole(bsp1->NbPoles(),mpnt);
 972       }
 973       first1 = bsp1->FirstParameter(); last1 = bsp1->LastParameter();
 974       pc1 = bsp1;
 975       done1 = Standard_True;
 976     }
 977     if (!bsp2.IsNull())
 978     {
 979       if(bsp2->Degree() == 1) bsp2->IncreaseDegree(2);
 980       if (reversed2) bsp2->SetPole(bsp2->NbPoles(),mpnt);
 981       else bsp2->SetPole(1,mpnt);
 982       first2 = bsp2->FirstParameter(); last2 = bsp2->LastParameter();
 983       pc2 = bsp2;
 984       done2 = Standard_True;
 985     }
 986   }
 987   else
 988   {
 989
 990     if (n1==n2)
 991     {
 992       if (pc1->IsKind(STANDARD_TYPE(Geom2d_Circle)) ||
 993           pc1->IsKind(STANDARD_TYPE(Geom2d_Ellipse)))
 994       {
 995         Standard_Real diff = M_PI - Abs(clast1-cfirst2)*0.5;
 996         first1 -= diff; last1 += diff;
 997         done1 = Standard_True;
 998       }
 999     }
1000     else
1001     {
1002
1003       // Determine domains for extremal points locating
1004       Standard_Real domfirst1 = first1, domlast1 = last1;
1005       if (pc1->IsKind(STANDARD_TYPE(Geom2d_BSplineCurve)) ||
1006           pc1->IsKind(STANDARD_TYPE(Geom2d_BezierCurve)))
1007       {
1008         domfirst1 = pc1->FirstParameter();
1009         domlast1  = pc1->LastParameter();
1010       }
1011       else if (pc1->IsKind(STANDARD_TYPE(Geom2d_Line)) ||
1012                pc1->IsKind(STANDARD_TYPE(Geom2d_Parabola)) ||
1013                pc1->IsKind(STANDARD_TYPE(Geom2d_Hyperbola)))
1014       {
1015         Standard_Real diff = domlast1 - domfirst1;
1016         if (reversed1) domfirst1 -= 10.*diff;
1017         else           domlast1 += 10.*diff;
1018       }
1019       else if (pc1->IsKind(STANDARD_TYPE(Geom2d_Circle)) ||
1020                pc1->IsKind(STANDARD_TYPE(Geom2d_Ellipse)))
1021       {
1022         domfirst1 = 0.; domlast1 = 2*M_PI;
1023       }
1024       Standard_Real domfirst2 = first2, domlast2 = last2;
1025       if (pc2->IsKind(STANDARD_TYPE(Geom2d_BSplineCurve)) ||
1026           pc2->IsKind(STANDARD_TYPE(Geom2d_BezierCurve)))
1027       {
1028         domfirst2 = pc2->FirstParameter();
1029         domlast2  = pc2->LastParameter();
1030       }
1031       else if (pc2->IsKind(STANDARD_TYPE(Geom2d_Line)) ||
1032                pc2->IsKind(STANDARD_TYPE(Geom2d_Parabola)) ||
1033                pc2->IsKind(STANDARD_TYPE(Geom2d_Hyperbola)))
1034       {
1035         Standard_Real diff = domlast2 - domfirst2;
1036         if (reversed2) domlast2 += 10.*diff;
1037         else           domfirst2 -= 10.*diff;
1038       }
1039       else if (pc2->IsKind(STANDARD_TYPE(Geom2d_Circle)) ||
1040                pc2->IsKind(STANDARD_TYPE(Geom2d_Ellipse)))
1041       {
1042         domfirst2 = 0.; domlast2 = 2*M_PI;
1043       }
1044
1045       Standard_Real ipar1 = clast1, ipar2 = cfirst2;
1046
1047       Geom2dInt_GInter Inter;
1048       Standard_Real tolint = ::Precision::PConfusion();
1049
1050       Geom2dAdaptor_Curve AC1(pc1), AC2(pc2);
1051
1052       // Try to find intersection points
1053       IntRes2d_Domain dom1(pc1->Value(domfirst1),domfirst1,tolint,
1054                            pc1->Value(domlast1),domlast1,tolint);
1055       IntRes2d_Domain dom2(pc2->Value(domfirst2),domfirst2,tolint,
1056                            pc2->Value(domlast2),domlast2,tolint);
1057       Inter.Perform( AC1, dom1, AC2, dom2, tolint, tolint );
1058       if (Inter.IsDone())
1059       {
1060         if (Inter.NbPoints() || Inter.NbSegments())
1061         {
1062           Standard_Integer i, index1 = 0, index2 = 0;
1063           Standard_Real pardist, pardist1=-1., pardist2=-1.;
1064           // iterate on intersection points
1065           IntRes2d_IntersectionPoint IP;
1066           for ( i=1; i<=Inter.NbPoints(); i++ )
1067           {
1068             IP = Inter.Point(i);
1069             // Adjust parameters on periodic curves
1070             Standard_Real u1 = AdjustOnPeriodic2d(pc1,reversed1,first1,last1,
1071                                                   IP.ParamOnFirst());
1072             Standard_Real u2 = AdjustOnPeriodic2d(pc2,!reversed2,first2,last2,
1073                                                   IP.ParamOnSecond());
1074             pardist = Abs(cfirst1-u1);
1075             if (pardist1>pardist || pardist1<0.)
1076             {
1077               index1 = i; pardist1 = pardist;
1078             }
1079               pardist = Abs(clast2-u2);
1080             if (pardist2>pardist || pardist2<0.)
1081             {
1082               index2 = i; pardist2 = pardist;
1083             }
1084           }
1085           Standard_Integer flag1 = 0, flag2 = 0;
1086           // iterate on intersection segments
1087           IntRes2d_IntersectionSegment IS;
1088           for ( i=1; i<=Inter.NbSegments(); i++ )
1089           {
1090             IS = Inter.Segment(i);
1091             for (Standard_Integer j=1; j<=2; j++)
1092             {
1093               if ((j==1 && IS.HasFirstPoint()) ||
1094                   (j==2 && IS.HasLastPoint()))
1095               {
1096                 if (j==1) IP = IS.FirstPoint();
1097                 else      IP = IS.LastPoint();
1098                 // Adjust parameters on periodic curves
1099                 Standard_Real u1 = AdjustOnPeriodic2d(pc1,reversed1,first1,last1,
1100                                                       IP.ParamOnFirst());
1101                 Standard_Real u2 = AdjustOnPeriodic2d(pc2,!reversed2,first2,last2,
1102                                                       IP.ParamOnSecond());
1103                 pardist = Abs(cfirst1-u1);
1104                 if (pardist1>pardist || pardist1<0.)
1105                 {
1106                   flag1 = j; index1 = i; pardist1 = pardist;
1107                 }
1108                 pardist = Abs(clast2-u2);
1109                 if (pardist2>pardist || pardist2<0.)
1110                 {
1111                   flag2 = j; index2 = i; pardist2 = pardist;
1112                 }
1113               }
1114             }
1115           }
1116           if (index1!=index2 || flag1!=flag2)
1117           {
1118             // take intersection closer to mean point
1119             gp_Pnt2d pt1, pt2;
1120             if (flag1==0) pt1 = Inter.Point(index1).Value();
1121             else
1122             {
1123               IS = Inter.Segment(index1);
1124               if (flag1==1) pt1 = IS.FirstPoint().Value();
1125               else          pt1 = IS.LastPoint().Value();
1126             }
1127             if (flag2==0) pt2 = Inter.Point(index2).Value();
1128             else
1129             {
1130               IS = Inter.Segment(index2);
1131               if (flag2==1) pt2 = IS.FirstPoint().Value();
1132               else          pt2 = IS.LastPoint().Value();
1133             }
1134             gp_Pnt2d mpnt((cpnt1.XY()+cpnt2.XY())*0.5);
1135             if (pt2.Distance(mpnt) < pt1.Distance(mpnt))
1136             {
1137               index1 = index2; flag1 = flag2;
1138             }
1139           }
1140           if (flag1==0) IP = Inter.Point(index1);
1141           else
1142           {
1143             IS = Inter.Segment(index1);
1144             if (flag1==1) IP = IS.FirstPoint();
1145             else          IP = IS.LastPoint();
1146           }
1147           // Ajust parameters on periodic curves
1148           Standard_Real u1 = AdjustOnPeriodic2d(pc1,reversed1,first1,last1,
1149                                                 IP.ParamOnFirst());
1150           Standard_Real u2 = AdjustOnPeriodic2d(pc2,!reversed2,first2,last2,
1151                                                 IP.ParamOnSecond());
1152           // Check points to satisfy distance criterium
1153           gp_Pnt2d p1 = pc1->Value(u1), p2 = pc2->Value(u2);
1154           if (p1.Distance(p2)<=gap &&
1155               Abs(cfirst1-u1) > ::Precision::PConfusion() &&
1156               Abs(clast2 -u2) > ::Precision::PConfusion() &&
1157               (((u1>first1) && (u1<last1)) || ((u2>first2) && (u2<last2)) ||
1158                (cpnt1.Distance(p1)<=gap) || (cpnt2.Distance(p2)<=gap)))
1159           {
1160             ipar1 = u1; ipar2 = u2;
1161             done1 = done2 = Standard_True;
1162           }
1163         }
1164       }
1165
1166       // Try to find closest points if nothing yet found
1167       if (!done1)
1168       {
1169         Geom2dAPI_ExtremaCurveCurve Extr(pc1,pc2,domfirst1,domlast1,domfirst2,domlast2);
1170         if (Extr.NbExtrema())
1171         {
1172           Standard_Real u1, u2;
1173           Extr.LowerDistanceParameters(u1,u2);
1174           // Ajust parameters on periodic curves
1175           u1 = AdjustOnPeriodic2d(pc1,reversed1,first1,last1,u1);
1176           u2 = AdjustOnPeriodic2d(pc2,!reversed2,first2,last2,u2);
1177           // Check points to satisfy distance criterium
1178           gp_Pnt2d p1 = pc1->Value(u1), p2 = pc2->Value(u2);
1179           if (p1.Distance(p2)<=gap &&
1180               Abs(cfirst1-u1) > ::Precision::PConfusion() &&
1181               Abs(clast2 -u2) > ::Precision::PConfusion() &&
1182               (((u1>first1) && (u1<last1)) || ((u2>first2) && (u2<last2)) ||
1183                (cpnt1.Distance(p1)<=gap) || (cpnt2.Distance(p2)<=gap)))
1184           {
1185             ipar1 = u1; ipar2 = u2;
1186             done1 = done2 = Standard_True;
1187           }
1188         }
1189       }
1190
1191       // Try to find projections if nothing yet found
1192       if (!done1)
1193       {
1194         Geom2dAPI_ProjectPointOnCurve Proj;
1195         gp_Pnt2d ipnt1 = cpnt1, ipnt2 = cpnt2;
1196         Standard_Real u1 = ipar1, u2 = ipar2;
1197         Proj.Init(cpnt2,pc1,domfirst1,domlast1);
1198         if (Proj.NbPoints())
1199         {
1200           Standard_Integer index = 1;
1201           Standard_Real dist, mindist=-1.;
1202           for (Standard_Integer i=1; i<=Proj.NbPoints(); i++)
1203           {
1204             dist = cpnt2.Distance(Proj.Point(i));
1205             if (mindist>dist || mindist<0.)
1206             {
1207               index = i; mindist = dist;
1208             }
1209             ipnt1 = Proj.Point(index);
1210             u1 = Proj.Parameter(index);
1211           }
1212         }
1213         Proj.Init(cpnt1,pc2,domfirst2,domlast2);
1214         if (Proj.NbPoints())
1215         {
1216           Standard_Integer index = 1;
1217           Standard_Real dist, mindist=-1.;
1218           for (Standard_Integer i=1; i<=Proj.NbPoints(); i++)
1219           {
1220             dist = cpnt1.Distance(Proj.Point(i));
1221             if (mindist>dist || mindist<0.)
1222             {
1223               index = i; mindist = dist;
1224             }
1225             ipnt2 = Proj.Point(index);
1226             u2 = Proj.Parameter(index);
1227           }
1228         }
1229         // Ajust parameters on periodic curves
1230         u1 = AdjustOnPeriodic2d(pc1,reversed1,first1,last1,u1);
1231         u2 = AdjustOnPeriodic2d(pc2,!reversed2,first2,last2,u2);
1232         // Process special case of projection
1233         if ((((reversed1 && u1>clast1)  || (!reversed1 && u1<clast1)) &&
1234              ((reversed2 && u2<cfirst2) || (!reversed2 && u2>cfirst2))) ||
1235             (((reversed1 && u1<clast1)  || (!reversed1 && u1>clast1)) &&
1236              ((reversed2 && u2>cfirst2) || (!reversed2 && u2<cfirst2))))
1237         {
1238           // both projections lie inside/outside initial domains
1239           // project mean point
1240           gp_Pnt2d mpnt((cpnt1.XY()+cpnt2.XY())*0.5);
1241           u1 = ipar1; u2 = ipar2;
1242           Proj.Init(mpnt,pc1,domfirst1,domlast1);
1243           if (Proj.NbPoints())
1244           {
1245             Standard_Integer index = 1;
1246             Standard_Real dist, mindist=-1.;
1247             for (Standard_Integer i=1; i<=Proj.NbPoints(); i++)
1248             {
1249               dist = mpnt.Distance(Proj.Point(i));
1250               if (mindist>dist || mindist<0.)
1251               {
1252                 index = i; mindist = dist;
1253               }
1254               ipnt1 = Proj.Point(index);
1255               u1 = Proj.Parameter(index);
1256             }
1257           }
1258           Proj.Init(mpnt,pc2,domfirst2,domlast2);
1259           if (Proj.NbPoints())
1260           {
1261             Standard_Integer index = 1;
1262             Standard_Real dist, mindist=-1.;
1263             for (Standard_Integer i=1; i<=Proj.NbPoints(); i++)
1264             {
1265               dist = mpnt.Distance(Proj.Point(i));
1266               if (mindist>dist || mindist<0.)
1267               {
1268                 index = i; mindist = dist;
1269               }
1270               ipnt2 = Proj.Point(index);
1271               u2 = Proj.Parameter(index);
1272             }
1273           }
1274         }
1275         else
1276         {
1277           if (cpnt1.Distance(ipnt2)<cpnt2.Distance(ipnt1)) u1 = ipar1;
1278           else                                             u2 = ipar2;
1279         }
1280         // Ajust parameters on periodic curves
1281         u1 = AdjustOnPeriodic2d(pc1,reversed1,first1,last1,u1);
1282         u2 = AdjustOnPeriodic2d(pc2,!reversed2,first2,last2,u2);
1283         // Check points to satisfy distance criterium
1284         gp_Pnt2d p1 = pc1->Value(u1), p2 = pc2->Value(u2);
1285         if (p1.Distance(p2)<=gap &&
1286             Abs(cfirst1-u1) > ::Precision::PConfusion() &&
1287             Abs(clast2 -u2) > ::Precision::PConfusion() &&
1288             (((u1>first1) && (u1<last1)) || ((u2>first2) && (u2<last2)) ||
1289              (cpnt1.Distance(p1)<=gap) || (cpnt2.Distance(p2)<=gap)))
1290         {
1291           ipar1 = u1; ipar2 = u2;
1292           done1 = done2 = Standard_True;
1293         }
1294       }
1295
1296       if (done1)
1297       {
1298
1299         if (ipar1<first1 || ipar1>last1 || ipar2<first2 || ipar2>last2)
1300         {
1301
1302           // Check whether new points lie inside the surface bounds
1303           Standard_Real umin, umax, vmin, vmax;
1304           myAnalyzer->Surface()->Surface()->Bounds(umin,umax,vmin,vmax);
1305           if (::Precision::IsInfinite(umin) || ::Precision::IsInfinite(umax) ||
1306               ::Precision::IsInfinite(vmin) || ::Precision::IsInfinite(vmax))
1307           {
1308             Standard_Real fumin, fumax, fvmin, fvmax;
1309             BRepTools::UVBounds(face,fumin,fumax,fvmin,fvmax);
1310             if (::Precision::IsInfinite(umin)) umin = fumin-preci;
1311             if (::Precision::IsInfinite(umax)) umax = fumax+preci;
1312             if (::Precision::IsInfinite(vmin)) vmin = fvmin-preci;
1313             if (::Precision::IsInfinite(vmax)) vmax = fvmax+preci;
1314           }
1315
1316           gp_Pnt2d ipnt, P1, P2;
1317           Standard_Real u, v;
1318           Standard_Boolean out;
1319           // iterate on curves
1320           for (Standard_Integer j=1; j<=2; j++)
1321           {
1322
1323             if (j==1)
1324             {
1325               if (ipar1>=first1 && ipar1<=last1) continue;
1326               ipnt = pc1->Value(ipar1);
1327             }
1328             else
1329             {
1330               if (ipar2>=first2 && ipar2<=last2) continue;
1331               ipnt = pc2->Value(ipar2);
1332             }
1333
1334             // iterate on bounding lines
1335             for (Standard_Integer k=1; k<=2; k++)
1336             {
1337
1338               u = ipnt.X(); v = ipnt.Y();
1339
1340               out = Standard_True;
1341               if (k==1)
1342               {
1343                 if (u<umin)
1344                 {
1345                   P1 = gp_Pnt2d(umin,vmin); P2 = gp_Pnt2d(umin,vmax);
1346                 }
1347                 else if (u>umax)
1348                 {
1349                   P1 = gp_Pnt2d(umax,vmin); P2 = gp_Pnt2d(umax,vmax);
1350                 }
1351                 else out = Standard_False;
1352               }
1353               else
1354               {
1355                 if (v<vmin)
1356                 {
1357                   P1 = gp_Pnt2d(umin,vmin); P2 = gp_Pnt2d(umax,vmin);
1358                 }
1359                 else if (v>vmax)
1360                 {
1361                   P1 = gp_Pnt2d(umin,vmax); P2 = gp_Pnt2d(umax,vmax);
1362                 }
1363                 else out = Standard_False;
1364               }
1365
1366               if (out)
1367               {
1368                 // Intersect pcurve with bounding line
1369                 Handle(Geom2d_Line) lin = new Geom2d_Line(P1,gp_Dir2d(gp_Vec2d(P1,P2)));
1370                 Geom2dAdaptor_Curve ACL(lin);
1371                 IntRes2d_Domain dlin(P1,0.,tolint,P2,P1.Distance(P2),tolint);
1372
1373                 Handle(Geom2d_Curve) pc;
1374                 Standard_Real fpar, lpar;
1375                 if (j==1)
1376                 {
1377                   if (cfirst1<ipar1)
1378                   {
1379                     fpar = cfirst1, lpar = ipar1;
1380                   }
1381                   else
1382                   {
1383                     fpar = ipar1, lpar = cfirst1;
1384                   }
1385                   pc = pc1;
1386                 }
1387                 else
1388                 {
1389                   if (clast2<ipar2)
1390                   {
1391                     fpar = clast2, lpar = ipar2;
1392                   }
1393                   else
1394                   {
1395                     fpar = ipar2, lpar = clast2;
1396                   }
1397                   pc = pc2;
1398                 }
1399                 Geom2dAdaptor_Curve ACC(pc);
1400                 IntRes2d_Domain domc(pc->Value(fpar),fpar,tolint,
1401                                      pc->Value(lpar),lpar,tolint);
1402
1403                 // Intersect line with the pcurve
1404                 Inter.Perform( ACL, dlin, ACC, domc, tolint, tolint );
1405                 if (Inter.IsDone())
1406                 {
1407                   if (Inter.NbPoints() || Inter.NbSegments())
1408                   {
1409                     Standard_Integer i, index = 1;
1410                     Standard_Real uu, dist, mindist=-1.;
1411                     // iterate on intersection points
1412                     for ( i=1; i<=Inter.NbPoints(); i++ )
1413                     {
1414                       // Adjust parameters on periodic curve
1415                       uu = AdjustOnPeriodic2d(pc,(j==1? reversed1 : !reversed2),
1416                                               fpar,lpar,Inter.Point(i).ParamOnSecond());
1417                       dist = Abs((j==1? cfirst1 : clast2)-uu);
1418                       if (mindist>dist || mindist<0.)
1419                       {
1420                         index = i; mindist = dist;
1421                       }
1422                     }
1423                     // iterate on intersection segments
1424                     Standard_Integer flag = 0;
1425                     IntRes2d_IntersectionPoint IP;
1426                     IntRes2d_IntersectionSegment IS;
1427                     for ( i=1; i<=Inter.NbSegments(); i++ )
1428                     {
1429                       IS = Inter.Segment(i);
1430                       for (Standard_Integer jj=1; jj<=2; jj++)
1431                       {
1432                         if ((jj==1 && IS.HasFirstPoint()) ||
1433                             (jj==2 && IS.HasLastPoint()))
1434                         {
1435                           if (jj==1) IP = IS.FirstPoint();
1436                           else       IP = IS.LastPoint();
1437                           // Adjust parameters on periodic curve
1438                           uu = AdjustOnPeriodic2d(pc,(jj==1? reversed1 : !reversed2),
1439                                                   fpar,lpar,IP.ParamOnSecond());
1440                           dist = Abs((jj==1? cfirst1 : clast2)-uu);
1441                           if (mindist>dist || mindist<0.)
1442                           {
1443                             flag = jj; index = i; mindist = dist;
1444                           }
1445                         }
1446                       }
1447                     }
1448                     if (flag==0) IP = Inter.Point(index);
1449                     else
1450                     {
1451                       IS = Inter.Segment(index);
1452                       if (flag==1) IP = IS.FirstPoint();
1453                       else         IP = IS.LastPoint();
1454                     }
1455                     // Ajust parameters on periodic curve
1456                     uu = AdjustOnPeriodic2d(pc,(j==1? reversed1 : !reversed2),
1457                                             fpar,lpar,IP.ParamOnSecond());
1458                     if (j==1 && Abs(cfirst1-uu) > ::Precision::PConfusion())
1459                     {
1460                       ipar1 = uu; ipnt = IP.Value();
1461                     }
1462                     if (j==2 && Abs(clast2-uu) > ::Precision::PConfusion())
1463                     {
1464                       ipar2 = uu; ipnt = IP.Value();
1465                     }
1466                   }
1467                 }
1468               }
1469             }
1470
1471             // Adjust if intersection lies inside old bounds
1472             if (j==1)
1473             {
1474               if (reversed1)
1475               {
1476                 if (ipar1>first1) ipar1 = first1;
1477               }
1478               else
1479               {
1480                 if (ipar1<last1) ipar1 = last1;
1481               }
1482             }
1483             else
1484             {
1485               if (reversed2)
1486               {
1487                 if (ipar2<last2) ipar2 = last2;
1488               }
1489               else
1490               {
1491                 if (ipar2>first2) ipar2 = first2;
1492               }
1493             }
1494           }
1495         }
1496       }
1497       try
1498       {
1499         OCC_CATCH_SIGNALS
1500         if (done1)
1501         {
1502           if (ipar1==clast1) done1 = Standard_False;
1503           else
1504           {
1505             // Set up new bounds for pcurve
1506             if (reversed1) first1 = ipar1; else last1 = ipar1;
1507             // Set new trim for old pcurve
1508             if (trimmed1) pc1 = new Geom2d_TrimmedCurve(pc1,first1,last1);
1509           }
1510         }
1511         if (done2)
1512         {
1513           if (ipar2==cfirst2) done2 = Standard_False;
1514           else
1515           {
1516             // Set up new bounds for pcurve
1517             if (reversed2) last2 = ipar2; else first2 = ipar2;
1518             // Set new trim for old pcurve
1519             if (trimmed2) pc2 = new Geom2d_TrimmedCurve(pc2,first2,last2);
1520           }
1521         }
1522       }
1523       catch (Standard_Failure)
1524       {
1525
1526 #ifdef SHAPEFIX_DEB
1527         cout << "Warning: ShapeFix_Wire_1::FixGap2d: Exception in TrimmedCurve2d  :"<<endl;
1528         Standard_Failure::Caught()->Print(cout); cout << endl;
1529 #endif
1530       }
1531     }
1532   }
1533
1534   if (done1 || done2)
1535   {
1536
1537     BRep_Builder B;
1538     ShapeBuild_Edge SBE;
1539     ShapeFix_ShapeTolerance SFST;
1540
1541     // Update vertices
1542     TopoDS_Vertex V1 = SAE.LastVertex(E1), V2 = SAE.FirstVertex(E2);
1543     TopoDS_Vertex nullV, newV1;
1544 //smh#8
1545     TopoDS_Shape emptyCopiedV2 = V2.EmptyCopied();
1546     TopoDS_Vertex newV2 = TopoDS::Vertex(emptyCopiedV2);
1547     SFST.SetTolerance(newV2,::Precision::Confusion());
1548     Context()->Replace(V2,newV2);
1549     if (V1.IsSame(V2))
1550 //smh#8
1551     {
1552       TopoDS_Shape tmpVertexRev = newV2.Oriented(TopAbs_REVERSED);
1553       newV1 = TopoDS::Vertex(tmpVertexRev);
1554     }
1555     else
1556     {
1557 //smh#8
1558       TopoDS_Shape emptyCopiedV1 = V1.EmptyCopied();
1559       newV1 = TopoDS::Vertex(emptyCopiedV1);
1560       SFST.SetTolerance(newV1,::Precision::Confusion());
1561       Context()->Replace(V1,newV1);
1562     }
1563
1564     if (done1)
1565     {
1566       // Update first edge
1567       TopoDS_Edge newE1 = SBE.CopyReplaceVertices(E1,nullV,newV1);
1568 //smh#8
1569       TopoDS_Shape tmpE1 = newE1.Oriented(TopAbs_FORWARD);
1570       B.UpdateEdge(TopoDS::Edge(tmpE1),pc1,face,0.);
1571       B.Range(TopoDS::Edge(tmpE1),face,first1,last1);
1572       SFST.SetTolerance(newE1,::Precision::Confusion(),TopAbs_EDGE);
1573       B.SameRange(newE1,Standard_False);
1574 //      B.SameParameter(newE1,Standard_False);
1575
1576       //To keep NM vertices belonging initial edges
1577       TopoDS_Iterator aItv(E1,Standard_False);
1578       for( ; aItv.More(); aItv.Next()) {
1579         if(aItv.Value().Orientation() == TopAbs_INTERNAL ||
1580            aItv.Value().Orientation() == TopAbs_EXTERNAL) {
1581           TopoDS_Vertex aOldV = TopoDS::Vertex(aItv.Value());
1582           TopoDS_Vertex anewV =  ShapeAnalysis_TransferParametersProj::CopyNMVertex(aOldV,newE1,E1);
1583           B.Add(newE1,anewV);
1584           Context()->Replace(aOldV,anewV);
1585         }
1586       }
1587
1588       Context()->Replace(E1,newE1);
1589       sbwd->Set(newE1,n1);
1590     }
1591
1592     if (done2)
1593     {
1594       // Update second edge
1595       TopoDS_Edge newE2 = SBE.CopyReplaceVertices(E2,newV2,nullV);
1596 //smh#8
1597       TopoDS_Shape tmpE2 = newE2.Oriented(TopAbs_FORWARD);
1598       B.UpdateEdge(TopoDS::Edge(tmpE2),pc2,face,0.);
1599       B.Range(TopoDS::Edge(tmpE2),face,first2,last2);
1600       SFST.SetTolerance(newE2,::Precision::Confusion(),TopAbs_EDGE);
1601       B.SameRange(newE2,Standard_False);
1602 //      B.SameParameter(newE2,Standard_False);
1603       //To keep NM vertices belonging initial edges
1604       TopoDS_Iterator aItv(E2,Standard_False);
1605       for( ; aItv.More(); aItv.Next()) {
1606         if(aItv.Value().Orientation() == TopAbs_INTERNAL ||
1607            aItv.Value().Orientation() == TopAbs_EXTERNAL) {
1608           TopoDS_Vertex aOldV = TopoDS::Vertex(aItv.Value());
1609           TopoDS_Vertex anewV =  ShapeAnalysis_TransferParametersProj::CopyNMVertex(aOldV,newE2,E2);
1610           B.Add(newE2,anewV);
1611           Context()->Replace(aOldV,anewV);
1612         }
1613       }
1614       Context()->Replace(E2,newE2);
1615       sbwd->Set(newE2,n2);
1616     }
1617
1618     myLastFixStatus |= ShapeExtend::EncodeStatus( ShapeExtend_DONE1 );
1619   }
1620   else
1621     if (convert)
1622       myLastFixStatus |= ShapeExtend::EncodeStatus( ShapeExtend_FAIL2 );
1623
1624   return (done1 || done2);
1625 }