src/GEOMImpl/GEOMImpl_IAdvancedOperations.cxx

   1 // Copyright (C) 2007-2012  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.
   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 <Standard_Stream.hxx>
  23
  24 #include "GEOMImpl_Types.hxx"
  25 #include "GEOMImpl_IAdvancedOperations.hxx"
  26 #include "GEOMImpl_IBasicOperations.hxx"
  27 #include "GEOMImpl_IBooleanOperations.hxx"
  28 #include "GEOMImpl_IShapesOperations.hxx"
  29 #include "GEOMImpl_ITransformOperations.hxx"
  30 #include "GEOMImpl_IBlocksOperations.hxx"
  31 #include "GEOMImpl_I3DPrimOperations.hxx"
  32 #include "GEOMImpl_ILocalOperations.hxx"
  33 #include "GEOMImpl_IHealingOperations.hxx"
  34
  35 #include "GEOMImpl_Gen.hxx"
  36
  37 #include <Basics_OCCTVersion.hxx>
  38
  39 #include <utilities.h>
  40 #include <OpUtil.hxx>
  41 #include <Utils_ExceptHandlers.hxx>
  42
  43 #include "GEOM_Function.hxx"
  44 #include "GEOM_PythonDump.hxx"
  45
  46 #include "GEOMImpl_PipeTShapeDriver.hxx"
  47 #include "GEOMImpl_IPipeTShape.hxx"
  48 #include <GEOMImpl_DividedDiskDriver.hxx>
  49 #include <GEOMImpl_IDividedDisk.hxx>
  50 // #include <GEOMImpl_DividedCylinderDriver.hxx>
  51 // #include <GEOMImpl_IDividedCylinder.hxx>
  52 /*@@ insert new functions before this line @@ do not remove this line @@ do not remove this line @@*/
  53
  54 #include <TopExp.hxx>
  55 #include <TopExp_Explorer.hxx>
  56 #include <TopoDS.hxx>
  57 #include <TopoDS_Vertex.hxx>
  58 #include <TopTools_IndexedMapOfShape.hxx>
  59
  60 #include <gp_Pnt.hxx>
  61 #include <gp_Vec.hxx>
  62 #include <gp_Ax3.hxx>
  63
  64 #include <BRepBuilderAPI_Transform.hxx>
  65 #include <BRep_Tool.hxx>
  66
  67 #include <cmath>
  68
  69 #include <TFunction_DriverTable.hxx>
  70 #include <TFunction_Driver.hxx>
  71 #include <TFunction_Logbook.hxx>
  72 #include <TDF_Tool.hxx>
  73 #include <Standard_Failure.hxx>
  74 #include <Standard_ErrorHandler.hxx> // CAREFUL ! position of this file is critic : see Lucien PIGNOLONI / OCC
  75
  76 #define HALF_LENGTH_MAIN_PIPE     "Main pipe half length" //"Tuyau principal - demi longueur"
  77 #define HALF_LENGTH_INCIDENT_PIPE "Incident pipe half length" //"Tuyau incident - demi longueur"
  78 #define CIRCULAR_QUARTER_PIPE     "Circular quarter of pipe" //"Circulaire - quart de tuyau"
  79 #define THICKNESS                 "Thickness" //"Epaisseur"
  80 #define FLANGE                    "Flange" // "Collerette"
  81 #define CHAMFER_OR_FILLET         "Chamfer or fillet" //"Chanfrein ou Raccord"
  82 #define JUNCTION_FACE_1           "Junction 1" //"Face de jonction 1"
  83 #define JUNCTION_FACE_2           "Junction 2" //"Face de jonction 2"
  84 #define JUNCTION_FACE_3           "Junction 3" //"Face de jonction 3"
  85
  86 //=============================================================================
  87 /*!
  88  *  Constructor
  89  */
  90 //=============================================================================
  91 GEOMImpl_IAdvancedOperations::GEOMImpl_IAdvancedOperations(GEOM_Engine* theEngine, int theDocID) :
  92   GEOM_IOperations(theEngine, theDocID)
  93 {
  94   MESSAGE("GEOMImpl_IAdvancedOperations::GEOMImpl_IAdvancedOperations");
  95   myBasicOperations     = new GEOMImpl_IBasicOperations(GetEngine(), GetDocID());
  96   myBooleanOperations   = new GEOMImpl_IBooleanOperations(GetEngine(), GetDocID());
  97   myShapesOperations    = new GEOMImpl_IShapesOperations(GetEngine(), GetDocID());
  98   myTransformOperations = new GEOMImpl_ITransformOperations(GetEngine(), GetDocID());
  99   myBlocksOperations    = new GEOMImpl_IBlocksOperations(GetEngine(), GetDocID());
 100   my3DPrimOperations    = new GEOMImpl_I3DPrimOperations(GetEngine(), GetDocID());
 101   myLocalOperations     = new GEOMImpl_ILocalOperations(GetEngine(), GetDocID());
 102   myHealingOperations   = new GEOMImpl_IHealingOperations(GetEngine(), GetDocID());
 103 }
 104
 105 //=============================================================================
 106 /*!
 107  *  Destructor
 108  */
 109 //=============================================================================
 110 GEOMImpl_IAdvancedOperations::~GEOMImpl_IAdvancedOperations()
 111 {
 112   MESSAGE("GEOMImpl_IAdvancedOperations::~GEOMImpl_IAdvancedOperations");
 113   delete myBasicOperations;
 114   delete myBooleanOperations;
 115   delete myShapesOperations;
 116   delete myTransformOperations;
 117   delete myBlocksOperations;
 118   delete my3DPrimOperations;
 119   delete myLocalOperations;
 120   delete myHealingOperations;
 121 }
 122
 123 //=============================================================================
 124 /*!
 125  *  SetPosition
 126  */
 127 //=============================================================================
 128 gp_Trsf GEOMImpl_IAdvancedOperations::GetPositionTrsf(double theL1, double theL2,
 129                                                       Handle(GEOM_Object) theP1,
 130                                                       Handle(GEOM_Object) theP2,
 131                                                       Handle(GEOM_Object) theP3)
 132 {
 133   // Old Local Coordinates System oldLCS
 134   gp_Pnt P0(0, 0, 0);
 135   gp_Pnt P1(-theL1, 0, 0);
 136   gp_Pnt P2(theL1, 0, 0);
 137   gp_Pnt P3(0, 0, theL2);
 138
 139   gp_Dir oldX(gp_Vec(P1, P2));
 140   gp_Dir oldZ(gp_Vec(P0, P3));
 141   gp_Ax3 oldLCS(P0, oldZ, oldX);
 142
 143   // New Local Coordinates System newLCS
 144   double LocX, LocY, LocZ;
 145   gp_Pnt newP1 = BRep_Tool::Pnt(TopoDS::Vertex(theP1->GetValue()));
 146   gp_Pnt newP2 = BRep_Tool::Pnt(TopoDS::Vertex(theP2->GetValue()));
 147   gp_Pnt newP3 = BRep_Tool::Pnt(TopoDS::Vertex(theP3->GetValue()));
 148   LocX = (newP1.X() + newP2.X()) / 2.;
 149   LocY = (newP1.Y() + newP2.Y()) / 2.;
 150   LocZ = (newP1.Z() + newP2.Z()) / 2.;
 151   gp_Pnt newO(LocX, LocY, LocZ);
 152
 153   gp_Dir newX(gp_Vec(newP1, newP2)); // P1P2 Vector
 154   gp_Dir newZ(gp_Vec(newO, newP3)); // OP3 Vector
 155   gp_Ax3 newLCS = gp_Ax3(newO, newZ, newX);
 156
 157   gp_Trsf aTrsf;
 158   aTrsf.SetDisplacement(oldLCS, newLCS);
 159
 160   return aTrsf;
 161 }
 162
 163 //=============================================================================
 164 /*!
 165  *  CheckCompatiblePosition
 166  *
 167  */
 168 //=============================================================================
 169 bool GEOMImpl_IAdvancedOperations::CheckCompatiblePosition(double& theL1, double& theL2,
 170                                                            Handle(GEOM_Object) theP1,
 171                                                            Handle(GEOM_Object) theP2,
 172                                                            Handle(GEOM_Object) theP3,
 173                                                            double theTolerance)
 174 {
 175   SetErrorCode(KO);
 176   gp_Pnt P1 = BRep_Tool::Pnt(TopoDS::Vertex(theP1->GetValue()));
 177   gp_Pnt P2 = BRep_Tool::Pnt(TopoDS::Vertex(theP2->GetValue()));
 178   gp_Pnt P3 = BRep_Tool::Pnt(TopoDS::Vertex(theP3->GetValue()));
 179
 180   double d12 = P1.Distance(P2);
 181   double d13 = P1.Distance(P3);
 182   double d23 = P2.Distance(P3);
 183   //    double d2 = newO.Distance(P3);
 184
 185   if (Abs(d12) <= Precision::Confusion()) {
 186     SetErrorCode("Junctions points P1 and P2 are identical");
 187     return false;
 188   }
 189   if (Abs(d13) <= Precision::Confusion()) {
 190     SetErrorCode("Junctions points P1 and P3 are identical");
 191     return false;
 192   }
 193   if (Abs(d23) <= Precision::Confusion()) {
 194     SetErrorCode("Junctions points P2 and P3 are identical");
 195     return false;
 196   }
 197
 198
 199   double newL1 = 0.5 * d12;
 200   double newL2 = sqrt(pow(d13,2)-pow(newL1,2));
 201   //
 202   // theL1*(1-theTolerance) <= newL1 <= theL1*(1+theTolerance)
 203   //
 204   if (fabs(newL1 - theL1) > Precision::Approximation()) {
 205     if ( (newL1 * (1 - theTolerance) -theL1 <= Precision::Approximation()) &&
 206          (newL1 * (1 + theTolerance) -theL1 >= Precision::Approximation()) ) {
 207       //            std::cerr << "theL1 = newL1" << std::endl;
 208       theL1 = newL1;
 209     } else {
 210       theL1 = -1;
 211       SetErrorCode("Dimension for main pipe (L1) is incompatible with new position");
 212       return false;
 213     }
 214   }
 215
 216   //
 217   // theL2*(1-theTolerance) <= newL2  <= theL2*(1+theTolerance)
 218   //
 219   if (fabs(newL2 - theL2) > Precision::Approximation()) {
 220     if ( (newL2 * (1 - theTolerance) -theL2 <= Precision::Approximation()) &&
 221          (newL2 * (1 + theTolerance) -theL2 >= Precision::Approximation()) ) {
 222       theL2 = newL2;
 223     } else {
 224       theL2 = -1;
 225       SetErrorCode("Dimension for incident pipe (L2) is incompatible with new position");
 226       return false;
 227     }
 228   }
 229
 230   SetErrorCode(OK);
 231   return true;
 232
 233 }
 234
 235 //=============================================================================
 236 /*!
 237  *  Generate the propagation groups of a Pipe T-Shape used for hexa mesh
 238  */
 239 //=============================================================================
 240 bool GEOMImpl_IAdvancedOperations::MakeGroups(Handle(GEOM_Object) theShape, int shapeType,
 241                                               double theR1, double theW1, double theL1,
 242                                               double theR2, double theW2, double theL2,
 243                                               Handle(TColStd_HSequenceOfTransient) theSeq,
 244                                               gp_Trsf aTrsf)
 245 {
 246   SetErrorCode(KO);
 247
 248   if (theShape.IsNull()) return false;
 249
 250   TopoDS_Shape aShape = theShape->GetValue();
 251   if (aShape.IsNull()) {
 252     SetErrorCode("Shape is not defined");
 253     return false;
 254   }
 255
 256   gp_Trsf aTrsfInv = aTrsf.Inverted();
 257
 258 //   int expectedGroups = 0;
 259 //   if (shapeType == TSHAPE_BASIC)
 260 //     if (Abs(theR2+theW2-theR1-theW1) <= Precision::Approximation())
 261 //       expectedGroups = 10;
 262 //     else
 263 //       expectedGroups = 11;
 264 //   else if (shapeType == TSHAPE_CHAMFER || shapeType == TSHAPE_FILLET)
 265 //     expectedGroups = 12;
 266
 267   double aR1Ext = theR1 + theW1;
 268   double aR2Ext = theR2 + theW2;
 269
 270   /////////////////////////
 271   //// Groups of Faces ////
 272   /////////////////////////
 273
 274   //
 275   // Comment the following lines when GetInPlace bug is solved
 276   // == BEGIN
 277   // Workaround of GetInPlace bug
 278   // Create a bounding box that fits the shape
 279   Handle(GEOM_Object) aBox = my3DPrimOperations->MakeBoxDXDYDZ(2*theL1, 2*aR1Ext, aR1Ext+theL2);
 280   aBox->GetLastFunction()->SetDescription("");
 281   myTransformOperations->TranslateDXDYDZ(aBox, -theL1, -aR1Ext, -aR1Ext);
 282   aBox->GetLastFunction()->SetDescription("");
 283   // Apply transformation to box
 284   BRepBuilderAPI_Transform aTransformationBox(aBox->GetValue(), aTrsf, Standard_False);
 285   TopoDS_Shape aBoxShapeTrsf = aTransformationBox.Shape();
 286   aBox->GetLastFunction()->SetValue(aBoxShapeTrsf);
 287
 288   // Get the shell of the box
 289   Handle(GEOM_Object) aShell = Handle(GEOM_Object)::DownCast
 290     (myShapesOperations->MakeExplode(aBox, TopAbs_SHELL, true)->Value(1));
 291   aBox->GetLastFunction()->SetDescription("");
 292   aShell->GetLastFunction()->SetDescription("");
 293   // Get the common shapes between shell and shape
 294   Handle(GEOM_Object) aCommonCompound = myBooleanOperations->MakeBoolean (theShape, aShell, 1); // MakeCommon
 295   if (aCommonCompound.IsNull()) {
 296     SetErrorCode(myBooleanOperations->GetErrorCode());
 297     return false;
 298   }
 299   aCommonCompound->GetLastFunction()->SetDescription("");
 300   // Explode the faces of common shapes => 3 faces
 301   Handle(TColStd_HSequenceOfTransient) aCommonFaces =
 302     myShapesOperations->MakeExplode(aCommonCompound, TopAbs_FACE, true);
 303   aCommonCompound->GetLastFunction()->SetDescription("");
 304   std::list<Handle(GEOM_Object)> aCompoundOfFacesList;
 305
 306   for (int i=0 ; i<= aCommonFaces->Length()-4 ; i+=4) {
 307     std::list<Handle(GEOM_Object)> aFacesList;
 308     for (int j = 1 ; j <= 4 ; j++) {
 309       Handle(GEOM_Object) aFace = Handle(GEOM_Object)::DownCast(aCommonFaces->Value(i+j)); // Junction faces
 310       if (!aFace.IsNull()) {
 311         aFace->GetLastFunction()->SetDescription("");
 312         aFacesList.push_back(aFace);
 313       }
 314     }
 315     Handle(GEOM_Object) aCompoundOfFaces = myShapesOperations->MakeCompound(aFacesList);
 316     if (!aCompoundOfFaces.IsNull()) {
 317       aCompoundOfFaces->GetLastFunction()->SetDescription("");
 318       aCompoundOfFacesList.push_back(aCompoundOfFaces);
 319     }
 320   }
 321
 322   if (aCompoundOfFacesList.size() == 3) {
 323     Handle(GEOM_Object) aPln1 = aCompoundOfFacesList.front();
 324     aCompoundOfFacesList.pop_front();
 325     Handle(GEOM_Object) aPln2 = aCompoundOfFacesList.front();
 326     aCompoundOfFacesList.pop_front();
 327     Handle(GEOM_Object) aPln3 = aCompoundOfFacesList.front();
 328     aCompoundOfFacesList.pop_front();
 329     // == END
 330     //
 331
 332
 333     //     Uncomment the following lines when GetInPlace bug is solved
 334     //     == BEGIN
 335 //     Handle(GEOM_Object) aP1 = myBasicOperations->MakePointXYZ(-theL1, 0, 0);
 336 //     Handle(GEOM_Object) aP2 = myBasicOperations->MakePointXYZ(-0, 0, theL2);
 337 //     Handle(GEOM_Object) aP3 = myBasicOperations->MakePointXYZ(theL1, 0, 0);
 338 //     aP1->GetLastFunction()->SetDescription("");
 339 //     aP2->GetLastFunction()->SetDescription("");
 340 //     aP3->GetLastFunction()->SetDescription("");
 341 //     Handle(GEOM_Object) aV1 = myBasicOperations->MakeVectorDXDYDZ(-1, 0, 0);
 342 //     Handle(GEOM_Object) aV2 = myBasicOperations->MakeVectorDXDYDZ(0, 0, 1);
 343 //     Handle(GEOM_Object) aV3 = myBasicOperations->MakeVectorDXDYDZ(1, 0, 0);
 344 //     aV1->GetLastFunction()->SetDescription("");
 345 //     aV2->GetLastFunction()->SetDescription("");
 346 //     aV3->GetLastFunction()->SetDescription("");
 347 //     Handle(GEOM_Object) aPln1 = myBasicOperations->MakePlanePntVec(aP1, aV1, 2*(aR1Ext+theL2));
 348 //     Handle(GEOM_Object) aPln2 = myBasicOperations->MakePlanePntVec(aP2, aV2, 2*(aR2Ext));
 349 //     Handle(GEOM_Object) aPln3 = myBasicOperations->MakePlanePntVec(aP3, aV3, 2*(aR1Ext+theL2));
 350 //     aPln1->GetLastFunction()->SetDescription("");
 351 //     aPln2->GetLastFunction()->SetDescription("");
 352 //     aPln3->GetLastFunction()->SetDescription("");
 353 //
 354 //     BRepBuilderAPI_Transform aTransformation1(aPln1->GetValue(), aTrsf, Standard_False);
 355 //     TopoDS_Shape aTrsf_Shape1 = aTransformation1.Shape();
 356 //     aPln1->GetLastFunction()->SetValue(aTrsf_Shape1);
 357 //     BRepBuilderAPI_Transform aTransformation2(aPln2->GetValue(), aTrsf, Standard_False);
 358 //     TopoDS_Shape aTrsf_Shape2 = aTransformation2.Shape();
 359 //     aPln2->GetLastFunction()->SetValue(aTrsf_Shape2);
 360 //     BRepBuilderAPI_Transform aTransformation3(aPln3->GetValue(), aTrsf, Standard_False);
 361 //     TopoDS_Shape aTrsf_Shape3 = aTransformation3.Shape();
 362 //     aPln3->GetLastFunction()->SetValue(aTrsf_Shape3);
 363     //     == END
 364     //
 365
 366     Handle(GEOM_Object) junctionFaces1 = myShapesOperations->GetInPlace(theShape, aPln1);
 367     if (junctionFaces1.IsNull())
 368       junctionFaces1 = myShapesOperations->GetShapesOnShapeAsCompound
 369         (aPln1, theShape, TopAbs_FACE,  GEOMAlgo_ST_ONIN);
 370     if (!junctionFaces1.IsNull()) {
 371       junctionFaces1->GetLastFunction()->SetDescription("");
 372       junctionFaces1->SetName("JUNCTION_FACE_1");
 373       theSeq->Append(junctionFaces1);
 374     }
 375     else {
 376       SetErrorCode("Junction face 1 not found");
 377       //        theSeq->Append(aPln1);
 378       //        return false;
 379     }
 380     Handle(GEOM_Object) junctionFaces2 = myShapesOperations->GetInPlace(theShape, aPln2);
 381     if (junctionFaces2.IsNull())
 382       junctionFaces2 = myShapesOperations->GetShapesOnShapeAsCompound
 383         (aPln2, theShape, TopAbs_FACE,  GEOMAlgo_ST_ONIN);
 384     if (!junctionFaces2.IsNull()) {
 385       junctionFaces2->GetLastFunction()->SetDescription("");
 386       junctionFaces2->SetName("JUNCTION_FACE_2");
 387       theSeq->Append(junctionFaces2);
 388     }
 389     else {
 390       SetErrorCode("Junction face 2 not found");
 391       //        theSeq->Append(aPln2);
 392       //        return false;
 393     }
 394     Handle(GEOM_Object) junctionFaces3 = myShapesOperations->GetInPlace(theShape, aPln3);
 395     if (junctionFaces3.IsNull())
 396       junctionFaces3 = myShapesOperations->GetShapesOnShapeAsCompound
 397         (aPln3, theShape, TopAbs_FACE,  GEOMAlgo_ST_ONIN);
 398     if (!junctionFaces3.IsNull()) {
 399       junctionFaces3->GetLastFunction()->SetDescription("");
 400       junctionFaces3->SetName("JUNCTION_FACE_3");
 401       theSeq->Append(junctionFaces3);
 402     }
 403     else {
 404       SetErrorCode("Junction face 3 not found");
 405       //        theSeq->Append(aPln3);
 406       //        return false;
 407     }
 408   // Comment the following lines when GetInPlace bug is solved
 409   // == BEGIN
 410   }
 411   //     == END
 412   /////////////////////////
 413   //// Groups of Edges ////
 414   /////////////////////////
 415   // Result of propagate
 416
 417   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
 418
 419   TCollection_AsciiString theDesc = aFunction->GetDescription();
 420   Handle(TColStd_HSequenceOfTransient) aSeqPropagate = myBlocksOperations->Propagate(theShape);
 421   if (aSeqPropagate.IsNull() || aSeqPropagate->Length() == 0) {
 422     SetErrorCode("Propagation groups not found");
 423     return false;
 424   }
 425   Standard_Integer nbEdges, aNbGroups = aSeqPropagate->Length();
 426   // Recover previous description to get rid of Propagate dump
 427   aFunction->SetDescription(theDesc);
 428
 429   bool addGroup;
 430   bool circularFoundAndAdded = false;
 431   bool circularFound10 = false;
 432   bool incidentPipeFound = false;
 433   bool mainPipeFound = false;
 434   bool mainPipeFoundAndAdded = false;
 435   bool radialFound =false;
 436   bool flangeFound = false;
 437   bool flangeFoundAndAdded = false;
 438   bool chamferOrFilletFound = false;
 439
 440   for (int i=1 ; i<= aNbGroups; i++) {
 441     addGroup = false;
 442
 443     Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(aSeqPropagate->Value(i));
 444     if(aGroup.IsNull())
 445       continue;
 446
 447     TopoDS_Shape aGroupShape = aGroup->GetValue();
 448     BRepBuilderAPI_Transform aTransformationShapeInv(aGroupShape, aTrsfInv, Standard_False);
 449     TopoDS_Shape aGroupShapeTrsfInv = aTransformationShapeInv.Shape();
 450
 451     TopTools_IndexedMapOfShape anEdgesMap;
 452     TopExp::MapShapes(aGroupShapeTrsfInv,TopAbs_EDGE, anEdgesMap);
 453     nbEdges = anEdgesMap.Extent();
 454
 455     if (shapeType == TSHAPE_BASIC) {
 456       if ((nbEdges == 21) || /*R1Ext = R2Ext*/(nbEdges == 17)){
 457         addGroup = true;
 458         aGroup->SetName("THICKNESS");
 459       }
 460       else if (nbEdges == 6) {
 461         if (!circularFoundAndAdded) {
 462           circularFoundAndAdded = true;
 463           addGroup = true;
 464           aGroup->SetName("CIRCULAR_QUARTER_PIPE");
 465         }
 466       }
 467       else if (nbEdges == 8) {
 468         incidentPipeFound = true;
 469         mainPipeFound = false;
 470         radialFound =false;
 471         flangeFound = false;
 472
 473         TopExp_Explorer Ex(aGroupShapeTrsfInv,TopAbs_VERTEX);
 474         while (Ex.More()) {
 475           gp_Pnt aP =  BRep_Tool::Pnt(TopoDS::Vertex(Ex.Current()));
 476           double x=aP.X(), y=aP.Y(), z=aP.Z();
 477
 478
 479           if ((Abs(x) > aR2Ext + Precision::Confusion()) ||
 480               (Abs(y) > aR2Ext + Precision::Confusion())) {
 481             incidentPipeFound = false;
 482           }
 483
 484           if ( z < -Precision::Confusion()) {
 485             // length of main pipe
 486             mainPipeFound = true;
 487             if (!mainPipeFoundAndAdded) {
 488               mainPipeFoundAndAdded = true;
 489               addGroup = true;
 490               aGroup->SetName("HALF_LENGTH_MAIN_PIPE");
 491             }
 492           }
 493
 494           else if (Abs(x) > (theL1-Precision::Confusion())) {
 495             // discretisation circulaire
 496             radialFound = true;
 497             if (!circularFoundAndAdded) {
 498               circularFoundAndAdded = true;
 499               addGroup = true;
 500               aGroup->SetName("CIRCULAR_QUARTER_PIPE");
 501             }
 502           }
 503           Ex.Next();
 504         }
 505         if (incidentPipeFound) {
 506           addGroup = true;
 507           aGroup->SetName("HALF_LENGTH_INCIDENT_PIPE");
 508         }
 509         if (!addGroup && (!incidentPipeFound &&
 510                           !radialFound &&
 511                           !mainPipeFound &&
 512                           !flangeFound)) {
 513           // Flange (collerette)
 514           flangeFound = true;
 515           addGroup = true;
 516           aGroup->SetName("FLANGE");
 517         }
 518       }
 519       else
 520         continue;
 521     }
 522     else if (shapeType == TSHAPE_CHAMFER || shapeType == TSHAPE_FILLET) {
 523       if (nbEdges == 25) {
 524         addGroup = true;
 525         aGroup->SetName("THICKNESS");
 526       }
 527       else if ((nbEdges == 10) || (nbEdges == 6)) {
 528         if (!circularFoundAndAdded) {
 529           addGroup = true;
 530           circularFoundAndAdded = true;
 531           aGroup->SetName("CIRCULAR_QUARTER_PIPE");
 532           if (nbEdges == 10) {
 533             circularFound10 = true;
 534           }
 535         }
 536         else if (!circularFound10 && nbEdges == 10) {
 537           circularFound10 = true;
 538           addGroup = true;
 539           aGroup->SetName("CIRCULAR_QUARTER_PIPE");
 540         }
 541       }
 542       else if (nbEdges == 8) {
 543         incidentPipeFound = true;
 544         mainPipeFound = true;
 545         flangeFound = false;
 546
 547         bool isNearZ0 = false;
 548         bool isBelowZ0 = false;
 549
 550         TopExp_Explorer Ex (aGroupShapeTrsfInv,TopAbs_VERTEX);
 551         while (Ex.More()) {
 552           gp_Pnt aP = BRep_Tool::Pnt(TopoDS::Vertex(Ex.Current()));
 553           double x=aP.X(), y=aP.Y(), z=aP.Z();
 554
 555           // tuy_princ_long_avant & tuy_princ_long_apres
 556           //bool isMain = (((z < Precision::Confusion()) || (x < Precision::Confusion())) &&
 557           //               ((y <= aR1Ext + Precision::Confusion()) ||
 558           //                (y <= -(aR1Ext + Precision::Confusion())) ||
 559           //                (y <= theR1 + Precision::Confusion()) ||
 560           //                (y == -(theR1 + Precision::Confusion()))));
 561           bool isMain = ((z < Precision::Confusion() || x < Precision::Confusion()) &&
 562                          (fabs(y) > theR1 - Precision::Confusion() ||
 563                           fabs(y) < Precision::Confusion()));
 564
 565           if (!isMain) {
 566             mainPipeFound = false;
 567           }
 568
 569           // collerette
 570           //if (z < Precision::Confusion() && !isMain) {
 571           //  flangeFound = true;
 572           //  if (!flangeFoundAndAdded) {
 573           //    flangeFoundAndAdded = true;
 574           //    addGroup = true;
 575           //    aGroup->SetName("FLANGE");
 576           //  }
 577           //}
 578           if (fabs(z) < Precision::Confusion()) isNearZ0 = true;
 579           if (z < - Precision::Confusion()) isBelowZ0 = true;
 580
 581           // tuyau incident
 582           if ((Abs(x) > aR2Ext + Precision::Confusion()) ||
 583               (Abs(y) > aR2Ext + Precision::Confusion())) {
 584             incidentPipeFound = false;
 585           }
 586           Ex.Next();
 587         }
 588         if (mainPipeFound) {
 589           addGroup = true;
 590           aGroup->SetName("HALF_LENGTH_MAIN_PIPE");
 591         }
 592         if (incidentPipeFound) {
 593           addGroup = true;
 594           aGroup->SetName("HALF_LENGTH_INCIDENT_PIPE");
 595         }
 596         if (isNearZ0 && !isBelowZ0) {
 597           flangeFound = true;
 598           if (!flangeFoundAndAdded) {
 599             flangeFoundAndAdded = true;
 600             addGroup = true;
 601             aGroup->SetName("FLANGE");
 602           }
 603         }
 604         if (!addGroup && (!incidentPipeFound &&
 605                           !mainPipeFound &&
 606                           !flangeFound &&
 607                           !chamferOrFilletFound)) {
 608           addGroup = true;
 609           chamferOrFilletFound = true;
 610           if (shapeType == TSHAPE_CHAMFER)
 611             aGroup->SetName("CHAMFER");
 612           else
 613             aGroup->SetName("FILLET");
 614         }
 615       }
 616       else
 617         continue;
 618     }
 619     // Add group to the list
 620     if (addGroup)
 621       theSeq->Append(aGroup);
 622   }
 623
 624   SetErrorCode(OK);
 625   return true;
 626 }
 627
 628 bool GEOMImpl_IAdvancedOperations::MakePipeTShapePartition(Handle(GEOM_Object) theShape,
 629                                                            double theR1, double theW1, double theL1,
 630                                                            double theR2, double theW2, double theL2,
 631                                                            double theH, double theW,
 632                                                            double theRF, bool isNormal)
 633 {
 634   SetErrorCode(KO);
 635
 636   // Build tools for partition operation:
 637   // 1 face and 2 planes
 638   // Face
 639   Handle(GEOM_Object) arete_intersect_int;
 640   Handle(GEOM_Object) wire_t, wire_t2, face_t, face_t2;
 641   Handle(GEOM_Object) chan_racc;
 642   Handle(GEOM_Object) vi1, vi2;
 643   Handle(GEOM_Object) Te3;
 644
 645   try {
 646 #if OCC_VERSION_LARGE > 0x06010000
 647     OCC_CATCH_SIGNALS;
 648 #endif
 649     Handle(GEOM_Object) Vector_Z = myBasicOperations->MakeVectorDXDYDZ(0, 0, 1);
 650     Vector_Z->GetLastFunction()->SetDescription("");
 651
 652     // Useful values
 653     double aSize = 2*(theL1 + theL2);
 654     double aR1Ext = theR1 + theW1;
 655     double aR2Ext = theR2 + theW2;
 656     double theVertCylinderRadius = aR2Ext + theW + theRF;
 657     double theHoriCylinderRadius = aR1Ext + theH + theRF;
 658
 659     // Common edges on internal cylinder
 660     Handle(GEOM_Object) box_i = my3DPrimOperations->MakeBoxDXDYDZ(theR2, theR2, theR1);
 661     box_i->GetLastFunction()->SetDescription("");
 662     box_i = myTransformOperations->TranslateDXDYDZ(box_i, -theR2, -theR2, 0);
 663     box_i->GetLastFunction()->SetDescription("");
 664
 665     Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
 666     TCollection_AsciiString theDesc = aFunction->GetDescription();
 667     Handle(TColStd_HSequenceOfTransient) edges_i =
 668       myShapesOperations->GetShapesOnBox(box_i, theShape, TopAbs_EDGE, GEOMAlgo_ST_IN);
 669     // Recover previous description to get rid of Propagate dump
 670     aFunction->SetDescription(theDesc);
 671     if (edges_i.IsNull() || edges_i->Length() == 0) {
 672       SetErrorCode("Internal edges not found");
 673       return false;
 674     }
 675     for (int i=1; i<=edges_i->Length();i++) {
 676       Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast(edges_i->Value(i));
 677       anObj->GetLastFunction()->SetDescription("");
 678     }
 679     arete_intersect_int = Handle(GEOM_Object)::DownCast(edges_i->Value(1));
 680
 681     // search for vertices located on both internal pipes
 682     aFunction = theShape->GetLastFunction();
 683     theDesc = aFunction->GetDescription();
 684     Handle(TColStd_HSequenceOfTransient) vertices_i =
 685       myShapesOperations->GetShapesOnBox(box_i, theShape, TopAbs_VERTEX, GEOMAlgo_ST_ONIN);
 686     // Recover previous description to get rid of Propagate dump
 687     aFunction->SetDescription(theDesc);
 688     if (vertices_i.IsNull() || vertices_i->Length() == 0) {
 689       SetErrorCode("Internal vertices not found");
 690       return false;
 691     }
 692
 693     for (int i = 1; i <= vertices_i->Length(); i++) {
 694       Handle(GEOM_Object) v = Handle(GEOM_Object)::DownCast(vertices_i->Value(i));
 695       v->GetLastFunction()->SetDescription("");
 696       TopoDS_Vertex aVertex = TopoDS::Vertex(v->GetValue());
 697       gp_Pnt aP = BRep_Tool::Pnt(aVertex);
 698 //       std::cout << "Coords: " << aP.X() << ", " << aP.Y() << ", " << aP.Z() << std::endl;
 699       if (Abs(aP.X()) <= Precision::Confusion()) {
 700         if (Abs(aP.Y()) - theR1 <= Precision::Confusion()) {
 701           vi1 = v;
 702         }
 703       } else if (Abs(aP.Y()) <= Precision::Confusion()) {
 704         if (Abs(aP.X()) - theR1 <= Precision::Confusion()) {
 705           vi2 = v;
 706         }
 707       }
 708     }
 709
 710     std::list<Handle(GEOM_Object)> theShapes;
 711
 712     if (isNormal) {
 713       Handle(GEOM_Object) ve1, ve2;
 714
 715       Handle(GEOM_Object) box_e = my3DPrimOperations->MakeBoxDXDYDZ(aR2Ext, aR2Ext, aR1Ext);
 716       box_e->GetLastFunction()->SetDescription("");
 717       box_e = myTransformOperations->TranslateDXDYDZ(box_e, -aR2Ext, -aR2Ext, 0);
 718       box_e->GetLastFunction()->SetDescription("");
 719       // Common edges on external cylinder
 720       aFunction = theShape->GetLastFunction();
 721       theDesc = aFunction->GetDescription();
 722       Handle(TColStd_HSequenceOfTransient) edges_e =
 723         myShapesOperations->GetShapesOnBox(box_e, theShape, TopAbs_EDGE, GEOMAlgo_ST_IN);
 724       // Recover previous description to get rid of Propagate dump
 725       aFunction->SetDescription(theDesc);
 726       if (edges_e.IsNull() || edges_e->Length() == 0) {
 727         SetErrorCode("External edges not found");
 728         return false;
 729       }
 730       for (int i=1; i<=edges_e->Length();i++) {
 731         Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast(edges_e->Value(i));
 732         anObj->GetLastFunction()->SetDescription("");
 733       }
 734
 735       // search for vertices located on both external pipes
 736       aFunction = theShape->GetLastFunction();
 737       theDesc = aFunction->GetDescription();
 738       Handle(TColStd_HSequenceOfTransient) vertices_e =
 739         myShapesOperations->GetShapesOnBox(box_e, theShape, TopAbs_VERTEX, GEOMAlgo_ST_ONIN);
 740       // Recover previous description to get rid of Propagate dump
 741       aFunction->SetDescription(theDesc);
 742       if (vertices_e.IsNull() || vertices_e->Length() == 0) {
 743         SetErrorCode("External vertices not found");
 744         return false;
 745       }
 746
 747       for (int i = 1; i <= vertices_e->Length(); i++) {
 748         Handle(GEOM_Object) v = Handle(GEOM_Object)::DownCast(vertices_e->Value(i));
 749         v->GetLastFunction()->SetDescription("");
 750         TopoDS_Vertex aVertex = TopoDS::Vertex(v->GetValue());
 751         gp_Pnt aP = BRep_Tool::Pnt(aVertex);
 752 //         std::cout << "Coords: " << aP.X() << ", " << aP.Y() << ", " << aP.Z() << std::endl;
 753         if (Abs(aP.X()) <= Precision::Confusion()) {
 754           if (Abs(aP.Y()) - theR2 > Precision::Confusion()) {
 755             ve1 = v;
 756           }
 757         } else if (Abs(aP.Y()) <= Precision::Confusion()) {
 758           if (Abs(aP.X()) - theR2 > Precision::Confusion()) {
 759             ve2 = v;
 760           }
 761         }
 762         if ( !ve1.IsNull() && !ve2.IsNull())
 763           break;
 764       }
 765       Handle(GEOM_Object) edge_e1, edge_e2;
 766
 767       edge_e1 = myBasicOperations->MakeLineTwoPnt(ve1, vi1);
 768       if (edge_e1.IsNull()) {
 769         SetErrorCode("Edge 1 could not be built");
 770         return false;
 771       }
 772
 773       edge_e2 = myBasicOperations->MakeLineTwoPnt(ve2, vi2);
 774       if (edge_e2.IsNull()) {
 775         SetErrorCode("Edge 2 could not be built");
 776         return false;
 777       }
 778
 779       edge_e1->GetLastFunction()->SetDescription("");
 780       edge_e2->GetLastFunction()->SetDescription("");
 781
 782       std::list<Handle(GEOM_Object)> edge_e_elist;
 783       edge_e_elist.push_back(arete_intersect_int);
 784       edge_e_elist.push_back(edge_e1);
 785       edge_e_elist.push_back(Handle(GEOM_Object)::DownCast(edges_e->Value(1)));
 786       edge_e_elist.push_back(edge_e2);
 787       wire_t = myShapesOperations->MakeWire(edge_e_elist, 1e-7);
 788       if (wire_t.IsNull()) {
 789         SetErrorCode("Impossible to build wire");
 790         return false;
 791       }
 792       wire_t->GetLastFunction()->SetDescription("");
 793       face_t = myShapesOperations->MakeFace(wire_t, false);
 794       if (face_t.IsNull()) {
 795         SetErrorCode("Impossible to build face");
 796         return false;
 797       }
 798       face_t->GetLastFunction()->SetDescription("");
 799
 800       theShapes.push_back(theShape);
 801       theShapes.push_back(vi1);
 802       theShapes.push_back(vi2);
 803       theShapes.push_back(ve1);
 804       theShapes.push_back(ve2);
 805       theShapes.push_back(edge_e1);
 806       theShapes.push_back(edge_e2);
 807       theShapes.push_back(wire_t);
 808       theShapes.push_back(face_t);
 809     }
 810     else {
 811       Handle(GEOM_Object) P1, P2, P3, P4, P5, P6;
 812       int idP1, idP2, idP3, idP4;
 813       int PZX, PZY;
 814       double ZX=0, ZY=0;
 815       std::vector<int> LX;
 816       std::vector<int> LY;
 817       Handle(GEOM_Object) box_e = my3DPrimOperations->MakeBoxDXDYDZ
 818         (theVertCylinderRadius, theVertCylinderRadius, theHoriCylinderRadius);
 819       box_e->GetLastFunction()->SetDescription("");
 820       box_e = myTransformOperations->TranslateDXDYDZ
 821         (box_e, -theVertCylinderRadius, -theVertCylinderRadius, 0);
 822       box_e->GetLastFunction()->SetDescription("");
 823
 824       aFunction = theShape->GetLastFunction();
 825       theDesc = aFunction->GetDescription();
 826       Handle(TColStd_HSequenceOfTransient) extremVertices =
 827         myShapesOperations->GetShapesOnBox(box_e, theShape, TopAbs_VERTEX, GEOMAlgo_ST_ONIN);
 828       // Recover previous description to get rid of Propagate dump
 829       aFunction->SetDescription(theDesc);
 830
 831       if (extremVertices.IsNull() || extremVertices->Length() == 0) {
 832         if (theRF == 0)
 833           SetErrorCode("Vertices on chamfer not found");
 834         else
 835           SetErrorCode("Vertices on fillet not found");
 836         return false;
 837       }
 838
 839       theShapes.push_back(theShape);
 840       theShapes.push_back(box_e);
 841       if (extremVertices->Length() != 6) {
 842         //           for (int i=1; i<=extremVertices->Length(); i++){
 843         //             theShapes.push_back(Handle(GEOM_Object)::DownCast(extremVertices->Value(i)));
 844         //           }
 845         //           Handle(GEOM_Object) aCompound = myShapesOperations->MakeCompound(theShapes);
 846         //           TopoDS_Shape aCompoundShape = aCompound->GetValue();
 847         //           theShape->GetLastFunction()->SetValue(aCompoundShape);
 848         SetErrorCode("Bad number of vertices on chamfer found");
 849         return false;
 850       }
 851
 852       for (int i=1; i<=extremVertices->Length(); i++){
 853         Handle(GEOM_Object) aV = Handle(GEOM_Object)::DownCast(extremVertices->Value(i));
 854         aV->GetLastFunction()->SetDescription("");
 855         gp_Pnt aP = BRep_Tool::Pnt(TopoDS::Vertex(aV->GetValue()));
 856
 857         if (Abs(aP.X()) <= Precision::Confusion()) {
 858           if (Abs(aP.Y()) - theR2 > Precision::Confusion()) {
 859             LX.push_back(i);
 860             if  (aP.Z()-ZX > Precision::Confusion()) {
 861               ZX = aP.Z();
 862               PZX = i;
 863             }
 864           }
 865         }
 866         else {
 867           if (Abs(aP.X()) - theR2 > Precision::Confusion()) {
 868             LY.push_back(i);
 869             if (aP.Z() - ZY > Precision::Confusion()) {
 870               ZY = aP.Z();
 871               PZY = i;
 872             }
 873           }
 874         }
 875       }
 876
 877       idP2 = PZX;
 878       idP4 = PZY;
 879       idP1 = LX.at(0);
 880       if (LX.at(0) == PZX)
 881         idP1 = LX.at(1);
 882       idP3 = LY.at(0);
 883       if (LY.at(0) == PZY)
 884         idP3 = LY.at(1);
 885
 886       P1 = Handle(GEOM_Object)::DownCast(extremVertices->Value(idP1));
 887       P2 = Handle(GEOM_Object)::DownCast(extremVertices->Value(idP2));
 888       P3 = Handle(GEOM_Object)::DownCast(extremVertices->Value(idP3));
 889       P4 = Handle(GEOM_Object)::DownCast(extremVertices->Value(idP4));
 890
 891       Handle(GEOM_Object) Cote_1 = myBasicOperations->MakeLineTwoPnt(P1, vi1);
 892       if (Cote_1.IsNull()) {
 893         SetErrorCode("Impossible to build edge in thickness");
 894         return false;
 895       }
 896       Cote_1->GetLastFunction()->SetDescription("");
 897
 898       Handle(GEOM_Object) Cote_2 = myBasicOperations->MakeLineTwoPnt(vi2, P3);
 899       if (Cote_2.IsNull()) {
 900         SetErrorCode("Impossible to build edge in thickness");
 901         return false;
 902       }
 903       Cote_2->GetLastFunction()->SetDescription("");
 904
 905       // edge_chan_princ = arete du chanfrein (ou raccord) sur le tuyau principal
 906       // edge_chan_inc = arete du chanfrein (ou raccord) sur le tuyau incident
 907       //         std::cerr << "Getting chamfer edge on main pipe" << std::endl;
 908       Handle(GEOM_Object) edge_chan_princ = myBlocksOperations->GetEdge(theShape, P1, P3);
 909       if (edge_chan_princ.IsNull()) {
 910         SetErrorCode("Impossible to find edge on main pipe");
 911         return false;
 912       }
 913       edge_chan_princ->GetLastFunction()->SetDescription("");
 914
 915       Handle(GEOM_Object) edge_chan_inc = myBlocksOperations->GetEdge(theShape, P2, P4);
 916       if (edge_chan_inc.IsNull()) {
 917         SetErrorCode("Impossible to find edge on incident pipe");
 918         return false;
 919       }
 920       edge_chan_inc->GetLastFunction()->SetDescription("");
 921
 922       std::list<Handle(GEOM_Object)> edgeList1;
 923       edgeList1.push_back(edge_chan_princ);
 924       edgeList1.push_back(Cote_1);
 925       edgeList1.push_back(arete_intersect_int);
 926       edgeList1.push_back(Cote_2);
 927
 928       //         std::cerr << "Creating wire 1" << std::endl;
 929       wire_t = myShapesOperations->MakeWire(edgeList1, 1e-7);
 930       if (wire_t.IsNull()) {
 931         SetErrorCode("Impossible to build wire");
 932         return false;
 933       }
 934       wire_t->GetLastFunction()->SetDescription("");
 935
 936       //         std::cerr << "Creating face 1" << std::endl;
 937       face_t = myShapesOperations->MakeFace(wire_t, false);
 938       if (face_t.IsNull()) {
 939         SetErrorCode("Impossible to build face");
 940         return false;
 941       }
 942       face_t->GetLastFunction()->SetDescription("");
 943       theShapes.push_back(face_t);
 944
 945       gp_Pnt aP2 = BRep_Tool::Pnt(TopoDS::Vertex(P2->GetValue()));
 946       gp_Pnt aP5 = BRep_Tool::Pnt(TopoDS::Vertex(vi1->GetValue()));
 947       double deltaZ = aP2.Z() - aP5.Z();
 948       //         std::cerr << "Creating new point from vi1 with deltaZ = " << deltaZ << std::endl;
 949       Handle(GEOM_Object) P5bis = myTransformOperations->TranslateDXDYDZCopy(vi1, 0, 0, deltaZ);
 950       if (P5bis.IsNull()) {
 951         SetErrorCode("Impossible to translate vertex");
 952         return false;
 953       }
 954       P5bis->GetLastFunction()->SetDescription("");
 955
 956       gp_Pnt aP4 = BRep_Tool::Pnt(TopoDS::Vertex(P4->GetValue()));
 957       gp_Pnt aP6 = BRep_Tool::Pnt(TopoDS::Vertex(vi2->GetValue()));
 958       deltaZ = aP4.Z() - aP6.Z();
 959       //         std::cerr << "Creating new point from vi2 with deltaZ = " << deltaZ << std::endl;
 960       Handle(GEOM_Object) P6bis = myTransformOperations->TranslateDXDYDZCopy(vi2, 0, 0, deltaZ);
 961       if (P6bis.IsNull()) {
 962         SetErrorCode("Impossible to translate vertex");
 963         return false;
 964       }
 965       P6bis->GetLastFunction()->SetDescription("");
 966
 967       //         std::cerr << "Creating new line 1 from 2 previous points" << std::endl;
 968       Handle(GEOM_Object) Cote_3 = myBasicOperations->MakeLineTwoPnt(P5bis, P2);
 969       if (Cote_3.IsNull()) {
 970         SetErrorCode("Impossible to build edge in thickness");
 971         return false;
 972       }
 973       Cote_3->GetLastFunction()->SetDescription("");
 974
 975       //         std::cerr << "Creating new line 2 from 2 previous points" << std::endl;
 976       Handle(GEOM_Object) Cote_4 = myBasicOperations->MakeLineTwoPnt(P6bis, P4);
 977       if (Cote_4.IsNull()) {
 978         SetErrorCode("Impossible to build edge in thickness");
 979         return false;
 980       }
 981       Cote_4->GetLastFunction()->SetDescription("");
 982
 983       //         std::cerr << "Creating new line 3 from 2 previous points" << std::endl;
 984       Handle(GEOM_Object) Cote_5 = myBasicOperations->MakeLineTwoPnt(P5bis, P6bis);
 985       if (Cote_4.IsNull()) {
 986         SetErrorCode("Impossible to build edge in thickness");
 987         return false;
 988       }
 989       Cote_5->GetLastFunction()->SetDescription("");
 990
 991       //std::list<Handle(GEOM_Object)> edgeList2;
 992       //edgeList2.push_back(edge_chan_inc);
 993       //edgeList2.push_back(Cote_3);
 994       //edgeList2.push_back(Cote_5);
 995       //edgeList2.push_back(Cote_4);
 996       //         std::cerr << "Creating wire 2" << std::endl;
 997       //wire_t2 = myShapesOperations->MakeWire(edgeList2, 1e-7);
 998       //if (wire_t2.IsNull()) {
 999       //  SetErrorCode("Impossible to build wire");
1000       //  return false;
1001       //}
1002       //wire_t2->GetLastFunction()->SetDescription("");
1003       //         std::cerr << "Creating face 2" << std::endl;
1004       //face_t2 = myShapesOperations->MakeFace(wire_t2, false);
1005
1006       // Mantis issue 0021682
1007       face_t2 = my3DPrimOperations->MakePrismVecH(edge_chan_inc, Cote_4, - (theR2 + theW2));
1008       //face_t2 = my3DPrimOperations->MakePrismVecH(edge_chan_inc, Cote_4, - 2.0*theR2);
1009       if (face_t2.IsNull()) {
1010         SetErrorCode("Impossible to build face");
1011         return false;
1012       }
1013       face_t2->GetLastFunction()->SetDescription("");
1014       theShapes.push_back(face_t2);
1015     }
1016
1017     // Planes
1018     Handle(GEOM_Object) aP0 = myBasicOperations->MakePointXYZ(0, 0, 0);
1019     Handle(GEOM_Object) aVZ = myBasicOperations->MakeVectorDXDYDZ(0, 0, 1);
1020     Handle(GEOM_Object) aVXZ = myBasicOperations->MakeVectorDXDYDZ(aR1Ext, 0, 0.5*(theL1+theVertCylinderRadius));
1021     Handle(GEOM_Object) aPlnOZ = myBasicOperations->MakePlanePntVec(aP0, aVZ, aSize);
1022     Handle(GEOM_Object) aPlnOXZ = myBasicOperations->MakePlanePntVec(aP0, aVXZ, aSize);
1023     aP0->GetLastFunction()->SetDescription("");
1024     aVZ->GetLastFunction()->SetDescription("");
1025     aVXZ->GetLastFunction()->SetDescription("");
1026     aPlnOZ->GetLastFunction()->SetDescription("");
1027     aPlnOXZ->GetLastFunction()->SetDescription("");
1028     theShapes.push_back(aPlnOZ);
1029     theShapes.push_back(aPlnOXZ);
1030
1031     // Partition
1032     Handle(TColStd_HSequenceOfTransient) partitionShapes = new TColStd_HSequenceOfTransient;
1033     Handle(TColStd_HSequenceOfTransient) theTools = new TColStd_HSequenceOfTransient;
1034     Handle(TColStd_HSequenceOfTransient) theKeepInside = new TColStd_HSequenceOfTransient;
1035     Handle(TColStd_HSequenceOfTransient) theRemoveInside = new TColStd_HSequenceOfTransient;
1036     Handle(TColStd_HArray1OfInteger) theMaterials;
1037
1038     partitionShapes->Append(theShape);
1039     theTools->Append(aPlnOZ);
1040     if (Abs(aR1Ext - aR2Ext) > Precision::Confusion())
1041       theTools->Append(aPlnOXZ);
1042     theTools->Append(face_t);
1043     if (!isNormal)
1044       theTools->Append(face_t2);
1045
1046     Te3 = myBooleanOperations->MakePartition
1047               (partitionShapes, theTools, theKeepInside, theRemoveInside,
1048               TopAbs_SOLID, false, theMaterials, 0, false);
1049     if (Te3.IsNull()) {
1050       SetErrorCode("Impossible to build partition of TShape");
1051       return false;
1052     }
1053     Te3->GetLastFunction()->SetDescription("");
1054
1055     // Last verification: result should be a block
1056     std::list<GEOMImpl_IBlocksOperations::BCError> errList;
1057     if (!myBlocksOperations->CheckCompoundOfBlocks(Te3,errList)) {
1058       SetErrorCode("TShape is not a compound of block");
1059       return false;
1060     }
1061
1062 //     // BEGIN Compound of created shapes - Only for debug purpose
1063 //     theShapes.clear();
1064 //     theShapes.push_back(theShape);
1065 //     theShapes.push_back(aPlnOZ);
1066 //     if (Abs(aR1Ext - aR2Ext) > Precision::Confusion() )
1067 //       theShapes.push_back(aPlnOXZ);
1068 //     theShapes.push_back(face_t);
1069 //     if (!isNormal)
1070 //       theShapes.push_back(face_t2);
1071 //
1072 //     Handle(GEOM_Object) aCompound = myShapesOperations->MakeCompound(theShapes);
1073 //     TopoDS_Shape aCompoundShape = aCompound->GetValue();
1074 //     theShape->GetLastFunction()->SetValue(aCompoundShape);
1075 //     // END Compound of created shapes - Only for debug purpose
1076
1077     TopoDS_Shape aShape = Te3->GetValue();
1078     theShape->GetLastFunction()->SetValue(aShape);
1079   } catch (Standard_Failure) {
1080     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1081     SetErrorCode(aFail->GetMessageString());
1082     return false;
1083   }
1084
1085   SetErrorCode(OK);
1086   return true;
1087 }
1088
1089 // Mirror and glue faces
1090 bool GEOMImpl_IAdvancedOperations::MakePipeTShapeMirrorAndGlue(Handle(GEOM_Object) theShape,
1091                                                                double theR1, double theW1, double theL1,
1092                                                                double theR2, double theW2, double theL2)
1093 {
1094   SetErrorCode(KO);
1095
1096   // Useful values
1097   double aSize = 2*(theL1 + theL2);
1098   double aR1Ext = theR1 + theW1;
1099
1100   // Planes
1101   Handle(GEOM_Object) aP0 = myBasicOperations->MakePointXYZ(0, 0, 0);
1102   aP0->GetLastFunction()->SetDescription("");
1103   Handle(GEOM_Object) aVX = myBasicOperations->MakeVectorDXDYDZ(1, 0, 0);
1104   Handle(GEOM_Object) aVY = myBasicOperations->MakeVectorDXDYDZ(0, 1, 0);
1105   aVX->GetLastFunction()->SetDescription("");
1106   aVY->GetLastFunction()->SetDescription("");
1107   Handle(GEOM_Object) aPlane_OX = myBasicOperations->MakePlanePntVec(aP0, aVX, 2*(aR1Ext + theL2));
1108   Handle(GEOM_Object) aPlane_OY = myBasicOperations->MakePlanePntVec(aP0, aVY, aSize);
1109   aPlane_OX->GetLastFunction()->SetDescription("");
1110   aPlane_OY->GetLastFunction()->SetDescription("");
1111
1112   Handle(GEOM_Object) Te4 = myTransformOperations->MirrorPlaneCopy(theShape, aPlane_OX);
1113   if (Te4.IsNull()) {
1114     SetErrorCode("Impossible to build mirror of quarter TShape");
1115     return false;
1116   }
1117
1118   Handle(GEOM_Object) Te5 = myTransformOperations->MirrorPlaneCopy(theShape, aPlane_OY);
1119   if (Te5.IsNull()) {
1120     SetErrorCode("Impossible to build mirror of half TShape");
1121     return false;
1122   }
1123
1124   Handle(GEOM_Object) Te6 = myTransformOperations->MirrorPlaneCopy(Te4, aPlane_OY);
1125   if (Te6.IsNull()) {
1126     SetErrorCode("Impossible to build mirror of half TShape");
1127     return false;
1128   }
1129
1130   std::list<Handle(GEOM_Object)> aShapesList;
1131   aShapesList.push_back(theShape);
1132   aShapesList.push_back(Te4);
1133   aShapesList.push_back(Te5);
1134   aShapesList.push_back(Te6);
1135   Handle(GEOM_Object) Te7 = myShapesOperations->MakeCompound(aShapesList);
1136   if (Te7.IsNull()) {
1137     SetErrorCode("Impossible to build compound");
1138     return false;
1139   }
1140
1141   Handle(GEOM_Object) Te8 = myShapesOperations->MakeGlueFaces(Te7, 1e-7, true);
1142   if (Te8.IsNull()) {
1143     SetErrorCode("Impossible to glue faces of TShape");
1144     return false;
1145   }
1146
1147   TopoDS_Shape aShape = Te8->GetValue();
1148
1149   theShape->GetLastFunction()->SetValue(aShape);
1150
1151   Te4->GetLastFunction()->SetDescription("");
1152   Te5->GetLastFunction()->SetDescription("");
1153   Te6->GetLastFunction()->SetDescription("");
1154   Te7->GetLastFunction()->SetDescription("");
1155   Te8->GetLastFunction()->SetDescription("");
1156
1157   SetErrorCode(OK);
1158   return true;
1159 }
1160
1161 //=============================================================================
1162 /*!
1163  *  MakePipeTShape
1164  *  \brief Create a T-shape object with specified caracteristics for the main and
1165  *         the incident pipes (radius, width, half-length).
1166  *         Center of the shape is (0,0,0). The main plane of the T-shape is XOY.
1167  *  \param theR1 Internal radius of main pipe
1168  *  \param theW1 Width of main pipe
1169  *  \param theL1 Half-length of main pipe
1170  *  \param theR2 Internal radius of incident pipe (R2 < R1)
1171  *  \param theW2 Width of incident pipe (R2+W2 < R1+W1)
1172  *  \param theL2 Half-length of incident pipe
1173  *  \param theHexMesh Boolean indicating if shape is prepared for hex mesh
1174  *  \return List of GEOM_Objects, containing the created shape and propagation groups.
1175  */
1176 //=============================================================================
1177 Handle(TColStd_HSequenceOfTransient)
1178   GEOMImpl_IAdvancedOperations::MakePipeTShape(double theR1, double theW1, double theL1,
1179                                                double theR2, double theW2, double theL2,
1180                                                bool theHexMesh)
1181 {
1182   MESSAGE("GEOMImpl_IAdvancedOperations::MakePipeTShape");
1183   SetErrorCode(KO);
1184   //Add a new object
1185   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_TSHAPE);
1186
1187   //Add a new shape function with parameters
1188   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_PipeTShapeDriver::GetID(), TSHAPE_BASIC);
1189   if (aFunction.IsNull()) return NULL;
1190
1191   //Check if the function is set correctly
1192   if (aFunction->GetDriverGUID() != GEOMImpl_PipeTShapeDriver::GetID()) return NULL;
1193
1194   GEOMImpl_IPipeTShape aData(aFunction);
1195
1196   aData.SetR1(theR1);
1197   aData.SetW1(theW1);
1198   aData.SetL1(theL1);
1199   aData.SetR2(theR2);
1200   aData.SetW2(theW2);
1201   aData.SetL2(theL2);
1202   aData.SetHexMesh(theHexMesh);
1203
1204   //Compute the resulting value
1205   try {
1206 #if OCC_VERSION_LARGE > 0x06010000
1207     OCC_CATCH_SIGNALS;
1208 #endif
1209     if (!GetSolver()->ComputeFunction(aFunction)) {
1210       SetErrorCode("TShape driver failed");
1211       return NULL;
1212     }
1213     if (theHexMesh) {
1214       if (!MakePipeTShapePartition(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
1215         return NULL;
1216       if (!MakePipeTShapeMirrorAndGlue(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
1217         return NULL;
1218     }
1219   } catch (Standard_Failure) {
1220     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1221     SetErrorCode(aFail->GetMessageString());
1222     return NULL;
1223   }
1224
1225
1226   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
1227   aSeq->Append(aShape);
1228
1229   if (theHexMesh) {
1230     /*
1231      * Get the groups: BEGIN
1232      */
1233     try {
1234       if (!MakeGroups(aShape, TSHAPE_BASIC, theR1, theW1, theL1, theR2, theW2, theL2, aSeq, gp_Trsf()))
1235         return NULL;
1236     }
1237     catch (Standard_Failure) {
1238       Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1239       SetErrorCode(aFail->GetMessageString());
1240       return NULL;
1241     }
1242
1243     TCollection_AsciiString aListRes, anEntry;
1244     // Iterate over the sequence aSeq
1245     Standard_Integer aNbGroups = aSeq->Length();
1246     Standard_Integer i = 2;
1247     for (; i <= aNbGroups; i++) {
1248       Handle(Standard_Transient) anItem = aSeq->Value(i);
1249       if (anItem.IsNull()) continue;
1250       Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(anItem);
1251       if (aGroup.IsNull()) continue;
1252       //Make a Python command
1253       TDF_Tool::Entry(aGroup->GetEntry(), anEntry);
1254       aListRes += anEntry + ", ";
1255     }
1256
1257     aListRes.Trunc(aListRes.Length() - 2);
1258
1259     //Make a Python command
1260     GEOM::TPythonDump(aFunction)
1261       << "[" << aShape << ", " << aListRes.ToCString() << "] = geompy.MakePipeTShape("
1262       << theR1 << ", " << theW1 << ", " << theL1 << ", " << theR2 << ", " << theW2 << ", " << theL2 << ", "
1263       << theHexMesh << ")";
1264   }
1265   /*
1266    * Get the groups: END
1267    */
1268   else {
1269     //Make a Python command
1270     GEOM::TPythonDump(aFunction)
1271       << "[" << aShape << "] = geompy.MakePipeTShape(" << theR1 << ", " << theW1 << ", "
1272       << theL1 << ", " << theR2 << ", " << theW2 << ", " << theL2 << ", " << theHexMesh << ")";
1273   }
1274
1275   SetErrorCode(OK);
1276
1277   return aSeq;
1278 }
1279
1280 //=============================================================================
1281 /*!
1282  *  MakePipeTShapeWithPosition
1283  *  Create a T-shape object with specified caracteristics for the main and
1284  *  the incident pipes (radius, width, half-length).
1285  *  The extremities of the main pipe are located on junctions points P1 and P2.
1286  *  The extremity of the incident pipe is located on junction point P3.
1287  *  \param theR1 Internal radius of main pipe
1288  *  \param theW1 Width of main pipe
1289  *  \param theL1 Half-length of main pipe
1290  *  \param theR2 Internal radius of incident pipe (R2 < R1)
1291  *  \param theW2 Width of incident pipe (R2+W2 < R1+W1)
1292  *  \param theL2 Half-length of incident pipe
1293  *  \param theHexMesh Boolean indicating if shape is prepared for hex mesh
1294  *  \param theP1 1st junction point of main pipe
1295  *  \param theP2 2nd junction point of main pipe
1296  *  \param theP3 Junction point of incident pipe
1297  *  \return List of GEOM_Objects, containing the created shape and propagation groups..
1298  */
1299 //=============================================================================
1300 Handle(TColStd_HSequenceOfTransient)
1301 GEOMImpl_IAdvancedOperations::MakePipeTShapeWithPosition(double theR1, double theW1, double theL1,
1302                                                          double theR2, double theW2, double theL2,
1303                                                          bool theHexMesh,
1304                                                          Handle(GEOM_Object) theP1,
1305                                                          Handle(GEOM_Object) theP2,
1306                                                          Handle(GEOM_Object) theP3)
1307 {
1308   SetErrorCode(KO);
1309   //Add a new object
1310   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_TSHAPE);
1311   /////////////////
1312   // TSHAPE CODE
1313   /////////////////
1314   //Add a new shape function with parameters
1315   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_PipeTShapeDriver::GetID(), TSHAPE_BASIC);
1316   if (aFunction.IsNull()) return NULL;
1317
1318   //Check if the function is set correctly
1319   if (aFunction->GetDriverGUID() != GEOMImpl_PipeTShapeDriver::GetID()) return NULL;
1320
1321   // Check new position
1322   if (!CheckCompatiblePosition(theL1, theL2, theP1, theP2, theP3, 0.01)) {
1323     return NULL;
1324   }
1325
1326   GEOMImpl_IPipeTShape aData(aFunction);
1327
1328   aData.SetR1(theR1);
1329   aData.SetW1(theW1);
1330   aData.SetL1(theL1);
1331   aData.SetR2(theR2);
1332   aData.SetW2(theW2);
1333   aData.SetL2(theL2);
1334   aData.SetHexMesh(theHexMesh);
1335
1336   //Compute the resulting value
1337   try {
1338 #if OCC_VERSION_LARGE > 0x06010000
1339     OCC_CATCH_SIGNALS;
1340 #endif
1341     if (!GetSolver()->ComputeFunction(aFunction)) {
1342       SetErrorCode("TShape driver failed");
1343       return NULL;
1344     }
1345   } catch (Standard_Failure) {
1346     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1347     SetErrorCode(aFail->GetMessageString());
1348     return NULL;
1349   }
1350
1351   if (theHexMesh) {
1352     if (!MakePipeTShapePartition(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
1353       return NULL;
1354     if (!MakePipeTShapeMirrorAndGlue(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
1355       return NULL;
1356   }
1357
1358   TopoDS_Shape Te = aShape->GetValue();
1359
1360   // Set Position
1361   gp_Trsf aTrsf = GetPositionTrsf(theL1, theL2, theP1, theP2, theP3);
1362   BRepBuilderAPI_Transform aTransformation(Te, aTrsf, Standard_False);
1363   TopoDS_Shape aTrsf_Shape = aTransformation.Shape();
1364   aFunction->SetValue(aTrsf_Shape);
1365   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
1366   aSeq->Append(aShape);
1367
1368   if (theHexMesh) {
1369     //
1370     // Get the groups: BEGIN
1371     //
1372     try {
1373       if (!MakeGroups(aShape,TSHAPE_BASIC, theR1, theW1, theL1, theR2, theW2, theL2, aSeq, aTrsf)) {
1374         return NULL;
1375       }
1376     }
1377     catch (Standard_Failure) {
1378       Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1379       SetErrorCode(aFail->GetMessageString());
1380       return NULL;
1381     }
1382
1383     TCollection_AsciiString aListRes, anEntry;
1384     // Iterate over the sequence aSeq
1385     Standard_Integer aNbGroups = aSeq->Length();
1386     Standard_Integer i = 2;
1387     for (; i <= aNbGroups; i++) {
1388       Handle(Standard_Transient) anItem = aSeq->Value(i);
1389       if (anItem.IsNull()) continue;
1390       Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(anItem);
1391       if (aGroup.IsNull()) continue;
1392       //Make a Python command
1393       TDF_Tool::Entry(aGroup->GetEntry(), anEntry);
1394       aListRes += anEntry + ", ";
1395     }
1396
1397     aListRes.Trunc(aListRes.Length() - 2);
1398
1399     //Make a Python command
1400     GEOM::TPythonDump(aFunction)
1401       << "[" << aShape << ", " << aListRes.ToCString() << "] = geompy.MakePipeTShape("
1402       << theR1 << ", " << theW1 << ", " << theL1 << ", " << theR2 << ", " << theW2 << ", " << theL2 << ", "
1403       << theHexMesh << ", " << theP1 << ", " << theP2 << ", " << theP3 << ")";
1404   }
1405   //
1406   // Get the groups: END
1407   //
1408
1409   else {
1410     //Make a Python command
1411     GEOM::TPythonDump(aFunction)
1412       << "[" << aShape << "] = geompy.MakePipeTShape(" << theR1 << ", " << theW1 << ", "
1413       << theL1 << ", " << theR2 << ", " << theW2 << ", " << theL2 << ", " << theHexMesh << ", " << theP1
1414       << ", " << theP2 << ", " << theP3 << ")";
1415   }
1416
1417   SetErrorCode(OK);
1418
1419   return aSeq;
1420 }
1421
1422 //=============================================================================
1423 /*!
1424  *  MakePipeTShapeChamfer
1425  *  Create a T-shape object with specified caracteristics for the main and
1426  *  the incident pipes (radius, width, half-length). A chamfer is created
1427  *  on the junction of the pipes.
1428  *  Center of the shape is (0,0,0). The main plane of the T-shape is XOY.
1429  *  \param theR1 Internal radius of main pipe
1430  *  \param theW1 Width of main pipe
1431  *  \param theL1 Half-length of main pipe
1432  *  \param theR2 Internal radius of incident pipe (R2 < R1)
1433  *  \param theW2 Width of incident pipe (R2+W2 < R1+W1)
1434  *  \param theL2 Half-length of incident pipe
1435  *  \param theH Height of chamfer.
1436  *  \param theW Width of chamfer.
1437  *  \param theHexMesh Boolean indicating if shape is prepared for hex mesh
1438  *  \return List of GEOM_Objects, containing the created shape and propagation groups.
1439  */
1440 //=============================================================================
1441 Handle(TColStd_HSequenceOfTransient)
1442 GEOMImpl_IAdvancedOperations::MakePipeTShapeChamfer(double theR1, double theW1, double theL1,
1443                                                     double theR2, double theW2, double theL2,
1444                                                     double theH, double theW,
1445                                                     bool theHexMesh)
1446 {
1447   SetErrorCode(KO);
1448   //Add a new object
1449   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_TSHAPE);
1450   //Add a new shape function with parameters
1451   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_PipeTShapeDriver::GetID(), TSHAPE_CHAMFER);
1452   if (aFunction.IsNull()) return NULL;
1453
1454   //Check if the function is set correctly
1455   if (aFunction->GetDriverGUID() != GEOMImpl_PipeTShapeDriver::GetID()) return NULL;
1456
1457   GEOMImpl_IPipeTShape aData(aFunction);
1458
1459   aData.SetR1(theR1);
1460   aData.SetW1(theW1);
1461   aData.SetL1(theL1);
1462   aData.SetR2(theR2);
1463   aData.SetW2(theW2);
1464   aData.SetL2(theL2);
1465   aData.SetH(theH);
1466   aData.SetW(theW);
1467   aData.SetHexMesh(theHexMesh);
1468
1469   //Compute the resulting value
1470   try {
1471 #if OCC_VERSION_LARGE > 0x06010000
1472     OCC_CATCH_SIGNALS;
1473 #endif
1474     if (!GetSolver()->ComputeFunction(aFunction)) {
1475       SetErrorCode("TShape driver failed");
1476       return NULL;
1477     }
1478   } catch (Standard_Failure) {
1479     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1480     SetErrorCode(aFail->GetMessageString());
1481     return NULL;
1482   }
1483
1484   // BEGIN of chamfer
1485   TopoDS_Shape aShapeShape = aShape->GetValue();
1486   TopTools_IndexedMapOfShape anEdgesIndices;
1487   TopExp::MapShapes(aShapeShape, anEdgesIndices);
1488   // Common edges on external cylinders
1489   Handle(GEOM_Object) box_e;
1490   if (theHexMesh) {
1491     box_e = my3DPrimOperations->MakeBoxDXDYDZ(theR2+theW2, theR2+theW2, theR1+theW1);
1492   }
1493   else {
1494     box_e = my3DPrimOperations->MakeBoxDXDYDZ(2*(theR2+theW2), 2*(theR2+theW2), theR1+theW1);
1495   }
1496   box_e->GetLastFunction()->SetDescription("");
1497   box_e = myTransformOperations->TranslateDXDYDZ(box_e, -(theR2+theW2), -(theR2+theW2), 0);
1498   box_e->GetLastFunction()->SetDescription("");
1499
1500   Handle(TColStd_HSequenceOfInteger) edges_e =
1501     myShapesOperations->GetShapesOnBoxIDs(box_e, aShape, TopAbs_EDGE, GEOMAlgo_ST_IN);
1502   box_e->GetLastFunction()->SetDescription("");
1503
1504   if (edges_e.IsNull() || edges_e->Length() == 0) {
1505     SetErrorCode("External edges not found");
1506     return NULL;
1507   }
1508   int nbEdgesInChamfer = 0;
1509   std::list<int> theEdges;
1510   for (int i=1; i<=edges_e->Length();i++) {
1511     int edgeID = edges_e->Value(i);
1512     TopoDS_Shape theEdge = anEdgesIndices.FindKey(edgeID);
1513     TopExp_Explorer Ex(theEdge,TopAbs_VERTEX);
1514     int iv=0;
1515     while (Ex.More()) {
1516       iv ++;
1517       gp_Pnt aPt = BRep_Tool::Pnt(TopoDS::Vertex(Ex.Current()));
1518       if (Abs(aPt.Z() - (theR1+theW1)) <= Precision::Confusion()) {
1519         nbEdgesInChamfer ++;
1520         theEdges.push_back(edgeID);
1521       }
1522       Ex.Next();
1523     }
1524     if (theHexMesh && nbEdgesInChamfer == 1)
1525       break;
1526   }
1527   Handle(GEOM_Object) aChamfer;
1528   try {
1529     aChamfer = myLocalOperations->MakeChamferEdges(aShape, theW, theH, theEdges);
1530   }
1531   catch (Standard_Failure) {
1532     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1533     SetErrorCode(aFail->GetMessageString());
1534     return NULL;
1535   }
1536   if (aChamfer.IsNull()) {
1537     SetErrorCode("Chamfer can not be computed on the given shape with the given parameters");
1538     return NULL;
1539   }
1540   aChamfer->GetLastFunction()->SetDescription("");
1541
1542   TopoDS_Shape aChamferShape = aChamfer->GetValue();
1543   aFunction->SetValue(aChamferShape);
1544   // END of chamfer
1545
1546   //   bool doMesh = false;
1547   if (theHexMesh) {
1548     //        doMesh = true;
1549     if (!MakePipeTShapePartition(aShape, theR1, theW1, theL1, theR2, theW2, theL2, theH, theW, 0, false)) {
1550       MESSAGE("PipeTShape partition failed");
1551       //            doMesh = false;
1552       return NULL;
1553     }
1554     if (!MakePipeTShapeMirrorAndGlue(aShape, theR1, theW1, theL1, theR2, theW2, theL2)) {
1555       MESSAGE("PipeTShape mirrors and glue failed");
1556       //          doMesh = false;
1557       return NULL;
1558     }
1559   }
1560
1561   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
1562   aSeq->Append(aShape);
1563
1564   //    if (doMesh) {
1565   if (theHexMesh) {
1566     //
1567     //         Get the groups: BEGIN
1568     //
1569     //if (!MakeGroups(aShape, TSHAPE_CHAMFER, theR1, theW1, theL1, theR2, theW2, theL2, aSeq, gp_Trsf())) {
1570     //  //Make a Python command
1571     //  GEOM::TPythonDump(aFunction)
1572     //    << "[" << aShape << "] = geompy.MakePipeTShapeChamfer(" << theR1 << ", " << theW1
1573     //    << ", " << theL1 << ", " << theR2 << ", " << theW2 << ", " << theL2 << ", " << theH << ", " << theW
1574     //    << ", " << theHexMesh << ")";
1575     //}
1576     //else {
1577     try {
1578       if (!MakeGroups(aShape, TSHAPE_CHAMFER, theR1, theW1, theL1, theR2, theW2, theL2, aSeq, gp_Trsf()))
1579         return NULL;
1580     }
1581     catch (Standard_Failure) {
1582       Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1583       SetErrorCode(aFail->GetMessageString());
1584       return NULL;
1585     }
1586
1587     TCollection_AsciiString aListRes, anEntry;
1588     // Iterate over the sequence aSeq
1589     Standard_Integer aNbGroups = aSeq->Length();
1590     Standard_Integer i = 2;
1591     for (; i <= aNbGroups; i++) {
1592       Handle(Standard_Transient) anItem = aSeq->Value(i);
1593       if (anItem.IsNull()) continue;
1594       Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(anItem);
1595       if (aGroup.IsNull()) continue;
1596       //Make a Python command
1597       TDF_Tool::Entry(aGroup->GetEntry(), anEntry);
1598       aListRes += anEntry + ", ";
1599     }
1600
1601     aListRes.Trunc(aListRes.Length() - 2);
1602
1603     //Make a Python command
1604     GEOM::TPythonDump(aFunction)
1605       << "[" << aShape << ", " << aListRes.ToCString()
1606       << "] = geompy.MakePipeTShapeChamfer(" << theR1 << ", " << theW1 << ", " << theL1 << ", " << theR2
1607       << ", " << theW2 << ", " << theL2 << ", " << theH << ", " << theW << ", " << theHexMesh << ")";
1608     //}
1609   }
1610   //
1611   //     Get the groups: END
1612   //
1613   else {
1614     //Make a Python command
1615     GEOM::TPythonDump(aFunction)
1616       << "[" << aShape << "] = geompy.MakePipeTShapeChamfer(" << theR1 << ", " << theW1
1617       << ", " << theL1 << ", " << theR2 << ", " << theW2 << ", " << theL2 << ", " << theH << ", " << theW
1618       << ", " << theHexMesh << ")";
1619   }
1620
1621   SetErrorCode(OK);
1622
1623   return aSeq;
1624 }
1625
1626 //=============================================================================
1627 /*!
1628  *  MakePipeTShapeChamferWithPosition
1629  *  Create a T-shape object with specified caracteristics for the main and
1630  *  the incident pipes (radius, width, half-length). A chamfer is created
1631  *  on the junction of the pipes.
1632  *  The extremities of the main pipe are located on junctions points P1 and P2.
1633  *  The extremity of the incident pipe is located on junction point P3.
1634  *  \param theR1 Internal radius of main pipe
1635  *  \param theW1 Width of main pipe
1636  *  \param theL1 Half-length of main pipe
1637  *  \param theR2 Internal radius of incident pipe (R2 < R1)
1638  *  \param theW2 Width of incident pipe (R2+W2 < R1+W1)
1639  *  \param theL2 Half-length of incident pipe
1640  *  \param theH Height of chamfer.
1641  *  \param theW Width of chamfer.
1642  *  \param theHexMesh Boolean indicating if shape is prepared for hex mesh
1643  *  \param theP1 1st junction point of main pipe
1644  *  \param theP2 2nd junction point of main pipe
1645  *  \param theP3 Junction point of incident pipe
1646  *  \return List of GEOM_Objects, containing the created shape and propagation groups.
1647  */
1648 //=============================================================================
1649 Handle(TColStd_HSequenceOfTransient)
1650 GEOMImpl_IAdvancedOperations::MakePipeTShapeChamferWithPosition(double theR1, double theW1, double theL1,
1651                                                                 double theR2, double theW2, double theL2,
1652                                                                 double theH, double theW,
1653                                                                 bool theHexMesh,
1654                                                                 Handle(GEOM_Object) theP1,
1655                                                                 Handle(GEOM_Object) theP2,
1656                                                                 Handle(GEOM_Object) theP3)
1657 {
1658   SetErrorCode(KO);
1659   //Add a new object
1660   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_TSHAPE);
1661   //Add a new shape function with parameters
1662   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_PipeTShapeDriver::GetID(), TSHAPE_CHAMFER);
1663   if (aFunction.IsNull()) return NULL;
1664
1665   //Check if the function is set correctly
1666   if (aFunction->GetDriverGUID() != GEOMImpl_PipeTShapeDriver::GetID()) return NULL;
1667
1668   // Check new position
1669   if (!CheckCompatiblePosition(theL1, theL2, theP1, theP2, theP3, 0.01)) {
1670     return NULL;
1671   }
1672
1673   GEOMImpl_IPipeTShape aData(aFunction);
1674
1675   aData.SetR1(theR1);
1676   aData.SetW1(theW1);
1677   aData.SetL1(theL1);
1678   aData.SetR2(theR2);
1679   aData.SetW2(theW2);
1680   aData.SetL2(theL2);
1681   aData.SetH(theH);
1682   aData.SetW(theW);
1683   aData.SetHexMesh(theHexMesh);
1684
1685   //Compute the resulting value
1686   try {
1687 #if OCC_VERSION_LARGE > 0x06010000
1688     OCC_CATCH_SIGNALS;
1689 #endif
1690     if (!GetSolver()->ComputeFunction(aFunction)) {
1691       SetErrorCode("TShape driver failed");
1692       return NULL;
1693     }
1694   } catch (Standard_Failure) {
1695     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1696     SetErrorCode(aFail->GetMessageString());
1697     return NULL;
1698   }
1699
1700   // BEGIN of chamfer
1701   TopoDS_Shape aShapeShape = aShape->GetValue();
1702   TopTools_IndexedMapOfShape anEdgesIndices;
1703   TopExp::MapShapes(aShapeShape, anEdgesIndices);
1704   // Common edges on external cylinders
1705   Handle(GEOM_Object) box_e;
1706   if (theHexMesh) {
1707     box_e = my3DPrimOperations->MakeBoxDXDYDZ(theR2+theW2, theR2+theW2, theR1+theW1);
1708   }
1709   else {
1710     box_e = my3DPrimOperations->MakeBoxDXDYDZ(2*(theR2+theW2), 2*(theR2+theW2), theR1+theW1);
1711   }
1712   box_e->GetLastFunction()->SetDescription("");
1713   box_e = myTransformOperations->TranslateDXDYDZ(box_e, -(theR2+theW2), -(theR2+theW2), 0);
1714   box_e->GetLastFunction()->SetDescription("");
1715
1716   Handle(TColStd_HSequenceOfInteger) edges_e =
1717     myShapesOperations->GetShapesOnBoxIDs(box_e, aShape, TopAbs_EDGE, GEOMAlgo_ST_IN);
1718   box_e->GetLastFunction()->SetDescription("");
1719
1720   if (edges_e.IsNull() || edges_e->Length() == 0) {
1721     SetErrorCode("External edges not found");
1722     return NULL;
1723   }
1724   int nbEdgesInChamfer = 0;
1725   std::list<int> theEdges;
1726   for (int i=1; i<=edges_e->Length();i++) {
1727     int edgeID = edges_e->Value(i);
1728     TopoDS_Shape theEdge = anEdgesIndices.FindKey(edgeID);
1729     TopExp_Explorer Ex(theEdge,TopAbs_VERTEX);
1730     while (Ex.More()) {
1731       gp_Pnt aPt = BRep_Tool::Pnt(TopoDS::Vertex(Ex.Current()));
1732       if (Abs(aPt.Z() - (theR1+theW1)) <= Precision::Confusion()) {
1733         nbEdgesInChamfer ++;
1734         theEdges.push_back(edgeID);
1735       }
1736       Ex.Next();
1737     }
1738     if (theHexMesh && nbEdgesInChamfer == 1)
1739       break;
1740   }
1741   Handle(GEOM_Object) aChamfer;
1742   try {
1743     aChamfer = myLocalOperations->MakeChamferEdges(aShape, theW, theH, theEdges);
1744   }
1745   catch (Standard_Failure) {
1746     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1747     SetErrorCode(aFail->GetMessageString());
1748     return NULL;
1749   }
1750   if (aChamfer.IsNull()) {
1751     SetErrorCode("Chamfer can not be computed on the given shape with the given parameters");
1752     return NULL;
1753   }
1754   aChamfer->GetLastFunction()->SetDescription("");
1755
1756   TopoDS_Shape aChamferShape = aChamfer->GetValue();
1757   aFunction->SetValue(aChamferShape);
1758   // END of chamfer
1759
1760   if (theHexMesh) {
1761     if (!MakePipeTShapePartition(aShape, theR1, theW1, theL1, theR2, theW2, theL2, theH, theW, 0, false))
1762       return NULL;
1763     if (!MakePipeTShapeMirrorAndGlue(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
1764       return NULL;
1765   }
1766
1767   TopoDS_Shape Te = aShape->GetValue();
1768
1769   // Set Position
1770   gp_Trsf aTrsf = GetPositionTrsf(theL1, theL2, theP1, theP2, theP3);
1771   BRepBuilderAPI_Transform aTransformation(Te, aTrsf, Standard_False);
1772   TopoDS_Shape aTrsf_Shape = aTransformation.Shape();
1773   aFunction->SetValue(aTrsf_Shape);
1774   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
1775   aSeq->Append(aShape);
1776   if (theHexMesh) {
1777     /*
1778      * Get the groups: BEGIN
1779      */
1780     try {
1781       if (!MakeGroups(aShape, TSHAPE_CHAMFER, theR1, theW1, theL1, theR2, theW2, theL2, aSeq, aTrsf))
1782         return NULL;
1783     }
1784     catch (Standard_Failure) {
1785       Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1786       SetErrorCode(aFail->GetMessageString());
1787       return NULL;
1788     }
1789
1790     TCollection_AsciiString aListRes, anEntry;
1791     // Iterate over the sequence aSeq
1792     Standard_Integer aNbGroups = aSeq->Length();
1793     Standard_Integer i = 2;
1794     for (; i <= aNbGroups; i++) {
1795       Handle(Standard_Transient) anItem = aSeq->Value(i);
1796       if (anItem.IsNull()) continue;
1797       Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(anItem);
1798       if (aGroup.IsNull()) continue;
1799       //Make a Python command
1800       TDF_Tool::Entry(aGroup->GetEntry(), anEntry);
1801       aListRes += anEntry + ", ";
1802     }
1803
1804     aListRes.Trunc(aListRes.Length() - 2);
1805
1806     //Make a Python command
1807     GEOM::TPythonDump(aFunction)
1808       << "[" << aShape << ", " << aListRes.ToCString()
1809       << "] = geompy.MakePipeTShapeChamfer(" << theR1 << ", " << theW1 << ", " << theL1 << ", " << theR2
1810       << ", " << theW2 << ", " << theL2 << ", " << theH << ", " << theW << ", " << theHexMesh << ", "
1811       << theP1 << ", " << theP2 << ", " << theP3 << ")";
1812   }
1813   /*
1814    * Get the groups: END
1815    */
1816   else {
1817     //Make a Python command
1818     GEOM::TPythonDump(aFunction)
1819       << "[" << aShape << "] = geompy.MakePipeTShapeChamfer(" << theR1 << ", " << theW1
1820       << ", " << theL1 << ", " << theR2 << ", " << theW2 << ", " << theL2 << ", " << theH << ", " << theW
1821       << ", " << theHexMesh << ", " << theP1 << ", " << theP2 << ", " << theP3 << ")";
1822   }
1823
1824   SetErrorCode(OK);
1825
1826   return aSeq;
1827 }
1828
1829 //=============================================================================
1830 /*!
1831  *  MakePipeTShapeFillet
1832  *  Create a T-shape object with specified caracteristics for the main and
1833  *  the incident pipes (radius, width, half-length). A fillet is created
1834  *  on the junction of the pipes.
1835  *  Center of the shape is (0,0,0). The main plane of the T-shape is XOY.
1836  *  \param theR1 Internal radius of main pipe
1837  *  \param theW1 Width of main pipe
1838  *  \param theL1 Half-length of main pipe
1839  *  \param theR2 Internal radius of incident pipe (R2 < R1)
1840  *  \param theW2 Width of incident pipe (R2+W2 < R1+W1)
1841  *  \param theL2 Half-length of incident pipe
1842  *  \param theRF Radius of curvature of fillet.
1843  *  \param theHexMesh Boolean indicating if shape is prepared for hex mesh
1844  *  \return List of GEOM_Objects, containing the created shape and propagation groups.
1845  */
1846 //=============================================================================
1847 Handle(TColStd_HSequenceOfTransient)
1848 GEOMImpl_IAdvancedOperations::MakePipeTShapeFillet(double theR1, double theW1, double theL1,
1849                                                    double theR2, double theW2, double theL2,
1850                                                    double theRF, bool theHexMesh)
1851 {
1852   SetErrorCode(KO);
1853   //Add a new object
1854   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_TSHAPE);
1855   //Add a new shape function with parameters
1856   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_PipeTShapeDriver::GetID(), TSHAPE_FILLET);
1857   if (aFunction.IsNull()) return NULL;
1858
1859   //Check if the function is set correctly
1860   if (aFunction->GetDriverGUID() != GEOMImpl_PipeTShapeDriver::GetID()) return NULL;
1861
1862   GEOMImpl_IPipeTShape aData(aFunction);
1863
1864   aData.SetR1(theR1);
1865   aData.SetW1(theW1);
1866   aData.SetL1(theL1);
1867   aData.SetR2(theR2);
1868   aData.SetW2(theW2);
1869   aData.SetL2(theL2);
1870   aData.SetRF(theRF);
1871   aData.SetHexMesh(theHexMesh);
1872
1873   //Compute the resulting value
1874   try {
1875 #if OCC_VERSION_LARGE > 0x06010000
1876     OCC_CATCH_SIGNALS;
1877 #endif
1878     if (!GetSolver()->ComputeFunction(aFunction)) {
1879       SetErrorCode("TShape driver failed");
1880       return NULL;
1881     }
1882   } catch (Standard_Failure) {
1883     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1884     SetErrorCode(aFail->GetMessageString());
1885     return NULL;
1886   }
1887
1888   // BEGIN of fillet
1889   TopoDS_Shape aShapeShape = aShape->GetValue();
1890   TopTools_IndexedMapOfShape anEdgesIndices;
1891   TopExp::MapShapes(aShapeShape, anEdgesIndices);
1892   // Common edges on external cylinders
1893   Handle(GEOM_Object) box_e;
1894   if (theHexMesh) {
1895     box_e = my3DPrimOperations->MakeBoxDXDYDZ(theR2+theW2, theR2+theW2, theR1+theW1);
1896   }
1897   else {
1898     box_e = my3DPrimOperations->MakeBoxDXDYDZ(2*(theR2+theW2), 2*(theR2+theW2), theR1+theW1);
1899   }
1900   box_e->GetLastFunction()->SetDescription("");
1901   box_e = myTransformOperations->TranslateDXDYDZ(box_e, -(theR2+theW2), -(theR2+theW2), 0);
1902   box_e->GetLastFunction()->SetDescription("");
1903
1904   Handle(TColStd_HSequenceOfInteger) edges_e =
1905     myShapesOperations->GetShapesOnBoxIDs(box_e, aShape, TopAbs_EDGE, GEOMAlgo_ST_IN);
1906   box_e->GetLastFunction()->SetDescription("");
1907
1908   if (edges_e.IsNull() || edges_e->Length() == 0) {
1909     SetErrorCode("External edges not found");
1910     return NULL;
1911   }
1912   int nbEdgesInFillet = 0;
1913   std::list<int> theEdges;
1914   for (int i=1; i<=edges_e->Length();i++) {
1915     int edgeID = edges_e->Value(i);
1916     TopoDS_Shape theEdge = anEdgesIndices.FindKey(edgeID);
1917     TopExp_Explorer Ex(theEdge,TopAbs_VERTEX);
1918     while (Ex.More()) {
1919       gp_Pnt aPt = BRep_Tool::Pnt(TopoDS::Vertex(Ex.Current()));
1920       if (Abs(aPt.Z() - (theR1+theW1)) <= Precision::Confusion()) {
1921         nbEdgesInFillet ++;
1922         theEdges.push_back(edgeID);
1923       }
1924       Ex.Next();
1925     }
1926     if (theHexMesh && nbEdgesInFillet == 1)
1927       break;
1928   }
1929
1930   Handle(GEOM_Object) aFillet;
1931   try {
1932     aFillet = myLocalOperations->MakeFilletEdges(aShape, theRF, theEdges);
1933   }
1934   catch (Standard_Failure) {
1935     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1936     SetErrorCode(aFail->GetMessageString());
1937     return NULL;
1938   }
1939   if (aFillet.IsNull()) {
1940     //SetErrorCode("Fillet can not be computed on the given shape with the given parameters");
1941     SetErrorCode(myLocalOperations->GetErrorCode());
1942     return NULL;
1943   }
1944   aFillet->GetLastFunction()->SetDescription("");
1945
1946   TopoDS_Shape aFilletShape = aFillet->GetValue();
1947   aFunction->SetValue(aFilletShape);
1948   // END of fillet
1949
1950 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - BEGIN (1)
1951 // the following block, when enabled, leads to partitioning problems
1952 #if 0
1953 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - END (1)
1954   // BEGIN: Limit tolerances (debug)
1955   Handle(GEOM_Object) aCorr1 = myHealingOperations->LimitTolerance(aShape, 1e-07);
1956   TopoDS_Shape aCorr1Shape = aCorr1->GetValue();
1957   aShape->GetLastFunction()->SetValue(aCorr1Shape);
1958   aCorr1->GetLastFunction()->SetDescription("");
1959   // END: Limit tolerances (debug)
1960 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - BEGIN (2)
1961 #endif
1962 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - END (2)
1963
1964   if (theHexMesh) {
1965     if (!MakePipeTShapePartition(aShape, theR1, theW1, theL1, theR2, theW2, theL2, 0, 0, theRF, false))
1966       return NULL;
1967     if (!MakePipeTShapeMirrorAndGlue(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
1968       return NULL;
1969   }
1970
1971   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
1972   aSeq->Append(aShape);
1973   if (theHexMesh) {
1974     /*
1975      * Get the groups: BEGIN
1976      */
1977     try {
1978       if (!MakeGroups(aShape, TSHAPE_FILLET, theR1, theW1, theL1, theR2, theW2, theL2, aSeq, gp_Trsf()))
1979         return NULL;
1980     }
1981     catch (Standard_Failure) {
1982       Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1983       SetErrorCode(aFail->GetMessageString());
1984       return NULL;
1985     }
1986
1987     TCollection_AsciiString aListRes, anEntry;
1988     // Iterate over the sequence aSeq
1989     Standard_Integer aNbGroups = aSeq->Length();
1990     Standard_Integer i = 2;
1991     for (; i <= aNbGroups; i++) {
1992       Handle(Standard_Transient) anItem = aSeq->Value(i);
1993       if (anItem.IsNull()) continue;
1994       Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(anItem);
1995       if (aGroup.IsNull()) continue;
1996       //Make a Python command
1997       TDF_Tool::Entry(aGroup->GetEntry(), anEntry);
1998       aListRes += anEntry + ", ";
1999     }
2000
2001     aListRes.Trunc(aListRes.Length() - 2);
2002
2003     //Make a Python command
2004     GEOM::TPythonDump(aFunction)
2005       << "[" << aShape << ", " << aListRes.ToCString()
2006       << "] = geompy.MakePipeTShapeFillet(" << theR1 << ", " << theW1 << ", " << theL1 << ", " << theR2
2007       << ", " << theW2 << ", " << theL2 << ", " << theRF << ", " << theHexMesh << ")";
2008   }
2009   /*
2010    * Get the groups: END
2011    */
2012   else {
2013     //Make a Python command
2014     GEOM::TPythonDump(aFunction)
2015       << "[" << aShape << "] = geompy.MakePipeTShapeFillet(" << theR1 << ", " << theW1
2016       << ", " << theL1 << ", " << theR2 << ", " << theW2 << ", " << theL2 << ", " << theRF << ", "
2017       << theHexMesh << ")";
2018   }
2019
2020   SetErrorCode(OK);
2021
2022   return aSeq;
2023 }
2024
2025 //=============================================================================
2026 /*!
2027  *  MakePipeTShapeFilletWithPosition
2028  *  \brief Create a T-shape object with specified caracteristics for the main and
2029  *  the incident pipes (radius, width, half-length). A fillet is created
2030  *  on the junction of the pipes.
2031  *  The extremities of the main pipe are located on junctions points P1 and P2.
2032  *  The extremity of the incident pipe is located on junction point P3.
2033  *  \param theR1 Internal radius of main pipe
2034  *  \param theW1 Width of main pipe
2035  *  \param theL1 Half-length of main pipe
2036  *  \param theR2 Internal radius of incident pipe (R2 < R1)
2037  *  \param theW2 Width of incident pipe (R2+W2 < R1+W1)
2038  *  \param theL2 Half-length of incident pipe
2039  *  \param theRF Radius of curvature of fillet
2040  *  \param theHexMesh Boolean indicating if shape is prepared for hex mesh
2041  *  \param theP1 1st junction point of main pipe
2042  *  \param theP2 2nd junction point of main pipe
2043  *  \param theP3 Junction point of incident pipe
2044  *  \return List of GEOM_Objects, containing the created shape and propagation groups.
2045  */
2046 //=============================================================================
2047 Handle(TColStd_HSequenceOfTransient)
2048 GEOMImpl_IAdvancedOperations::MakePipeTShapeFilletWithPosition(double theR1, double theW1, double theL1,
2049                                                                double theR2, double theW2, double theL2,
2050                                                                double theRF, bool theHexMesh,
2051                                                                Handle(GEOM_Object) theP1,
2052                                                                Handle(GEOM_Object) theP2,
2053                                                                Handle(GEOM_Object) theP3)
2054 {
2055   SetErrorCode(KO);
2056   //Add a new object
2057   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_TSHAPE);
2058   //Add a new shape function with parameters
2059   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_PipeTShapeDriver::GetID(), TSHAPE_FILLET);
2060   if (aFunction.IsNull()) return NULL;
2061
2062   //Check if the function is set correctly
2063   if (aFunction->GetDriverGUID() != GEOMImpl_PipeTShapeDriver::GetID()) return NULL;
2064
2065   // Check new position
2066   if (!CheckCompatiblePosition(theL1, theL2, theP1, theP2, theP3, 0.01)) {
2067     return NULL;
2068   }
2069
2070   GEOMImpl_IPipeTShape aData(aFunction);
2071
2072   aData.SetR1(theR1);
2073   aData.SetW1(theW1);
2074   aData.SetL1(theL1);
2075   aData.SetR2(theR2);
2076   aData.SetW2(theW2);
2077   aData.SetL2(theL2);
2078   aData.SetRF(theRF);
2079   aData.SetHexMesh(theHexMesh);
2080
2081   //Compute the resulting value
2082   try {
2083 #if OCC_VERSION_LARGE > 0x06010000
2084     OCC_CATCH_SIGNALS;
2085 #endif
2086     if (!GetSolver()->ComputeFunction(aFunction)) {
2087       SetErrorCode("TShape driver failed");
2088       return NULL;
2089     }
2090   } catch (Standard_Failure) {
2091     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2092     SetErrorCode(aFail->GetMessageString());
2093     return NULL;
2094   }
2095
2096   // BEGIN of fillet
2097   TopoDS_Shape aShapeShape = aShape->GetValue();
2098   TopTools_IndexedMapOfShape anEdgesIndices;
2099   TopExp::MapShapes(aShapeShape, anEdgesIndices);
2100   // Common edges on external cylinders
2101   Handle(GEOM_Object) box_e;
2102   if (theHexMesh) {
2103     box_e = my3DPrimOperations->MakeBoxDXDYDZ(theR2+theW2, theR2+theW2, theR1+theW1);
2104   }
2105   else {
2106     box_e = my3DPrimOperations->MakeBoxDXDYDZ(2*(theR2+theW2), 2*(theR2+theW2), theR1+theW1);
2107   }
2108   box_e->GetLastFunction()->SetDescription("");
2109   box_e = myTransformOperations->TranslateDXDYDZ(box_e, -(theR2+theW2), -(theR2+theW2), 0);
2110   box_e->GetLastFunction()->SetDescription("");
2111
2112   Handle(TColStd_HSequenceOfInteger) edges_e =
2113     myShapesOperations->GetShapesOnBoxIDs(box_e, aShape, TopAbs_EDGE, GEOMAlgo_ST_IN);
2114   box_e->GetLastFunction()->SetDescription("");
2115
2116   if (edges_e.IsNull() || edges_e->Length() == 0) {
2117     SetErrorCode("External edges not found");
2118     return NULL;
2119   }
2120   int nbEdgesInFillet = 0;
2121   std::list<int> theEdges;
2122   for (int i=1; i<=edges_e->Length();i++) {
2123     int edgeID = edges_e->Value(i);
2124     TopoDS_Shape theEdge = anEdgesIndices.FindKey(edgeID);
2125     TopExp_Explorer Ex(theEdge,TopAbs_VERTEX);
2126     while (Ex.More()) {
2127       gp_Pnt aPt = BRep_Tool::Pnt(TopoDS::Vertex(Ex.Current()));
2128       if (Abs(aPt.Z() - (theR1+theW1)) <= Precision::Confusion()) {
2129         nbEdgesInFillet ++;
2130         theEdges.push_back(edgeID);
2131       }
2132       Ex.Next();
2133     }
2134     if (theHexMesh && nbEdgesInFillet == 1)
2135       break;
2136   }
2137
2138   Handle(GEOM_Object) aFillet;
2139   try {
2140     aFillet = myLocalOperations->MakeFilletEdges(aShape, theRF, theEdges);
2141   }
2142   catch (Standard_Failure) {
2143     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2144     SetErrorCode(aFail->GetMessageString());
2145     return NULL;
2146   }
2147   if (aFillet.IsNull()) {
2148     SetErrorCode("Fillet can not be computed on the given shape with the given parameters");
2149     return NULL;
2150   }
2151   aFillet->GetLastFunction()->SetDescription("");
2152
2153   TopoDS_Shape aFilletShape = aFillet->GetValue();
2154   aFunction->SetValue(aFilletShape);
2155   // END of fillet
2156
2157 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - BEGIN (3)
2158 // the following block, when enabled, leads to partitioning problems
2159 #if 0
2160 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - END (3)
2161   // BEGIN: Limit tolerances (debug)
2162   Handle(GEOM_Object) aCorr1 = myHealingOperations->LimitTolerance(aShape, 1e-07);
2163   TopoDS_Shape aCorr1Shape = aCorr1->GetValue();
2164   aShape->GetLastFunction()->SetValue(aCorr1Shape);
2165   aCorr1->GetLastFunction()->SetDescription("");
2166   // END: Limit tolerances (debug)
2167 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - BEGIN (4)
2168 #endif
2169 // VSR: debug issues 0021568 and 0021550 (15/05/2012) - END (4)
2170
2171   if (theHexMesh) {
2172     if (!MakePipeTShapePartition(aShape, theR1, theW1, theL1, theR2, theW2, theL2, 0, 0, theRF, false))
2173       return NULL;
2174     if (!MakePipeTShapeMirrorAndGlue(aShape, theR1, theW1, theL1, theR2, theW2, theL2))
2175       return NULL;
2176   }
2177
2178   TopoDS_Shape Te = aShape->GetValue();
2179
2180   // Set Position
2181   gp_Trsf aTrsf = GetPositionTrsf(theL1, theL2, theP1, theP2, theP3);
2182   BRepBuilderAPI_Transform aTransformation(Te, aTrsf, Standard_False);
2183   TopoDS_Shape aTrsf_Shape = aTransformation.Shape();
2184   aFunction->SetValue(aTrsf_Shape);
2185   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
2186   aSeq->Append(aShape);
2187   if (theHexMesh) {
2188     /*
2189      * Get the groups: BEGIN
2190      */
2191     try {
2192       if (!MakeGroups(aShape, TSHAPE_FILLET, theR1, theW1, theL1, theR2, theW2, theL2, aSeq, aTrsf))
2193         return NULL;
2194     }
2195     catch (Standard_Failure) {
2196       Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2197       SetErrorCode(aFail->GetMessageString());
2198       return NULL;
2199     }
2200
2201     TCollection_AsciiString aListRes, anEntry;
2202     // Iterate over the sequence aSeq
2203     Standard_Integer aNbGroups = aSeq->Length();
2204     Standard_Integer i = 2;
2205     for (; i <= aNbGroups; i++) {
2206       Handle(Standard_Transient) anItem = aSeq->Value(i);
2207       if (anItem.IsNull()) continue;
2208       Handle(GEOM_Object) aGroup = Handle(GEOM_Object)::DownCast(anItem);
2209       if (aGroup.IsNull()) continue;
2210       //Make a Python command
2211       TDF_Tool::Entry(aGroup->GetEntry(), anEntry);
2212       aListRes += anEntry + ", ";
2213     }
2214
2215     aListRes.Trunc(aListRes.Length() - 2);
2216
2217     //Make a Python command
2218     GEOM::TPythonDump(aFunction)
2219       << "[" << aShape << ", " << aListRes.ToCString()
2220       << "] = geompy.MakePipeTShapeFillet(" << theR1 << ", " << theW1 << ", " << theL1 << ", " << theR2
2221       << ", " << theW2 << ", " << theL2 << ", " << theRF << ", " << theHexMesh << ", " << theP1 << ", "
2222       << theP2 << ", " << theP3 << ")";
2223   }
2224   /*
2225    * Get the groups: END
2226    */
2227   else {
2228     //Make a Python command
2229     GEOM::TPythonDump(aFunction)
2230       << "[" << aShape << "] = geompy.MakePipeTShapeFillet(" << theR1 << ", " << theW1
2231       << ", " << theL1 << ", " << theR2 << ", " << theW2 << ", " << theL2 << ", " << theRF << ", "
2232       << theHexMesh << ", " << theP1 << ", " << theP2 << ", " << theP3 << ")";
2233   }
2234
2235   SetErrorCode(OK);
2236
2237   return aSeq;
2238 }
2239
2240 //=============================================================================
2241 /*!
2242  *  This function allows to create a disk already divided into blocks. It can be
2243  *  used to create divided pipes for later meshing in hexaedra.
2244  *  \param theR Radius of the disk
2245  *  \param theRatio Relative size of the central square diagonal against the disk diameter
2246  *  \param theOrientation Plane on which the disk will be built
2247  *  \param thePattern The division pattern of the disk (hexagon or square in the center)
2248  *  \return New GEOM_Object, containing the created shape.
2249  */
2250 //=============================================================================
2251 Handle(GEOM_Object) GEOMImpl_IAdvancedOperations::MakeDividedDisk (double theR, double theRatio,
2252                                                                    int theOrientation, int thePattern)
2253 {
2254   SetErrorCode(KO);
2255
2256   //Add a new object
2257   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_DIVIDEDDISK);
2258
2259   //Add a new shape function with parameters
2260   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_DividedDiskDriver::GetID(), DIVIDEDDISK_R_RATIO);
2261   if (aFunction.IsNull()) return NULL;
2262
2263   //Check if the function is set correctly
2264   if (aFunction->GetDriverGUID() != GEOMImpl_DividedDiskDriver::GetID()) return NULL;
2265
2266   GEOMImpl_IDividedDisk aData (aFunction);
2267
2268   aData.SetR(theR);
2269   aData.SetRatio(theRatio);
2270   aData.SetOrientation(theOrientation);
2271   aData.SetType(thePattern);
2272
2273   //Compute the resulting value
2274   try {
2275 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
2276     OCC_CATCH_SIGNALS;
2277 #endif
2278     if (!GetSolver()->ComputeFunction(aFunction)) {
2279       SetErrorCode("DividedDisk driver failed");
2280       return NULL;
2281     }
2282   }
2283   catch (Standard_Failure) {
2284     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2285     SetErrorCode(aFail->GetMessageString());
2286     return NULL;
2287   }
2288
2289   std::string aPatternStr;
2290
2291   switch(thePattern)
2292   {
2293     case 0:
2294       aPatternStr = "GEOM.SQUARE";
2295       break;
2296     case 1:
2297       aPatternStr = "GEOM.HEXAGON";
2298       break;
2299   }
2300
2301   //Make a Python command
2302   GEOM::TPythonDump(aFunction) << aShape << " = geompy.MakeDividedDisk(" << theR << ", " << theOrientation << ", " << aPatternStr.c_str() << ")";
2303
2304   SetErrorCode(OK);
2305
2306   return aShape;
2307 }
2308
2309 //=============================================================================
2310 /*!
2311  *  This function allows to create a disk already divided into blocks. It can be
2312  *  used to create divided pipes for later meshing in hexaedra.
2313  *  \param theR Radius of the disk
2314  *  \param theRatio Relative size of the central square diagonal against the disk diameter
2315  *  \return New GEOM_Object, containing the created shape.
2316  */
2317 //=============================================================================
2318 Handle(GEOM_Object) GEOMImpl_IAdvancedOperations::MakeDividedDiskPntVecR (Handle(GEOM_Object) thePnt,
2319                                                                           Handle(GEOM_Object) theVec,
2320                                                                           double theR,
2321                                                                           double theRatio,
2322                                                                           int    thePattern)
2323 {
2324   SetErrorCode(KO);
2325
2326   //Add a new object
2327   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_DIVIDEDDISK);
2328
2329   //Add a new shape function with parameters
2330   Handle(GEOM_Function) aFunction = aShape->AddFunction(GEOMImpl_DividedDiskDriver::GetID(), DIVIDEDDISK_R_VECTOR_PNT);
2331   if (aFunction.IsNull()) return NULL;
2332
2333   //Check if the function is set correctly
2334   if (aFunction->GetDriverGUID() != GEOMImpl_DividedDiskDriver::GetID()) return NULL;
2335
2336   GEOMImpl_IDividedDisk aData (aFunction);
2337
2338   Handle(GEOM_Function) aRefPnt = thePnt->GetLastFunction();
2339   Handle(GEOM_Function) aRefVec = theVec->GetLastFunction();
2340
2341   if (aRefPnt.IsNull() || aRefVec.IsNull()) return NULL;
2342
2343   aData.SetCenter(aRefPnt);
2344   aData.SetVector(aRefVec);
2345
2346   aData.SetR(theR);
2347   aData.SetRatio(theRatio);
2348   aData.SetType(thePattern);
2349
2350   //Compute the resulting value
2351   try {
2352 #if (OCC_VERSION_MAJOR << 16 | OCC_VERSION_MINOR << 8 | OCC_VERSION_MAINTENANCE) > 0x060100
2353     OCC_CATCH_SIGNALS;
2354 #endif
2355     if (!GetSolver()->ComputeFunction(aFunction)) {
2356       SetErrorCode("DividedDisk driver failed");
2357       return NULL;
2358     }
2359   }
2360   catch (Standard_Failure) {
2361     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2362     SetErrorCode(aFail->GetMessageString());
2363     return NULL;
2364   }
2365
2366   std::string aPatternStr;
2367
2368   switch(thePattern)
2369   {
2370     case 0:
2371       aPatternStr = "GEOM.SQUARE";
2372       break;
2373     case 1:
2374       aPatternStr = "GEOM.HEXAGON";
2375       break;
2376   }
2377
2378
2379   //Make a Python command
2380   GEOM::TPythonDump(aFunction) << aShape << " = geompy.MakeDividedDiskPntVecR(" << thePnt << ", " << theVec << ", " << theR << ", " << aPatternStr.c_str() << ")";
2381
2382   SetErrorCode(OK);
2383
2384   return aShape;
2385 }
2386
2387 //=============================================================================
2388 /*!
2389  *  Builds a cylinder prepared for hexa meshes
2390  *  \param theR Radius of the cylinder
2391  *  \param theH Height of the cylinder
2392  *  \return New GEOM_Object, containing the created shape.
2393  */
2394 //=============================================================================
2395 Handle(GEOM_Object) GEOMImpl_IAdvancedOperations::MakeDividedCylinder (double theR,
2396                                                                        double theH,
2397                                                                        int    thePattern)
2398 {
2399   SetErrorCode(KO);
2400
2401   //Add a new object
2402   Handle(GEOM_Object) aShape = GetEngine()->AddObject(GetDocID(), GEOM_DIVIDEDCYLINDER);
2403
2404   Handle(GEOM_Object) aBaseShape = MakeDividedDisk(theR, 67.0, 1, thePattern);
2405   aBaseShape->GetLastFunction()->SetDescription("");   // Erase dump of MakeDividedDisk
2406
2407   aShape = my3DPrimOperations->MakePrismDXDYDZ(aBaseShape,0.0,0.0,theH, -1.0);
2408
2409   Handle(GEOM_Function) aFunction =  aShape->GetLastFunction();
2410   aFunction->SetDescription("");   // Erase dump of MakePrismDXDYDZ
2411   aShape->SetType(GEOM_DIVIDEDCYLINDER);
2412
2413   std::string aPatternStr;
2414
2415   switch(thePattern)
2416   {
2417     case 0:
2418       aPatternStr = "GEOM.SQUARE";
2419       break;
2420     case 1:
2421       aPatternStr = "GEOM.HEXAGON";
2422       break;
2423   }
2424
2425   //Make a Python command
2426   GEOM::TPythonDump(aFunction) << aShape << " = geompy.MakeDividedCylinder(" << theR << ", " << theH << ", " << aPatternStr.c_str() << ")";
2427
2428   SetErrorCode(OK);
2429
2430   return aShape;
2431 }
2432 /*@@ insert new functions before this line @@ do not remove this line @@ do not remove this line @@*/