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