src/AdvancedEngine/GEOMImpl_IAdvancedOperations.cxx

   1 // Copyright (C) 2007-2014  CEA/DEN, EDF R&D, OPEN CASCADE
   2 //
   3 // This library is free software; you can redistribute it and/or
   4 // modify it under the terms of the GNU Lesser General Public
   5 // License as published by the Free Software Foundation; either
   6 // version 2.1 of the License, or (at your option) any later version.
   7 //
   8 // This library is distributed in the hope that it will be useful,
   9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11 // Lesser General Public License for more details.
  12 //
  13 // You should have received a copy of the GNU Lesser General Public
  14 // License along with this library; if not, write to the Free Software
  15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19 //  File   : GEOMImpl_IAdvancedOperations.cxx
  20 //  Author : Vadim SANDLER, Open CASCADE S.A.S. (vadim.sandler@opencascade.com)
  21
  22 #include "GEOMImpl_IAdvancedOperations.hxx"
  23
  24 #include <Basics_OCCTVersion.hxx>
  25
  26 #include <utilities.h>
  27 #include <OpUtil.hxx>
  28 #include <Utils_ExceptHandlers.hxx>
  29
  30 #include "GEOM_Function.hxx"
  31 #include "GEOM_PythonDump.hxx"
  32 #include "GEOMUtils.hxx"
  33 #include "GEOMAlgo_Splitter.hxx"
  34 #include "GEOMAlgo_FinderShapeOn1.hxx"
  35
  36 #include "GEOMImpl_Gen.hxx"
  37 #include "GEOMImpl_Types.hxx"
  38
  39 #include "GEOMImpl_IBasicOperations.hxx"
  40 #include "GEOMImpl_IBooleanOperations.hxx"
  41 #include "GEOMImpl_IShapesOperations.hxx"
  42 #include "GEOMImpl_ITransformOperations.hxx"
  43 #include "GEOMImpl_IBlocksOperations.hxx"
  44 #include "GEOMImpl_I3DPrimOperations.hxx"
  45 #include "GEOMImpl_ILocalOperations.hxx"
  46 #include "GEOMImpl_IHealingOperations.hxx"
  47 #include "GEOMImpl_IGroupOperations.hxx"
  48
  49 #include "GEOMImpl_GlueDriver.hxx"
  50 #include "GEOMImpl_PipeTShapeDriver.hxx"
  51 #include "GEOMImpl_IPipeTShape.hxx"
  52 #include "GEOMImpl_DividedDiskDriver.hxx"
  53 #include "GEOMImpl_IDividedDisk.hxx"
  54 #include <GEOMImpl_SmoothingSurfaceDriver.hxx>
  55 #include <GEOMImpl_ISmoothingSurface.hxx>
  56 /*@@ insert new functions before this line @@ do not remove this line @@ do not remove this line @@*/
  57
  58 #include <TDF_Tool.hxx>
  59 #include <TFunction_DriverTable.hxx>
  60 #include <TFunction_Driver.hxx>
  61 #include <TFunction_Logbook.hxx>
  62 #include <TNaming_CopyShape.hxx>
  63
  64 #include <TopExp.hxx>
  65 #include <TopExp_Explorer.hxx>
  66 #include <TopoDS.hxx>
  67 #include <TopoDS_Vertex.hxx>
  68 #include <TopTools_IndexedMapOfShape.hxx>
  69 #include <TopTools_ListIteratorOfListOfShape.hxx>
  70 #include <TColStd_IndexedDataMapOfTransientTransient.hxx>
  71
  72 #include <BRep_Builder.hxx>
  73 #include <BRep_Tool.hxx>
  74
  75 #include <BRepAdaptor_Surface.hxx>
  76 #include <BRepAlgoAPI_Cut.hxx>
  77 #include <BRepAlgoAPI_Fuse.hxx>
  78 #include <BRepBuilderAPI_MakeFace.hxx>
  79 #include <BRepBuilderAPI_MakeVertex.hxx>
  80 #include <BRepBuilderAPI_Transform.hxx>
  81 #include <BRepPrimAPI_MakeCone.hxx>
  82 #include <BRepPrimAPI_MakeCylinder.hxx>
  83
  84 #include <gp_Ax3.hxx>
  85 #include <gp_Pln.hxx>
  86 #include <gp_Pnt.hxx>
  87 #include <gp_Vec.hxx>
  88 #include <GC_MakeConicalSurface.hxx>
  89 #include <Geom_CylindricalSurface.hxx>
  90
  91 #include <ShapeAnalysis_Edge.hxx>
  92
  93 #include <cmath>
  94
  95 #include "AdvancedEngine_Types.hxx"
  96
  97 #include <Standard_Stream.hxx>
  98 #include <Standard_Failure.hxx>
  99 #include <StdFail_NotDone.hxx>
 100 #include <Standard_ErrorHandler.hxx> // CAREFUL ! position of this file is critic : see Lucien PIGNOLONI / OCC
 101
 102 #define HALF_LENGTH_MAIN_PIPE     "Main pipe half length" //"Tuyau principal - demi longueur"
 103 #define HALF_LENGTH_INCIDENT_PIPE "Incident pipe half length" //"Tuyau incident - demi longueur"
 104 #define CIRCULAR_QUARTER_PIPE     "Circular quarter of pipe" //"Circulaire - quart de tuyau"
 105 #define THICKNESS                 "Thickness" //"Epaisseur"
 106 #define FLANGE                    "Flange" // "Collerette"
 107 #define CHAMFER_OR_FILLET         "Chamfer or fillet" //"Chanfrein ou Raccord"
 108 #define JUNCTION_FACE_1           "Junction 1" //"Face de jonction 1"
 109 #define JUNCTION_FACE_2           "Junction 2" //"Face de jonction 2"
 110 #define JUNCTION_FACE_3           "Junction 3" //"Face de jonction 3"
 111
 112 #define FIND_GROUPS_BY_POINTS 1
 113
 114 //=============================================================================
 115 /*!
 116  *  Constructor
 117  */
 118 //=============================================================================
 119 GEOMImpl_IAdvancedOperations::GEOMImpl_IAdvancedOperations(GEOM_Engine* theEngine, int theDocID) :
 120   GEOM_IOperations(theEngine, theDocID)
 121 {
 122   MESSAGE("GEOMImpl_IAdvancedOperations::GEOMImpl_IAdvancedOperations");
 123   myBasicOperations     = new GEOMImpl_IBasicOperations(GetEngine(), GetDocID());
 124   myBooleanOperations   = new GEOMImpl_IBooleanOperations(GetEngine(), GetDocID());
 125   myShapesOperations    = new GEOMImpl_IShapesOperations(GetEngine(), GetDocID());
 126   myTransformOperations = new GEOMImpl_ITransformOperations(GetEngine(), GetDocID());
 127   myBlocksOperations    = new GEOMImpl_IBlocksOperations(GetEngine(), GetDocID());
 128   my3DPrimOperations    = new GEOMImpl_I3DPrimOperations(GetEngine(), GetDocID());
 129   myLocalOperations     = new GEOMImpl_ILocalOperations(GetEngine(), GetDocID());
 130   myHealingOperations   = new GEOMImpl_IHealingOperations(GetEngine(), GetDocID());
 131   myGroupOperations     = new GEOMImpl_IGroupOperations(GetEngine(), GetDocID());
 132 }
 133
 134 //=============================================================================
 135 /*!
 136  *  Destructor
 137  */
 138 //=============================================================================
 139 GEOMImpl_IAdvancedOperations::~GEOMImpl_IAdvancedOperations()
 140 {
 141   MESSAGE("GEOMImpl_IAdvancedOperations::~GEOMImpl_IAdvancedOperations");
 142   delete myBasicOperations;
 143   delete myBooleanOperations;
 144   delete myShapesOperations;
 145   delete myTransformOperations;
 146   delete myBlocksOperations;
 147   delete my3DPrimOperations;
 148   delete myLocalOperations;
 149   delete myHealingOperations;
 150   delete myGroupOperations;
 151 }
 152
 153 //=============================================================================
 154 /*!
 155  *  SetPosition
 156  */
 157 //=============================================================================
 158 gp_Trsf GEOMImpl_IAdvancedOperations::GetPositionTrsf(double theL1, double theL2,
 159                                                       Handle(GEOM_Object) theP1,
 160                                                       Handle(GEOM_Object) theP2,
 161                                                       Handle(GEOM_Object) theP3)
 162 {
 163   // Old Local Coordinates System oldLCS
 164   gp_Pnt P0(0, 0, 0);
 165   gp_Pnt P1(-theL1, 0, 0);
 166   gp_Pnt P2(theL1, 0, 0);
 167   gp_Pnt P3(0, 0, theL2);
 168
 169   gp_Dir oldX(gp_Vec(P1, P2));
 170   gp_Dir oldZ(gp_Vec(P0, P3));
 171   gp_Ax3 oldLCS(P0, oldZ, oldX);
 172
 173   // New Local Coordinates System newLCS
 174   double LocX, LocY, LocZ;
 175   gp_Pnt newP1 = BRep_Tool::Pnt(TopoDS::Vertex(theP1->GetValue()));
 176   gp_Pnt newP2 = BRep_Tool::Pnt(TopoDS::Vertex(theP2->GetValue()));
 177   gp_Pnt newP3 = BRep_Tool::Pnt(TopoDS::Vertex(theP3->GetValue()));
 178   LocX = (newP1.X() + newP2.X()) / 2.;
 179   LocY = (newP1.Y() + newP2.Y()) / 2.;
 180   LocZ = (newP1.Z() + newP2.Z()) / 2.;
 181   gp_Pnt newO(LocX, LocY, LocZ);
 182
 183   gp_Dir newX(gp_Vec(newP1, newP2)); // P1P2 Vector
 184   gp_Dir newZ(gp_Vec(newO, newP3)); // OP3 Vector
 185   gp_Ax3 newLCS = gp_Ax3(newO, newZ, newX);
 186
 187   gp_Trsf aTrsf;
 188   aTrsf.SetDisplacement(oldLCS, newLCS);
 189
 190   return aTrsf;
 191 }
 192
 193 //=============================================================================
 194 /*!
 195  *  CheckCompatiblePosition
 196  *
 197  */
 198 //=============================================================================
 199 bool GEOMImpl_IAdvancedOperations::CheckCompatiblePosition(double& theL1, double& theL2,
 200                                                            Handle(GEOM_Object) theP1,
 201                                                            Handle(GEOM_Object) theP2,
 202                                                            Handle(GEOM_Object) theP3,
 203                                                            double theTolerance)
 204 {
 205   SetErrorCode(KO);
 206   gp_Pnt P1 = BRep_Tool::Pnt(TopoDS::Vertex(theP1->GetValue()));
 207   gp_Pnt P2 = BRep_Tool::Pnt(TopoDS::Vertex(theP2->GetValue()));
 208   gp_Pnt P3 = BRep_Tool::Pnt(TopoDS::Vertex(theP3->GetValue()));
 209
 210   double d12 = P1.Distance(P2);
 211   double d13 = P1.Distance(P3);
 212   double d23 = P2.Distance(P3);
 213   //    double d2 = newO.Distance(P3);
 214
 215   if (Abs(d12) <= Precision::Confusion()) {
 216     SetErrorCode("Junctions points P1 and P2 are identical");
 217     return false;
 218   }
 219   if (Abs(d13) <= Precision::Confusion()) {
 220     SetErrorCode("Junctions points P1 and P3 are identical");
 221     return false;
 222   }
 223   if (Abs(d23) <= Precision::Confusion()) {
 224     SetErrorCode("Junctions points P2 and P3 are identical");
 225     return false;
 226   }
 227
 228
 229   double newL1 = 0.5 * d12;
 230   double newL2 = sqrt(pow(d13,2)-pow(newL1,2));
 231   //
 232   // theL1*(1-theTolerance) <= newL1 <= theL1*(1+theTolerance)
 233   //
 234   if (fabs(newL1 - theL1) > Precision::Approximation()) {
 235     if ( (newL1 * (1 - theTolerance) -theL1 <= Precision::Approximation()) &&
 236          (newL1 * (1 + theTolerance) -theL1 >= Precision::Approximation()) ) {
 237       //            std::cerr << "theL1 = newL1" << std::endl;
 238       theL1 = newL1;
 239     } else {
 240       theL1 = -1;
 241       SetErrorCode("Dimension for main pipe (L1) is incompatible with new position");
 242       return false;
 243     }
 244   }
 245
 246   //
 247   // theL2*(1-theTolerance) <= newL2  <= theL2*(1+theTolerance)
 248   //
 249   if (fabs(newL2 - theL2) > Precision::Approximation()) {
 250     if ( (newL2 * (1 - theTolerance) -theL2 <= Precision::Approximation()) &&
 251          (newL2 * (1 + theTolerance) -theL2 >= Precision::Approximation()) ) {
 252       theL2 = newL2;
 253     } else {
 254       theL2 = -1;
 255       SetErrorCode("Dimension for incident pipe (L2) is incompatible with new position");
 256       return false;
 257     }
 258   }
 259
 260   SetErrorCode(OK);
 261   return true;
 262
 263 }
 264
 265 //=============================================================================
 266 /*!
 267  *  Generate the propagation groups of a Pipe T-Shape used for hexa mesh
 268  */
 269 //=============================================================================
 270 bool GEOMImpl_IAdvancedOperations::MakeGroups(Handle(GEOM_Object) theShape, int shapeType,
 271                                               double theR1, double theW1, double theL1,
 272                                               double theR2, double theW2, double theL2,
 273                                               double theH, double theW, double theRF,
 274                                               Handle(TColStd_HSequenceOfTransient) theSeq,
 275                                               gp_Trsf aTrsf)
 276 {
 277   SetErrorCode(KO);
 278
 279   if (theShape.IsNull()) return false;
 280
 281   TopoDS_Shape aShape = theShape->GetValue();
 282   if (aShape.IsNull()) {
 283     SetErrorCode("Shape is not defined");
 284     return false;
 285   }
 286
 287   gp_Trsf aTrsfInv = aTrsf.Inverted();
 288
 289 //   int expectedGroups = 0;
 290 //   if (shapeType == TSHAPE_BASIC)
 291 //     if (Abs(theR2+theW2-theR1-theW1) <= Precision::Approximation())
 292 //       expectedGroups = 10;
 293 //     else
 294 //       expectedGroups = 11;
 295 //   else if (shapeType == TSHAPE_CHAMFER || shapeType == TSHAPE_FILLET)
 296 //     expectedGroups = 12;
 297
 298   double aR1Ext = theR1 + theW1;
 299   double aR2Ext = theR2 + theW2;
 300
 301   /////////////////////////
 302   //// Groups of Faces ////
 303   /////////////////////////
 304
 305   //
 306   // Comment the following lines when GetInPlace bug is solved
 307   // == BEGIN
 308   // Workaround of GetInPlace bug
 309   // Create a bounding box that fits the shape
 310   Handle(GEOM_Object) aBox = my3DPrimOperations->MakeBoxDXDYDZ(2*theL1, 2*aR1Ext, aR1Ext+theL2);
 311   aBox->GetLastFunction()->SetDescription("");
 312   myTransformOperations->TranslateDXDYDZ(aBox, -theL1, -aR1Ext, -aR1Ext);
 313   aBox->GetLastFunction()->SetDescription("");
 314   // Apply transformation to box
 315   BRepBuilderAPI_Transform aTransformationBox(aBox->GetValue(), aTrsf, Standard_False);
 316   TopoDS_Shape aBoxShapeTrsf = aTransformationBox.Shape();
 317   aBox->GetLastFunction()->SetValue(aBoxShapeTrsf);
 318
 319   // Get the shell of the box
 320   Handle(GEOM_Object) aShell = Handle(GEOM_Object)::DownCast
 321     (myShapesOperations->MakeExplode(aBox, TopAbs_SHELL, true)->Value(1));
 322   aBox->GetLastFunction()->SetDescription("");
 323   aShell->GetLastFunction()->SetDescription("");
 324   // Get the common shapes between shell and shape
 325   Handle(GEOM_Object) aCommonCompound = myBooleanOperations->MakeBoolean
 326                             (theShape, aShell, 1, Standard_False); // MakeCommon
 327   if (aCommonCompound.IsNull()) {
 328     SetErrorCode(myBooleanOperations->GetErrorCode());
 329     return false;
 330   }
 331   aCommonCompound->GetLastFunction()->SetDescription("");
 332   // Explode the faces of common shapes => 3 faces
 333   Handle(TColStd_HSequenceOfTransient) aCommonFaces =
 334     myShapesOperations->MakeExplode(aCommonCompound, TopAbs_FACE, true);
 335   aCommonCompound->GetLastFunction()->SetDescription("");
 336   std::list<Handle(GEOM_Object)> aCompoundOfFacesList;
 337
 338   for (int i=0 ; i<= aCommonFaces->Length()-4 ; i+=4) {
 339     std::list<Handle(GEOM_Object)> aFacesList;
 340     for (int j = 1 ; j <= 4 ; j++) {
 341       Handle(GEOM_Object) aFace = Handle(GEOM_Object)::DownCast(aCommonFaces->Value(i+j)); // Junction faces
 342       if (!aFace.IsNull()) {
 343         aFace->GetLastFunction()->SetDescription("");
 344         aFacesList.push_back(aFace);
 345       }
 346     }
 347     Handle(GEOM_Object) aCompoundOfFaces = myShapesOperations->MakeCompound(aFacesList);
 348     if (!aCompoundOfFaces.IsNull()) {
 349       aCompoundOfFaces->GetLastFunction()->SetDescription("");
 350       aCompoundOfFacesList.push_back(aCompoundOfFaces);
 351     }
 352   }
 353
 354   if (aCompoundOfFacesList.size() == 3) {
 355     Handle(GEOM_Object) aPln1 = aCompoundOfFacesList.front();
 356     aCompoundOfFacesList.pop_front();
 357     Handle(GEOM_Object) aPln2 = aCompoundOfFacesList.front();
 358     aCompoundOfFacesList.pop_front();
 359     Handle(GEOM_Object) aPln3 = aCompoundOfFacesList.front();
 360     aCompoundOfFacesList.pop_front();
 361     // == END
 362     //
 363
 364
 365     //     Uncomment the following lines when GetInPlace bug is solved
 366     //     == BEGIN
 367 //     Handle(GEOM_Object) aP1 = myBasicOperations->MakePointXYZ(-theL1, 0, 0);
 368 //     Handle(GEOM_Object) aP2 = myBasicOperations->MakePointXYZ(-0, 0, theL2);
 369 //     Handle(GEOM_Object) aP3 = myBasicOperations->MakePointXYZ(theL1, 0, 0);
 370 //     aP1->GetLastFunction()->SetDescription("");
 371 //     aP2->GetLastFunction()->SetDescription("");
 372 //     aP3->GetLastFunction()->SetDescription("");
 373 //     Handle(GEOM_Object) aV1 = myBasicOperations->MakeVectorDXDYDZ(-1, 0, 0);
 374 //     Handle(GEOM_Object) aV2 = myBasicOperations->MakeVectorDXDYDZ(0, 0, 1);
 375 //     Handle(GEOM_Object) aV3 = myBasicOperations->MakeVectorDXDYDZ(1, 0, 0);
 376 //     aV1->GetLastFunction()->SetDescription("");
 377 //     aV2->GetLastFunction()->SetDescription("");
 378 //     aV3->GetLastFunction()->SetDescription("");
 379 //     Handle(GEOM_Object) aPln1 = myBasicOperations->MakePlanePntVec(aP1, aV1, 2*(aR1Ext+theL2));
 380 //     Handle(GEOM_Object) aPln2 = myBasicOperations->MakePlanePntVec(aP2, aV2, 2*(aR2Ext));
 381 //     Handle(GEOM_Object) aPln3 = myBasicOperations->MakePlanePntVec(aP3, aV3, 2*(aR1Ext+theL2));
 382 //     aPln1->GetLastFunction()->SetDescription("");
 383 //     aPln2->GetLastFunction()->SetDescription("");
 384 //     aPln3->GetLastFunction()->SetDescription("");
 385 //
 386 //     BRepBuilderAPI_Transform aTransformation1(aPln1->GetValue(), aTrsf, Standard_False);
 387 //     TopoDS_Shape aTrsf_Shape1 = aTransformation1.Shape();
 388 //     aPln1->GetLastFunction()->SetValue(aTrsf_Shape1);
 389 //     BRepBuilderAPI_Transform aTransformation2(aPln2->GetValue(), aTrsf, Standard_False);
 390 //     TopoDS_Shape aTrsf_Shape2 = aTransformation2.Shape();
 391 //     aPln2->GetLastFunction()->SetValue(aTrsf_Shape2);
 392 //     BRepBuilderAPI_Transform aTransformation3(aPln3->GetValue(), aTrsf, Standard_False);
 393 //     TopoDS_Shape aTrsf_Shape3 = aTransformation3.Shape();
 394 //     aPln3->GetLastFunction()->SetValue(aTrsf_Shape3);
 395     //     == END
 396     //
 397
 398     Handle(GEOM_Object) junctionFaces1 = myShapesOperations->GetInPlace(theShape, aPln1);
 399     if (junctionFaces1.IsNull())
 400       junctionFaces1 = myShapesOperations->GetShapesOnShapeAsCompound
 401         (aPln1, theShape, TopAbs_FACE,  GEOMAlgo_ST_ONIN);
 402     if (!junctionFaces1.IsNull()) {
 403       junctionFaces1->GetLastFunction()->SetDescription("");
 404       junctionFaces1->SetName("JUNCTION_FACE_1");
 405       theSeq->Append(junctionFaces1);
 406     }
 407     else {
 408       SetErrorCode("Junction face 1 not found");
 409       //        theSeq->Append(aPln1);
 410       //        return false;
 411     }
 412     Handle(GEOM_Object) junctionFaces2 = myShapesOperations->GetInPlace(theShape, aPln2);
 413     if (junctionFaces2.IsNull())
 414       junctionFaces2 = myShapesOperations->GetShapesOnShapeAsCompound
 415         (aPln2, theShape, TopAbs_FACE,  GEOMAlgo_ST_ONIN);
 416     if (!junctionFaces2.IsNull()) {
 417       junctionFaces2->GetLastFunction()->SetDescription("");
 418       junctionFaces2->SetName("JUNCTION_FACE_2");
 419       theSeq->Append(junctionFaces2);
 420     }
 421     else {
 422       SetErrorCode("Junction face 2 not found");
 423       //        theSeq->Append(aPln2);
 424       //        return false;
 425     }
 426     Handle(GEOM_Object) junctionFaces3 = myShapesOperations->GetInPlace(theShape, aPln3);
 427     if (junctionFaces3.IsNull())
 428       junctionFaces3 = myShapesOperations->GetShapesOnShapeAsCompound
 429         (aPln3, theShape, TopAbs_FACE,  GEOMAlgo_ST_ONIN);
 430     if (!junctionFaces3.IsNull()) {
 431       junctionFaces3->GetLastFunction()->SetDescription("");
 432       junctionFaces3->SetName("JUNCTION_FACE_3");
 433       theSeq->Append(junctionFaces3);
 434     }
 435     else {
 436       SetErrorCode("Junction face 3 not found");
 437       //        theSeq->Append(aPln3);
 438       //        return false;
 439     }
 440   // Comment the following lines when GetInPlace bug is solved
 441   // == BEGIN
 442   }
 443   //     == END
 444
 445   /////////////////////////
 446   //// Groups of Edges ////
 447   /////////////////////////
 448   // Result of propagate
 449
 450   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
 451
 452   TCollection_AsciiString theDesc = aFunction->GetDescription();
 453   Handle(TColStd_HSequenceOfTransient) aSeqPropagate = myBlocksOperations->Propagate(theShape);
 454   if (aSeqPropagate.IsNull() || aSeqPropagate->Length() == 0) {
 455     SetErrorCode("Propagation groups not found");
 456     return false;
 457   }
 458   Standard_Integer aNbGroups = aSeqPropagate->Length();
 459   // Recover previous description to get rid of Propagate dump
 460   aFunction->SetDescription(theDesc);
 461
 462 #ifdef FIND_GROUPS_BY_POINTS
 463   // BEGIN: new groups search
 464
 465   //              W2  R2
 466   //            .----.-----.----.
 467   //           e|    |  |  |    |
 468   //            |    |  |  |    |
 469   //            .    |  |  |    .
 470   //         g / ''..|  |  |..'' \
 471   //       f  /      '''''''      \
 472   //  .---.--'..     |  |  |     ..'--.---.
 473   //  |a    \   '''...........'''   /     |
 474   //  |-------\------'  |  '------/-------.
 475   //  |         \       |       /         |
 476   // c|           \     |     /           |
 477   //  |    R1       \   |   /             |
 478   //  |               \ | /               |
 479   //  ._________________|_________________.
 480   //  |       L1        |                 |
 481   //  |                 |                 |
 482   //  |                 |                 |
 483   // b|                 |                 |
 484   //  |                 |                 |
 485   //  |-----------------|-----------------|
 486   //  |    W1           |                 |
 487   //  '-----------------'-----------------'
 488   //          d
 489
 490   // "Thickness" group (a)
 491   gp_Pnt aPntA (-theL1, 0, theR1 + theW1/2.);
 492   aPntA.Transform(aTrsf);
 493   BRepBuilderAPI_MakeVertex mkVertexA (aPntA);
 494   TopoDS_Vertex aVertA = TopoDS::Vertex(mkVertexA.Shape());
 495   TopoDS_Shape anEdgeA = GEOMUtils::GetEdgeNearPoint(aShape, aVertA);
 496
 497   // "Circular quarter of pipe" group (b)
 498   gp_Pnt aPntB (-theL1, -aR1Ext * sin(M_PI/4.), -aR1Ext * sin(M_PI/4.));
 499   aPntB.Transform(aTrsf);
 500   BRepBuilderAPI_MakeVertex mkVertexB (aPntB);
 501   TopoDS_Vertex aVertB = TopoDS::Vertex(mkVertexB.Shape());
 502   TopoDS_Shape anEdgeB = GEOMUtils::GetEdgeNearPoint(aShape, aVertB);
 503
 504   // "Circular quarter of pipe" group (c)
 505   gp_Pnt aPntC (-theL1, -aR1Ext * sin(M_PI/4.), aR1Ext * sin(M_PI/4.));
 506   aPntC.Transform(aTrsf);
 507   BRepBuilderAPI_MakeVertex mkVertexC (aPntC);
 508   TopoDS_Vertex aVertC = TopoDS::Vertex(mkVertexC.Shape());
 509   TopoDS_Shape anEdgeC = GEOMUtils::GetEdgeNearPoint(aShape, aVertC);
 510
 511   // "Main pipe half length" group (d)
 512   gp_Pnt aPntD (-theL1/2., 0, -aR1Ext);
 513   aPntD.Transform(aTrsf);
 514   BRepBuilderAPI_MakeVertex mkVertexD (aPntD);
 515   TopoDS_Vertex aVertD = TopoDS::Vertex(mkVertexD.Shape());
 516   TopoDS_Shape anEdgeD = GEOMUtils::GetEdgeNearPoint(aShape, aVertD);
 517
 518   // "Incident pipe half length" group (e)
 519   double aTol10 = Precision::Confusion() * 10.;
 520   gp_Pnt aPntE (-aR2Ext, 0, theL2 - aTol10);
 521   aPntE.Transform(aTrsf);
 522   BRepBuilderAPI_MakeVertex mkVertexE (aPntE);
 523   TopoDS_Vertex aVertE = TopoDS::Vertex(mkVertexE.Shape());
 524   TopoDS_Shape anEdgeE = GEOMUtils::GetEdgeNearPoint(aShape, aVertE);
 525
 526   // "Flange" group (f)
 527   double aFx = - aR2Ext - aTol10;
 528   if (shapeType == TSHAPE_CHAMFER)
 529     aFx -= theW;
 530   else if (shapeType == TSHAPE_FILLET)
 531     aFx -= theRF;
 532   gp_Pnt aPntF (aFx, 0, aR1Ext);
 533   aPntF.Transform(aTrsf);
 534   BRepBuilderAPI_MakeVertex mkVertexF (aPntF);
 535   TopoDS_Vertex aVertF = TopoDS::Vertex(mkVertexF.Shape());
 536   TopoDS_Shape anEdgeF = GEOMUtils::GetEdgeNearPoint(aShape, aVertF);
 537
 538   // "Chamfer or Fillet" group (g)
 539   TopoDS_Shape anEdgeG;
 540   if (shapeType == TSHAPE_CHAMFER) {
 541     gp_Pnt aPntG (-aR2Ext - theW/2., 0, aR1Ext + theH/2.);
 542     aPntG.Transform(aTrsf);
 543     BRepBuilderAPI_MakeVertex mkVertexG (aPntG);
 544     TopoDS_Vertex aVertG = TopoDS::Vertex(mkVertexG.Shape());
 545     anEdgeG = GEOMUtils::GetEdgeNearPoint(aShape, aVertG);
 546   }
 547   else if (shapeType == TSHAPE_FILLET) {
 548     gp_Pnt aPntG (-aR2Ext - theRF/2., 0, aR1Ext + theRF/2.);
 549     aPntG.Transform(aTrsf);
 550     BRepBuilderAPI_MakeVertex mkVertexG (aPntG);
 551     TopoDS_Vertex aVertG = TopoDS::Vertex(mkVertexG.Shape());
 552     anEdgeG = GEOMUtils::GetEdgeNearPoint(aShape, aVertG);
 553   }
 554
 555   for (int i = 1 ; i <= aNbGroups; i++) {
 556     Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(aSeqPropagate->Value(i));
 557     if (aGroup.IsNull())
 558       continue;
 559
 560     TopoDS_Shape aGroupShape = aGroup->GetValue();
 561     TopTools_IndexedMapOfShape anEdgesMap;
 562     TopExp::MapShapes(aGroupShape, TopAbs_EDGE, anEdgesMap);
 563
 564     if (anEdgesMap.Contains(anEdgeA)) { // a
 565       aGroup->SetName("THICKNESS");
 566       theSeq->Append(aGroup);
 567     }
 568     else if (anEdgesMap.Contains(anEdgeB)) { // b
 569       aGroup->SetName("CIRCULAR_QUARTER_PIPE");
 570       theSeq->Append(aGroup);
 571     }
 572     else if (anEdgesMap.Contains(anEdgeC)) { // c
 573       aGroup->SetName("CIRCULAR_QUARTER_PIPE");
 574       theSeq->Append(aGroup);
 575     }
 576     else if (anEdgesMap.Contains(anEdgeD)) { // d
 577       aGroup->SetName("HALF_LENGTH_MAIN_PIPE");
 578       theSeq->Append(aGroup);
 579     }
 580     else if (anEdgesMap.Contains(anEdgeE)) { // e
 581       aGroup->SetName("HALF_LENGTH_INCIDENT_PIPE");
 582       theSeq->Append(aGroup);
 583     }
 584     else if (anEdgesMap.Contains(anEdgeF)) { // f
 585       aGroup->SetName("FLANGE");
 586       theSeq->Append(aGroup);
 587     }
 588     else if (shapeType == TSHAPE_CHAMFER) { // g
 589       if (anEdgesMap.Contains(anEdgeG)) {
 590         aGroup->SetName("CHAMFER");
 591         theSeq->Append(aGroup);
 592       }
 593     }
 594     else if (shapeType == TSHAPE_FILLET) { // g
 595       if (anEdgesMap.Contains(anEdgeG)) {
 596         aGroup->SetName("FILLET");
 597         theSeq->Append(aGroup);
 598       }
 599     }
 600     else {
 601     }
 602   }
 603   // END: new groups search
 604 #else
 605   bool addGroup;
 606   bool circularFoundAndAdded = false;
 607   bool circularFound10 = false;
 608   bool incidentPipeFound = false;
 609   bool mainPipeFound = false;
 610   bool mainPipeFoundAndAdded = false;
 611   bool radialFound =false;
 612   bool flangeFound = false;
 613   bool flangeFoundAndAdded = false;
 614   bool chamferOrFilletFound = false;
 615
 616   for (int i = 1 ; i <= aNbGroups; i++) {
 617     addGroup = false;
 618
 619     Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(aSeqPropagate->Value(i));
 620     if (aGroup.IsNull())
 621       continue;
 622
 623     TopoDS_Shape aGroupShape = aGroup->GetValue();
 624     BRepBuilderAPI_Transform aTransformationShapeInv (aGroupShape, aTrsfInv, Standard_False);
 625     TopoDS_Shape aGroupShapeTrsfInv = aTransformationShapeInv.Shape();
 626
 627     TopTools_IndexedMapOfShape anEdgesMap;
 628     TopExp::MapShapes(aGroupShapeTrsfInv,TopAbs_EDGE, anEdgesMap);
 629     Standard_Integer nbEdges = anEdgesMap.Extent();
 630
 631     if (shapeType == TSHAPE_BASIC) {
 632       if ((nbEdges >= 21) || /*R1Ext = R2Ext*/(nbEdges == 17)) { // 17, 17+8*{1,2,3}, 21, 21+8*{1,2,3}
 633         addGroup = true;
 634         aGroup->SetName("THICKNESS");
 635       }
 636       else if (nbEdges == 6) {
 637         if (!circularFoundAndAdded) {
 638           circularFoundAndAdded = true;
 639           addGroup = true;
 640           aGroup->SetName("CIRCULAR_QUARTER_PIPE");
 641         }
 642       }
 643       else if (nbEdges == 8) {
 644         incidentPipeFound = true;
 645         mainPipeFound = false;
 646         radialFound = false;
 647         flangeFound = false;
 648
 649         TopExp_Explorer Ex(aGroupShapeTrsfInv,TopAbs_VERTEX);
 650         while (Ex.More()) {
 651           gp_Pnt aP =  BRep_Tool::Pnt(TopoDS::Vertex(Ex.Current()));
 652           double x=aP.X(), y=aP.Y(), z=aP.Z();
 653
 654
 655           if ((Abs(x) > aR2Ext + Precision::Confusion()) ||
 656               (Abs(y) > aR2Ext + Precision::Confusion())) {
 657             incidentPipeFound = false;
 658           }
 659
 660           if ( z < -Precision::Confusion()) {
 661             // length of main pipe
 662             mainPipeFound = true;
 663             if (!mainPipeFoundAndAdded) {
 664               mainPipeFoundAndAdded = true;
 665               addGroup = true;
 666               aGroup->SetName("HALF_LENGTH_MAIN_PIPE");
 667             }
 668           }
 669
 670           else if (Abs(x) > (theL1-Precision::Confusion())) {
 671             // discretisation circulaire
 672             radialFound = true;
 673             if (!circularFoundAndAdded) {
 674               circularFoundAndAdded = true;
 675               addGroup = true;
 676               aGroup->SetName("CIRCULAR_QUARTER_PIPE");
 677             }
 678           }
 679           Ex.Next();
 680         }
 681         if (incidentPipeFound) {
 682           addGroup = true;
 683           aGroup->SetName("HALF_LENGTH_INCIDENT_PIPE");
 684         }
 685         if (!addGroup && (!incidentPipeFound &&
 686                           !radialFound &&
 687                           !mainPipeFound &&
 688                           !flangeFound)) {
 689           // Flange (collerette)
 690           flangeFound = true;
 691           addGroup = true;
 692           aGroup->SetName("FLANGE");
 693         }
 694       }
 695       else
 696         continue;
 697     }
 698     else if (shapeType == TSHAPE_CHAMFER || shapeType == TSHAPE_FILLET) {
 699       if (nbEdges >= 25) { // 25, 25+8, 25+16, 25+24
 700         addGroup = true;
 701         aGroup->SetName("THICKNESS");
 702       }
 703       else if ((nbEdges == 10) || (nbEdges == 6)) {
 704         if (!circularFoundAndAdded) {
 705           addGroup = true;
 706           circularFoundAndAdded = true;
 707           aGroup->SetName("CIRCULAR_QUARTER_PIPE");
 708           if (nbEdges == 10) {
 709             circularFound10 = true;
 710           }
 711         }
 712         else if (!circularFound10 && nbEdges == 10) {
 713           circularFound10 = true;
 714           addGroup = true;
 715           aGroup->SetName("CIRCULAR_QUARTER_PIPE");
 716         }
 717       }
 718       else if (nbEdges == 8) {
 719         incidentPipeFound = true;
 720         mainPipeFound = true;
 721         flangeFound = false;
 722
 723         bool isNearZ0 = false;
 724         bool isBelowZ0 = false;
 725
 726         TopExp_Explorer Ex (aGroupShapeTrsfInv,TopAbs_VERTEX);
 727         while (Ex.More()) {
 728           gp_Pnt aP = BRep_Tool::Pnt(TopoDS::Vertex(Ex.Current()));
 729           double x=aP.X(), y=aP.Y(), z=aP.Z();
 730
 731           // tuy_princ_long_avant & tuy_princ_long_apres
 732           //bool isMain = (((z < Precision::Confusion()) || (x < Precision::Confusion())) &&
 733           //               ((y <= aR1Ext + Precision::Confusion()) ||
 734           //                (y <= -(aR1Ext + Precision::Confusion())) ||
 735           //                (y <= theR1 + Precision::Confusion()) ||
 736           //                (y == -(theR1 + Precision::Confusion()))));
 737           bool isMain = ((z < Precision::Confusion() || x < Precision::Confusion()) &&
 738                          (fabs(y) > theR1 - Precision::Confusion() ||
 739                           fabs(y) < Precision::Confusion()));
 740
 741           if (!isMain) {
 742             mainPipeFound = false;
 743           }
 744
 745           // collerette
 746           //if (z < Precision::Confusion() && !isMain) {
 747           //  flangeFound = true;
 748           //  if (!flangeFoundAndAdded) {
 749           //    flangeFoundAndAdded = true;
 750           //    addGroup = true;
 751           //    aGroup->SetName("FLANGE");
 752           //  }
 753           //}
 754           if (fabs(z) < Precision::Confusion()) isNearZ0 = true;
 755           if (z < - Precision::Confusion()) isBelowZ0 = true;
 756
 757           // tuyau incident
 758           if ((Abs(x) > aR2Ext + Precision::Confusion()) ||
 759               (Abs(y) > aR2Ext + Precision::Confusion())) {
 760             incidentPipeFound = false;
 761           }
 762           Ex.Next();
 763         }
 764         if (mainPipeFound) {
 765           addGroup = true;
 766           aGroup->SetName("HALF_LENGTH_MAIN_PIPE");
 767         }
 768         if (incidentPipeFound) {
 769           addGroup = true;
 770           aGroup->SetName("HALF_LENGTH_INCIDENT_PIPE");
 771         }
 772         if (isNearZ0 && !isBelowZ0) {
 773           flangeFound = true;
 774           if (!flangeFoundAndAdded) {
 775             flangeFoundAndAdded = true;
 776             addGroup = true;
 777             aGroup->SetName("FLANGE");
 778           }
 779         }
 780         if (!addGroup && (!incidentPipeFound &&
 781                           !mainPipeFound &&
 782                           !flangeFound &&
 783                           !chamferOrFilletFound)) {
 784           addGroup = true;
 785           chamferOrFilletFound = true;
 786           if (shapeType == TSHAPE_CHAMFER)
 787             aGroup->SetName("CHAMFER");
 788           else
 789             aGroup->SetName("FILLET");
 790         }
 791       }
 792       else
 793         continue;
 794     }
 795     // Add group to the list
 796     if (addGroup)
 797       theSeq->Append(aGroup);
 798   }
 799 #endif
 800
 801   SetErrorCode(OK);
 802   return true;
 803 }
 804
 805 //=============================================================================
 806 /*!
 807  *  Return faces that are laying on surface.
 808  */
 809 //=============================================================================
 810 bool GEOMImpl_IAdvancedOperations::GetFacesOnSurf
 811                      (const TopoDS_Shape &theShape,
 812                       const Handle_Geom_Surface& theSurface,
 813                       const Standard_Real theTolerance,
 814                       TopTools_ListOfShape &theFaces)
 815 {
 816   GEOMAlgo_FinderShapeOn1 aFinder;
 817
 818   aFinder.SetShape(theShape);
 819   aFinder.SetTolerance(theTolerance);
 820   aFinder.SetSurface(theSurface);
 821   aFinder.SetShapeType(TopAbs_FACE);
 822   aFinder.SetState(GEOMAlgo_ST_ON);
 823
 824   // Sets the minimal number of inner points for the faces that do not have own
 825   // inner points at all (for e.g. rectangular planar faces have just 2 triangles).
 826   // Default value=3
 827   aFinder.SetNbPntsMin(3);
 828   // Sets the maximal number of inner points for edges or faces.
 829   // It is usefull for the cases when this number is very big (e.g =2000) to improve
 830   // the performance. If this value =0, all inner points will be taken into account.
 831   // Default value=0
 832   aFinder.SetNbPntsMax(100);
 833   aFinder.Perform();
 834
 835   // Interprete results
 836   Standard_Integer iErr = aFinder.ErrorStatus();
 837   // the detailed description of error codes is in GEOMAlgo_FinderShapeOn1.cxx
 838   if (iErr) {
 839     MESSAGE(" iErr : " << iErr);
 840     TCollection_AsciiString aMsg (" iErr : ");
 841     aMsg += TCollection_AsciiString(iErr);
 842     SetErrorCode(aMsg);
 843     return false;
 844   }
 845   Standard_Integer iWrn = aFinder.WarningStatus();
 846   // the detailed description of warning codes is in GEOMAlgo_FinderShapeOn1.cxx
 847   if (iWrn) {
 848     MESSAGE(" *** iWrn : " << iWrn);
 849   }
 850
 851   const TopTools_ListOfShape &aListRes = aFinder.Shapes(); // the result
 852   TopTools_ListIteratorOfListOfShape anIter (aListRes);
 853
 854   for (; anIter.More(); anIter.Next()) {
 855     theFaces.Append(anIter.Value());
 856   }
 857
 858   return true;
 859 }
 860
 861 //=============================================================================
 862 /*!
 863  *  Creates and returns conical face.
 864  */
 865 //=============================================================================
 866 TopoDS_Shape GEOMImpl_IAdvancedOperations::MakeConicalFace
 867                                   (const gp_Ax2 &theAxis,
 868                                    const double theRadius,
 869                                    const double theRadiusThin,
 870                                    const double theHeight,
 871                                    const gp_Trsf &theTrsf)
 872 {
 873   BRepPrimAPI_MakeCone aMkCone (theAxis, theRadius, theRadiusThin, theHeight);
 874   TopoDS_Shape aResult;
 875
 876   aMkCone.Build();
 877   if (aMkCone.IsDone()) {
 878     TopExp_Explorer anExp(aMkCone.Shape(), TopAbs_FACE);
 879
 880     for (; anExp.More(); anExp.Next()) {
 881       TopoDS_Face aFace = TopoDS::Face(anExp.Current());
 882
 883       if (aFace.IsNull() == Standard_False) {
 884         BRepAdaptor_Surface anAdaptor(aFace, Standard_False);
 885
 886         if (anAdaptor.GetType() == GeomAbs_Cone) {
 887           // This is a conical face. Transform and return it.
 888           BRepBuilderAPI_Transform aTransf(aFace, theTrsf, Standard_False);
 889
 890           aResult = aTransf.Shape();
 891           break;
 892         }
 893       }
 894     }
 895   }
 896
 897   return aResult;
 898 }
 899
 900 //=============================================================================
 901 /*!
 902  *  Generate the internal group of a Pipe T-Shape
 903  */
 904 //=============================================================================
 905 bool GEOMImpl_IAdvancedOperations::MakeInternalGroup
 906                       (const Handle(GEOM_Object) &theShape,
 907                        const double theR1, const double theLen1,
 908                        const double theR2, const double theLen2,
 909                        const double theRL, double theTransLenL,
 910                        const double theRR, double theTransLenR,
 911                        const double theRI, double theTransLenI,
 912                        const Handle(TColStd_HSequenceOfTransient) &theSeq,
 913                        const gp_Trsf &theTrsf)
 914 {
 915   SetErrorCode(KO);
 916
 917   if (theShape.IsNull()) {
 918     return false;
 919   }
 920
 921   TopoDS_Shape aShape = theShape->GetValue();
 922
 923   if (aShape.IsNull()) {
 924     SetErrorCode("Shape is not defined");
 925     return false;
 926   }
 927
 928   // Compute tolerance
 929   Standard_Real aMaxTol = -RealLast();
 930   TopExp_Explorer anExp(aShape, TopAbs_VERTEX);
 931
 932   for (; anExp.More(); anExp.Next()) {
 933     TopoDS_Vertex aVertex = TopoDS::Vertex(anExp.Current());
 934
 935     if (aVertex.IsNull() == Standard_False) {
 936       const Standard_Real aTol = BRep_Tool::Tolerance(aVertex);
 937
 938       if (aTol > aMaxTol) {
 939         aMaxTol = aTol;
 940       }
 941     }
 942   }
 943
 944   // Construct internal surfaces.
 945   Standard_Integer i = 0;
 946   const Standard_Integer aMaxNbSurf = 5;
 947   Handle(Geom_Surface) aSurface[aMaxNbSurf];
 948   TopTools_ListOfShape aConicalFaces;
 949   Standard_Real aTolConf = Precision::Confusion();
 950
 951   // 1. Construct the internal surface of main pipe.
 952   gp_Ax2 anAxis1 (gp::Origin(), gp::DX(), gp::DZ());
 953   gp_Ax2 anAxis2 (gp::Origin(), gp::DZ(), gp::DX());
 954
 955   aSurface[i++] = new Geom_CylindricalSurface(anAxis1, theR1);
 956
 957   // 2. Construct the internal surface of incident pipe.
 958   aSurface[i++] = new Geom_CylindricalSurface(anAxis2, theR2);
 959
 960   // 3. Construct the internal surface of left reduction pipe.
 961   if (theRL > aTolConf) {
 962     aSurface[i++] = new Geom_CylindricalSurface(anAxis1, theRL);
 963
 964     if (theTransLenL > aTolConf) {
 965       // 3.1. Construct the internal surface of left transition pipe.
 966       gp_Pnt aPLeft (-theLen1, 0., 0.);
 967       gp_Ax2 anAxisLeft (aPLeft, -gp::DX(), gp::DZ());
 968       TopoDS_Shape aConeLeft =
 969         MakeConicalFace(anAxisLeft, theR1, theRL, theTransLenL, theTrsf);
 970
 971       if (aConeLeft.IsNull() == Standard_False) {
 972         aConicalFaces.Append(aConeLeft);
 973       }
 974     }
 975   }
 976
 977   // 4. Construct the internal surface of right reduction pipe.
 978   if (theRR > aTolConf) {
 979     // There is no need to construct another cylinder of the same radius. Skip it.
 980     if (Abs(theRR - theRL) > aTolConf) {
 981       aSurface[i++] = new Geom_CylindricalSurface(anAxis1, theRR);
 982     }
 983
 984     if (theTransLenL > aTolConf) {
 985       // 4.1. Construct the internal surface of right transition pipe.
 986       gp_Pnt aPRight (theLen1, 0., 0.);
 987       gp_Ax2 anAxisRight (aPRight, gp::DX(), gp::DZ());
 988       TopoDS_Shape aConeRight =
 989         MakeConicalFace(anAxisRight, theR1, theRR, theTransLenR, theTrsf);
 990
 991       if (aConeRight.IsNull() == Standard_False) {
 992         aConicalFaces.Append(aConeRight);
 993       }
 994     }
 995   }
 996
 997   // 5. Construct the internal surface of incident reduction pipe.
 998   if (theRI > aTolConf) {
 999     aSurface[i++] = new Geom_CylindricalSurface(anAxis2, theRI);
1000
1001     if (theTransLenI > aTolConf) {
1002       // 5.1. Construct the internal surface of incident transition pipe.
1003       gp_Pnt aPInci (0., 0., theLen2);
1004       gp_Ax2 anAxisInci (aPInci, gp::DZ(), gp::DX());
1005       TopoDS_Shape aConeInci =
1006         MakeConicalFace(anAxisInci, theR2, theRI, theTransLenI, theTrsf);
1007
1008       if (aConeInci.IsNull() == Standard_False) {
1009         aConicalFaces.Append(aConeInci);
1010       }
1011     }
1012   }
1013
1014   // Get faces that are laying on cylindrical surfaces.
1015   TopTools_ListOfShape aFaces;
1016   gp_Trsf anInvTrsf = theTrsf.Inverted();
1017
1018   for (i = 0; i < aMaxNbSurf; i++) {
1019     if (aSurface[i].IsNull()) {
1020       break;
1021     }
1022
1023     aSurface[i]->Transform(theTrsf);
1024
1025     TopTools_ListOfShape aLocalFaces;
1026
1027     if (!GetFacesOnSurf(aShape, aSurface[i], aMaxTol, aLocalFaces)) {
1028       return false;
1029     }
1030
1031     if (i < 2) {
1032       // Check if the result contains outer cylinders.
1033       // It is required for main and incident pipes.
1034       TopTools_ListIteratorOfListOfShape anIter(aLocalFaces);
1035
1036       while (anIter.More()) {
1037         TopExp_Explorer anExp(anIter.Value(), TopAbs_VERTEX);
1038         Standard_Boolean isInside = Standard_False;
1039
1040         // Get a vertex from this shape
1041         if (anExp.More()) {
1042           TopoDS_Vertex aVtx = TopoDS::Vertex(anExp.Current());
1043
1044           if (aVtx.IsNull() == Standard_False) {
1045             gp_Pnt aPnt = BRep_Tool::Pnt(aVtx);
1046
1047             aPnt.Transform(anInvTrsf);
1048
1049             if (i == 0) {
1050               // Check if the point is inside the main pipe.
1051               isInside = (Abs(aPnt.X()) <= theLen1);
1052             } else { // i == 1
1053               // Check if the point is inside the incident pipe.
1054               isInside = (aPnt.Z() <= theLen2);
1055             }
1056           }
1057         }
1058
1059         if (isInside) {
1060           // Keep this face.
1061           anIter.Next();
1062         } else {
1063           // Remove this face.
1064           aLocalFaces.Remove(anIter);
1065         }
1066       }
1067     }
1068
1069     aFaces.Append(aLocalFaces);
1070   }
1071
1072   // Get faces that are laying on conical faces.
1073   if (aConicalFaces.IsEmpty() == Standard_False) {
1074     Handle(GEOM_Object) aCone =
1075       GetEngine()->AddObject(GetDocID(), GEOM_TSHAPE);
1076     Handle(GEOM_Function) aFunction =
1077       aCone->AddFunction(GEOMImpl_PipeTShapeDriver::GetID(), TSHAPE_BASIC);
1078     TopTools_ListIteratorOfListOfShape aFIter(aConicalFaces);
1079     Handle(GEOM_Object) aConeFromShape;
1080
1081     for (; aFIter.More(); aFIter.Next()) {
1082       aFunction->SetValue(aFIter.Value());
1083       aConeFromShape = myShapesOperations->GetInPlace(theShape, aCone);
1084
1085       if (aConeFromShape.IsNull() == Standard_False) {
1086         aConeFromShape->GetLastFunction()->SetDescription("");
1087         TopoDS_Shape aConeFaces = aConeFromShape->GetValue();
1088         TopExp_Explorer anExp(aConeFaces, TopAbs_FACE);
1089
1090         for (; anExp.More(); anExp.Next()) {
1091           TopoDS_Face aConeFace = TopoDS::Face(anExp.Current());
1092
1093           if (aConeFace.IsNull() == Standard_False) {
1094             aFaces.Append(aConeFace);
1095           }
1096         }
1097       }
1098     }
1099   }
1100
1101   // Create a group of internal faces.
1102   if (aFaces.IsEmpty() == Standard_False) {
1103     Handle(GEOM_Object) aGroup = myGroupOperations->CreateGroup(theShape, TopAbs_FACE);
1104
1105     if (aGroup.IsNull() == Standard_False) {
1106       aGroup->GetLastFunction()->SetDescription("");
1107       aGroup->SetName("INTERNAL_FACES");
1108
1109       TopTools_IndexedMapOfShape anIndices;
1110       Handle(TColStd_HSequenceOfInteger) aSeqIDs = new TColStd_HSequenceOfInteger;
1111
1112       TopExp::MapShapes(aShape, anIndices);
1113
1114       TopTools_ListIteratorOfListOfShape anIter(aFaces);
1115
1116       for (; anIter.More(); anIter.Next()) {
1117         const TopoDS_Shape &aFace = anIter.Value();
1118         const Standard_Integer anIndex = anIndices.FindIndex(aFace);
1119
1120         if (anIndex > 0) {
1121           aSeqIDs->Append(anIndex);
1122         }
1123       }
1124
1125       myGroupOperations->UnionIDs(aGroup, aSeqIDs);
1126       aGroup->GetLastFunction()->SetDescription("");
1127       theSeq->Append(aGroup);
1128     }
1129   }
1130
1131   SetErrorCode(OK);
1132
1133   return true;
1134 }
1135
1136 bool GEOMImpl_IAdvancedOperations::MakePipeTShapePartition(Handle(GEOM_Object) theShape,
1137                                                            double theR1, double theW1, double theL1,
1138                                                            double theR2, double theW2, double theL2,
1139                                                            double theH, double theW,
1140                                                            double theRF, bool isNormal)
1141 {
1142   SetErrorCode(KO);
1143
1144   // Build tools for partition operation:
1145   // 1 face and 2 planes
1146   // Face
1147   Handle(GEOM_Object) arete_intersect_int, arete_intersect_ext;
1148   Handle(GEOM_Object) wire_t, wire_t2, face_t, face_t2;
1149   Handle(GEOM_Object) chan_racc;
1150   Handle(GEOM_Object) vi1, vi2;
1151   Handle(GEOM_Object) Te3;
1152
1153   try {
1154     OCC_CATCH_SIGNALS;
1155     Handle(GEOM_Object) Vector_Z = myBasicOperations->MakeVectorDXDYDZ(0, 0, 1);
1156     Vector_Z->GetLastFunction()->SetDescription("");
1157
1158     // Useful values
1159     double aSize = 2*(theL1 + theL2);
1160     double aR1Ext = theR1 + theW1;
1161     double aR2Ext = theR2 + theW2;
1162     double theVertCylinderRadius = aR2Ext + theW + theRF;
1163     double theHoriCylinderRadius = aR1Ext + theH + theRF;
1164
1165     // Common edges on internal cylinder
1166     Handle(GEOM_Object) box_i = my3DPrimOperations->MakeBoxDXDYDZ(theR2, theR2, theR1);
1167     box_i->GetLastFunction()->SetDescription("");
1168     box_i = myTransformOperations->TranslateDXDYDZ(box_i, -theR2, -theR2, 0);
1169     box_i->GetLastFunction()->SetDescription("");
1170
1171     Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
1172     TCollection_AsciiString theDesc = aFunction->GetDescription();
1173     Handle(TColStd_HSequenceOfTransient) edges_i =
1174       myShapesOperations->GetShapesOnBox(box_i, theShape, TopAbs_EDGE, GEOMAlgo_ST_IN);
1175     // Recover previous description to get rid of Propagate dump
1176     aFunction->SetDescription(theDesc);
1177     if (edges_i.IsNull() || edges_i->Length() == 0) {
1178       SetErrorCode("Internal edges not found");
1179       return false;
1180     }
1181     for (int i=1; i<=edges_i->Length();i++) {
1182       Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast(edges_i->Value(i));
1183       anObj->GetLastFunction()->SetDescription("");
1184     }
1185     arete_intersect_int = Handle(GEOM_Object)::DownCast(edges_i->Value(1));
1186
1187     // search for vertices located on both internal pipes
1188     aFunction = theShape->GetLastFunction();
1189     theDesc = aFunction->GetDescription();
1190     Handle(TColStd_HSequenceOfTransient) vertices_i =
1191       myShapesOperations->GetShapesOnBox(box_i, theShape, TopAbs_VERTEX, GEOMAlgo_ST_ONIN);
1192     // Recover previous description to get rid of Propagate dump
1193     aFunction->SetDescription(theDesc);
1194     if (vertices_i.IsNull() || vertices_i->Length() == 0) {
1195       SetErrorCode("Internal vertices not found");
1196       return false;
1197     }
1198
1199     double d1min = theR2+theW2, d2min=theR2+theW2;
1200     for (int i = 1; i <= vertices_i->Length(); i++) {
1201       Handle(GEOM_Object) v = Handle(GEOM_Object)::DownCast(vertices_i->Value(i));
1202       v->GetLastFunction()->SetDescription("");
1203       TopoDS_Vertex aVertex = TopoDS::Vertex(v->GetValue());
1204       gp_Pnt aP = BRep_Tool::Pnt(aVertex);
1205       if (Abs(aP.X()) <= Precision::Confusion()) {
1206         if (Abs(aP.Y()) < d1min) {
1207           vi1 = v;
1208           d1min = Abs(aP.Y());
1209         }
1210       } else if (Abs(aP.Y()) <= Precision::Confusion()) {
1211         if (Abs(aP.X()) < d2min) {
1212           vi2 = v;
1213           d2min = Abs(aP.X());
1214         }
1215       }
1216     }
1217     if (vi1.IsNull() || vi2.IsNull()) {
1218       SetErrorCode("Cannot find internal intersection vertices");
1219       return false;
1220     }
1221
1222     std::list<Handle(GEOM_Object)> theShapes;
1223
1224     if (isNormal) {
1225       Handle(GEOM_Object) ve1, ve2;
1226       TopoDS_Vertex vertex1, vertex2;
1227
1228       Handle(GEOM_Object) box_e = my3DPrimOperations->MakeBoxDXDYDZ(aR2Ext, aR2Ext, aR1Ext);
1229       box_e->GetLastFunction()->SetDescription("");
1230       box_e = myTransformOperations->TranslateDXDYDZ(box_e, -aR2Ext, -aR2Ext, 0);
1231       box_e->GetLastFunction()->SetDescription("");
1232
1233       // search for vertices located on both external pipes
1234       aFunction = theShape->GetLastFunction();
1235       theDesc = aFunction->GetDescription();
1236       Handle(TColStd_HSequenceOfTransient) vertices_e =
1237         myShapesOperations->GetShapesOnBox(box_e, theShape, TopAbs_VERTEX, GEOMAlgo_ST_ONIN);
1238       // Recover previous description to get rid of Propagate dump
1239       aFunction->SetDescription(theDesc);
1240       if (vertices_e.IsNull() || vertices_e->Length() == 0) {
1241         SetErrorCode("External vertices not found");
1242         return false;
1243       }
1244
1245       double d1max = 0, d2max = 0;
1246       for (int i = 1; i <= vertices_e->Length(); i++) {
1247         Handle(GEOM_Object) v = Handle(GEOM_Object)::DownCast(vertices_e->Value(i));
1248         v->GetLastFunction()->SetDescription("");
1249         TopoDS_Vertex aVertex = TopoDS::Vertex(v->GetValue());
1250         gp_Pnt aP = BRep_Tool::Pnt(aVertex);
1251         if (Abs(aP.X()) <= Precision::Confusion()) {
1252           if (Abs(aP.Y()) > d1max) {
1253             ve1 = v;
1254             vertex1 = aVertex;
1255             d1max = Abs(aP.Y());
1256           }
1257         } else if (Abs(aP.Y()) <= Precision::Confusion()) {
1258           if (Abs(aP.X()) > d2max) {
1259             ve2 = v;
1260             vertex2 = aVertex;
1261             d2max = Abs(aP.X());
1262           }
1263         }
1264       }
1265       if (ve1.IsNull() || ve2.IsNull()) {
1266         SetErrorCode("Cannot find external intersection vertices");
1267         return false;
1268       }
1269       Handle(GEOM_Object) edge_e1, edge_e2;
1270
1271       // Common edges on external cylinder
1272       aFunction = theShape->GetLastFunction();
1273       theDesc = aFunction->GetDescription();
1274       Handle(TColStd_HSequenceOfTransient) edges_e =
1275         myShapesOperations->GetShapesOnBox(box_e, theShape, TopAbs_EDGE, GEOMAlgo_ST_IN);
1276       // Recover previous description to get rid of Propagate dump
1277       aFunction->SetDescription(theDesc);
1278       if (edges_e.IsNull() || edges_e->Length() == 0) {
1279         SetErrorCode("External edges not found");
1280         return false;
1281       }
1282       ShapeAnalysis_Edge sae;
1283       for (int i=1; i<=edges_e->Length();i++) {
1284         Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast(edges_e->Value(i));
1285         anObj->GetLastFunction()->SetDescription("");
1286         TopoDS_Edge anEdge = TopoDS::Edge(anObj->GetValue());
1287         if ( !anEdge.IsNull() &&
1288              (sae.FirstVertex(anEdge).IsSame(vertex1) || sae.LastVertex(anEdge).IsSame(vertex1)) &&
1289              (sae.FirstVertex(anEdge).IsSame(vertex2) || sae.LastVertex(anEdge).IsSame(vertex2))) {
1290           arete_intersect_ext = anObj;
1291         }
1292       }
1293
1294       edge_e1 = myBasicOperations->MakeLineTwoPnt(ve1, vi1);
1295       if (edge_e1.IsNull()) {
1296         SetErrorCode("Edge 1 could not be built");
1297         return false;
1298       }
1299
1300       edge_e2 = myBasicOperations->MakeLineTwoPnt(ve2, vi2);
1301       if (edge_e2.IsNull()) {
1302         SetErrorCode("Edge 2 could not be built");
1303         return false;
1304       }
1305
1306       edge_e1->GetLastFunction()->SetDescription("");
1307       edge_e2->GetLastFunction()->SetDescription("");
1308
1309       std::list<Handle(GEOM_Object)> edge_e_elist;
1310       edge_e_elist.push_back(arete_intersect_int);
1311       edge_e_elist.push_back(edge_e1);
1312       edge_e_elist.push_back(arete_intersect_ext);
1313       edge_e_elist.push_back(edge_e2);
1314       wire_t = myShapesOperations->MakeWire(edge_e_elist, 1e-7);
1315       if (wire_t.IsNull()) {
1316         SetErrorCode("Impossible to build wire");
1317         return false;
1318       }
1319       wire_t->GetLastFunction()->SetDescription("");
1320       face_t = myShapesOperations->MakeFace(wire_t, false);
1321       if (face_t.IsNull()) {
1322         SetErrorCode("Impossible to build face");
1323         return false;
1324       }
1325       face_t->GetLastFunction()->SetDescription("");
1326
1327       theShapes.push_back(theShape);
1328       theShapes.push_back(vi1);
1329       theShapes.push_back(vi2);
1330       theShapes.push_back(ve1);
1331       theShapes.push_back(ve2);
1332       theShapes.push_back(edge_e1);
1333       theShapes.push_back(edge_e2);
1334       theShapes.push_back(wire_t);
1335       theShapes.push_back(face_t);
1336     }
1337     else {
1338       Handle(GEOM_Object) P1, P2, P3, P4, P5, P6;
1339       int idP1, idP2, idP3, idP4;
1340       int PZX, PZY;
1341       double ZX=0, ZY=0;
1342       std::vector<int> LX;
1343       std::vector<int> LY;
1344       Handle(GEOM_Object) box_e = my3DPrimOperations->MakeBoxDXDYDZ
1345         (theVertCylinderRadius, theVertCylinderRadius, theHoriCylinderRadius);
1346       box_e->GetLastFunction()->SetDescription("");
1347       box_e = myTransformOperations->TranslateDXDYDZ
1348         (box_e, -theVertCylinderRadius, -theVertCylinderRadius, 0);
1349       box_e->GetLastFunction()->SetDescription("");
1350
1351       aFunction = theShape->GetLastFunction();
1352       theDesc = aFunction->GetDescription();
1353       Handle(TColStd_HSequenceOfTransient) extremVertices =
1354         myShapesOperations->GetShapesOnBox(box_e, theShape, TopAbs_VERTEX, GEOMAlgo_ST_ONIN);
1355       // Recover previous description to get rid of Propagate dump
1356       aFunction->SetDescription(theDesc);
1357
1358       if (extremVertices.IsNull() || extremVertices->Length() == 0) {
1359         if (theRF == 0)
1360           SetErrorCode("Vertices on chamfer not found");
1361         else
1362           SetErrorCode("Vertices on fillet not found");
1363         return false;
1364       }
1365
1366       theShapes.push_back(theShape);
1367       theShapes.push_back(box_e);
1368       if (extremVertices->Length() != 6) {
1369         //           for (int i=1; i<=extremVertices->Length(); i++){
1370         //             theShapes.push_back(Handle(GEOM_Object)::DownCast(extremVertices->Value(i)));
1371         //           }
1372         //           Handle(GEOM_Object) aCompound = myShapesOperations->MakeCompound(theShapes);
1373         //           TopoDS_Shape aCompoundShape = aCompound->GetValue();
1374         //           theShape->GetLastFunction()->SetValue(aCompoundShape);
1375         SetErrorCode("Bad number of vertices on chamfer found");
1376         return false;
1377       }
1378
1379       for (int i=1; i<=extremVertices->Length(); i++){
1380         Handle(GEOM_Object) aV = Handle(GEOM_Object)::DownCast(extremVertices->Value(i));
1381         aV->GetLastFunction()->SetDescription("");
1382         gp_Pnt aP = BRep_Tool::Pnt(TopoDS::Vertex(aV->GetValue()));
1383
1384         if (Abs(aP.X()) <= Precision::Confusion()) {
1385           if (Abs(aP.Y()) - theR2 > Precision::Confusion()) {
1386             LX.push_back(i);
1387             if  (aP.Z()-ZX > Precision::Confusion()) {
1388               ZX = aP.Z();
1389               PZX = i;
1390             }
1391           }
1392         }
1393         else {
1394           if (Abs(aP.X()) - theR2 > Precision::Confusion()) {
1395             LY.push_back(i);
1396             if (aP.Z() - ZY > Precision::Confusion()) {
1397               ZY = aP.Z();
1398               PZY = i;
1399             }
1400           }
1401         }
1402       }
1403
1404       idP2 = PZX;
1405       idP4 = PZY;
1406       idP1 = LX.at(0);
1407       if (LX.at(0) == PZX)
1408         idP1 = LX.at(1);
1409       idP3 = LY.at(0);
1410       if (LY.at(0) == PZY)
1411         idP3 = LY.at(1);
1412
1413       P1 = Handle(GEOM_Object)::DownCast(extremVertices->Value(idP1));
1414       P2 = Handle(GEOM_Object)::DownCast(extremVertices->Value(idP2));
1415       P3 = Handle(GEOM_Object)::DownCast(extremVertices->Value(idP3));
1416       P4 = Handle(GEOM_Object)::DownCast(extremVertices->Value(idP4));
1417
1418       Handle(GEOM_Object) Cote_1 = myBasicOperations->MakeLineTwoPnt(P1, vi1);
1419       if (Cote_1.IsNull()) {
1420         SetErrorCode("Impossible to build edge in thickness");
1421         return false;
1422       }
1423       Cote_1->GetLastFunction()->SetDescription("");
1424
1425       Handle(GEOM_Object) Cote_2 = myBasicOperations->MakeLineTwoPnt(vi2, P3);
1426       if (Cote_2.IsNull()) {
1427         SetErrorCode("Impossible to build edge in thickness");
1428         return false;
1429       }
1430       Cote_2->GetLastFunction()->SetDescription("");
1431
1432       // edge_chan_princ = arete du chanfrein (ou raccord) sur le tuyau principal
1433       // edge_chan_inc = arete du chanfrein (ou raccord) sur le tuyau incident
1434       //         std::cerr << "Getting chamfer edge on main pipe" << std::endl;
1435       Handle(GEOM_Object) edge_chan_princ = myBlocksOperations->GetEdge(theShape, P1, P3);
1436       if (edge_chan_princ.IsNull()) {
1437         SetErrorCode("Impossible to find edge on main pipe");
1438         return false;
1439       }
1440       edge_chan_princ->GetLastFunction()->SetDescription("");
1441
1442       Handle(GEOM_Object) edge_chan_inc = myBlocksOperations->GetEdge(theShape, P2, P4);
1443       if (edge_chan_inc.IsNull()) {
1444         SetErrorCode("Impossible to find edge on incident pipe");
1445         return false;
1446       }
1447       edge_chan_inc->GetLastFunction()->SetDescription("");
1448
1449       std::list<Handle(GEOM_Object)> edgeList1;
1450       edgeList1.push_back(edge_chan_princ);
1451       edgeList1.push_back(Cote_1);
1452       edgeList1.push_back(arete_intersect_int);
1453       edgeList1.push_back(Cote_2);
1454
1455       //         std::cerr << "Creating wire 1" << std::endl;
1456       wire_t = myShapesOperations->MakeWire(edgeList1, 1e-7);
1457       if (wire_t.IsNull()) {
1458         SetErrorCode("Impossible to build wire");
1459         return false;
1460       }
1461       wire_t->GetLastFunction()->SetDescription("");
1462
1463       //         std::cerr << "Creating face 1" << std::endl;
1464       face_t = myShapesOperations->MakeFace(wire_t, false);
1465       if (face_t.IsNull()) {
1466         SetErrorCode("Impossible to build face");
1467         return false;
1468       }
1469       face_t->GetLastFunction()->SetDescription("");
1470       theShapes.push_back(face_t);
1471
1472       gp_Pnt aP2 = BRep_Tool::Pnt(TopoDS::Vertex(P2->GetValue()));
1473       gp_Pnt aP5 = BRep_Tool::Pnt(TopoDS::Vertex(vi1->GetValue()));
1474       double deltaZ = aP2.Z() - aP5.Z();
1475       //         std::cerr << "Creating new point from vi1 with deltaZ = " << deltaZ << std::endl;
1476       Handle(GEOM_Object) P5bis = myTransformOperations->TranslateDXDYDZCopy(vi1, 0, 0, deltaZ);
1477       if (P5bis.IsNull()) {
1478         SetErrorCode("Impossible to translate vertex");
1479         return false;
1480       }
1481       P5bis->GetLastFunction()->SetDescription("");
1482
1483       gp_Pnt aP4 = BRep_Tool::Pnt(TopoDS::Vertex(P4->GetValue()));
1484       gp_Pnt aP6 = BRep_Tool::Pnt(TopoDS::Vertex(vi2->GetValue()));
1485       deltaZ = aP4.Z() - aP6.Z();
1486       //         std::cerr << "Creating new point from vi2 with deltaZ = " << deltaZ << std::endl;
1487       Handle(GEOM_Object) P6bis = myTransformOperations->TranslateDXDYDZCopy(vi2, 0, 0, deltaZ);
1488       if (P6bis.IsNull()) {
1489         SetErrorCode("Impossible to translate vertex");
1490         return false;
1491       }
1492       P6bis->GetLastFunction()->SetDescription("");
1493
1494       //         std::cerr << "Creating new line 1 from 2 previous points" << std::endl;
1495       Handle(GEOM_Object) Cote_3 = myBasicOperations->MakeLineTwoPnt(P5bis, P2);
1496       if (Cote_3.IsNull()) {
1497         SetErrorCode("Impossible to build edge in thickness");
1498         return false;
1499       }
1500       Cote_3->GetLastFunction()->SetDescription("");
1501
1502       //         std::cerr << "Creating new line 2 from 2 previous points" << std::endl;
1503       Handle(GEOM_Object) Cote_4 = myBasicOperations->MakeLineTwoPnt(P6bis, P4);
1504       if (Cote_4.IsNull()) {
1505         SetErrorCode("Impossible to build edge in thickness");
1506         return false;
1507       }
1508       Cote_4->GetLastFunction()->SetDescription("");
1509
1510       //         std::cerr << "Creating new line 3 from 2 previous points" << std::endl;
1511       Handle(GEOM_Object) Cote_5 = myBasicOperations->MakeLineTwoPnt(P5bis, P6bis);
1512       if (Cote_4.IsNull()) {
1513         SetErrorCode("Impossible to build edge in thickness");
1514         return false;
1515       }
1516       Cote_5->GetLastFunction()->SetDescription("");
1517
1518       //std::list<Handle(GEOM_Object)> edgeList2;
1519       //edgeList2.push_back(edge_chan_inc);
1520       //edgeList2.push_back(Cote_3);
1521       //edgeList2.push_back(Cote_5);
1522       //edgeList2.push_back(Cote_4);
1523       //         std::cerr << "Creating wire 2" << std::endl;
1524       //wire_t2 = myShapesOperations->MakeWire(edgeList2, 1e-7);
1525       //if (wire_t2.IsNull()) {
1526       //  SetErrorCode("Impossible to build wire");
1527       //  return false;
1528       //}
1529       //wire_t2->GetLastFunction()->SetDescription("");
1530       //         std::cerr << "Creating face 2" << std::endl;
1531       //face_t2 = myShapesOperations->MakeFace(wire_t2, false);
1532
1533       // Mantis issue 0021682
1534       face_t2 = my3DPrimOperations->MakePrismVecH(edge_chan_inc, Cote_4, - (theR2 + theW2));
1535       //face_t2 = my3DPrimOperations->MakePrismVecH(edge_chan_inc, Cote_4, - 2.0*theR2);
1536       if (face_t2.IsNull()) {
1537         SetErrorCode("Impossible to build face");
1538         return false;
1539       }
1540       face_t2->GetLastFunction()->SetDescription("");
1541       theShapes.push_back(face_t2);
1542     }
1543
1544     // Planes
1545     Handle(GEOM_Object) aP0 = myBasicOperations->MakePointXYZ(0, 0, 0);
1546     Handle(GEOM_Object) aVZ = myBasicOperations->MakeVectorDXDYDZ(0, 0, 1);
1547     Handle(GEOM_Object) aVXZ = myBasicOperations->MakeVectorDXDYDZ(aR1Ext, 0, 0.5*(theL1+theVertCylinderRadius));
1548     Handle(GEOM_Object) aPlnOZ = myBasicOperations->MakePlanePntVec(aP0, aVZ, aSize);
1549     Handle(GEOM_Object) aPlnOXZ = myBasicOperations->MakePlanePntVec(aP0, aVXZ, aSize);
1550     aP0->GetLastFunction()->SetDescription("");
1551     aVZ->GetLastFunction()->SetDescription("");
1552     aVXZ->GetLastFunction()->SetDescription("");
1553     aPlnOZ->GetLastFunction()->SetDescription("");
1554     aPlnOXZ->GetLastFunction()->SetDescription("");
1555     theShapes.push_back(aPlnOZ);
1556     theShapes.push_back(aPlnOXZ);
1557
1558     // Partition
1559     Handle(TColStd_HSequenceOfTransient) partitionShapes = new TColStd_HSequenceOfTransient;
1560     Handle(TColStd_HSequenceOfTransient) theTools = new TColStd_HSequenceOfTransient;
1561     Handle(TColStd_HSequenceOfTransient) theKeepInside = new TColStd_HSequenceOfTransient;
1562     Handle(TColStd_HSequenceOfTransient) theRemoveInside = new TColStd_HSequenceOfTransient;
1563     Handle(TColStd_HArray1OfInteger) theMaterials;
1564
1565     partitionShapes->Append(theShape);
1566     theTools->Append(aPlnOZ);
1567     if (Abs(aR1Ext - aR2Ext) > Precision::Confusion())
1568       theTools->Append(aPlnOXZ);
1569     theTools->Append(face_t);
1570     if (!isNormal)
1571       theTools->Append(face_t2);
1572
1573     Te3 = myBooleanOperations->MakePartition
1574               (partitionShapes, theTools, theKeepInside, theRemoveInside,
1575               TopAbs_SOLID, false, theMaterials, 0, false, Standard_False);
1576     if (Te3.IsNull()) {
1577       SetErrorCode("Impossible to build partition of TShape");
1578       return false;
1579     }
1580     Te3->GetLastFunction()->SetDescription("");
1581
1582     // Last verification: result should be a block
1583     std::list<GEOMImpl_IBlocksOperations::BCError> errList;
1584     if (!myBlocksOperations->CheckCompoundOfBlocks(Te3,errList)) {
1585       SetErrorCode("TShape is not a compound of block");
1586       return false;
1587     }
1588
1589 //     // BEGIN Compound of created shapes - Only for debug purpose
1590 //     theShapes.clear();
1591 //     theShapes.push_back(theShape);
1592 //     theShapes.push_back(aPlnOZ);
1593 //     if (Abs(aR1Ext - aR2Ext) > Precision::Confusion() )
1594 //       theShapes.push_back(aPlnOXZ);
1595 //     theShapes.push_back(face_t);
1596 //     if (!isNormal)
1597 //       theShapes.push_back(face_t2);
1598 //
1599 //     Handle(GEOM_Object) aCompound = myShapesOperations->MakeCompound(theShapes);
1600 //     TopoDS_Shape aCompoundShape = aCompound->GetValue();
1601 //     theShape->GetLastFunction()->SetValue(aCompoundShape);
1602 //     // END Compound of created shapes - Only for debug purpose
1603
1604     TopoDS_Shape aShape = Te3->GetValue();
1605     theShape->GetLastFunction()->SetValue(aShape);
1606   }
1607   catch (Standard_Failure) {
1608     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1609     SetErrorCode(aFail->GetMessageString());
1610     return false;
1611   }
1612
1613   SetErrorCode(OK);
1614   return true;
1615 }
1616
1617 // Mirror and glue faces
1618 bool GEOMImpl_IAdvancedOperations::MakePipeTShapeMirrorAndGlue(Handle(GEOM_Object) theShape,
1619                                                                double theR1, double theW1, double theL1,
1620                                                                double theR2, double theW2, double theL2)
1621 {
1622   SetErrorCode(KO);
1623
1624   // Useful values
1625   double aSize = 2*(theL1 + theL2);
1626   double aR1Ext = theR1 + theW1;
1627
1628   // Planes
1629   Handle(GEOM_Object) aP0 = myBasicOperations->MakePointXYZ(0, 0, 0);
1630   aP0->GetLastFunction()->SetDescription("");
1631   Handle(GEOM_Object) aVX = myBasicOperations->MakeVectorDXDYDZ(1, 0, 0);
1632   Handle(GEOM_Object) aVY = myBasicOperations->MakeVectorDXDYDZ(0, 1, 0);
1633   aVX->GetLastFunction()->SetDescription("");
1634   aVY->GetLastFunction()->SetDescription("");
1635   Handle(GEOM_Object) aPlane_OX = myBasicOperations->MakePlanePntVec(aP0, aVX, 2*(aR1Ext + theL2));
1636   Handle(GEOM_Object) aPlane_OY = myBasicOperations->MakePlanePntVec(aP0, aVY, aSize);
1637   aPlane_OX->GetLastFunction()->SetDescription("");
1638   aPlane_OY->GetLastFunction()->SetDescription("");
1639
1640   Handle(GEOM_Object) Te4 = myTransformOperations->MirrorPlaneCopy(theShape, aPlane_OX);
1641   if (Te4.IsNull()) {
1642     SetErrorCode("Impossible to build mirror of quarter TShape");
1643     return false;
1644   }
1645
1646   Handle(GEOM_Object) Te5 = myTransformOperations->MirrorPlaneCopy(theShape, aPlane_OY);
1647   if (Te5.IsNull()) {
1648     SetErrorCode("Impossible to build mirror of half TShape");
1649     return false;
1650   }
1651
1652   Handle(GEOM_Object) Te6 = myTransformOperations->MirrorPlaneCopy(Te4, aPlane_OY);
1653   if (Te6.IsNull()) {
1654     SetErrorCode("Impossible to build mirror of half TShape");
1655     return false;
1656   }
1657
1658   std::list<Handle(GEOM_Object)> aShapesList;
1659   aShapesList.push_back(theShape);
1660   aShapesList.push_back(Te4);
1661   aShapesList.push_back(Te5);
1662   aShapesList.push_back(Te6);
1663   Handle(GEOM_Object) Te7 = myShapesOperations->MakeCompound(aShapesList);
1664   if (Te7.IsNull()) {
1665     SetErrorCode("Impossible to build compound");
1666     return false;
1667   }
1668
1669   // Copy source shape
1670   TopoDS_Shape aShapeCopy;
1671   TColStd_IndexedDataMapOfTransientTransient aMapTShapes;
1672   TNaming_CopyShape::CopyTool(Te7->GetValue(), aMapTShapes, aShapeCopy);
1673
1674   Handle(GEOM_Object) Te8 = myShapesOperations->MakeGlueFaces(Te7, 1e-7, true);
1675   if (Te8.IsNull()) {
1676     SetErrorCode("Impossible to glue faces of TShape");
1677     return false;
1678   }
1679
1680   TopoDS_Shape aShape = Te8->GetValue();
1681   BRepCheck_Analyzer anAna (aShape, Standard_True);
1682
1683   if (!anAna.IsValid()) {
1684     // Try to do gluing with the tolerance equal to maximal
1685     // tolerance of vertices of the source shape.
1686     Standard_Real aTolMax = -RealLast();
1687
1688     for (TopExp_Explorer ExV (aShapeCopy, TopAbs_VERTEX); ExV.More(); ExV.Next()) {
1689       TopoDS_Vertex aVertex = TopoDS::Vertex(ExV.Current());
1690       Standard_Real aTol = BRep_Tool::Tolerance(aVertex);
1691
1692       if (aTol > aTolMax) {
1693         aTolMax = aTol;
1694       }
1695     }
1696
1697     // Perform gluing
1698     Te7->GetLastFunction()->SetValue(aShapeCopy);
1699     Te8 = myShapesOperations->MakeGlueFaces(Te7, aTolMax, true);
1700
1701     if (Te8.IsNull()) {
1702       SetErrorCode("Impossible to glue faces of TShape");
1703       return false;
1704     }
1705
1706     aShape = Te8->GetValue();
1707   }
1708
1709
1710   theShape->GetLastFunction()->SetValue(aShape);
1711
1712   Te4->GetLastFunction()->SetDescription("");
1713   Te5->GetLastFunction()->SetDescription("");
1714   Te6->GetLastFunction()->SetDescription("");
1715   Te7->GetLastFunction()->SetDescription("");
1716   Te8->GetLastFunction()->SetDescription("");
1717
1718   SetErrorCode(OK);
1719   return true;
1720 }
1721
1722 //=======================================================================
1723 //function : MakePipeTShapeThicknessReduction
1724 //purpose  : Static method. Add thiskness reduction elements at the three
1725 //                          open ends of the T-Shape.
1726 //=======================================================================
1727 TopoDS_Shape GEOMImpl_IAdvancedOperations::MakePipeTShapeThicknessReduction
1728                                   (TopoDS_Shape theShape,
1729                                    double r1, double w1, double l1,
1730                                    double r2, double w2, double l2,
1731                                    double rL, double wL, double ltransL, double lthinL,
1732                                    double rR, double wR, double ltransR, double lthinR,
1733                                    double rI, double wI, double ltransI, double lthinI,
1734                                    bool fuseReductions)
1735 {
1736   // Add thickness reduction elements
1737   // at the three extremities: Left, Right and Incident
1738   //
1739   // ---------------------.
1740   //   W                   \
1741   // ---------------------. \
1742   //   ^                   \ '-----------------.
1743   //   |R                   \        Wthin     |
1744   //   |                     '-----------------'
1745   //   v                             Rthin
1746   // --.--.--.--.--.--.--.--.--.--.--.--.--.--.--
1747   //                     Ltrans    Lthin
1748
1749   TopoDS_Shape aResult = theShape;
1750   double aTol = Precision::Confusion();
1751
1752   gp_Vec aVX = gp::DX(), aVZ = gp::DZ();
1753
1754   // Left reduction (rL, wL, ltransL, lthinL)
1755   if (rL > aTol && wL > aTol && ltransL > aTol) {
1756     gp_Pnt aPLeft (-l1, 0, 0);
1757     gp_Ax2 anAxesLeft (aPLeft, -aVX, aVZ);
1758     TopoDS_Shape aReductionLeft = GEOMImpl_IAdvancedOperations::MakeThicknessReduction
1759       (anAxesLeft, r1, w1, rL, wL, ltransL, lthinL, fuseReductions);
1760
1761     if (fuseReductions) {
1762       BRepAlgoAPI_Fuse fuseL (aResult, aReductionLeft);
1763       if (!fuseL.IsDone())
1764         StdFail_NotDone::Raise("Cannot fuse Te with left reduction");
1765       aResult = fuseL.Shape();
1766     }
1767     else {
1768       BRep_Builder B;
1769       TopoDS_Compound C;
1770       B.MakeCompound(C);
1771       B.Add(C, aResult);
1772       B.Add(C, aReductionLeft);
1773       aResult = C;
1774     }
1775   }
1776
1777   // Right reduction
1778   if (rR > aTol && wR > aTol && ltransR > aTol) {
1779     gp_Pnt aPRight (l1, 0, 0);
1780     gp_Ax2 anAxesRight (aPRight, aVX, aVZ);
1781     TopoDS_Shape aReductionRight = GEOMImpl_IAdvancedOperations::MakeThicknessReduction
1782       (anAxesRight, r1, w1, rR, wR, ltransR, lthinR, fuseReductions);
1783
1784     if (fuseReductions) {
1785       BRepAlgoAPI_Fuse fuseR (aResult, aReductionRight);
1786       if (!fuseR.IsDone())
1787         StdFail_NotDone::Raise("Cannot fuse Te with right reduction");
1788       aResult = fuseR.Shape();
1789     }
1790     else {
1791       BRep_Builder B;
1792       TopoDS_Compound C;
1793       B.MakeCompound(C);
1794       B.Add(C, aResult);
1795       B.Add(C, aReductionRight);
1796       aResult = C;
1797     }
1798   }
1799
1800   // Incident reduction
1801   if (rI > aTol && wI > aTol && ltransI > aTol) {
1802     gp_Pnt aPInci (0, 0, l2);
1803     gp_Ax2 anAxesInci (aPInci, aVZ, aVX);
1804     TopoDS_Shape aReductionInci = GEOMImpl_IAdvancedOperations::MakeThicknessReduction
1805       (anAxesInci, r2, w2, rI, wI, ltransI, lthinI, fuseReductions);
1806
1807     if (fuseReductions) {
1808       BRepAlgoAPI_Fuse fuseInci (aResult, aReductionInci);
1809       if (!fuseInci.IsDone())
1810         StdFail_NotDone::Raise("Cannot fuse Te with incident reduction");
1811       aResult = fuseInci.Shape();
1812     }
1813     else {
1814       BRep_Builder B;
1815       TopoDS_Compound C;
1816       B.MakeCompound(C);
1817       B.Add(C, aResult);
1818       B.Add(C, aReductionInci);
1819       aResult = C;
1820     }
1821   }
1822
1823   // Get rid of extra compounds
1824   TopTools_ListOfShape listShapeRes;
1825   GEOMUtils::AddSimpleShapes(aResult, listShapeRes);
1826   aResult = listShapeRes.First(); // useful for the case "fuseReductions == true"
1827
1828   if (!fuseReductions && listShapeRes.Extent() > 1) {
1829     // Simplify T-Shape compound (get rid of sub-compounds) and glue duplicated faces
1830     BRep_Builder B;
1831     TopoDS_Compound C;
1832     B.MakeCompound(C);
1833
1834     TopTools_ListIteratorOfListOfShape itSub (listShapeRes);
1835     for (; itSub.More(); itSub.Next())
1836       B.Add(C, itSub.Value());
1837
1838     // GlueFaces
1839     aResult = GEOMImpl_GlueDriver::GlueFaces(C, Precision::Confusion(), Standard_True);
1840   }
1841
1842   return aResult;
1843 }
1844
1845 //=======================================================================
1846 //function : MakeThicknessReduction
1847 //purpose  : Static method. Create one thickness reduction element.
1848 //=======================================================================
1849 TopoDS_Shape GEOMImpl_IAdvancedOperations::MakeThicknessReduction (gp_Ax2 theAxes,
1850                                                                    const double R, const double W,
1851                                                                    const double Rthin, const double Wthin,
1852                                                                    const double Ltrans, const double Lthin,
1853                                                                    bool fuse)
1854 {
1855   double aTol = Precision::Confusion();
1856   if (Rthin < aTol || Wthin < aTol || Ltrans < aTol) {
1857     StdFail_NotDone::Raise("Cannot build thickness reduction: too small values");
1858   }
1859   bool isThinPart = (Lthin > aTol);
1860
1861   //     .
1862   //   W |\
1863   //     . \
1864   //   ^  \ '-----------------.
1865   //   |R  \|                 | Wthin
1866   //   |    '-----------------'
1867   //   v                        Rthin
1868   // --.--.--.--.--.--.--.--.--.--.--.--.--> theAxes.Direction()
1869   //     Ltrans     Lthin
1870
1871   double RExt = R + W;
1872   double RthinExt = Rthin + Wthin;
1873
1874   gp_Dir aNormal = theAxes.Direction();
1875   gp_Dir anXDir  = theAxes.XDirection();
1876   gp_Pnt aPntCyl (theAxes.Location().XYZ() + aNormal.XYZ()*Ltrans);
1877   gp_Ax2 anAxesCyl (aPntCyl, aNormal, anXDir);
1878
1879   // Build the transition part
1880   BRepPrimAPI_MakeCone ConeExt (theAxes, RExt, RthinExt, Ltrans);
1881   BRepPrimAPI_MakeCone ConeInt (theAxes, R, Rthin, Ltrans);
1882   ConeExt.Build();
1883   ConeInt.Build();
1884   if (!ConeExt.IsDone() || !ConeInt.IsDone())
1885     StdFail_NotDone::Raise("Cannot build cones of thickness reduction");
1886   BRepAlgoAPI_Cut cut1 (ConeExt.Shape(), ConeInt.Shape());
1887   if (!cut1.IsDone())
1888     StdFail_NotDone::Raise("Coudn't build transition part of thickness reduction");
1889   TopoDS_Shape aReduction = cut1.Shape();
1890
1891   // Build the thin part, if required
1892   TopoDS_Shape aThinPart;
1893   if (isThinPart) {
1894     BRepPrimAPI_MakeCylinder CExt (anAxesCyl, RthinExt, Lthin);
1895     BRepPrimAPI_MakeCylinder CInt (anAxesCyl, Rthin, Lthin);
1896     CExt.Build();
1897     CInt.Build();
1898     if (!CExt.IsDone() || !CInt.IsDone())
1899       StdFail_NotDone::Raise("Cannot build cylinders of thickness reduction");
1900     BRepAlgoAPI_Cut cut2 (CExt.Shape(), CInt.Shape());
1901     if (!cut2.IsDone())
1902       StdFail_NotDone::Raise("Coudn't build thin part of thickness reduction");
1903     aThinPart = cut2.Shape();
1904   }
1905
1906   // Join parts
1907   if (fuse) {
1908     if (isThinPart) {
1909       BRepAlgoAPI_Fuse fuse1 (aReduction, aThinPart);
1910       if (!fuse1.IsDone())
1911         StdFail_NotDone::Raise("Cannot fuse parts of thickness reduction");
1912       aReduction = fuse1.Shape();
1913     }
1914   }
1915   else {
1916     // Partition the reduction on blocks
1917     gp_Ax3 anAxesPln1 (aPntCyl, theAxes.XDirection(), aNormal);
1918     gp_Ax3 anAxesPln2 (aPntCyl, theAxes.YDirection(), aNormal);
1919     gp_Pln aPln1 (anAxesPln1);
1920     gp_Pln aPln2 (anAxesPln2);
1921     double aSize = Ltrans + Lthin + R + Rthin + Wthin; // to guarantee enough size in all directions
1922     TopoDS_Shape aTool1 = BRepBuilderAPI_MakeFace(aPln1, -aSize, +aSize, -aSize, +aSize).Shape();
1923     TopoDS_Shape aTool2 = BRepBuilderAPI_MakeFace(aPln2, -aSize, +aSize, -aSize, +aSize).Shape();
1924
1925     GEOMAlgo_Splitter PS;
1926     PS.AddArgument(aReduction);
1927     if (isThinPart)
1928       PS.AddArgument(aThinPart);
1929     PS.AddTool(aTool1);
1930     PS.AddTool(aTool2);
1931     PS.SetLimit(TopAbs_SOLID);
1932     PS.Perform();
1933
1934     aReduction = PS.Shape();
1935   }
1936
1937   return aReduction;
1938 }
1939
1940 //=============================================================================
1941 /*!
1942  *  MakePipeTShape
1943  *  \brief Create a T-shape object with specified caracteristics for the main and
1944  *         the incident pipes (radius, width, half-length).
1945  *         Center of the shape is (0,0,0). The main plane of the T-shape is XOY.
1946  *  \param theR1 Internal radius of main pipe
1947  *  \param theW1 Width of main pipe
1948  *  \param theL1 Half-length of main pipe
1949  *  \param theR2 Internal radius of incident pipe (R2 < R1)
1950  *  \param theW2 Width of incident pipe (R2+W2 < R1+W1)
1951  *  \param theL2 Half-length of incident pipe
1952  *  \param theHexMesh Boolean indicating if shape is prepared for hex mesh
1953  *  \return List of GEOM_Objects, containing the created shape and propagation groups.
1954  */
1955 //=============================================================================
1956 Handle(TColStd_HSequenceOfTransient)
1957   GEOMImpl_IAdvancedOperations::MakePipeTShape(double theR1, double theW1, double theL1,
1958                                                double theR2, double theW2, double theL2,
1959                                                double theRL, double theWL, double theLtransL, double theLthinL,
1960                                                double theRR, double theWR, double theLtransR, double theLthinR,
1961                                                double theRI, double theWI, double theLtransI, double theLthinI,
1962                                                bool theHexMesh)
1963 {
1964   MESSAGE("GEOMImpl_IAdvancedOperations::MakePipeTShape");
1965   SetErrorCode(KO);
1966   //Add a new object
1967   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_TSHAPE);
1968
1969   //Add a new shape function with parameters
1970   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_PipeTShapeDriver::GetID(), TSHAPE_BASIC);
1971   if (aFunction.IsNull()) return NULL;
1972
1973   //Check if the function is set correctly
1974   if (aFunction->GetDriverGUID() != GEOMImpl_PipeTShapeDriver::GetID()) return NULL;
1975
1976   GEOMImpl_IPipeTShape aData (aFunction);
1977
1978   aData.SetR1(theR1);
1979   aData.SetW1(theW1);
1980   aData.SetL1(theL1);
1981   aData.SetR2(theR2);
1982   aData.SetW2(theW2);
1983   aData.SetL2(theL2);
1984   aData.SetHexMesh(theHexMesh);
1985
1986   bool isTRL = (theRL + theWL + theLtransL + theLthinL) > Precision::Confusion();
1987   bool isTRR = (theRR + theWR + theLtransR + theLthinR) > Precision::Confusion();
1988   bool isTRI = (theRI + theWI + theLtransI + theLthinI) > Precision::Confusion();
1989
1990   //Compute the resulting value
1991   try {
1992     OCC_CATCH_SIGNALS;
1993     if (!GetSolver()->ComputeFunction(aFunction)) {
1994       SetErrorCode("TShape driver failed");
1995       return NULL;
1996     }
1997
1998     if (theHexMesh) {
1999       if (!MakePipeTShapePartition(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
2000         return NULL;
2001       if (!MakePipeTShapeMirrorAndGlue(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
2002         return NULL;
2003     }
2004
2005     if (isTRL || isTRR || isTRI) {
2006       // Add thickness reduction elements
2007       // at the three extremities: Left, Right and Incident
2008       TopoDS_Shape aResShape =
2009         MakePipeTShapeThicknessReduction(aShape->GetValue(), theR1, theW1, theL1, theR2, theW2, theL2,
2010                                          theRL, theWL, theLtransL, theLthinL,
2011                                          theRR, theWR, theLtransR, theLthinR,
2012                                          theRI, theWI, theLtransI, theLthinI,
2013                                          !theHexMesh);
2014       aFunction->SetValue(aResShape);
2015     }
2016   }
2017   catch (Standard_Failure) {
2018     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2019     SetErrorCode(aFail->GetMessageString());
2020     return NULL;
2021   }
2022
2023   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
2024   aSeq->Append(aShape);
2025
2026   try {
2027     if (theHexMesh) {
2028       // Get the groups
2029       if (!MakeGroups(aShape, TSHAPE_BASIC, theR1, theW1, theL1, theR2, theW2, theL2,
2030                       0., 0., 0., aSeq, gp_Trsf()))
2031         return NULL;
2032     }
2033
2034     // Get internal group.
2035     if (!MakeInternalGroup(aShape, theR1, theL1, theR2, theL2, theRL, theLtransL,
2036                            theRR, theLtransR, theRI, theLtransI,
2037                            aSeq, gp_Trsf())) {
2038       return NULL;
2039     }
2040   }
2041   catch (Standard_Failure) {
2042     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2043     SetErrorCode(aFail->GetMessageString());
2044     return NULL;
2045   }
2046
2047   //Make a Python command
2048   TCollection_AsciiString anEntry, aListRes("[");
2049   // Iterate over the sequence aSeq
2050   Standard_Integer aNbGroups = aSeq->Length();
2051   Standard_Integer i = 1;
2052   for (; i <= aNbGroups; i++) {
2053     Handle(Standard_Transient) anItem = aSeq->Value(i);
2054     if (anItem.IsNull()) continue;
2055     Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(anItem);
2056     if (aGroup.IsNull()) continue;
2057     //Make a Python command
2058     TDF_Tool::Entry(aGroup->GetEntry(), anEntry);
2059     aListRes += anEntry + ", ";
2060   }
2061   aListRes.Trunc(aListRes.Length() - 2);
2062
2063   GEOM::TPythonDump pd (aFunction);
2064
2065   pd << aListRes.ToCString() << "] = geompy.MakePipeTShape("
2066      << theR1 << ", " << theW1 << ", " << theL1 << ", "
2067      << theR2 << ", " << theW2 << ", " << theL2 << ", "
2068      << theHexMesh;
2069
2070   // thickness reduction
2071   if (isTRL)
2072     pd << ", theRL=" << theRL << ", theWL=" << theWL
2073        << ", theLtransL=" << theLtransL << ", theLthinL=" << theLthinL;
2074   if (isTRR)
2075     pd << ", theRR=" << theRR << ", theWR=" << theWR
2076        << ", theLtransR=" << theLtransR << ", theLthinR=" << theLthinR;
2077   if (isTRI)
2078     pd << ", theRI=" << theRI << ", theWI=" << theWI
2079        << ", theLtransI=" << theLtransI << ", theLthinI=" << theLthinI;
2080
2081   pd << ")";
2082
2083   SetErrorCode(OK);
2084
2085   return aSeq;
2086 }
2087
2088 //=============================================================================
2089 /*!
2090  *  MakePipeTShapeWithPosition
2091  *  Create a T-shape object with specified caracteristics for the main and
2092  *  the incident pipes (radius, width, half-length).
2093  *  The extremities of the main pipe are located on junctions points P1 and P2.
2094  *  The extremity of the incident pipe is located on junction point P3.
2095  *  \param theR1 Internal radius of main pipe
2096  *  \param theW1 Width of main pipe
2097  *  \param theL1 Half-length of main pipe
2098  *  \param theR2 Internal radius of incident pipe (R2 < R1)
2099  *  \param theW2 Width of incident pipe (R2+W2 < R1+W1)
2100  *  \param theL2 Half-length of incident pipe
2101  *  \param theHexMesh Boolean indicating if shape is prepared for hex mesh
2102  *  \param theP1 1st junction point of main pipe
2103  *  \param theP2 2nd junction point of main pipe
2104  *  \param theP3 Junction point of incident pipe
2105  *  \return List of GEOM_Objects, containing the created shape and propagation groups..
2106  */
2107 //=============================================================================
2108 Handle(TColStd_HSequenceOfTransient)
2109 GEOMImpl_IAdvancedOperations::MakePipeTShapeWithPosition
2110                              (double theR1, double theW1, double theL1,
2111                               double theR2, double theW2, double theL2,
2112                               double theRL, double theWL, double theLtransL, double theLthinL,
2113                               double theRR, double theWR, double theLtransR, double theLthinR,
2114                               double theRI, double theWI, double theLtransI, double theLthinI,
2115                               bool theHexMesh,
2116                               Handle(GEOM_Object) theP1,
2117                               Handle(GEOM_Object) theP2,
2118                               Handle(GEOM_Object) theP3)
2119 {
2120   SetErrorCode(KO);
2121   //Add a new object
2122   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_TSHAPE);
2123   /////////////////
2124   // TSHAPE CODE
2125   /////////////////
2126   //Add a new shape function with parameters
2127   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_PipeTShapeDriver::GetID(), TSHAPE_BASIC);
2128   if (aFunction.IsNull()) return NULL;
2129
2130   //Check if the function is set correctly
2131   if (aFunction->GetDriverGUID() != GEOMImpl_PipeTShapeDriver::GetID()) return NULL;
2132
2133   // Check new position
2134   if (!CheckCompatiblePosition(theL1, theL2, theP1, theP2, theP3, 0.01)) {
2135     return NULL;
2136   }
2137
2138   GEOMImpl_IPipeTShape aData(aFunction);
2139
2140   aData.SetR1(theR1);
2141   aData.SetW1(theW1);
2142   aData.SetL1(theL1);
2143   aData.SetR2(theR2);
2144   aData.SetW2(theW2);
2145   aData.SetL2(theL2);
2146   aData.SetHexMesh(theHexMesh);
2147
2148   bool isTRL = (theRL + theWL + theLtransL + theLthinL) > Precision::Confusion();
2149   bool isTRR = (theRR + theWR + theLtransR + theLthinR) > Precision::Confusion();
2150   bool isTRI = (theRI + theWI + theLtransI + theLthinI) > Precision::Confusion();
2151
2152   //Compute the resulting value
2153   try {
2154     OCC_CATCH_SIGNALS;
2155     if (!GetSolver()->ComputeFunction(aFunction)) {
2156       SetErrorCode("TShape driver failed");
2157       return NULL;
2158     }
2159
2160     if (theHexMesh) {
2161       if (!MakePipeTShapePartition(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
2162         return NULL;
2163       if (!MakePipeTShapeMirrorAndGlue(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
2164         return NULL;
2165     }
2166
2167     if (isTRL || isTRR || isTRI) {
2168       // Add thickness reduction elements
2169       // at the three extremities: Left, Right and Incident
2170       TopoDS_Shape aResShape =
2171         MakePipeTShapeThicknessReduction(aShape->GetValue(), theR1, theW1, theL1, theR2, theW2, theL2,
2172                                          theRL, theWL, theLtransL, theLthinL,
2173                                          theRR, theWR, theLtransR, theLthinR,
2174                                          theRI, theWI, theLtransI, theLthinI,
2175                                          !theHexMesh);
2176       aFunction->SetValue(aResShape);
2177     }
2178   }
2179   catch (Standard_Failure) {
2180     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2181     SetErrorCode(aFail->GetMessageString());
2182     return NULL;
2183   }
2184
2185   TopoDS_Shape Te = aShape->GetValue();
2186
2187   // Set Position
2188   gp_Trsf aTrsf = GetPositionTrsf(theL1, theL2, theP1, theP2, theP3);
2189   BRepBuilderAPI_Transform aTransformation(Te, aTrsf, Standard_False);
2190   TopoDS_Shape aTrsf_Shape = aTransformation.Shape();
2191   aFunction->SetValue(aTrsf_Shape);
2192
2193   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
2194   aSeq->Append(aShape);
2195
2196   try {
2197     if (theHexMesh) {
2198       // Get the groups
2199       if (!MakeGroups(aShape,TSHAPE_BASIC, theR1, theW1, theL1, theR2, theW2, theL2,
2200                       0., 0., 0., aSeq, aTrsf)) {
2201         return NULL;
2202       }
2203     }
2204
2205     // Get internal group.
2206     if (!MakeInternalGroup(aShape, theR1, theL1, theR2, theL2, theRL, theLtransL,
2207                            theRR, theLtransR, theRI, theLtransI,
2208                            aSeq, aTrsf)) {
2209       return NULL;
2210     }
2211   }
2212   catch (Standard_Failure) {
2213     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2214     SetErrorCode(aFail->GetMessageString());
2215     return NULL;
2216   }
2217
2218   //Make a Python command
2219   TCollection_AsciiString anEntry, aListRes("[");
2220   // Iterate over the sequence aSeq
2221   Standard_Integer aNbGroups = aSeq->Length();
2222   Standard_Integer i = 1;
2223   for (; i <= aNbGroups; i++) {
2224     Handle(Standard_Transient) anItem = aSeq->Value(i);
2225     if (anItem.IsNull()) continue;
2226     Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(anItem);
2227     if (aGroup.IsNull()) continue;
2228     //Make a Python command
2229     TDF_Tool::Entry(aGroup->GetEntry(), anEntry);
2230     aListRes += anEntry + ", ";
2231   }
2232   aListRes.Trunc(aListRes.Length() - 2);
2233
2234   GEOM::TPythonDump pd (aFunction);
2235
2236   pd << aListRes.ToCString() << "] = geompy.MakePipeTShape("
2237      << theR1 << ", " << theW1 << ", " << theL1 << ", "
2238      << theR2 << ", " << theW2 << ", " << theL2 << ", "
2239      << theHexMesh << ", " << theP1 << ", " << theP2 << ", " << theP3;
2240
2241   // thickness reduction
2242   if (isTRL)
2243     pd << ", theRL=" << theRL << ", theWL=" << theWL
2244        << ", theLtransL=" << theLtransL << ", theLthinL=" << theLthinL;
2245   if (isTRR)
2246     pd << ", theRR=" << theRR << ", theWR=" << theWR
2247        << ", theLtransR=" << theLtransR << ", theLthinR=" << theLthinR;
2248   if (isTRI)
2249     pd << ", theRI=" << theRI << ", theWI=" << theWI
2250        << ", theLtransI=" << theLtransI << ", theLthinI=" << theLthinI;
2251
2252   pd << ")";
2253
2254   SetErrorCode(OK);
2255
2256   return aSeq;
2257 }
2258
2259 //=============================================================================
2260 /*!
2261  *  MakePipeTShapeChamfer
2262  *  Create a T-shape object with specified caracteristics for the main and
2263  *  the incident pipes (radius, width, half-length). A chamfer is created
2264  *  on the junction of the pipes.
2265  *  Center of the shape is (0,0,0). The main plane of the T-shape is XOY.
2266  *  \param theR1 Internal radius of main pipe
2267  *  \param theW1 Width of main pipe
2268  *  \param theL1 Half-length of main pipe
2269  *  \param theR2 Internal radius of incident pipe (R2 < R1)
2270  *  \param theW2 Width of incident pipe (R2+W2 < R1+W1)
2271  *  \param theL2 Half-length of incident pipe
2272  *  \param theH Height of chamfer.
2273  *  \param theW Width of chamfer.
2274  *  \param theHexMesh Boolean indicating if shape is prepared for hex mesh
2275  *  \return List of GEOM_Objects, containing the created shape and propagation groups.
2276  */
2277 //=============================================================================
2278 Handle(TColStd_HSequenceOfTransient)
2279 GEOMImpl_IAdvancedOperations::MakePipeTShapeChamfer
2280                              (double theR1, double theW1, double theL1,
2281                               double theR2, double theW2, double theL2,
2282                               double theRL, double theWL, double theLtransL, double theLthinL,
2283                               double theRR, double theWR, double theLtransR, double theLthinR,
2284                               double theRI, double theWI, double theLtransI, double theLthinI,
2285                               double theH, double theW,
2286                               bool theHexMesh)
2287 {
2288   SetErrorCode(KO);
2289   //Add a new object
2290   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_TSHAPE);
2291   //Add a new shape function with parameters
2292   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_PipeTShapeDriver::GetID(), TSHAPE_CHAMFER);
2293   if (aFunction.IsNull()) return NULL;
2294
2295   //Check if the function is set correctly
2296   if (aFunction->GetDriverGUID() != GEOMImpl_PipeTShapeDriver::GetID()) return NULL;
2297
2298   GEOMImpl_IPipeTShape aData(aFunction);
2299
2300   aData.SetR1(theR1);
2301   aData.SetW1(theW1);
2302   aData.SetL1(theL1);
2303   aData.SetR2(theR2);
2304   aData.SetW2(theW2);
2305   aData.SetL2(theL2);
2306   aData.SetH(theH);
2307   aData.SetW(theW);
2308   aData.SetHexMesh(theHexMesh);
2309
2310   bool isTRL = (theRL + theWL + theLtransL + theLthinL) > Precision::Confusion();
2311   bool isTRR = (theRR + theWR + theLtransR + theLthinR) > Precision::Confusion();
2312   bool isTRI = (theRI + theWI + theLtransI + theLthinI) > Precision::Confusion();
2313
2314   //Compute the resulting value
2315   try {
2316     OCC_CATCH_SIGNALS;
2317     if (!GetSolver()->ComputeFunction(aFunction)) {
2318       SetErrorCode("TShape driver failed");
2319       return NULL;
2320     }
2321   }
2322   catch (Standard_Failure) {
2323     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2324     SetErrorCode(aFail->GetMessageString());
2325     return NULL;
2326   }
2327
2328   // BEGIN of chamfer
2329   TopoDS_Shape aShapeShape = aShape->GetValue();
2330   TopTools_IndexedMapOfShape anEdgesIndices;
2331   TopExp::MapShapes(aShapeShape, anEdgesIndices);
2332   // Common edges on external cylinders
2333   Handle(GEOM_Object) box_e;
2334   if (theHexMesh) {
2335     box_e = my3DPrimOperations->MakeBoxDXDYDZ(theR2+theW2, theR2+theW2, theR1+theW1);
2336   }
2337   else {
2338     box_e = my3DPrimOperations->MakeBoxDXDYDZ(2*(theR2+theW2), 2*(theR2+theW2), theR1+theW1);
2339   }
2340   box_e->GetLastFunction()->SetDescription("");
2341   box_e = myTransformOperations->TranslateDXDYDZ(box_e, -(theR2+theW2), -(theR2+theW2), 0);
2342   box_e->GetLastFunction()->SetDescription("");
2343
2344   Handle(TColStd_HSequenceOfInteger) edges_e =
2345     myShapesOperations->GetShapesOnBoxIDs(box_e, aShape, TopAbs_EDGE, GEOMAlgo_ST_IN);
2346   box_e->GetLastFunction()->SetDescription("");
2347
2348   if (edges_e.IsNull() || edges_e->Length() == 0) {
2349     SetErrorCode("External edges not found");
2350     return NULL;
2351   }
2352   int nbEdgesInChamfer = 0;
2353   std::list<int> theEdges;
2354   for (int i=1; i<=edges_e->Length();i++) {
2355     int edgeID = edges_e->Value(i);
2356     TopoDS_Shape theEdge = anEdgesIndices.FindKey(edgeID);
2357     TopExp_Explorer Ex(theEdge,TopAbs_VERTEX);
2358     int iv=0;
2359     while (Ex.More()) {
2360       iv ++;
2361       gp_Pnt aPt = BRep_Tool::Pnt(TopoDS::Vertex(Ex.Current()));
2362       if (Abs(aPt.Z() - (theR1+theW1)) <= Precision::Confusion()) {
2363         nbEdgesInChamfer ++;
2364         theEdges.push_back(edgeID);
2365       }
2366       Ex.Next();
2367     }
2368     if (theHexMesh && nbEdgesInChamfer == 1)
2369       break;
2370   }
2371   Handle(GEOM_Object) aChamfer;
2372   try {
2373     aChamfer = myLocalOperations->MakeChamferEdges(aShape, theW, theH, theEdges);
2374   }
2375   catch (Standard_Failure) {
2376     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2377     SetErrorCode(aFail->GetMessageString());
2378     return NULL;
2379   }
2380   if (aChamfer.IsNull()) {
2381     SetErrorCode("Chamfer can not be computed on the given shape with the given parameters");
2382     return NULL;
2383   }
2384   aChamfer->GetLastFunction()->SetDescription("");
2385
2386   TopoDS_Shape aChamferShape = aChamfer->GetValue();
2387   aFunction->SetValue(aChamferShape);
2388   // END of chamfer
2389
2390   if (theHexMesh) {
2391     if (!MakePipeTShapePartition(aShape, theR1, theW1, theL1, theR2, theW2, theL2, theH, theW, 0, false))
2392       return NULL;
2393     if (!MakePipeTShapeMirrorAndGlue(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
2394       return NULL;
2395   }
2396
2397   // Add thickness reduction elements
2398   // at the three extremities: Left, Right and Incident
2399   try {
2400     OCC_CATCH_SIGNALS;
2401     if (isTRL || isTRR || isTRI) {
2402       TopoDS_Shape aResShape =
2403         MakePipeTShapeThicknessReduction(aShape->GetValue(), theR1, theW1, theL1, theR2, theW2, theL2,
2404                                          theRL, theWL, theLtransL, theLthinL,
2405                                          theRR, theWR, theLtransR, theLthinR,
2406                                          theRI, theWI, theLtransI, theLthinI,
2407                                          !theHexMesh);
2408       aFunction->SetValue(aResShape);
2409     }
2410   }
2411   catch (Standard_Failure) {
2412     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2413     SetErrorCode(aFail->GetMessageString());
2414     return NULL;
2415   }
2416
2417   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
2418   aSeq->Append(aShape);
2419
2420   try {
2421     if (theHexMesh) {
2422       // Get the groups
2423       if (!MakeGroups(aShape, TSHAPE_CHAMFER, theR1, theW1, theL1, theR2, theW2, theL2,
2424                       theH, theW, 0., aSeq, gp_Trsf()))
2425         return NULL;
2426     }
2427
2428     // Get internal group.
2429     if (!MakeInternalGroup(aShape, theR1, theL1, theR2, theL2, theRL, theLtransL,
2430                            theRR, theLtransR, theRI, theLtransI,
2431                            aSeq, gp_Trsf())) {
2432       return NULL;
2433     }
2434   }
2435   catch (Standard_Failure) {
2436     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2437     SetErrorCode(aFail->GetMessageString());
2438     return NULL;
2439   }
2440
2441   //Make a Python command
2442   TCollection_AsciiString anEntry, aListRes("[");
2443   // Iterate over the sequence aSeq
2444   Standard_Integer aNbGroups = aSeq->Length();
2445   Standard_Integer i = 1;
2446   for (; i <= aNbGroups; i++) {
2447     Handle(Standard_Transient) anItem = aSeq->Value(i);
2448     if (anItem.IsNull()) continue;
2449     Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(anItem);
2450     if (aGroup.IsNull()) continue;
2451     //Make a Python command
2452     TDF_Tool::Entry(aGroup->GetEntry(), anEntry);
2453     aListRes += anEntry + ", ";
2454   }
2455   aListRes.Trunc(aListRes.Length() - 2);
2456
2457   GEOM::TPythonDump pd (aFunction);
2458
2459   pd << aListRes.ToCString() << "] = geompy.MakePipeTShapeChamfer("
2460      << theR1 << ", " << theW1 << ", " << theL1 << ", "
2461      << theR2 << ", " << theW2 << ", " << theL2 << ", "
2462      << theH << ", " << theW << ", " << theHexMesh;
2463
2464   // thickness reduction
2465   if (isTRL)
2466     pd << ", theRL=" << theRL << ", theWL=" << theWL
2467        << ", theLtransL=" << theLtransL << ", theLthinL=" << theLthinL;
2468   if (isTRR)
2469     pd << ", theRR=" << theRR << ", theWR=" << theWR
2470        << ", theLtransR=" << theLtransR << ", theLthinR=" << theLthinR;
2471   if (isTRI)
2472     pd << ", theRI=" << theRI << ", theWI=" << theWI
2473        << ", theLtransI=" << theLtransI << ", theLthinI=" << theLthinI;
2474
2475   pd << ")";
2476
2477   SetErrorCode(OK);
2478
2479   return aSeq;
2480 }
2481
2482 //=============================================================================
2483 /*!
2484  *  MakePipeTShapeChamferWithPosition
2485  *  Create a T-shape object with specified caracteristics for the main and
2486  *  the incident pipes (radius, width, half-length). A chamfer is created
2487  *  on the junction of the pipes.
2488  *  The extremities of the main pipe are located on junctions points P1 and P2.
2489  *  The extremity of the incident pipe is located on junction point P3.
2490  *  \param theR1 Internal radius of main pipe
2491  *  \param theW1 Width of main pipe
2492  *  \param theL1 Half-length of main pipe
2493  *  \param theR2 Internal radius of incident pipe (R2 < R1)
2494  *  \param theW2 Width of incident pipe (R2+W2 < R1+W1)
2495  *  \param theL2 Half-length of incident pipe
2496  *  \param theH Height of chamfer.
2497  *  \param theW Width of chamfer.
2498  *  \param theHexMesh Boolean indicating if shape is prepared for hex mesh
2499  *  \param theP1 1st junction point of main pipe
2500  *  \param theP2 2nd junction point of main pipe
2501  *  \param theP3 Junction point of incident pipe
2502  *  \return List of GEOM_Objects, containing the created shape and propagation groups.
2503  */
2504 //=============================================================================
2505 Handle(TColStd_HSequenceOfTransient)
2506 GEOMImpl_IAdvancedOperations::MakePipeTShapeChamferWithPosition
2507                              (double theR1, double theW1, double theL1,
2508                               double theR2, double theW2, double theL2,
2509                               double theRL, double theWL, double theLtransL, double theLthinL,
2510                               double theRR, double theWR, double theLtransR, double theLthinR,
2511                               double theRI, double theWI, double theLtransI, double theLthinI,
2512                               double theH, double theW,
2513                               bool theHexMesh,
2514                               Handle(GEOM_Object) theP1,
2515                               Handle(GEOM_Object) theP2,
2516                               Handle(GEOM_Object) theP3)
2517 {
2518   SetErrorCode(KO);
2519   //Add a new object
2520   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_TSHAPE);
2521   //Add a new shape function with parameters
2522   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_PipeTShapeDriver::GetID(), TSHAPE_CHAMFER);
2523   if (aFunction.IsNull()) return NULL;
2524
2525   //Check if the function is set correctly
2526   if (aFunction->GetDriverGUID() != GEOMImpl_PipeTShapeDriver::GetID()) return NULL;
2527
2528   // Check new position
2529   if (!CheckCompatiblePosition(theL1, theL2, theP1, theP2, theP3, 0.01)) {
2530     return NULL;
2531   }
2532
2533   GEOMImpl_IPipeTShape aData(aFunction);
2534
2535   aData.SetR1(theR1);
2536   aData.SetW1(theW1);
2537   aData.SetL1(theL1);
2538   aData.SetR2(theR2);
2539   aData.SetW2(theW2);
2540   aData.SetL2(theL2);
2541   aData.SetH(theH);
2542   aData.SetW(theW);
2543   aData.SetHexMesh(theHexMesh);
2544
2545   bool isTRL = (theRL + theWL + theLtransL + theLthinL) > Precision::Confusion();
2546   bool isTRR = (theRR + theWR + theLtransR + theLthinR) > Precision::Confusion();
2547   bool isTRI = (theRI + theWI + theLtransI + theLthinI) > Precision::Confusion();
2548
2549   //Compute the resulting value
2550   try {
2551     OCC_CATCH_SIGNALS;
2552     if (!GetSolver()->ComputeFunction(aFunction)) {
2553       SetErrorCode("TShape driver failed");
2554       return NULL;
2555     }
2556   }
2557   catch (Standard_Failure) {
2558     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2559     SetErrorCode(aFail->GetMessageString());
2560     return NULL;
2561   }
2562
2563   // BEGIN of chamfer
2564   TopoDS_Shape aShapeShape = aShape->GetValue();
2565   TopTools_IndexedMapOfShape anEdgesIndices;
2566   TopExp::MapShapes(aShapeShape, anEdgesIndices);
2567   // Common edges on external cylinders
2568   Handle(GEOM_Object) box_e;
2569   if (theHexMesh) {
2570     box_e = my3DPrimOperations->MakeBoxDXDYDZ(theR2+theW2, theR2+theW2, theR1+theW1);
2571   }
2572   else {
2573     box_e = my3DPrimOperations->MakeBoxDXDYDZ(2*(theR2+theW2), 2*(theR2+theW2), theR1+theW1);
2574   }
2575   box_e->GetLastFunction()->SetDescription("");
2576   box_e = myTransformOperations->TranslateDXDYDZ(box_e, -(theR2+theW2), -(theR2+theW2), 0);
2577   box_e->GetLastFunction()->SetDescription("");
2578
2579   Handle(TColStd_HSequenceOfInteger) edges_e =
2580     myShapesOperations->GetShapesOnBoxIDs(box_e, aShape, TopAbs_EDGE, GEOMAlgo_ST_IN);
2581   box_e->GetLastFunction()->SetDescription("");
2582
2583   if (edges_e.IsNull() || edges_e->Length() == 0) {
2584     SetErrorCode("External edges not found");
2585     return NULL;
2586   }
2587   int nbEdgesInChamfer = 0;
2588   std::list<int> theEdges;
2589   for (int i=1; i<=edges_e->Length();i++) {
2590     int edgeID = edges_e->Value(i);
2591     TopoDS_Shape theEdge = anEdgesIndices.FindKey(edgeID);
2592     TopExp_Explorer Ex(theEdge,TopAbs_VERTEX);
2593     while (Ex.More()) {
2594       gp_Pnt aPt = BRep_Tool::Pnt(TopoDS::Vertex(Ex.Current()));
2595       if (Abs(aPt.Z() - (theR1+theW1)) <= Precision::Confusion()) {
2596         nbEdgesInChamfer ++;
2597         theEdges.push_back(edgeID);
2598       }
2599       Ex.Next();
2600     }
2601     if (theHexMesh && nbEdgesInChamfer == 1)
2602       break;
2603   }
2604   Handle(GEOM_Object) aChamfer;
2605   try {
2606     aChamfer = myLocalOperations->MakeChamferEdges(aShape, theW, theH, theEdges);
2607   }
2608   catch (Standard_Failure) {
2609     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2610     SetErrorCode(aFail->GetMessageString());
2611     return NULL;
2612   }
2613   if (aChamfer.IsNull()) {
2614     SetErrorCode("Chamfer can not be computed on the given shape with the given parameters");
2615     return NULL;
2616   }
2617   aChamfer->GetLastFunction()->SetDescription("");
2618
2619   TopoDS_Shape aChamferShape = aChamfer->GetValue();
2620   aFunction->SetValue(aChamferShape);
2621   // END of chamfer
2622
2623   if (theHexMesh) {
2624     if (!MakePipeTShapePartition(aShape, theR1, theW1, theL1, theR2, theW2, theL2, theH, theW, 0, false))
2625       return NULL;
2626     if (!MakePipeTShapeMirrorAndGlue(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
2627       return NULL;
2628   }
2629
2630   // Add thickness reduction elements
2631   // at the three extremities: Left, Right and Incident
2632   try {
2633     OCC_CATCH_SIGNALS;
2634     if (isTRL || isTRR || isTRI) {
2635       TopoDS_Shape aResShape =
2636         MakePipeTShapeThicknessReduction(aShape->GetValue(), theR1, theW1, theL1, theR2, theW2, theL2,
2637                                          theRL, theWL, theLtransL, theLthinL,
2638                                          theRR, theWR, theLtransR, theLthinR,
2639                                          theRI, theWI, theLtransI, theLthinI,
2640                                          !theHexMesh);
2641       aFunction->SetValue(aResShape);
2642     }
2643   }
2644   catch (Standard_Failure) {
2645     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2646     SetErrorCode(aFail->GetMessageString());
2647     return NULL;
2648   }
2649
2650   // Set Position
2651   gp_Trsf aTrsf = GetPositionTrsf(theL1, theL2, theP1, theP2, theP3);
2652   BRepBuilderAPI_Transform aTransformation (aShape->GetValue(), aTrsf, Standard_False);
2653   TopoDS_Shape aTrsf_Shape = aTransformation.Shape();
2654   aFunction->SetValue(aTrsf_Shape);
2655
2656   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
2657   aSeq->Append(aShape);
2658
2659   try {
2660     if (theHexMesh) {
2661       // Get the groups
2662       if (!MakeGroups(aShape, TSHAPE_CHAMFER, theR1, theW1, theL1, theR2, theW2, theL2,
2663                       theH, theW, 0., aSeq, aTrsf))
2664         return NULL;
2665     }
2666
2667     // Get internal group.
2668     if (!MakeInternalGroup(aShape, theR1, theL1, theR2, theL2, theRL, theLtransL,
2669                            theRR, theLtransR, theRI, theLtransI,
2670                            aSeq, aTrsf)) {
2671       return NULL;
2672     }
2673   }
2674   catch (Standard_Failure) {
2675     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2676     SetErrorCode(aFail->GetMessageString());
2677     return NULL;
2678   }
2679
2680   //Make a Python command
2681   TCollection_AsciiString anEntry, aListRes("[");
2682   // Iterate over the sequence aSeq
2683   Standard_Integer aNbGroups = aSeq->Length();
2684   Standard_Integer i = 1;
2685   for (; i <= aNbGroups; i++) {
2686     Handle(Standard_Transient) anItem = aSeq->Value(i);
2687     if (anItem.IsNull()) continue;
2688     Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(anItem);
2689     if (aGroup.IsNull()) continue;
2690     //Make a Python command
2691     TDF_Tool::Entry(aGroup->GetEntry(), anEntry);
2692     aListRes += anEntry + ", ";
2693   }
2694   aListRes.Trunc(aListRes.Length() - 2);
2695
2696   GEOM::TPythonDump pd (aFunction);
2697
2698   pd << aListRes.ToCString() << "] = geompy.MakePipeTShapeChamfer("
2699      << theR1 << ", " << theW1 << ", " << theL1 << ", "
2700      << theR2 << ", " << theW2 << ", " << theL2 << ", "
2701      << theH << ", " << theW << ", " << theHexMesh << ", "
2702      << theP1 << ", " << theP2 << ", " << theP3;
2703
2704   // thickness reduction
2705   if (isTRL)
2706     pd << ", theRL=" << theRL << ", theWL=" << theWL
2707        << ", theLtransL=" << theLtransL << ", theLthinL=" << theLthinL;
2708   if (isTRR)
2709     pd << ", theRR=" << theRR << ", theWR=" << theWR
2710        << ", theLtransR=" << theLtransR << ", theLthinR=" << theLthinR;
2711   if (isTRI)
2712     pd << ", theRI=" << theRI << ", theWI=" << theWI
2713        << ", theLtransI=" << theLtransI << ", theLthinI=" << theLthinI;
2714
2715   pd << ")";
2716
2717   SetErrorCode(OK);
2718
2719   return aSeq;
2720 }
2721
2722 //=============================================================================
2723 /*!
2724  *  MakePipeTShapeFillet
2725  *  Create a T-shape object with specified caracteristics for the main and
2726  *  the incident pipes (radius, width, half-length). A fillet is created
2727  *  on the junction of the pipes.
2728  *  Center of the shape is (0,0,0). The main plane of the T-shape is XOY.
2729  *  \param theR1 Internal radius of main pipe
2730  *  \param theW1 Width of main pipe
2731  *  \param theL1 Half-length of main pipe
2732  *  \param theR2 Internal radius of incident pipe (R2 < R1)
2733  *  \param theW2 Width of incident pipe (R2+W2 < R1+W1)
2734  *  \param theL2 Half-length of incident pipe
2735  *  \param theRF Radius of curvature of fillet.
2736  *  \param theHexMesh Boolean indicating if shape is prepared for hex mesh
2737  *  \return List of GEOM_Objects, containing the created shape and propagation groups.
2738  */
2739 //=============================================================================
2740 Handle(TColStd_HSequenceOfTransient)
2741 GEOMImpl_IAdvancedOperations::MakePipeTShapeFillet
2742                              (double theR1, double theW1, double theL1,
2743                               double theR2, double theW2, double theL2,
2744                               double theRL, double theWL, double theLtransL, double theLthinL,
2745                               double theRR, double theWR, double theLtransR, double theLthinR,
2746                               double theRI, double theWI, double theLtransI, double theLthinI,
2747                               double theRF, bool theHexMesh)
2748 {
2749   SetErrorCode(KO);
2750   //Add a new object
2751   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_TSHAPE);
2752   //Add a new shape function with parameters
2753   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_PipeTShapeDriver::GetID(), TSHAPE_FILLET);
2754   if (aFunction.IsNull()) return NULL;
2755
2756   //Check if the function is set correctly
2757   if (aFunction->GetDriverGUID() != GEOMImpl_PipeTShapeDriver::GetID()) return NULL;
2758
2759   GEOMImpl_IPipeTShape aData(aFunction);
2760
2761   aData.SetR1(theR1);
2762   aData.SetW1(theW1);
2763   aData.SetL1(theL1);
2764   aData.SetR2(theR2);
2765   aData.SetW2(theW2);
2766   aData.SetL2(theL2);
2767   aData.SetRF(theRF);
2768   aData.SetHexMesh(theHexMesh);
2769
2770   bool isTRL = (theRL + theWL + theLtransL + theLthinL) > Precision::Confusion();
2771   bool isTRR = (theRR + theWR + theLtransR + theLthinR) > Precision::Confusion();
2772   bool isTRI = (theRI + theWI + theLtransI + theLthinI) > Precision::Confusion();
2773
2774   //Compute the resulting value
2775   try {
2776     OCC_CATCH_SIGNALS;
2777     if (!GetSolver()->ComputeFunction(aFunction)) {
2778       SetErrorCode("TShape driver failed");
2779       return NULL;
2780     }
2781   }
2782   catch (Standard_Failure) {
2783     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2784     SetErrorCode(aFail->GetMessageString());
2785     return NULL;
2786   }
2787
2788   // BEGIN of fillet
2789   TopoDS_Shape aShapeShape = aShape->GetValue();
2790   TopTools_IndexedMapOfShape anEdgesIndices;
2791   TopExp::MapShapes(aShapeShape, anEdgesIndices);
2792   // Common edges on external cylinders
2793   Handle(GEOM_Object) box_e;
2794   if (theHexMesh) {
2795     box_e = my3DPrimOperations->MakeBoxDXDYDZ(theR2+theW2, theR2+theW2, theR1+theW1);
2796   }
2797   else {
2798     box_e = my3DPrimOperations->MakeBoxDXDYDZ(2*(theR2+theW2), 2*(theR2+theW2), theR1+theW1);
2799   }
2800   box_e->GetLastFunction()->SetDescription("");
2801   box_e = myTransformOperations->TranslateDXDYDZ(box_e, -(theR2+theW2), -(theR2+theW2), 0);
2802   box_e->GetLastFunction()->SetDescription("");
2803
2804   Handle(TColStd_HSequenceOfInteger) edges_e =
2805     myShapesOperations->GetShapesOnBoxIDs(box_e, aShape, TopAbs_EDGE, GEOMAlgo_ST_IN);
2806   box_e->GetLastFunction()->SetDescription("");
2807
2808   if (edges_e.IsNull() || edges_e->Length() == 0) {
2809     SetErrorCode("External edges not found");
2810     return NULL;
2811   }
2812   int nbEdgesInFillet = 0;
2813   std::list<int> theEdges;
2814   for (int i=1; i<=edges_e->Length();i++) {
2815     int edgeID = edges_e->Value(i);
2816     TopoDS_Shape theEdge = anEdgesIndices.FindKey(edgeID);
2817     TopExp_Explorer Ex(theEdge,TopAbs_VERTEX);
2818     while (Ex.More()) {
2819       gp_Pnt aPt = BRep_Tool::Pnt(TopoDS::Vertex(Ex.Current()));
2820       if (Abs(aPt.Z() - (theR1+theW1)) <= Precision::Confusion()) {
2821         nbEdgesInFillet ++;
2822         theEdges.push_back(edgeID);
2823       }
2824       Ex.Next();
2825     }
2826     if (theHexMesh && nbEdgesInFillet == 1)
2827       break;
2828   }
2829
2830   Handle(GEOM_Object) aFillet;
2831   try {
2832     aFillet = myLocalOperations->MakeFilletEdges(aShape, theRF, theEdges);
2833   }
2834   catch (Standard_Failure) {
2835     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2836     SetErrorCode(aFail->GetMessageString());
2837     return NULL;
2838   }
2839   if (aFillet.IsNull()) {
2840     //SetErrorCode("Fillet can not be computed on the given shape with the given parameters");
2841     SetErrorCode(myLocalOperations->GetErrorCode());
2842     return NULL;
2843   }
2844   aFillet->GetLastFunction()->SetDescription("");
2845
2846   TopoDS_Shape aFilletShape = aFillet->GetValue();
2847   aFunction->SetValue(aFilletShape);
2848   // END of fillet
2849
2850 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - BEGIN (1)
2851 // the following block, when enabled, leads to partitioning problems
2852 #if 0
2853 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - END (1)
2854   // BEGIN: Limit tolerances (debug)
2855   Handle(GEOM_Object) aCorr1 = myHealingOperations->LimitTolerance(aShape, 1e-07);
2856   TopoDS_Shape aCorr1Shape = aCorr1->GetValue();
2857   aShape->GetLastFunction()->SetValue(aCorr1Shape);
2858   aCorr1->GetLastFunction()->SetDescription("");
2859   // END: Limit tolerances (debug)
2860 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - BEGIN (2)
2861 #endif
2862 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - END (2)
2863
2864   if (theHexMesh) {
2865     if (!MakePipeTShapePartition(aShape, theR1, theW1, theL1, theR2, theW2, theL2, 0, 0, theRF, false))
2866       return NULL;
2867     if (!MakePipeTShapeMirrorAndGlue(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
2868       return NULL;
2869   }
2870
2871   // Add thickness reduction elements
2872   // at the three extremities: Left, Right and Incident
2873   try {
2874     OCC_CATCH_SIGNALS;
2875     if (isTRL || isTRR || isTRI) {
2876       TopoDS_Shape aResShape =
2877         MakePipeTShapeThicknessReduction(aShape->GetValue(), theR1, theW1, theL1, theR2, theW2, theL2,
2878                                          theRL, theWL, theLtransL, theLthinL,
2879                                          theRR, theWR, theLtransR, theLthinR,
2880                                          theRI, theWI, theLtransI, theLthinI,
2881                                          !theHexMesh);
2882       aFunction->SetValue(aResShape);
2883     }
2884   }
2885   catch (Standard_Failure) {
2886     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2887     SetErrorCode(aFail->GetMessageString());
2888     return NULL;
2889   }
2890
2891   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
2892   aSeq->Append(aShape);
2893
2894   try {
2895     if (theHexMesh) {
2896       // Get the groups
2897       if (!MakeGroups(aShape, TSHAPE_FILLET, theR1, theW1, theL1, theR2, theW2, theL2,
2898                       0., 0., theRF, aSeq, gp_Trsf()))
2899         return NULL;
2900     }
2901
2902     // Get internal group.
2903     if (!MakeInternalGroup(aShape, theR1, theL1, theR2, theL2, theRL, theLtransL,
2904                            theRR, theLtransR, theRI, theLtransI,
2905                            aSeq, gp_Trsf())) {
2906       return NULL;
2907     }
2908   }
2909   catch (Standard_Failure) {
2910     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2911     SetErrorCode(aFail->GetMessageString());
2912     return NULL;
2913   }
2914
2915   //Make a Python command
2916   TCollection_AsciiString anEntry, aListRes("[");
2917   // Iterate over the sequence aSeq
2918   Standard_Integer aNbGroups = aSeq->Length();
2919   Standard_Integer i = 1;
2920   for (; i <= aNbGroups; i++) {
2921     Handle(Standard_Transient) anItem = aSeq->Value(i);
2922     if (anItem.IsNull()) continue;
2923     Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(anItem);
2924     if (aGroup.IsNull()) continue;
2925     //Make a Python command
2926     TDF_Tool::Entry(aGroup->GetEntry(), anEntry);
2927     aListRes += anEntry + ", ";
2928   }
2929   aListRes.Trunc(aListRes.Length() - 2);
2930
2931   GEOM::TPythonDump pd (aFunction);
2932
2933   pd << aListRes.ToCString() << "] = geompy.MakePipeTShapeFillet("
2934      << theR1 << ", " << theW1 << ", " << theL1 << ", "
2935      << theR2 << ", " << theW2 << ", " << theL2 << ", "
2936      << theRF << ", " << theHexMesh;
2937
2938   // thickness reduction
2939   if (isTRL)
2940     pd << ", theRL=" << theRL << ", theWL=" << theWL
2941        << ", theLtransL=" << theLtransL << ", theLthinL=" << theLthinL;
2942   if (isTRR)
2943     pd << ", theRR=" << theRR << ", theWR=" << theWR
2944        << ", theLtransR=" << theLtransR << ", theLthinR=" << theLthinR;
2945   if (isTRI)
2946     pd << ", theRI=" << theRI << ", theWI=" << theWI
2947        << ", theLtransI=" << theLtransI << ", theLthinI=" << theLthinI;
2948
2949   pd << ")";
2950
2951   SetErrorCode(OK);
2952
2953   return aSeq;
2954 }
2955
2956 //=============================================================================
2957 /*!
2958  *  MakePipeTShapeFilletWithPosition
2959  *  \brief Create a T-shape object with specified caracteristics for the main and
2960  *  the incident pipes (radius, width, half-length). A fillet is created
2961  *  on the junction of the pipes.
2962  *  The extremities of the main pipe are located on junctions points P1 and P2.
2963  *  The extremity of the incident pipe is located on junction point P3.
2964  *  \param theR1 Internal radius of main pipe
2965  *  \param theW1 Width of main pipe
2966  *  \param theL1 Half-length of main pipe
2967  *  \param theR2 Internal radius of incident pipe (R2 < R1)
2968  *  \param theW2 Width of incident pipe (R2+W2 < R1+W1)
2969  *  \param theL2 Half-length of incident pipe
2970  *  \param theRF Radius of curvature of fillet
2971  *  \param theHexMesh Boolean indicating if shape is prepared for hex mesh
2972  *  \param theP1 1st junction point of main pipe
2973  *  \param theP2 2nd junction point of main pipe
2974  *  \param theP3 Junction point of incident pipe
2975  *  \return List of GEOM_Objects, containing the created shape and propagation groups.
2976  */
2977 //=============================================================================
2978 Handle(TColStd_HSequenceOfTransient)
2979 GEOMImpl_IAdvancedOperations::MakePipeTShapeFilletWithPosition
2980                              (double theR1, double theW1, double theL1,
2981                               double theR2, double theW2, double theL2,
2982                               double theRL, double theWL, double theLtransL, double theLthinL,
2983                               double theRR, double theWR, double theLtransR, double theLthinR,
2984                               double theRI, double theWI, double theLtransI, double theLthinI,
2985                               double theRF, bool theHexMesh,
2986                               Handle(GEOM_Object) theP1,
2987                               Handle(GEOM_Object) theP2,
2988                               Handle(GEOM_Object) theP3)
2989 {
2990   SetErrorCode(KO);
2991   //Add a new object
2992   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_TSHAPE);
2993   //Add a new shape function with parameters
2994   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_PipeTShapeDriver::GetID(), TSHAPE_FILLET);
2995   if (aFunction.IsNull()) return NULL;
2996
2997   //Check if the function is set correctly
2998   if (aFunction->GetDriverGUID() != GEOMImpl_PipeTShapeDriver::GetID()) return NULL;
2999
3000   // Check new position
3001   if (!CheckCompatiblePosition(theL1, theL2, theP1, theP2, theP3, 0.01)) {
3002     return NULL;
3003   }
3004
3005   GEOMImpl_IPipeTShape aData(aFunction);
3006
3007   aData.SetR1(theR1);
3008   aData.SetW1(theW1);
3009   aData.SetL1(theL1);
3010   aData.SetR2(theR2);
3011   aData.SetW2(theW2);
3012   aData.SetL2(theL2);
3013   aData.SetRF(theRF);
3014   aData.SetHexMesh(theHexMesh);
3015
3016   bool isTRL = (theRL + theWL + theLtransL + theLthinL) > Precision::Confusion();
3017   bool isTRR = (theRR + theWR + theLtransR + theLthinR) > Precision::Confusion();
3018   bool isTRI = (theRI + theWI + theLtransI + theLthinI) > Precision::Confusion();
3019
3020   //Compute the resulting value
3021   try {
3022     OCC_CATCH_SIGNALS;
3023     if (!GetSolver()->ComputeFunction(aFunction)) {
3024       SetErrorCode("TShape driver failed");
3025       return NULL;
3026     }
3027   }
3028   catch (Standard_Failure) {
3029     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
3030     SetErrorCode(aFail->GetMessageString());
3031     return NULL;
3032   }
3033
3034   // BEGIN of fillet
3035   TopoDS_Shape aShapeShape = aShape->GetValue();
3036   TopTools_IndexedMapOfShape anEdgesIndices;
3037   TopExp::MapShapes(aShapeShape, anEdgesIndices);
3038   // Common edges on external cylinders
3039   Handle(GEOM_Object) box_e;
3040   if (theHexMesh) {
3041     box_e = my3DPrimOperations->MakeBoxDXDYDZ(theR2+theW2, theR2+theW2, theR1+theW1);
3042   }
3043   else {
3044     box_e = my3DPrimOperations->MakeBoxDXDYDZ(2*(theR2+theW2), 2*(theR2+theW2), theR1+theW1);
3045   }
3046   box_e->GetLastFunction()->SetDescription("");
3047   box_e = myTransformOperations->TranslateDXDYDZ(box_e, -(theR2+theW2), -(theR2+theW2), 0);
3048   box_e->GetLastFunction()->SetDescription("");
3049
3050   Handle(TColStd_HSequenceOfInteger) edges_e =
3051     myShapesOperations->GetShapesOnBoxIDs(box_e, aShape, TopAbs_EDGE, GEOMAlgo_ST_IN);
3052   box_e->GetLastFunction()->SetDescription("");
3053
3054   if (edges_e.IsNull() || edges_e->Length() == 0) {
3055     SetErrorCode("External edges not found");
3056     return NULL;
3057   }
3058   int nbEdgesInFillet = 0;
3059   std::list<int> theEdges;
3060   for (int i=1; i<=edges_e->Length();i++) {
3061     int edgeID = edges_e->Value(i);
3062     TopoDS_Shape theEdge = anEdgesIndices.FindKey(edgeID);
3063     TopExp_Explorer Ex(theEdge,TopAbs_VERTEX);
3064     while (Ex.More()) {
3065       gp_Pnt aPt = BRep_Tool::Pnt(TopoDS::Vertex(Ex.Current()));
3066       if (Abs(aPt.Z() - (theR1+theW1)) <= Precision::Confusion()) {
3067         nbEdgesInFillet ++;
3068         theEdges.push_back(edgeID);
3069       }
3070       Ex.Next();
3071     }
3072     if (theHexMesh && nbEdgesInFillet == 1)
3073       break;
3074   }
3075
3076   Handle(GEOM_Object) aFillet;
3077   try {
3078     aFillet = myLocalOperations->MakeFilletEdges(aShape, theRF, theEdges);
3079   }
3080   catch (Standard_Failure) {
3081     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
3082     SetErrorCode(aFail->GetMessageString());
3083     return NULL;
3084   }
3085   if (aFillet.IsNull()) {
3086     SetErrorCode("Fillet can not be computed on the given shape with the given parameters");
3087     return NULL;
3088   }
3089   aFillet->GetLastFunction()->SetDescription("");
3090
3091   TopoDS_Shape aFilletShape = aFillet->GetValue();
3092   aFunction->SetValue(aFilletShape);
3093   // END of fillet
3094
3095 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - BEGIN (3)
3096 // the following block, when enabled, leads to partitioning problems
3097 #if 0
3098 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - END (3)
3099   // BEGIN: Limit tolerances (debug)
3100   Handle(GEOM_Object) aCorr1 = myHealingOperations->LimitTolerance(aShape, 1e-07);
3101   TopoDS_Shape aCorr1Shape = aCorr1->GetValue();
3102   aShape->GetLastFunction()->SetValue(aCorr1Shape);
3103   aCorr1->GetLastFunction()->SetDescription("");
3104   // END: Limit tolerances (debug)
3105 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - BEGIN (4)
3106 #endif
3107 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - END (4)
3108
3109   if (theHexMesh) {
3110     if (!MakePipeTShapePartition(aShape, theR1, theW1, theL1, theR2, theW2, theL2, 0, 0, theRF, false))
3111       return NULL;
3112     if (!MakePipeTShapeMirrorAndGlue(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
3113       return NULL;
3114   }
3115
3116   // Add thickness reduction elements
3117   // at the three extremities: Left, Right and Incident
3118   try {
3119     OCC_CATCH_SIGNALS;
3120     if (isTRL || isTRR || isTRI) {
3121       TopoDS_Shape aResShape =
3122         MakePipeTShapeThicknessReduction(aShape->GetValue(), theR1, theW1, theL1, theR2, theW2, theL2,
3123                                          theRL, theWL, theLtransL, theLthinL,
3124                                          theRR, theWR, theLtransR, theLthinR,
3125                                          theRI, theWI, theLtransI, theLthinI,
3126                                          !theHexMesh);
3127       aFunction->SetValue(aResShape);
3128     }
3129   }
3130   catch (Standard_Failure) {
3131     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
3132     SetErrorCode(aFail->GetMessageString());
3133     return NULL;
3134   }
3135
3136   // Set Position
3137   gp_Trsf aTrsf = GetPositionTrsf(theL1, theL2, theP1, theP2, theP3);
3138   BRepBuilderAPI_Transform aTransformation (aShape->GetValue(), aTrsf, Standard_False);
3139   TopoDS_Shape aTrsf_Shape = aTransformation.Shape();
3140   aFunction->SetValue(aTrsf_Shape);
3141
3142   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
3143   aSeq->Append(aShape);
3144
3145   try {
3146     if (theHexMesh) {
3147       // Get the groups
3148       if (!MakeGroups(aShape, TSHAPE_FILLET, theR1, theW1, theL1, theR2, theW2, theL2,
3149                       0., 0., theRF, aSeq, aTrsf))
3150         return NULL;
3151     }
3152
3153     // Get internal group.
3154     if (!MakeInternalGroup(aShape, theR1, theL1, theR2, theL2, theRL, theLtransL,
3155                            theRR, theLtransR, theRI, theLtransI,
3156                            aSeq, aTrsf)) {
3157       return NULL;
3158     }
3159   }
3160   catch (Standard_Failure) {
3161     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
3162     SetErrorCode(aFail->GetMessageString());
3163     return NULL;
3164   }
3165
3166   //Make a Python command
3167   TCollection_AsciiString anEntry, aListRes("[");
3168   // Iterate over the sequence aSeq
3169   Standard_Integer aNbGroups = aSeq->Length();
3170   Standard_Integer i = 1;
3171   for (; i <= aNbGroups; i++) {
3172     Handle(Standard_Transient) anItem = aSeq->Value(i);
3173     if (anItem.IsNull()) continue;
3174     Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(anItem);
3175     if (aGroup.IsNull()) continue;
3176     //Make a Python command
3177     TDF_Tool::Entry(aGroup->GetEntry(), anEntry);
3178     aListRes += anEntry + ", ";
3179   }
3180   aListRes.Trunc(aListRes.Length() - 2);
3181
3182   GEOM::TPythonDump pd (aFunction);
3183
3184   pd << aListRes.ToCString() << "] = geompy.MakePipeTShapeFillet("
3185      << theR1 << ", " << theW1 << ", " << theL1 << ", "
3186      << theR2  << ", " << theW2 << ", " << theL2 << ", "
3187      << theRF << ", " << theHexMesh << ", "
3188      << theP1 << ", " << theP2 << ", " << theP3;
3189
3190   // thickness reduction
3191   if (isTRL)
3192     pd << ", theRL=" << theRL << ", theWL=" << theWL
3193        << ", theLtransL=" << theLtransL << ", theLthinL=" << theLthinL;
3194   if (isTRR)
3195     pd << ", theRR=" << theRR << ", theWR=" << theWR
3196        << ", theLtransR=" << theLtransR << ", theLthinR=" << theLthinR;
3197   if (isTRI)
3198     pd << ", theRI=" << theRI << ", theWI=" << theWI
3199        << ", theLtransI=" << theLtransI << ", theLthinI=" << theLthinI;
3200
3201   pd << ")";
3202
3203   SetErrorCode(OK);
3204
3205   return aSeq;
3206 }
3207
3208 //=============================================================================
3209 /*!
3210  *  This function allows to create a disk already divided into blocks. It can be
3211  *  used to create divided pipes for later meshing in hexaedra.
3212  *  \param theR Radius of the disk
3213  *  \param theRatio Relative size of the central square diagonal against the disk diameter
3214  *  \param theOrientation Plane on which the disk will be built
3215  *  \param thePattern The division pattern of the disk (hexagon or square in the center)
3216  *  \return New GEOM_Object, containing the created shape.
3217  */
3218 //=============================================================================
3219 Handle(GEOM_Object) GEOMImpl_IAdvancedOperations::MakeDividedDisk (double theR, double theRatio,
3220                                                                    int theOrientation, int thePattern)
3221 {
3222   SetErrorCode(KO);
3223
3224   if (theOrientation != 1 &&
3225       theOrientation != 2 &&
3226       theOrientation != 3)
3227   {
3228     SetErrorCode("theOrientation must be 1(=OXY), 2(=OYZ) or 3(=OZX)");
3229     return NULL;
3230   }
3231   //Add a new object
3232   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_DIVIDEDDISK);
3233
3234   //Add a new shape function with parameters
3235   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_DividedDiskDriver::GetID(), DIVIDEDDISK_R_RATIO);
3236   if (aFunction.IsNull()) return NULL;
3237
3238   //Check if the function is set correctly
3239   if (aFunction->GetDriverGUID() != GEOMImpl_DividedDiskDriver::GetID()) return NULL;
3240
3241   GEOMImpl_IDividedDisk aData (aFunction);
3242
3243   aData.SetR(theR);
3244   aData.SetRatio(theRatio);
3245   aData.SetOrientation(theOrientation);
3246   aData.SetType(thePattern);
3247
3248   //Compute the resulting value
3249   try {
3250 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
3251     OCC_CATCH_SIGNALS;
3252 #endif
3253     if (!GetSolver()->ComputeFunction(aFunction)) {
3254       SetErrorCode("DividedDisk driver failed");
3255       return NULL;
3256     }
3257   }
3258   catch (Standard_Failure) {
3259     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
3260     SetErrorCode(aFail->GetMessageString());
3261     return NULL;
3262   }
3263
3264   std::string aPatternStr;
3265
3266   switch(thePattern)
3267   {
3268     case 0:
3269       aPatternStr = "GEOM.SQUARE";
3270       break;
3271     case 1:
3272       aPatternStr = "GEOM.HEXAGON";
3273       break;
3274   }
3275
3276   //Make a Python command
3277   GEOM::TPythonDump(aFunction) << aShape << " = geompy.MakeDividedDisk(" << theR << ", " << theOrientation << ", " << aPatternStr.c_str() << ")";
3278
3279   SetErrorCode(OK);
3280
3281   return aShape;
3282 }
3283
3284 //=============================================================================
3285 /*!
3286  *  This function allows to create a disk already divided into blocks. It can be
3287  *  used to create divided pipes for later meshing in hexaedra.
3288  *  \param theR Radius of the disk
3289  *  \param theRatio Relative size of the central square diagonal against the disk diameter
3290  *  \return New GEOM_Object, containing the created shape.
3291  */
3292 //=============================================================================
3293 Handle(GEOM_Object) GEOMImpl_IAdvancedOperations::MakeDividedDiskPntVecR (Handle(GEOM_Object) thePnt,
3294                                                                           Handle(GEOM_Object) theVec,
3295                                                                           double theR,
3296                                                                           double theRatio,
3297                                                                           int    thePattern)
3298 {
3299   SetErrorCode(KO);
3300
3301   //Add a new object
3302   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_DIVIDEDDISK);
3303
3304   //Add a new shape function with parameters
3305   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_DividedDiskDriver::GetID(), DIVIDEDDISK_R_VECTOR_PNT);
3306   if (aFunction.IsNull()) return NULL;
3307
3308   //Check if the function is set correctly
3309   if (aFunction->GetDriverGUID() != GEOMImpl_DividedDiskDriver::GetID()) return NULL;
3310
3311   GEOMImpl_IDividedDisk aData (aFunction);
3312
3313   Handle(GEOM_Function) aRefPnt = thePnt->GetLastFunction();
3314   Handle(GEOM_Function) aRefVec = theVec->GetLastFunction();
3315
3316   if (aRefPnt.IsNull() || aRefVec.IsNull()) return NULL;
3317
3318   aData.SetCenter(aRefPnt);
3319   aData.SetVector(aRefVec);
3320
3321   aData.SetR(theR);
3322   aData.SetRatio(theRatio);
3323   aData.SetType(thePattern);
3324
3325   //Compute the resulting value
3326   try {
3327 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
3328     OCC_CATCH_SIGNALS;
3329 #endif
3330     if (!GetSolver()->ComputeFunction(aFunction)) {
3331       SetErrorCode("DividedDisk driver failed");
3332       return NULL;
3333     }
3334   }
3335   catch (Standard_Failure) {
3336     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
3337     SetErrorCode(aFail->GetMessageString());
3338     return NULL;
3339   }
3340
3341   std::string aPatternStr;
3342
3343   switch(thePattern)
3344   {
3345     case 0:
3346       aPatternStr = "GEOM.SQUARE";
3347       break;
3348     case 1:
3349       aPatternStr = "GEOM.HEXAGON";
3350       break;
3351   }
3352
3353
3354   //Make a Python command
3355   GEOM::TPythonDump(aFunction) << aShape << " = geompy.MakeDividedDiskPntVecR(" << thePnt << ", " << theVec << ", " << theR << ", " << aPatternStr.c_str() << ")";
3356
3357   SetErrorCode(OK);
3358
3359   return aShape;
3360 }
3361
3362 //=============================================================================
3363 /*!
3364  *  Builds a cylinder prepared for hexa meshes
3365  *  \param theR Radius of the cylinder
3366  *  \param theH Height of the cylinder
3367  *  \return New GEOM_Object, containing the created shape.
3368  */
3369 //=============================================================================
3370 Handle(GEOM_Object) GEOMImpl_IAdvancedOperations::MakeDividedCylinder (double theR,
3371                                                                        double theH,
3372                                                                        int    thePattern)
3373 {
3374   SetErrorCode(KO);
3375
3376   //Add a new object
3377   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_DIVIDEDCYLINDER);
3378
3379   Handle(GEOM_Object) aBaseShape = MakeDividedDisk(theR, 67.0, 1, thePattern);
3380   aBaseShape->GetLastFunction()->SetDescription("");   // Erase dump of MakeDividedDisk
3381
3382   aShape = my3DPrimOperations->MakePrismDXDYDZ(aBaseShape,0.0,0.0,theH, -1.0);
3383
3384   Handle(GEOM_Function) aFunction =  aShape->GetLastFunction();
3385   aFunction->SetDescription("");   // Erase dump of MakePrismDXDYDZ
3386   aShape->SetType(GEOM_DIVIDEDCYLINDER);
3387
3388   std::string aPatternStr;
3389
3390   switch(thePattern)
3391   {
3392     case 0:
3393       aPatternStr = "GEOM.SQUARE";
3394       break;
3395     case 1:
3396       aPatternStr = "GEOM.HEXAGON";
3397       break;
3398   }
3399
3400   //Make a Python command
3401   GEOM::TPythonDump(aFunction) << aShape << " = geompy.MakeDividedCylinder(" << theR << ", " << theH << ", " << aPatternStr.c_str() << ")";
3402
3403   SetErrorCode(OK);
3404
3405   return aShape;
3406 }
3407 //=============================================================================
3408 /*!
3409  *  Create a smoothing surface from a set of points
3410  *  \param thelPoints list of points or compounds of points
3411  *  \param theNbMax maximum number of Bezier pieces in the resulting surface.
3412  *  \param theDegMax maximum degree of the resulting BSpline surface
3413  *  \param theDMax specifies maximum value of the GeomPlate_PlateG0Criterion criterion.
3414  *  \return New GEOM_Object, containing the created shape.
3415  */
3416 //=============================================================================
3417 Handle(GEOM_Object) GEOMImpl_IAdvancedOperations::MakeSmoothingSurface (std::list<Handle(GEOM_Object)> thelPoints,
3418                                                                         int                            theNbMax,
3419                                                                         int                            theDegMax,
3420                                                                         double                         theDMax)
3421 {
3422   SetErrorCode(KO);
3423
3424   //Add a new object
3425   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_SMOOTHINGSURFACE);
3426
3427   //Add a new shape function with parameters
3428   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_SmoothingSurfaceDriver::GetID(), SMOOTHINGSURFACE_LPOINTS);
3429   if (aFunction.IsNull()) return NULL;
3430
3431   //Check if the function is set correctly
3432   if (aFunction->GetDriverGUID() != GEOMImpl_SmoothingSurfaceDriver::GetID()) return NULL;
3433
3434   GEOMImpl_ISmoothingSurface aData (aFunction);
3435
3436   int aLen = thelPoints.size();
3437   aData.SetLength(aLen);
3438   int ind = 1;
3439   std::list<Handle(GEOM_Object)>::iterator it = thelPoints.begin();
3440   for (; it != thelPoints.end(); it++, ind++) {
3441     Handle(GEOM_Function) aRefObj = (*it)->GetLastFunction();
3442     if (aRefObj.IsNull()) {
3443       SetErrorCode("NULL point or compound for bSplineFaceShape");
3444       return NULL;
3445     }
3446     aData.SetPntOrComp(ind, aRefObj);
3447   }
3448
3449   aData.SetNbMax(theNbMax);
3450   aData.SetDegMax(theDegMax);
3451   aData.SetDMax(theDMax);
3452
3453   //Compute the resulting value
3454   try {
3455 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
3456     OCC_CATCH_SIGNALS;
3457 #endif
3458     if (!GetSolver()->ComputeFunction(aFunction)) {
3459       SetErrorCode("SmoothingSurface driver failed");
3460       return NULL;
3461     }
3462   }
3463   catch (Standard_Failure) {
3464     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
3465     SetErrorCode(aFail->GetMessageString());
3466     return NULL;
3467   }
3468
3469   //Make a Python command
3470   GEOM::TPythonDump pd (aFunction);
3471   pd << aShape << " = geompy.MakeSmoothingSurface([";
3472   it = thelPoints.begin();
3473   pd << (*it++);
3474   while (it != thelPoints.end()) {
3475     pd << ", " << (*it++);
3476   }
3477   pd << "], "
3478      << theNbMax << ", "
3479      << theDegMax << ", "
3480      << theDMax <<")";
3481
3482   SetErrorCode(OK);
3483
3484   return aShape;
3485 }
3486 /*@@ insert new functions before this line @@ do not remove this line @@ do not remove this line @@*/