src/IntCurveSurface/IntCurveSurface_Polyhedron.gxx

   1 // File:        IntCurveSurface_Polyhedron.gxx
   2 // Created:     Wed Feb  3 11:33:09 1993
   3 // Author:      Laurent BUCHARD
   4 //              <lbr@topsn3>
   5
   6 #include <gp_Pnt.hxx>
   7 #include <gp_Vec.hxx>
   8 #include <gp_Dir.hxx>
   9 #include <gp_Lin.hxx>
  10 #include <TColgp_Array2OfPnt.hxx>
  11 #include <TColStd_Array2OfReal.hxx>
  12
  13 #include <Bnd_Array1OfBox.hxx>
  14
  15 #include <Standard_ConstructionError.hxx>
  16
  17 //#if  defined (WNT) || !defined (DEB)
  18 #include <stdio.h>
  19 //#endif
  20
  21 #define CHECKBOUNDS 0
  22
  23
  24 //-----------------------------------------------------
  25 #define LONGUEUR_MINI_EDGE_TRIANGLE 1e-15
  26
  27
  28 //=======================================================================
  29 //function : IntCurveSurface_Polyhedron
  30 //purpose  :
  31 //=======================================================================
  32 IntCurveSurface_Polyhedron::IntCurveSurface_Polyhedron (const ThePSurface&     Surface,
  33                                                         const Standard_Integer nbdU,
  34                                                         const Standard_Integer nbdV,
  35                                                         const Standard_Real    u1,
  36                                                         const Standard_Real    v1,
  37                                                         const Standard_Real    u2,
  38                                                         const Standard_Real    v2)
  39      : nbdeltaU((nbdU<3)? 3 : nbdU),
  40        nbdeltaV((nbdV<3)? 3 : nbdV),
  41        TheDeflection(Epsilon(100.)),
  42        C_MyPnts(NULL),C_MyU(NULL),C_MyV(NULL),C_MyIsOnBounds(NULL)
  43 {
  44   Standard_Integer t = (nbdeltaU+1)*(nbdeltaV+1)+1;
  45   gp_Pnt           *CMyPnts = new gp_Pnt[t];         C_MyPnts = (void *)CMyPnts;
  46   Standard_Real    *CMyU    = new Standard_Real[t];  C_MyU    = (void *)CMyU;
  47   Standard_Real    *CMyV    = new Standard_Real[t];  C_MyV    = (void *)CMyV;
  48
  49 //  Modified by Sergey KHROMOV - Fri Dec  7 12:03:46 2001 Begin
  50   Standard_Boolean *CMyIsOnBounds = new Standard_Boolean[t];
  51
  52   C_MyIsOnBounds = (void *)CMyIsOnBounds;
  53 //  Modified by Sergey KHROMOV - Fri Dec  7 12:03:47 2001 End
  54   Init(Surface,u1,v1,u2,v2);
  55 }
  56
  57 //=======================================================================
  58 //function : IntCurveSurface_Polyhedron
  59 //purpose  :
  60 //=======================================================================
  61 IntCurveSurface_Polyhedron::IntCurveSurface_Polyhedron (const ThePSurface&     Surface,
  62                                                         const TColStd_Array1OfReal&    Upars,
  63                                                         const TColStd_Array1OfReal&    Vpars)
  64     : nbdeltaU(Upars.Length()-1),
  65       nbdeltaV(Vpars.Length()-1),
  66       TheDeflection(Epsilon(100.)),
  67       C_MyPnts(NULL),C_MyU(NULL),C_MyV(NULL),C_MyIsOnBounds(NULL)
  68 {
  69   Standard_Integer t = (nbdeltaU+1)*(nbdeltaV+1)+1;
  70   gp_Pnt           *CMyPnts = new gp_Pnt[t];         C_MyPnts = (void *)CMyPnts;
  71   Standard_Real    *CMyU    = new Standard_Real[t];  C_MyU    = (void *)CMyU;
  72   Standard_Real    *CMyV    = new Standard_Real[t];  C_MyV    = (void *)CMyV;
  73
  74 //  Modified by Sergey KHROMOV - Fri Dec  7 12:03:46 2001 Begin
  75   Standard_Boolean *CMyIsOnBounds = new Standard_Boolean[t];
  76
  77   C_MyIsOnBounds = (void *)CMyIsOnBounds;
  78 //  Modified by Sergey KHROMOV - Fri Dec  7 12:03:47 2001 End
  79   Init(Surface, Upars, Vpars);
  80 }
  81
  82
  83 void IntCurveSurface_Polyhedron::Destroy() {
  84   //-- printf("\n IntCurveSurface_Polyhedron::Destroy\n");
  85   gp_Pnt *CMyPnts     = (gp_Pnt *)C_MyPnts;       if(C_MyPnts) delete [] CMyPnts;
  86   Standard_Real *CMyU = (Standard_Real *)C_MyU;   if(C_MyU)    delete [] CMyU;
  87   Standard_Real *CMyV = (Standard_Real *)C_MyV;   if(C_MyV)    delete [] CMyV;
  88
  89 //  Modified by Sergey KHROMOV - Fri Dec  7 12:03:46 2001 Begin
  90   Standard_Boolean *CMyIsOnBounds = (Standard_Boolean *)C_MyIsOnBounds;
  91
  92   if(C_MyIsOnBounds)    delete [] CMyIsOnBounds;
  93 //  Modified by Sergey KHROMOV - Fri Dec  7 12:03:47 2001 End
  94
  95   C_MyPnts=C_MyU=C_MyV=C_MyIsOnBounds=NULL;
  96 }
  97
  98 //=======================================================================
  99 //function : Init
 100 //purpose  :
 101 //=======================================================================
 102 void IntCurveSurface_Polyhedron::Init(const ThePSurface&     Surface,
 103                                       const Standard_Real U0,
 104                                       const Standard_Real V0,
 105                                       const Standard_Real U1,
 106                                       const Standard_Real V1) {
 107   const Standard_Integer DebugDump = 0;
 108   Standard_Integer i1,i2;
 109   Standard_Real    U,V;
 110   Standard_Real    U1mU0sNbdeltaU = (U1-U0)/(Standard_Real)nbdeltaU;
 111   Standard_Real    V1mV0sNbdeltaV = (V1-V0)/(Standard_Real)nbdeltaV;
 112   gp_Pnt TP;
 113   Standard_Integer Index=1;
 114   //-- --------------------------------------------------
 115   //-- Index varie de 1 -> (nbdu+1)*(nbdv+1)
 116   //-- V est la colonne
 117   //-- U est la ligne
 118   //-- --------------------------------------------------
 119   gp_Pnt *CMyPnts     = (gp_Pnt *)C_MyPnts;
 120   Standard_Real *CMyU = (Standard_Real *)C_MyU;
 121   Standard_Real *CMyV = (Standard_Real *)C_MyV;
 122   Standard_Boolean *CMyIsOnBounds = (Standard_Boolean *)C_MyIsOnBounds;
 123
 124   for (i1 = 0, U = U0; i1 <= nbdeltaU; i1++, U+= U1mU0sNbdeltaU) {
 125     for (i2 = 0, V = V0; i2 <= nbdeltaV; i2++, V+= V1mV0sNbdeltaV ) {
 126       ThePSurfaceTool::D0(Surface,U,V,TP);
 127       //-- Point(TP,i1, i2,U,V);
 128       CMyPnts[Index] = TP;
 129       CMyU[Index]    = U;
 130       CMyV[Index]    = V;
 131 //  Modified by Sergey KHROMOV - Fri Dec  7 12:07:51 2001
 132       CMyIsOnBounds[Index] = (i1 == 0 || i1 == nbdeltaU ||
 133                               i2 == 0 || i2 == nbdeltaV);
 134 //  Modified by Sergey KHROMOV - Fri Dec  7 12:07:52 2001
 135       TheBnd.Add(TP);
 136       Index++;
 137     }
 138   }
 139   //-- Calcul de la deflection Triangle <-> Point milieu
 140   Standard_Real tol=0.0; Standard_Integer nbtriangles = NbTriangles();
 141   for (i1=1; i1<=nbtriangles; i1++) {
 142     Standard_Real tol1 = DeflectionOnTriangle(Surface,i1);
 143     if(tol1>tol) tol=tol1;
 144   }
 145   //-- Calcul de la deflection Bord <-> Point milieu
 146
 147
 148   DeflectionOverEstimation(tol*1.2);
 149   FillBounding();
 150
 151 //  Modified by Sergey KHROMOV - Fri Dec  7 11:23:33 2001 Begin
 152   Standard_Real aDeflection;
 153
 154   TheBorderDeflection = RealFirst();
 155
 156 // Compute the deflection on the lower bound (U-isoline) of the surface.
 157   aDeflection = ComputeBorderDeflection(Surface, U0, V0, V1, Standard_True);
 158
 159   if (aDeflection > TheBorderDeflection)
 160     TheBorderDeflection = aDeflection;
 161
 162 // Compute the deflection on the upper bound (U-isoline) of the surface.
 163   aDeflection = ComputeBorderDeflection(Surface, U1, V0, V1, Standard_True);
 164
 165   if (aDeflection > TheBorderDeflection)
 166     TheBorderDeflection = aDeflection;
 167
 168 // Compute the deflection on the lower bound (V-isoline) of the surface.
 169   aDeflection = ComputeBorderDeflection(Surface, V0, U0, U1, Standard_False);
 170
 171   if (aDeflection > TheBorderDeflection)
 172     TheBorderDeflection = aDeflection;
 173
 174 // Compute the deflection on the upper bound (V-isoline) of the surface.
 175   aDeflection = ComputeBorderDeflection(Surface, V1, U0, U1, Standard_False);
 176
 177   if (aDeflection > TheBorderDeflection)
 178     TheBorderDeflection = aDeflection;
 179
 180 //  Modified by Sergey KHROMOV - Fri Dec  7 11:23:34 2001 End
 181
 182   if(DebugDump) {
 183     Dump();
 184   }
 185 }
 186 //=======================================================================
 187 //function : Init
 188 //purpose  :
 189 //=======================================================================
 190 void IntCurveSurface_Polyhedron::Init(const ThePSurface&     Surface,
 191                                       const TColStd_Array1OfReal& Upars,
 192                                       const TColStd_Array1OfReal& Vpars) {
 193   const Standard_Integer DebugDump = 0;
 194   Standard_Integer i1,i2;
 195   Standard_Real    U,V;
 196   gp_Pnt TP;
 197   Standard_Integer Index=1;
 198   //-- --------------------------------------------------
 199   //-- Index varie de 1 -> (nbdu+1)*(nbdv+1)
 200   //-- V est la colonne
 201   //-- U est la ligne
 202   //-- --------------------------------------------------
 203   gp_Pnt *CMyPnts     = (gp_Pnt *)C_MyPnts;
 204   Standard_Real *CMyU = (Standard_Real *)C_MyU;
 205   Standard_Real *CMyV = (Standard_Real *)C_MyV;
 206   Standard_Boolean *CMyIsOnBounds = (Standard_Boolean *)C_MyIsOnBounds;
 207   Standard_Integer i0 = Upars.Lower(), j0 = Vpars.Lower();
 208
 209   for (i1 = 0; i1 <= nbdeltaU; i1++) {
 210     U = Upars(i1+i0);
 211     for (i2 = 0; i2 <= nbdeltaV; i2++) {
 212       V = Vpars(i2+j0);
 213       ThePSurfaceTool::D0(Surface,U,V,TP);
 214       //-- Point(TP,i1, i2,U,V);
 215       CMyPnts[Index] = TP;
 216       CMyU[Index]    = U;
 217       CMyV[Index]    = V;
 218 //  Modified by Sergey KHROMOV - Fri Dec  7 12:07:51 2001
 219       CMyIsOnBounds[Index] = (i1 == 0 || i1 == nbdeltaU ||
 220                               i2 == 0 || i2 == nbdeltaV);
 221 //  Modified by Sergey KHROMOV - Fri Dec  7 12:07:52 2001
 222       TheBnd.Add(TP);
 223       Index++;
 224     }
 225   }
 226   //-- Calcul de la deflection Triangle <-> Point milieu
 227   Standard_Real tol=0.0; Standard_Integer nbtriangles = NbTriangles();
 228   for (i1=1; i1<=nbtriangles; i1++) {
 229     Standard_Real tol1 = DeflectionOnTriangle(Surface,i1);
 230     if(tol1>tol) tol=tol1;
 231   }
 232   //-- Calcul de la deflection Bord <-> Point milieu
 233
 234
 235   DeflectionOverEstimation(tol*1.2);
 236   FillBounding();
 237
 238 //  Modified by Sergey KHROMOV - Fri Dec  7 11:23:33 2001 Begin
 239   Standard_Real aDeflection;
 240
 241   TheBorderDeflection = RealFirst();
 242   Standard_Real U0 = Upars(i0);
 243   Standard_Real V0 = Vpars(j0);
 244   Standard_Real U1 = Upars(Upars.Upper());
 245   Standard_Real V1 = Vpars(Vpars.Upper());
 246
 247 // Compute the deflection on the lower bound (U-isoline) of the surface.
 248   aDeflection = ComputeBorderDeflection(Surface, U0, V0, V1, Standard_True);
 249
 250   if (aDeflection > TheBorderDeflection)
 251     TheBorderDeflection = aDeflection;
 252
 253 // Compute the deflection on the upper bound (U-isoline) of the surface.
 254   aDeflection = ComputeBorderDeflection(Surface, U1, V0, V1, Standard_True);
 255
 256   if (aDeflection > TheBorderDeflection)
 257     TheBorderDeflection = aDeflection;
 258
 259 // Compute the deflection on the lower bound (V-isoline) of the surface.
 260   aDeflection = ComputeBorderDeflection(Surface, V0, U0, U1, Standard_False);
 261
 262   if (aDeflection > TheBorderDeflection)
 263     TheBorderDeflection = aDeflection;
 264
 265 // Compute the deflection on the upper bound (V-isoline) of the surface.
 266   aDeflection = ComputeBorderDeflection(Surface, V1, U0, U1, Standard_False);
 267
 268   if (aDeflection > TheBorderDeflection)
 269     TheBorderDeflection = aDeflection;
 270
 271 //  Modified by Sergey KHROMOV - Fri Dec  7 11:23:34 2001 End
 272
 273   if(DebugDump) {
 274     Dump();
 275   }
 276 }
 277 //=======================================================================
 278 //function : DeflectionOnTriangle
 279 //purpose  :
 280 //=======================================================================
 281 Standard_Real IntCurveSurface_Polyhedron::DeflectionOnTriangle (const ThePSurface& Surface,
 282                                                                 const Standard_Integer Triang) const
 283 {
 284   Standard_Integer i1,i2,i3;
 285   Triangle(Triang,i1,i2,i3);
 286   //-- Calcul de l equation du plan
 287   Standard_Real u1,v1,u2,v2,u3,v3;
 288   gp_Pnt P1,P2,P3;
 289   P1 = Point(i1,u1,v1);
 290   P2 = Point(i2,u2,v2);
 291   P3 = Point(i3,u3,v3);
 292   if(P1.SquareDistance(P2)<=LONGUEUR_MINI_EDGE_TRIANGLE) return(0);
 293   if(P1.SquareDistance(P3)<=LONGUEUR_MINI_EDGE_TRIANGLE) return(0);
 294   if(P2.SquareDistance(P3)<=LONGUEUR_MINI_EDGE_TRIANGLE) return(0);
 295   gp_XYZ XYZ1=P2.XYZ()-P1.XYZ();
 296   gp_XYZ XYZ2=P3.XYZ()-P2.XYZ();
 297   gp_XYZ XYZ3=P1.XYZ()-P3.XYZ();
 298   gp_Vec NormalVector((XYZ1^XYZ2)+(XYZ2^XYZ3)+(XYZ3^XYZ1));
 299   NormalVector.Normalize();
 300   //-- Standard_Real PolarDistance = NormalVector * P1.XYZ();
 301   //-- Calcul du point u,v  au centre du triangle
 302   Standard_Real u = (u1+u2+u3)/3.0;
 303   Standard_Real v = (v1+v2+v3)/3.0;
 304   gp_Pnt P =  ThePSurfaceTool::Value(Surface,u,v);
 305   gp_Vec P1P(P1,P);
 306   return(Abs(P1P.Dot(NormalVector)));
 307 }
 308 //=======================================================================
 309 //function : Parameters
 310 //purpose  :
 311 //=======================================================================
 312 void IntCurveSurface_Polyhedron::Parameters( const Standard_Integer Index
 313                                             ,Standard_Real &U
 314                                             ,Standard_Real &V) const
 315 {
 316 #if CHECKBOUNDS
 317   if(Index<0 || Index>((nbdeltaU+1)*(nbdeltaV+1))) {
 318     printf("\n Erreur IntCurveSurface_Polyhedron::Parameters\n");
 319   }
 320 #endif
 321   Standard_Real *CMyU = (Standard_Real *)C_MyU;
 322   U = CMyU[Index];
 323   Standard_Real *CMyV = (Standard_Real *)C_MyV;
 324   V = CMyV[Index];
 325 }
 326 //=======================================================================
 327 //function : DeflectionOverEstimation
 328 //purpose  : Set
 329 //=======================================================================
 330 void IntCurveSurface_Polyhedron::DeflectionOverEstimation(const Standard_Real flec)
 331 {
 332   if(flec<0.0001) {
 333     TheDeflection=0.0001;
 334     TheBnd.Enlarge(0.0001);
 335   }
 336   else {
 337     TheDeflection=flec;
 338     TheBnd.Enlarge(flec);
 339   }
 340 }
 341 //=======================================================================
 342 //function : DeflectionOverEstimation
 343 //purpose  : Get
 344 //=======================================================================
 345 Standard_Real IntCurveSurface_Polyhedron::DeflectionOverEstimation() const
 346 {
 347   return TheDeflection;
 348 }
 349 //=======================================================================
 350 //function : Bounding
 351 //purpose  :
 352 //=======================================================================
 353 const Bnd_Box& IntCurveSurface_Polyhedron::Bounding() const
 354 {
 355   return TheBnd;
 356 }
 357 //=======================================================================
 358 //function : FillBounding
 359 //purpose  :
 360 //=======================================================================
 361 void IntCurveSurface_Polyhedron::FillBounding()
 362 {
 363   TheComponentsBnd=new Bnd_HArray1OfBox(1, NbTriangles());
 364   Bnd_Box Boite;
 365   Standard_Integer np1, np2, np3;
 366   Standard_Integer nbtriangles = NbTriangles();
 367   for (Standard_Integer iTri=1; iTri<=nbtriangles; iTri++) {
 368     Triangle(iTri, np1, np2, np3);
 369     gp_Pnt p1(Point(np1));
 370     gp_Pnt p2(Point(np2));
 371     gp_Pnt p3(Point(np3));
 372     Boite.SetVoid();
 373     if(p1.SquareDistance(p2)>LONGUEUR_MINI_EDGE_TRIANGLE) {
 374       if(p1.SquareDistance(p3)>LONGUEUR_MINI_EDGE_TRIANGLE) {
 375         if(p2.SquareDistance(p3)>LONGUEUR_MINI_EDGE_TRIANGLE) {
 376           Boite.Add(p1);
 377           Boite.Add(p2);
 378           Boite.Add(p3);
 379           Boite.Enlarge(TheDeflection);
 380         }
 381       }
 382     }
 383     Boite.Enlarge(TheDeflection);
 384     TheComponentsBnd->SetValue(iTri,Boite);
 385   }
 386 }
 387 //=======================================================================
 388 //function : ComponentsBounding
 389 //purpose  :
 390 //=======================================================================
 391 const Handle(Bnd_HArray1OfBox)&
 392       IntCurveSurface_Polyhedron::ComponentsBounding() const
 393 {
 394  return TheComponentsBnd;
 395 }
 396 //=======================================================================
 397 //function : NbTriangles
 398 //purpose  :
 399 //=======================================================================
 400 Standard_Integer IntCurveSurface_Polyhedron::NbTriangles  () const
 401 {
 402   return nbdeltaU*nbdeltaV*2;
 403 }
 404 //=======================================================================
 405 //function : NbPoints
 406 //purpose  :
 407 //=======================================================================
 408 Standard_Integer IntCurveSurface_Polyhedron::NbPoints () const
 409 {
 410   return (nbdeltaU+1)*(nbdeltaV+1);
 411 }
 412 //=======================================================================
 413 //function : TriConnex
 414 //purpose  :
 415 //=======================================================================
 416 Standard_Integer IntCurveSurface_Polyhedron::TriConnex
 417   (const Standard_Integer Triang,
 418    const Standard_Integer Pivot,
 419    const Standard_Integer Pedge,
 420    Standard_Integer&      TriCon,
 421    Standard_Integer&      OtherP) const
 422 {
 423 #ifdef DEB
 424   Standard_Integer nbdeltaUp1 = nbdeltaU + 1;
 425   Standard_Integer nbdeltaUm2 = nbdeltaU + nbdeltaU;
 426 #endif
 427   Standard_Integer Pivotm1    = Pivot-1;
 428   Standard_Integer nbdeltaVp1 = nbdeltaV+1;
 429   Standard_Integer nbdeltaVm2 = nbdeltaV + nbdeltaV;
 430
 431 // Pivot position in the MaTriangle :
 432   Standard_Integer ligP = Pivotm1/nbdeltaVp1;
 433   Standard_Integer colP = Pivotm1 - ligP * nbdeltaVp1;
 434
 435 // Point sur Edge position in the MaTriangle and edge typ :
 436 #ifndef DEB
 437   Standard_Integer ligE =0, colE =0, typE =0;
 438 #else
 439   Standard_Integer ligE, colE, typE;
 440 #endif
 441   if (Pedge!=0) {
 442     ligE= (Pedge-1)/nbdeltaVp1;
 443     colE= (Pedge-1) - (ligE * nbdeltaVp1);
 444   // Horizontal
 445     if      (ligP==ligE) typE=1;
 446   // Vertical
 447     else if (colP==colE) typE=2;
 448   // Oblique
 449     else                 typE=3;
 450   }
 451   else {
 452     typE=0;
 453   }
 454
 455 // Triangle position General case :
 456 #ifndef DEB
 457   Standard_Integer linT =0, colT =0;
 458   Standard_Integer linO =0, colO =0;
 459   Standard_Integer t =0, tt =0;
 460 #else
 461   Standard_Integer linT, colT;
 462   Standard_Integer linO, colO;
 463   Standard_Integer t,tt;
 464 #endif
 465   if (Triang!=0) {
 466     t = (Triang-1)/(nbdeltaVm2);
 467     tt= (Triang-1)-t*nbdeltaVm2;
 468     linT= 1+t;
 469     colT= 1+tt;
 470     if (typE==0) {
 471       if (ligP==linT) {
 472         ligE=ligP-1;
 473         colE=colP-1;
 474         typE=3;
 475       }
 476       else {
 477         if (colT==ligP+ligP) {
 478           ligE=ligP;
 479           colE=colP-1;
 480           typE=1;
 481         }
 482         else {
 483           ligE=ligP+1;
 484           colE=colP+1;
 485           typE=3;
 486         }
 487       }
 488     }
 489     switch (typE) {
 490     case 1:  // Horizontal
 491       if (linT==ligP) {
 492         linT++;
 493         linO=ligP+1;
 494         colO=(colP>colE)? colP : colE; //--colO=Max(colP, colE);
 495       }
 496       else {
 497         linT--;
 498         linO=ligP-1;
 499         colO=(colP<colE)? colP : colE;  //--colO=Min(colP, colE);
 500       }
 501       break;
 502     case 2:  // Vertical
 503       if (colT==(colP+colP)) {
 504         colT++;
 505         linO=(ligP>ligE)? ligP : ligE;  //--linO=Max(ligP, ligE);
 506         colO=colP+1;;
 507       }
 508       else {
 509         colT--;
 510         linO=(ligP<ligE)? ligP : ligE;  //--linO=Min(ligP, ligE);
 511         colO=colP-1;;
 512       }
 513       break;
 514     case 3:  // Oblique
 515       if ((colT&1)==0) {
 516         colT--;
 517         linO=(ligP>ligE)? ligP : ligE;  //--linO=Max(ligP, ligE);
 518         colO=(colP<colE)? colP : colE;  //--colO=Min(colP, colE);
 519       }
 520       else {
 521         colT++;
 522         linO=(ligP<ligE)? ligP : ligE;  //--linO=Min(ligP, ligE);
 523         colO=(colP>colE)? colP : colE;  //--colO=Max(colP, colE);
 524       }
 525       break;
 526     }
 527   }
 528   else {
 529     // Unknown Triangle position :
 530     if (Pedge==0) {
 531       // Unknown edge :
 532       linT=(1>ligP)? 1 : ligP;      //--linT=Max(1, ligP);
 533       colT=(1>(colP+colP))? 1 : (colP+colP);       //--colT=Max(1, colP+colP);
 534       if (ligP==0) linO=ligP+1;
 535       else         linO=ligP-1;
 536       colO=colP;
 537     }
 538     else {
 539       // Known edge We take the left or down connectivity :
 540       switch (typE) {
 541       case 1:  // Horizontal
 542         linT=ligP+1;
 543         colT=(colP>colE)? colP : colE;  //--colT=Max(colP,colE);
 544         colT+=colT;
 545         linO=ligP+1;
 546         colO=(colP>colE)? colP : colE;  //--colO=Max(colP,colE);
 547         break;
 548       case 2:  // Vertical
 549         linT=(ligP>ligE)? ligP : ligE;  //--linT=Max(ligP, ligE);
 550         colT=colP+colP;
 551         linO=(ligP<ligE)? ligP : ligE;  //--linO=Min(ligP, ligE);
 552         colO=colP-1;
 553         break;
 554       case 3:  // Oblique
 555         linT=(ligP>ligE)? ligP : ligE;  //--linT=Max(ligP, ligE);
 556         colT=colP+colE;
 557         linO=(ligP>ligE)? ligP : ligE;  //--linO=Max(ligP, ligE);
 558         colO=(colP<colE)? colP : colE;  //--colO=Min(colP, colE);
 559         break;
 560       }
 561     }
 562   }
 563
 564   TriCon=(linT-1)*nbdeltaVm2 + colT;
 565
 566   if (linT<1) {
 567     linO=0;
 568     colO=colP+colP-colE;
 569     if (colO<0) {colO=0;linO=1;}
 570     else if (colO>nbdeltaV) {colO=nbdeltaV;linO=1;}
 571     TriCon=0;
 572   }
 573   else if (linT>nbdeltaU) {
 574     linO=nbdeltaU;
 575     colO=colP+colP-colE;
 576     if (colO<0) {colO=0;linO=nbdeltaU-1;}
 577     else if (colO>nbdeltaV) {colO=nbdeltaV;linO=nbdeltaU-1;}
 578     TriCon=0;
 579   }
 580
 581   if (colT<1) {
 582     colO=0;
 583     linO=ligP+ligP-ligE;
 584     if (linO<0) {linO=0;colO=1;}
 585     else if (linO>nbdeltaU) {linO=nbdeltaU;colO=1;}
 586     TriCon=0;
 587   }
 588   else if (colT>nbdeltaV) {
 589     colO=nbdeltaV;
 590     linO=ligP+ligP-ligE;
 591     if (linO<0) {linO=0;colO=nbdeltaV-1;}
 592     else if (linO>nbdeltaU) {linO=nbdeltaU;colO=nbdeltaV-1;}
 593     TriCon=0;
 594   }
 595
 596   OtherP=linO*nbdeltaVp1 + colO+1;
 597
 598   return TriCon;
 599 }
 600
 601
 602 //=======================================================================
 603 //function : PlaneEquation
 604 //purpose  :
 605 //=======================================================================
 606 void IntCurveSurface_Polyhedron::PlaneEquation (const Standard_Integer Triang,
 607                                         gp_XYZ&        NormalVector,
 608                                         Standard_Real& PolarDistance)  const
 609 {
 610   Standard_Integer i1,i2,i3;
 611   Triangle(Triang,i1,i2,i3);
 612
 613   //--  gp_XYZ v1=Point(i2).XYZ()-Point(i1).XYZ();
 614   //--  gp_XYZ v2=Point(i3).XYZ()-Point(i2).XYZ();
 615   //--  gp_XYZ v3=Point(i1).XYZ()-Point(i3).XYZ();
 616
 617   gp_XYZ Pointi1(Point(i1).XYZ());
 618   gp_XYZ Pointi2(Point(i2).XYZ());
 619   gp_XYZ Pointi3(Point(i3).XYZ());
 620
 621
 622   gp_XYZ v1= Pointi2 - Pointi1;
 623   gp_XYZ v2= Pointi3 - Pointi2;
 624   gp_XYZ v3= Pointi1 - Pointi3;
 625
 626   if(v1.SquareModulus()<=LONGUEUR_MINI_EDGE_TRIANGLE) { NormalVector.SetCoord(1.0,0.0,0.0); return; }
 627   if(v2.SquareModulus()<=LONGUEUR_MINI_EDGE_TRIANGLE) { NormalVector.SetCoord(1.0,0.0,0.0); return; }
 628   if(v3.SquareModulus()<=LONGUEUR_MINI_EDGE_TRIANGLE) { NormalVector.SetCoord(1.0,0.0,0.0); return; }
 629
 630   NormalVector= (v1^v2)+(v2^v3)+(v3^v1);
 631   NormalVector.Normalize();
 632   PolarDistance = NormalVector * Point(i1).XYZ();
 633 }
 634 //=======================================================================
 635 //function : Contain
 636 //purpose  :
 637 //=======================================================================
 638 Standard_Boolean IntCurveSurface_Polyhedron::Contain (const Standard_Integer Triang,
 639                                                       const gp_Pnt& ThePnt) const
 640 {
 641   Standard_Integer i1,i2,i3;
 642   Triangle(Triang,i1,i2,i3);
 643   gp_XYZ Pointi1(Point(i1).XYZ());
 644   gp_XYZ Pointi2(Point(i2).XYZ());
 645   gp_XYZ Pointi3(Point(i3).XYZ());
 646
 647   gp_XYZ v1=(Pointi2-Pointi1)^(ThePnt.XYZ()-Pointi1);
 648   gp_XYZ v2=(Pointi3-Pointi2)^(ThePnt.XYZ()-Pointi2);
 649   gp_XYZ v3=(Pointi1-Pointi3)^(ThePnt.XYZ()-Pointi3);
 650   if (v1*v2 >= 0. && v2*v3 >= 0. && v3*v1>=0.)
 651     return Standard_True;
 652   else
 653     return Standard_False;
 654 }
 655 //=======================================================================
 656 //function : Dump
 657 //purpose  :
 658 //=======================================================================
 659 void IntCurveSurface_Polyhedron::Dump() const
 660 {
 661
 662 }
 663 //=======================================================================
 664 //function : Size
 665 //purpose  :
 666 //=======================================================================
 667 void IntCurveSurface_Polyhedron::Size (Standard_Integer& nbdu,
 668                               Standard_Integer& nbdv) const
 669 {
 670   nbdu=nbdeltaU;
 671   nbdv=nbdeltaV;
 672 }
 673 //=======================================================================
 674 //function : Triangle
 675 //purpose  :
 676 //=======================================================================
 677 void IntCurveSurface_Polyhedron::Triangle (const Standard_Integer Index,
 678                                   Standard_Integer& P1,
 679                                   Standard_Integer& P2,
 680                                   Standard_Integer& P3) const
 681 {
 682   Standard_Integer line=1+((Index-1)/(nbdeltaV*2));
 683   Standard_Integer colon=1+((Index-1)%(nbdeltaV*2));
 684   Standard_Integer colpnt=(colon+1)/2;
 685
 686 // General formula = (line-1)*(nbdeltaV+1)+colpnt
 687
 688 //  Position of P1 = MesXYZ(line,colpnt);
 689   P1= (line-1)*(nbdeltaV+1) + colpnt;
 690
 691 //  Position of P2= MesXYZ(line+1,colpnt+((colon-1)%2));
 692   P2= line*(nbdeltaV+1) + colpnt+((colon-1)%2);
 693
 694 //  Position of P3= MesXYZ(line+(colon%2),colpnt+1);
 695   P3= (line-1+(colon%2))*(nbdeltaV+1) + colpnt + 1;
 696 }
 697 //=======================================================================
 698 //function : Point
 699 //=======================================================================
 700 const gp_Pnt& IntCurveSurface_Polyhedron::Point(const Standard_Integer Index
 701                                                 ,Standard_Real& U
 702                                                 ,Standard_Real& V) const
 703 {
 704 #if CHECKBOUNDS
 705   if(Index<0 || Index>((nbdeltaU+1)*(nbdeltaV+1))) {
 706     printf("\n Erreur IntCurveSurface_Polyhedron::Parameters\n");
 707   }
 708 #endif
 709   gp_Pnt *CMyPnts     = (gp_Pnt *)C_MyPnts;
 710   Standard_Real *CMyU = (Standard_Real *)C_MyU;
 711   Standard_Real *CMyV = (Standard_Real *)C_MyV;
 712   U=CMyU[Index];
 713   V=CMyV[Index];
 714   return CMyPnts[Index];
 715 }
 716 //=======================================================================
 717 //function : Point
 718 //=======================================================================
 719 const gp_Pnt& IntCurveSurface_Polyhedron::Point(const Standard_Integer Index) const {
 720 #if CHECKBOUNDS
 721   if(Index<0 || Index>((nbdeltaU+1)*(nbdeltaV+1))) {
 722     printf("\n Erreur IntCurveSurface_Polyhedron::Parameters\n");
 723   }
 724 #endif
 725
 726   gp_Pnt *CMyPnts     = (gp_Pnt *)C_MyPnts;
 727   return CMyPnts[Index];
 728 }
 729 //=======================================================================
 730 //function : Point
 731 //=======================================================================
 732 //void IntCurveSurface_Polyhedron::Point (const gp_Pnt& p,
 733 //                                      const Standard_Integer lig,
 734 //                                      const Standard_Integer col,
 735 //                                      const Standard_Real u,
 736 //                                      const Standard_Real v)
 737 void IntCurveSurface_Polyhedron::Point (const gp_Pnt& ,
 738                                         const Standard_Integer ,
 739                                         const Standard_Integer ,
 740                                         const Standard_Real ,
 741                                         const Standard_Real )
 742 {
 743   printf("\n IntCurveSurface_Polyhedron::Point : Ne dois pas etre appelle\n");
 744 }
 745 //=======================================================================
 746 //function : Point
 747 //=======================================================================
 748 void IntCurveSurface_Polyhedron::Point (const Standard_Integer Index,gp_Pnt& P) const
 749 {
 750 #if CHECKBOUNDS
 751   if(Index<0 || Index>((nbdeltaU+1)*(nbdeltaV+1))) {
 752     printf("\n Erreur IntCurveSurface_Polyhedron::Parameters\n");
 753   }
 754 #endif
 755
 756   gp_Pnt *CMyPnts     = (gp_Pnt *)C_MyPnts;
 757   P = CMyPnts[Index];
 758 }
 759
 760 //  Modified by Sergey KHROMOV - Fri Dec  7 10:12:47 2001 Begin
 761
 762 //=======================================================================
 763 //function : IsOnBound
 764 //purpose  : This method returns true if the edge based on points with
 765 //           indices Index1 and Index2 represents a boundary edge.
 766 //=======================================================================
 767
 768 Standard_Boolean IntCurveSurface_Polyhedron::IsOnBound
 769                                  (const Standard_Integer Index1,
 770                                   const Standard_Integer Index2) const
 771 {
 772 #if CHECKBOUNDS
 773   if(Index1<0 || Index1>((nbdeltaU+1)*(nbdeltaV+1))) {
 774     printf("\n Erreur IntCurveSurface_Polyhedron::IsOnBound\n");
 775   }
 776   if(Index2<0 || Index2>((nbdeltaU+1)*(nbdeltaV+1))) {
 777     printf("\n Erreur IntCurveSurface_Polyhedron::IsOnBound\n");
 778   }
 779 #endif
 780   Standard_Boolean *CMyIsOnBounds = (Standard_Boolean *)C_MyIsOnBounds;
 781   Standard_Integer  aDiff         = Abs(Index1 - Index2);
 782   Standard_Integer  i;
 783
 784 // Check if points are neighbour ones.
 785   if (aDiff != 1 && aDiff != nbdeltaV + 1)
 786     return Standard_False;
 787
 788   for (i = 0; i <= nbdeltaU; i++) {
 789     if ((Index1 == 1 + i*(nbdeltaV + 1)) && (Index2 == Index1 - 1))
 790       return Standard_False;
 791
 792     if ((Index1 == (1 + i)*(nbdeltaV + 1)) && (Index2 == Index1 + 1))
 793       return Standard_False;
 794   }
 795
 796   return (CMyIsOnBounds[Index1] && CMyIsOnBounds[Index2]);
 797 }
 798
 799 //=======================================================================
 800 //function : ComputeBorderDeflection
 801 //purpose  : This method computes and returns a deflection of isoline
 802 //           of given parameter on Surface.
 803 //=======================================================================
 804
 805 Standard_Real IntCurveSurface_Polyhedron::ComputeBorderDeflection
 806                               (const ThePSurface      &Surface,
 807                                const Standard_Real     Parameter,
 808                                const Standard_Real     PMin,
 809                                const Standard_Real     PMax,
 810                                const Standard_Boolean  isUIso) const
 811 {
 812   Standard_Integer aNbSamples;
 813   Standard_Integer i;
 814
 815   if (isUIso)
 816     aNbSamples = nbdeltaV;
 817   else
 818     aNbSamples = nbdeltaU;
 819
 820   Standard_Real aDelta      = (PMax - PMin)/aNbSamples;
 821   Standard_Real aPar        = PMin;
 822   Standard_Real aDeflection = RealFirst();
 823   gp_XYZ        aP1;
 824   gp_XYZ        aP2;
 825   gp_XYZ        aPMid;
 826   gp_XYZ        aPParMid;
 827
 828   for (i = 0; i <= aNbSamples; i++, aPar+= aDelta) {
 829     if (isUIso) {
 830       aP1      = ThePSurfaceTool::Value(Surface, Parameter, aPar).XYZ();
 831       aP2      = ThePSurfaceTool::Value(Surface, Parameter, aPar + aDelta).XYZ();
 832       aPParMid = ThePSurfaceTool::Value(Surface, Parameter, aPar + aDelta/2.).XYZ();
 833     } else {
 834       aP1      = ThePSurfaceTool::Value(Surface, aPar,             Parameter).XYZ();
 835       aP2      = ThePSurfaceTool::Value(Surface, aPar + aDelta,    Parameter).XYZ();
 836       aPParMid = ThePSurfaceTool::Value(Surface, aPar + aDelta/2., Parameter).XYZ();
 837     }
 838     aPMid = (aP2 + aP1)/2.;
 839
 840     Standard_Real aDist = (aPMid - aPParMid).Modulus();
 841
 842     if (aDist > aDeflection)
 843       aDeflection = aDist;
 844   }
 845
 846   return aDeflection;
 847 }
 848
 849 //  Modified by Sergey KHROMOV - Fri Dec  7 11:21:52 2001 End