src/GEOMImpl/GEOMImpl_IBlocksOperations.cxx

   1 using namespace std;
   2
   3 #include "GEOMImpl_IBlocksOperations.hxx"
   4
   5 #include "GEOMImpl_Types.hxx"
   6
   7 #include "GEOMImpl_BlockDriver.hxx"
   8 #include "GEOMImpl_IBlocks.hxx"
   9 #include "GEOMImpl_IBlockTrsf.hxx"
  10 #include "GEOMImpl_CopyDriver.hxx"
  11 #include "GEOMImpl_Block6Explorer.hxx"
  12
  13 #include "GEOM_Function.hxx"
  14 #include "GEOM_PythonDump.hxx"
  15
  16 #include "GEOMAlgo_GlueAnalyser.hxx"
  17 #include "GEOMAlgo_CoupleOfShapes.hxx"
  18 #include "GEOMAlgo_ListOfCoupleOfShapes.hxx"
  19 #include "GEOMAlgo_ListIteratorOfListOfCoupleOfShapes.hxx"
  20 #include "BlockFix_CheckTool.hxx"
  21
  22 #include "utilities.h"
  23 #include "OpUtil.hxx"
  24 #include "Utils_ExceptHandlers.hxx"
  25
  26 #include <TFunction_DriverTable.hxx>
  27 #include <TFunction_Driver.hxx>
  28 #include <TFunction_Logbook.hxx>
  29 #include <TDataStd_Integer.hxx>
  30 #include <TDF_Tool.hxx>
  31
  32 #include <BRep_Tool.hxx>
  33 #include <BRep_Builder.hxx>
  34 #include <BRepTools.hxx>
  35 #include <BRepTools_WireExplorer.hxx>
  36 #include <BRepGProp.hxx>
  37 #include <BRepBndLib.hxx>
  38 #include <BRepAdaptor_Surface.hxx>
  39 #include <BRepClass_FaceClassifier.hxx>
  40 #include <BRepClass3d_SolidClassifier.hxx>
  41 #include <BRepExtrema_DistShapeShape.hxx>
  42
  43 #include <TopAbs.hxx>
  44 #include <TopoDS.hxx>
  45 #include <TopoDS_Edge.hxx>
  46 #include <TopoDS_Vertex.hxx>
  47 #include <TopoDS_Compound.hxx>
  48 #include <TopoDS_Iterator.hxx>
  49 #include <TopExp.hxx>
  50 #include <TopExp_Explorer.hxx>
  51 #include <TopTools_MapOfShape.hxx>
  52 #include <TopTools_Array1OfShape.hxx>
  53 #include <TopTools_IndexedMapOfShape.hxx>
  54 #include <TopTools_DataMapOfShapeInteger.hxx>
  55 #include <TopTools_ListIteratorOfListOfShape.hxx>
  56 #include <TopTools_IndexedDataMapOfShapeListOfShape.hxx>
  57 #include <TopTools_DataMapIteratorOfDataMapOfShapeInteger.hxx>
  58
  59 #include <Bnd_Box.hxx>
  60 #include <GProp_GProps.hxx>
  61
  62 #include <Geom_Surface.hxx>
  63 #include <ShapeAnalysis_Surface.hxx>
  64
  65 #include <TColStd_MapOfInteger.hxx>
  66 #include <TColStd_Array1OfReal.hxx>
  67 #include <TColStd_Array1OfInteger.hxx>
  68 #include <TColStd_Array2OfInteger.hxx>
  69
  70 //#include <OSD_Timer.hxx>
  71
  72 #include <Precision.hxx>
  73
  74 #include <Standard_ErrorHandler.hxx> // CAREFUL ! position of this file is critic : see Lucien PIGNOLONI / OCC
  75
  76 //=============================================================================
  77 /*!
  78  *   constructor:
  79  */
  80 //=============================================================================
  81 GEOMImpl_IBlocksOperations::GEOMImpl_IBlocksOperations (GEOM_Engine* theEngine, int theDocID)
  82 : GEOM_IOperations(theEngine, theDocID)
  83 {
  84   MESSAGE("GEOMImpl_IBlocksOperations::GEOMImpl_IBlocksOperations");
  85 }
  86
  87 //=============================================================================
  88 /*!
  89  *  destructor
  90  */
  91 //=============================================================================
  92 GEOMImpl_IBlocksOperations::~GEOMImpl_IBlocksOperations()
  93 {
  94   MESSAGE("GEOMImpl_IBlocksOperations::~GEOMImpl_IBlocksOperations");
  95 }
  96
  97
  98 //=============================================================================
  99 /*!
 100  *  MakeQuad
 101  */
 102 //=============================================================================
 103 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::MakeQuad
 104                      (Handle(GEOM_Object) theEdge1, Handle(GEOM_Object) theEdge2,
 105                       Handle(GEOM_Object) theEdge3, Handle(GEOM_Object) theEdge4)
 106 {
 107   SetErrorCode(KO);
 108
 109   if (theEdge1.IsNull() || theEdge2.IsNull() ||
 110       theEdge3.IsNull() || theEdge4.IsNull()) return NULL;
 111
 112   //Add a new Face object
 113   Handle(GEOM_Object) aFace = GetEngine()->AddObject(GetDocID(), GEOM_FACE);
 114
 115   //Add a new Face function
 116   Handle(GEOM_Function) aFunction =
 117     aFace->AddFunction(GEOMImpl_BlockDriver::GetID(), BLOCK_FACE_FOUR_EDGES);
 118
 119   //Check if the function is set correctly
 120   if (aFunction->GetDriverGUID() != GEOMImpl_BlockDriver::GetID()) return NULL;
 121
 122   GEOMImpl_IBlocks aPI (aFunction);
 123
 124   Handle(GEOM_Function) aRef1 = theEdge1->GetLastFunction();
 125   Handle(GEOM_Function) aRef2 = theEdge2->GetLastFunction();
 126   Handle(GEOM_Function) aRef3 = theEdge3->GetLastFunction();
 127   Handle(GEOM_Function) aRef4 = theEdge4->GetLastFunction();
 128   if (aRef1.IsNull() || aRef2.IsNull() ||
 129       aRef3.IsNull() || aRef4.IsNull()) return NULL;
 130
 131   Handle(TColStd_HSequenceOfTransient) aShapesSeq = new TColStd_HSequenceOfTransient;
 132   aShapesSeq->Append(aRef1);
 133   aShapesSeq->Append(aRef2);
 134   aShapesSeq->Append(aRef3);
 135   aShapesSeq->Append(aRef4);
 136
 137   aPI.SetShapes(aShapesSeq);
 138
 139   //Compute the Face value
 140   try {
 141     if (!GetSolver()->ComputeFunction(aFunction)) {
 142       SetErrorCode("Block driver failed to compute a face");
 143       return NULL;
 144     }
 145   }
 146   catch (Standard_Failure) {
 147     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 148     SetErrorCode(aFail->GetMessageString());
 149     return NULL;
 150   }
 151
 152   //Make a Python command
 153   GEOM::TPythonDump(aFunction) << aFace << " = geompy.MakeQuad("
 154     << theEdge1 << ", " << theEdge2 << ", " << theEdge3 << ", " << theEdge4 << ")";
 155
 156   SetErrorCode(OK);
 157   return aFace;
 158 }
 159
 160 //=============================================================================
 161 /*!
 162  *  MakeQuad2Edges
 163  */
 164 //=============================================================================
 165 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::MakeQuad2Edges
 166                      (Handle(GEOM_Object) theEdge1, Handle(GEOM_Object) theEdge2)
 167 {
 168   SetErrorCode(KO);
 169
 170   if (theEdge1.IsNull() || theEdge2.IsNull()) return NULL;
 171
 172   //Add a new Face object
 173   Handle(GEOM_Object) aFace = GetEngine()->AddObject(GetDocID(), GEOM_FACE);
 174
 175   //Add a new Face function
 176   Handle(GEOM_Function) aFunction =
 177     aFace->AddFunction(GEOMImpl_BlockDriver::GetID(), BLOCK_FACE_TWO_EDGES);
 178
 179   //Check if the function is set correctly
 180   if (aFunction->GetDriverGUID() != GEOMImpl_BlockDriver::GetID()) return NULL;
 181
 182   GEOMImpl_IBlocks aPI (aFunction);
 183
 184   Handle(GEOM_Function) aRef1 = theEdge1->GetLastFunction();
 185   Handle(GEOM_Function) aRef2 = theEdge2->GetLastFunction();
 186   if (aRef1.IsNull() || aRef2.IsNull()) return NULL;
 187
 188   Handle(TColStd_HSequenceOfTransient) aShapesSeq = new TColStd_HSequenceOfTransient;
 189   aShapesSeq->Append(aRef1);
 190   aShapesSeq->Append(aRef2);
 191
 192   aPI.SetShapes(aShapesSeq);
 193
 194   //Compute the Face value
 195   try {
 196     if (!GetSolver()->ComputeFunction(aFunction)) {
 197       SetErrorCode("Block driver failed to compute a face");
 198       return NULL;
 199     }
 200   }
 201   catch (Standard_Failure) {
 202     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 203     SetErrorCode(aFail->GetMessageString());
 204     return NULL;
 205   }
 206
 207   //Make a Python command
 208   GEOM::TPythonDump(aFunction) << aFace << " = geompy.MakeQuad2Edges("
 209                                << theEdge1 << ", " << theEdge2 << ")";
 210
 211   SetErrorCode(OK);
 212   return aFace;
 213 }
 214
 215 //=============================================================================
 216 /*!
 217  *  MakeQuad4Vertices
 218  */
 219 //=============================================================================
 220 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::MakeQuad4Vertices
 221                      (Handle(GEOM_Object) thePnt1, Handle(GEOM_Object) thePnt2,
 222                       Handle(GEOM_Object) thePnt3, Handle(GEOM_Object) thePnt4)
 223 {
 224   SetErrorCode(KO);
 225
 226   if (thePnt1.IsNull() || thePnt2.IsNull() ||
 227       thePnt3.IsNull() || thePnt4.IsNull()) return NULL;
 228
 229   //Add a new Face object
 230   Handle(GEOM_Object) aFace = GetEngine()->AddObject(GetDocID(), GEOM_FACE);
 231
 232   //Add a new Face function
 233   Handle(GEOM_Function) aFunction =
 234     aFace->AddFunction(GEOMImpl_BlockDriver::GetID(), BLOCK_FACE_FOUR_PNT);
 235
 236   //Check if the function is set correctly
 237   if (aFunction->GetDriverGUID() != GEOMImpl_BlockDriver::GetID()) return NULL;
 238
 239   GEOMImpl_IBlocks aPI (aFunction);
 240
 241   Handle(GEOM_Function) aRef1 = thePnt1->GetLastFunction();
 242   Handle(GEOM_Function) aRef2 = thePnt2->GetLastFunction();
 243   Handle(GEOM_Function) aRef3 = thePnt3->GetLastFunction();
 244   Handle(GEOM_Function) aRef4 = thePnt4->GetLastFunction();
 245   if (aRef1.IsNull() || aRef2.IsNull() ||
 246       aRef3.IsNull() || aRef4.IsNull()) return NULL;
 247
 248   Handle(TColStd_HSequenceOfTransient) aShapesSeq = new TColStd_HSequenceOfTransient;
 249   aShapesSeq->Append(aRef1);
 250   aShapesSeq->Append(aRef2);
 251   aShapesSeq->Append(aRef3);
 252   aShapesSeq->Append(aRef4);
 253
 254   aPI.SetShapes(aShapesSeq);
 255
 256   //Compute the Face value
 257   try {
 258     if (!GetSolver()->ComputeFunction(aFunction)) {
 259       SetErrorCode("Block driver failed to compute a face");
 260       return NULL;
 261     }
 262   }
 263   catch (Standard_Failure) {
 264     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 265     SetErrorCode(aFail->GetMessageString());
 266     return NULL;
 267   }
 268
 269   //Make a Python command
 270   GEOM::TPythonDump(aFunction) << aFace << " = geompy.MakeQuad4Vertices("
 271     << thePnt1 << ", " << thePnt2 << ", " << thePnt3 << ", " << thePnt4 << ")";
 272
 273   SetErrorCode(OK);
 274   return aFace;
 275 }
 276
 277 //=============================================================================
 278 /*!
 279  *  MakeHexa
 280  */
 281 //=============================================================================
 282 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::MakeHexa
 283                      (Handle(GEOM_Object) theFace1, Handle(GEOM_Object) theFace2,
 284                       Handle(GEOM_Object) theFace3, Handle(GEOM_Object) theFace4,
 285                       Handle(GEOM_Object) theFace5, Handle(GEOM_Object) theFace6)
 286 {
 287   SetErrorCode(KO);
 288
 289   if (theFace1.IsNull() || theFace2.IsNull() ||
 290       theFace3.IsNull() || theFace4.IsNull() ||
 291       theFace5.IsNull() || theFace6.IsNull()) return NULL;
 292
 293   //Add a new Solid object
 294   Handle(GEOM_Object) aBlock = GetEngine()->AddObject(GetDocID(), GEOM_BLOCK);
 295
 296   //Add a new Block function
 297   Handle(GEOM_Function) aFunction =
 298     aBlock->AddFunction(GEOMImpl_BlockDriver::GetID(), BLOCK_SIX_FACES);
 299
 300   //Check if the function is set correctly
 301   if (aFunction->GetDriverGUID() != GEOMImpl_BlockDriver::GetID()) return NULL;
 302
 303   GEOMImpl_IBlocks aPI (aFunction);
 304
 305   Handle(GEOM_Function) aRef1 = theFace1->GetLastFunction();
 306   Handle(GEOM_Function) aRef2 = theFace2->GetLastFunction();
 307   Handle(GEOM_Function) aRef3 = theFace3->GetLastFunction();
 308   Handle(GEOM_Function) aRef4 = theFace4->GetLastFunction();
 309   Handle(GEOM_Function) aRef5 = theFace5->GetLastFunction();
 310   Handle(GEOM_Function) aRef6 = theFace6->GetLastFunction();
 311   if (aRef1.IsNull() || aRef2.IsNull() ||
 312       aRef3.IsNull() || aRef4.IsNull() ||
 313       aRef5.IsNull() || aRef6.IsNull()) return NULL;
 314
 315   Handle(TColStd_HSequenceOfTransient) aShapesSeq = new TColStd_HSequenceOfTransient;
 316   aShapesSeq->Append(aRef1);
 317   aShapesSeq->Append(aRef2);
 318   aShapesSeq->Append(aRef3);
 319   aShapesSeq->Append(aRef4);
 320   aShapesSeq->Append(aRef5);
 321   aShapesSeq->Append(aRef6);
 322
 323   aPI.SetShapes(aShapesSeq);
 324
 325   //Compute the Block value
 326   try {
 327     if (!GetSolver()->ComputeFunction(aFunction)) {
 328       SetErrorCode("Block driver failed to compute a block");
 329       return NULL;
 330     }
 331   }
 332   catch (Standard_Failure) {
 333     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 334     SetErrorCode(aFail->GetMessageString());
 335     return NULL;
 336   }
 337
 338   //Make a Python command
 339   GEOM::TPythonDump(aFunction) << aBlock << " = geompy.MakeHexa("
 340     << theFace1 << ", " << theFace2 << ", " << theFace3 << ", "
 341       << theFace4 << ", " << theFace5 << ", " << theFace6 << ")";
 342
 343   SetErrorCode(OK);
 344   return aBlock;
 345 }
 346
 347 //=============================================================================
 348 /*!
 349  *  MakeHexa2Faces
 350  */
 351 //=============================================================================
 352 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::MakeHexa2Faces
 353                    (Handle(GEOM_Object) theFace1, Handle(GEOM_Object) theFace2)
 354 {
 355   SetErrorCode(KO);
 356
 357   if (theFace1.IsNull() || theFace2.IsNull()) return NULL;
 358
 359   //Add a new Solid object
 360   Handle(GEOM_Object) aBlock = GetEngine()->AddObject(GetDocID(), GEOM_BLOCK);
 361
 362   //Add a new Block function
 363   Handle(GEOM_Function) aFunction =
 364     aBlock->AddFunction(GEOMImpl_BlockDriver::GetID(), BLOCK_TWO_FACES);
 365
 366   //Check if the function is set correctly
 367   if (aFunction->GetDriverGUID() != GEOMImpl_BlockDriver::GetID()) return NULL;
 368
 369   GEOMImpl_IBlocks aPI (aFunction);
 370
 371   Handle(GEOM_Function) aRef1 = theFace1->GetLastFunction();
 372   Handle(GEOM_Function) aRef2 = theFace2->GetLastFunction();
 373   if (aRef1.IsNull() || aRef2.IsNull()) return NULL;
 374
 375   Handle(TColStd_HSequenceOfTransient) aShapesSeq = new TColStd_HSequenceOfTransient;
 376   aShapesSeq->Append(aRef1);
 377   aShapesSeq->Append(aRef2);
 378
 379   aPI.SetShapes(aShapesSeq);
 380
 381   //Compute the Block value
 382   try {
 383     if (!GetSolver()->ComputeFunction(aFunction)) {
 384       SetErrorCode("Block driver failed to compute a block");
 385       return NULL;
 386     }
 387   }
 388   catch (Standard_Failure) {
 389     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 390     SetErrorCode(aFail->GetMessageString());
 391     return NULL;
 392   }
 393
 394   //Make a Python command
 395   GEOM::TPythonDump(aFunction) << aBlock << " = geompy.MakeHexa2Faces("
 396                                << theFace1 << ", " << theFace2 << ")";
 397
 398   SetErrorCode(OK);
 399   return aBlock;
 400 }
 401
 402 //=============================================================================
 403 /*!
 404  *  MakeBlockCompound
 405  */
 406 //=============================================================================
 407 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::MakeBlockCompound
 408                                               (Handle(GEOM_Object) theCompound)
 409 {
 410   SetErrorCode(KO);
 411
 412   if (theCompound.IsNull()) return NULL;
 413
 414   //Add a new object
 415   Handle(GEOM_Object) aBlockComp = GetEngine()->AddObject(GetDocID(), GEOM_COMPOUND);
 416
 417   //Add a new BlocksComp function
 418   Handle(GEOM_Function) aFunction =
 419     aBlockComp->AddFunction(GEOMImpl_BlockDriver::GetID(), BLOCK_COMPOUND_GLUE);
 420
 421   //Check if the function is set correctly
 422   if (aFunction->GetDriverGUID() != GEOMImpl_BlockDriver::GetID()) return NULL;
 423
 424   GEOMImpl_IBlocks aPI (aFunction);
 425
 426   Handle(GEOM_Function) aRef1 = theCompound->GetLastFunction();
 427   if (aRef1.IsNull()) return NULL;
 428
 429   Handle(TColStd_HSequenceOfTransient) aShapesSeq = new TColStd_HSequenceOfTransient;
 430   aShapesSeq->Append(aRef1);
 431
 432   aPI.SetShapes(aShapesSeq);
 433
 434   //Compute the Blocks Compound value
 435   try {
 436     if (!GetSolver()->ComputeFunction(aFunction)) {
 437       SetErrorCode("Block driver failed to compute a blocks compound");
 438       return NULL;
 439     }
 440   }
 441   catch (Standard_Failure) {
 442     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 443     SetErrorCode(aFail->GetMessageString());
 444     return NULL;
 445   }
 446
 447   //Make a Python command
 448   GEOM::TPythonDump(aFunction) << aBlockComp
 449     << " = geompy.MakeBlockCompound(" << theCompound << ")";
 450
 451   SetErrorCode(OK);
 452   return aBlockComp;
 453 }
 454
 455 //=============================================================================
 456 /*!
 457  *  GetEdge
 458  */
 459 //=============================================================================
 460 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::GetPoint
 461                                                (Handle(GEOM_Object) theShape,
 462                                                 const Standard_Real theX,
 463                                                 const Standard_Real theY,
 464                                                 const Standard_Real theZ,
 465                                                 const Standard_Real theEpsilon)
 466 {
 467   SetErrorCode(KO);
 468
 469   //New Point object
 470   Handle(GEOM_Object) aResult;
 471
 472   // Arguments
 473   if (theShape.IsNull()) return NULL;
 474
 475   TopoDS_Shape aBlockOrComp = theShape->GetValue();
 476   if (aBlockOrComp.IsNull()) {
 477     SetErrorCode("Block or compound is null");
 478     return NULL;
 479   }
 480   if (aBlockOrComp.ShapeType() != TopAbs_SOLID &&
 481       aBlockOrComp.ShapeType() != TopAbs_COMPOUND &&
 482       aBlockOrComp.ShapeType() != TopAbs_COMPSOLID) {
 483     SetErrorCode("Shape is neither a block, nor a compound of blocks");
 484     return NULL;
 485   }
 486
 487   //Compute the Vertex value
 488   gp_Pnt P (theX, theY, theZ);
 489   Standard_Real eps = Max(theEpsilon, Precision::Confusion());
 490
 491   TopoDS_Shape V;
 492   Standard_Integer isFound = 0;
 493   TopTools_MapOfShape mapShape;
 494   TopExp_Explorer exp (aBlockOrComp, TopAbs_VERTEX);
 495
 496   for (; exp.More(); exp.Next()) {
 497     if (mapShape.Add(exp.Current())) {
 498       TopoDS_Vertex aVi = TopoDS::Vertex(exp.Current());
 499       gp_Pnt aPi = BRep_Tool::Pnt(aVi);
 500       if (aPi.Distance(P) < eps) {
 501         V = aVi;
 502         isFound++;
 503       }
 504     }
 505   }
 506
 507   if (isFound == 0) {
 508     SetErrorCode("Vertex has not been found");
 509     return NULL;
 510   } else if (isFound > 1) {
 511     SetErrorCode("Multiple vertices found by the given coordinates and epsilon");
 512     return NULL;
 513   } else {
 514     TopTools_IndexedMapOfShape anIndices;
 515     TopExp::MapShapes(aBlockOrComp, anIndices);
 516     Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
 517     anArray->SetValue(1, anIndices.FindIndex(V));
 518     aResult = GetEngine()->AddSubShape(theShape, anArray);
 519   }
 520
 521   //The GetPoint() doesn't change object so no new function is required.
 522   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
 523   TCollection_AsciiString anOldDescr = aFunction->GetDescription();
 524
 525   //Make a Python command
 526   GEOM::TPythonDump(aFunction) << anOldDescr.ToCString() << "\n\t"
 527     << aResult << " = geompy.GetPoint(" << theShape << ", "
 528       << theX << ", " << theY << ", " << theZ << ", " << theEpsilon << ")";
 529
 530   SetErrorCode(OK);
 531   return aResult;
 532 }
 533
 534 //=============================================================================
 535 /*!
 536  *  GetEdge
 537  */
 538 //=============================================================================
 539 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::GetEdge
 540                                                 (Handle(GEOM_Object) theShape,
 541                                                  Handle(GEOM_Object) thePoint1,
 542                                                  Handle(GEOM_Object) thePoint2)
 543 {
 544   SetErrorCode(KO);
 545
 546   //New Edge object
 547   Handle(GEOM_Object) aResult;
 548
 549   // Arguments
 550   if (theShape.IsNull() || thePoint1.IsNull() || thePoint2.IsNull()) return NULL;
 551
 552   TopoDS_Shape aBlockOrComp = theShape->GetValue();
 553   if (aBlockOrComp.IsNull()) {
 554     SetErrorCode("Block or compound is null");
 555     return NULL;
 556   }
 557   if (aBlockOrComp.ShapeType() != TopAbs_SOLID &&
 558       aBlockOrComp.ShapeType() != TopAbs_COMPOUND &&
 559       aBlockOrComp.ShapeType() != TopAbs_COMPSOLID) {
 560     SetErrorCode("Shape is neither a block, nor a compound of blocks");
 561     return NULL;
 562   }
 563
 564   TopoDS_Shape anArg1 = thePoint1->GetValue();
 565   TopoDS_Shape anArg2 = thePoint2->GetValue();
 566   if (anArg1.IsNull() || anArg2.IsNull()) {
 567     SetErrorCode("Null shape is given as argument");
 568     return NULL;
 569   }
 570   if (anArg1.ShapeType() != TopAbs_VERTEX ||
 571       anArg2.ShapeType() != TopAbs_VERTEX) {
 572     SetErrorCode("Element for edge identification is not a vertex");
 573     return NULL;
 574   }
 575
 576   //Compute the Edge value
 577   try {
 578     TopTools_IndexedDataMapOfShapeListOfShape MVE;
 579     GEOMImpl_Block6Explorer::MapShapesAndAncestors
 580       (aBlockOrComp, TopAbs_VERTEX, TopAbs_EDGE, MVE);
 581
 582     TopoDS_Shape V1,V2;
 583     Standard_Integer ish, ext = MVE.Extent();
 584
 585     if (MVE.Contains(anArg1)) {
 586       V1 = anArg1;
 587     } else {
 588       for (ish = 1; ish <= ext; ish++) {
 589         TopoDS_Shape aShi = MVE.FindKey(ish);
 590         if (BRepTools::Compare(TopoDS::Vertex(anArg1), TopoDS::Vertex(aShi))) {
 591           V1 = aShi;
 592           break;
 593         }
 594       }
 595     }
 596
 597     if (MVE.Contains(anArg2)) {
 598       V2 = anArg2;
 599     } else {
 600       for (ish = 1; ish <= ext; ish++) {
 601         TopoDS_Shape aShi = MVE.FindKey(ish);
 602         if (BRepTools::Compare(TopoDS::Vertex(anArg2), TopoDS::Vertex(aShi))) {
 603           V2 = aShi;
 604           break;
 605         }
 606       }
 607     }
 608
 609     if (V1.IsNull() || V2.IsNull()) {
 610       SetErrorCode("The given vertex does not belong to the shape");
 611       return NULL;
 612     }
 613
 614     TopoDS_Shape anEdge;
 615     Standard_Integer isFound =
 616       GEOMImpl_Block6Explorer::FindEdge(anEdge, V1, V2, MVE, Standard_True);
 617     if (isFound == 0) {
 618       SetErrorCode("The given vertices do not belong to one edge of the given shape");
 619       return NULL;
 620     } else if (isFound > 1) {
 621       SetErrorCode("Multiple edges found by the given vertices of the shape");
 622       return NULL;
 623     } else {
 624       TopTools_IndexedMapOfShape anIndices;
 625       TopExp::MapShapes(aBlockOrComp, anIndices);
 626       Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
 627       anArray->SetValue(1, anIndices.FindIndex(anEdge));
 628       aResult = GetEngine()->AddSubShape(theShape, anArray);
 629     }
 630   } catch (Standard_Failure) {
 631     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 632     SetErrorCode(aFail->GetMessageString());
 633     return NULL;
 634   }
 635
 636   Handle(GEOM_Function) aFunction = aResult->GetLastFunction();
 637
 638   //Make a Python command
 639   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetEdge("
 640     << theShape << ", " << thePoint1 << ", " << thePoint2 << ")";
 641
 642   SetErrorCode(OK);
 643   return aResult;
 644 }
 645
 646 //=============================================================================
 647 /*!
 648  *  GetEdgeNearPoint
 649  */
 650 //=============================================================================
 651 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::GetEdgeNearPoint
 652                                                 (Handle(GEOM_Object) theShape,
 653                                                  Handle(GEOM_Object) thePoint)
 654 {
 655   SetErrorCode(KO);
 656
 657   //New object
 658   Handle(GEOM_Object) aResult;
 659
 660   // Arguments
 661   if (theShape.IsNull() || thePoint.IsNull()) return NULL;
 662
 663   TopoDS_Shape aBlockOrComp = theShape->GetValue();
 664   if (aBlockOrComp.IsNull()) {
 665     SetErrorCode("Block or compound is null");
 666     return NULL;
 667   }
 668   if (aBlockOrComp.ShapeType() != TopAbs_SOLID &&
 669       aBlockOrComp.ShapeType() != TopAbs_COMPOUND &&
 670       aBlockOrComp.ShapeType() != TopAbs_COMPSOLID) {
 671     SetErrorCode("Shape is neither a block, nor a compound of blocks");
 672     return NULL;
 673   }
 674
 675   TopoDS_Shape anArg = thePoint->GetValue();
 676   if (anArg.IsNull()) {
 677     SetErrorCode("Null shape is given as argument");
 678     return NULL;
 679   }
 680   if (anArg.ShapeType() != TopAbs_VERTEX) {
 681     SetErrorCode("Element for edge identification is not a vertex");
 682     return NULL;
 683   }
 684
 685   //Compute the Edge value
 686   try {
 687     TopoDS_Shape aShape;
 688
 689     TopoDS_Vertex aVert = TopoDS::Vertex(anArg);
 690
 691     // 1. Explode blocks on edges
 692     TopTools_MapOfShape mapShape;
 693     Standard_Integer nbEdges = 0;
 694     TopExp_Explorer exp (aBlockOrComp, TopAbs_EDGE);
 695     for (; exp.More(); exp.Next()) {
 696       if (mapShape.Add(exp.Current())) {
 697         nbEdges++;
 698       }
 699     }
 700
 701     mapShape.Clear();
 702     Standard_Integer ind = 1;
 703     TopTools_Array1OfShape anEdges (1, nbEdges);
 704     TColStd_Array1OfReal aDistances (1, nbEdges);
 705     for (exp.Init(aBlockOrComp, TopAbs_EDGE); exp.More(); exp.Next()) {
 706       if (mapShape.Add(exp.Current())) {
 707         TopoDS_Shape anEdge = exp.Current();
 708         anEdges(ind) = anEdge;
 709
 710         // 2. Classify the point relatively each edge
 711         BRepExtrema_DistShapeShape aDistTool (aVert, anEdges(ind));
 712         if (!aDistTool.IsDone()) {
 713           SetErrorCode("Can not find a distance from the given point to one of edges");
 714           return NULL;
 715         }
 716         aDistances(ind) = aDistTool.Value();
 717         ind++;
 718       }
 719     }
 720
 721     // 3. Define edge, having minimum distance to the point
 722     Standard_Real nearest = RealLast(), nbFound = 0;
 723     Standard_Real prec = Precision::Confusion();
 724     for (ind = 1; ind <= nbEdges; ind++) {
 725       if (Abs(aDistances(ind) - nearest) < prec) {
 726         nbFound++;
 727       } else if (aDistances(ind) < nearest) {
 728         nearest = aDistances(ind);
 729         aShape = anEdges(ind);
 730         nbFound = 1;
 731       } else {
 732       }
 733     }
 734     if (nbFound > 1) {
 735       SetErrorCode("Multiple edges near the given point are found");
 736       return NULL;
 737     } else if (nbFound == 0) {
 738       SetErrorCode("There are no edges near the given point");
 739       return NULL;
 740     } else {
 741       TopTools_IndexedMapOfShape anIndices;
 742       TopExp::MapShapes(aBlockOrComp, anIndices);
 743       Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
 744       anArray->SetValue(1, anIndices.FindIndex(aShape));
 745       aResult = GetEngine()->AddSubShape(theShape, anArray);
 746     }
 747   }
 748   catch (Standard_Failure) {
 749     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 750     SetErrorCode(aFail->GetMessageString());
 751     return NULL;
 752   }
 753
 754   Handle(GEOM_Function) aFunction = aResult->GetLastFunction();
 755
 756   //Make a Python command
 757   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetEdgeNearPoint("
 758                                << theShape << ", " << thePoint << ")";
 759
 760   SetErrorCode(OK);
 761   return aResult;
 762 }
 763
 764 //=============================================================================
 765 /*!
 766  *  GetFaceByPoints
 767  */
 768 //=============================================================================
 769 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::GetFaceByPoints
 770                                                 (Handle(GEOM_Object) theShape,
 771                                                  Handle(GEOM_Object) thePoint1,
 772                                                  Handle(GEOM_Object) thePoint2,
 773                                                  Handle(GEOM_Object) thePoint3,
 774                                                  Handle(GEOM_Object) thePoint4)
 775 {
 776   SetErrorCode(KO);
 777
 778   //New object
 779   Handle(GEOM_Object) aResult;
 780
 781   // Arguments
 782   if (theShape.IsNull() ||
 783       thePoint1.IsNull() || thePoint2.IsNull() ||
 784       thePoint3.IsNull() || thePoint4.IsNull()) return NULL;
 785
 786   TopoDS_Shape aBlockOrComp = theShape->GetValue();
 787   if (aBlockOrComp.IsNull()) {
 788     SetErrorCode("Block or compound is null");
 789     return NULL;
 790   }
 791   if (aBlockOrComp.ShapeType() != TopAbs_SOLID &&
 792       aBlockOrComp.ShapeType() != TopAbs_COMPOUND &&
 793       aBlockOrComp.ShapeType() != TopAbs_COMPSOLID) {
 794     SetErrorCode("Shape is neither a block, nor a compound of blocks");
 795     return NULL;
 796   }
 797
 798   TopoDS_Shape anArg1 = thePoint1->GetValue();
 799   TopoDS_Shape anArg2 = thePoint2->GetValue();
 800   TopoDS_Shape anArg3 = thePoint3->GetValue();
 801   TopoDS_Shape anArg4 = thePoint4->GetValue();
 802   if (anArg1.IsNull() || anArg2.IsNull() ||
 803       anArg3.IsNull() || anArg4.IsNull()) {
 804     SetErrorCode("Null shape is given as argument");
 805     return NULL;
 806   }
 807   if (anArg1.ShapeType() != TopAbs_VERTEX ||
 808       anArg2.ShapeType() != TopAbs_VERTEX ||
 809       anArg3.ShapeType() != TopAbs_VERTEX ||
 810       anArg4.ShapeType() != TopAbs_VERTEX) {
 811     SetErrorCode("Element for face identification is not a vertex");
 812     return NULL;
 813   }
 814
 815   //Compute the Face value
 816   try {
 817     TopoDS_Shape aShape;
 818
 819     TopTools_IndexedDataMapOfShapeListOfShape MVF;
 820     GEOMImpl_Block6Explorer::MapShapesAndAncestors(aBlockOrComp, TopAbs_VERTEX, TopAbs_FACE, MVF);
 821
 822     TopoDS_Shape V1,V2,V3,V4;
 823     Standard_Integer ish, ext = MVF.Extent();
 824
 825     if (MVF.Contains(anArg1)) {
 826       V1 = anArg1;
 827     } else {
 828       for (ish = 1; ish <= ext; ish++) {
 829         TopoDS_Shape aShi = MVF.FindKey(ish);
 830         if (BRepTools::Compare(TopoDS::Vertex(anArg1), TopoDS::Vertex(aShi))) {
 831           V1 = aShi;
 832           break;
 833         }
 834       }
 835     }
 836
 837     if (MVF.Contains(anArg2)) {
 838       V2 = anArg2;
 839     } else {
 840       for (ish = 1; ish <= ext; ish++) {
 841         TopoDS_Shape aShi = MVF.FindKey(ish);
 842         if (BRepTools::Compare(TopoDS::Vertex(anArg2), TopoDS::Vertex(aShi))) {
 843           V2 = aShi;
 844           break;
 845         }
 846       }
 847     }
 848
 849     if (MVF.Contains(anArg3)) {
 850       V3 = anArg3;
 851     } else {
 852       for (ish = 1; ish <= ext; ish++) {
 853         TopoDS_Shape aShi = MVF.FindKey(ish);
 854         if (BRepTools::Compare(TopoDS::Vertex(anArg3), TopoDS::Vertex(aShi))) {
 855           V3 = aShi;
 856           break;
 857         }
 858       }
 859     }
 860
 861     if (MVF.Contains(anArg4)) {
 862       V4 = anArg4;
 863     } else {
 864       for (ish = 1; ish <= ext; ish++) {
 865         TopoDS_Shape aShi = MVF.FindKey(ish);
 866         if (BRepTools::Compare(TopoDS::Vertex(anArg4), TopoDS::Vertex(aShi))) {
 867           V4 = aShi;
 868           break;
 869         }
 870       }
 871     }
 872
 873     if (V1.IsNull() || V2.IsNull() || V3.IsNull() || V4.IsNull()) {
 874       SetErrorCode("The given vertex does not belong to the shape");
 875       return NULL;
 876     }
 877
 878     Standard_Integer isFound =
 879       GEOMImpl_Block6Explorer::FindFace(aShape, V1, V2, V3, V4, MVF, Standard_True);
 880     if (isFound == 0) {
 881       SetErrorCode("The given vertices do not belong to one face of the given shape");
 882       return NULL;
 883     } else if (isFound > 1) {
 884       SetErrorCode("The given vertices belong to several faces of the given shape");
 885       return NULL;
 886     } else {
 887       TopTools_IndexedMapOfShape anIndices;
 888       TopExp::MapShapes(aBlockOrComp, anIndices);
 889       Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
 890       anArray->SetValue(1, anIndices.FindIndex(aShape));
 891       aResult = GetEngine()->AddSubShape(theShape, anArray);
 892     }
 893   }
 894   catch (Standard_Failure) {
 895     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
 896     SetErrorCode(aFail->GetMessageString());
 897     return NULL;
 898   }
 899
 900   Handle(GEOM_Function) aFunction = aResult->GetLastFunction();
 901
 902   //Make a Python command
 903   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetFaceByPoints("
 904     << theShape << ", " << thePoint1 << ", " << thePoint2
 905       << ", " << thePoint3 << ", " << thePoint4 << ")";
 906
 907   SetErrorCode(OK);
 908   return aResult;
 909 }
 910
 911 //=============================================================================
 912 /*!
 913  *  GetFaceByEdges
 914  */
 915 //=============================================================================
 916 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::GetFaceByEdges
 917                                                 (Handle(GEOM_Object) theShape,
 918                                                  Handle(GEOM_Object) theEdge1,
 919                                                  Handle(GEOM_Object) theEdge2)
 920 {
 921   SetErrorCode(KO);
 922
 923   //New object
 924   Handle(GEOM_Object) aResult;
 925
 926   // Arguments
 927   if (theShape.IsNull() || theEdge1.IsNull() || theEdge2.IsNull()) return NULL;
 928
 929   TopoDS_Shape aBlockOrComp = theShape->GetValue();
 930   if (aBlockOrComp.IsNull()) {
 931     SetErrorCode("Block or compound is null");
 932     return NULL;
 933   }
 934   if (aBlockOrComp.ShapeType() != TopAbs_SOLID &&
 935       aBlockOrComp.ShapeType() != TopAbs_COMPOUND &&
 936       aBlockOrComp.ShapeType() != TopAbs_COMPSOLID) {
 937     SetErrorCode("Shape is neither a block, nor a compound of blocks");
 938     return NULL;
 939   }
 940
 941   TopoDS_Shape anArg1 = theEdge1->GetValue();
 942   TopoDS_Shape anArg2 = theEdge2->GetValue();
 943   if (anArg1.IsNull() || anArg2.IsNull()) {
 944     SetErrorCode("Null shape is given as argument");
 945     return NULL;
 946   }
 947   if (anArg1.ShapeType() != TopAbs_EDGE ||
 948       anArg2.ShapeType() != TopAbs_EDGE) {
 949     SetErrorCode("Element for face identification is not an edge");
 950     return NULL;
 951   }
 952
 953   //Compute the Face value
 954   try {
 955     TopoDS_Shape aShape;
 956
 957     TopTools_IndexedDataMapOfShapeListOfShape MEF;
 958     GEOMImpl_Block6Explorer::MapShapesAndAncestors(aBlockOrComp, TopAbs_EDGE, TopAbs_FACE, MEF);
 959
 960     TopoDS_Shape E1,E2;
 961     Standard_Integer ish, ext = MEF.Extent();
 962
 963     if (MEF.Contains(anArg1)) {
 964       E1 = anArg1;
 965     } else {
 966       for (ish = 1; ish <= ext; ish++) {
 967         TopoDS_Shape aShi = MEF.FindKey(ish);
 968         if (GEOMImpl_Block6Explorer::IsSimilarEdges(anArg1, aShi)) {
 969           E1 = aShi;
 970         }
 971       }
 972     }
 973
 974     if (MEF.Contains(anArg2)) {
 975       E2 = anArg2;
 976     } else {
 977       for (ish = 1; ish <= ext; ish++) {
 978         TopoDS_Shape aShi = MEF.FindKey(ish);
 979         if (GEOMImpl_Block6Explorer::IsSimilarEdges(anArg2, aShi)) {
 980           E2 = aShi;
 981         }
 982       }
 983     }
 984
 985     if (E1.IsNull() || E2.IsNull()) {
 986       SetErrorCode("The given edge does not belong to the shape");
 987       return NULL;
 988     }
 989
 990     const TopTools_ListOfShape& aFacesOfE1 = MEF.FindFromKey(E1);
 991     const TopTools_ListOfShape& aFacesOfE2 = MEF.FindFromKey(E2);
 992
 993     Standard_Integer isFound = 0;
 994     TopTools_ListIteratorOfListOfShape anIterF1 (aFacesOfE1);
 995     for (; anIterF1.More(); anIterF1.Next()) {
 996
 997       TopTools_ListIteratorOfListOfShape anIterF2 (aFacesOfE2);
 998       for (; anIterF2.More(); anIterF2.Next()) {
 999
1000         if (anIterF1.Value().IsSame(anIterF2.Value())) {
1001           isFound++;
1002
1003           // Store the face, defined by two edges
1004           aShape = anIterF1.Value();
1005         }
1006       }
1007     }
1008     if (isFound == 0) {
1009       SetErrorCode("The given edges do not belong to one face of the given shape");
1010       return NULL;
1011     } else if (isFound > 1) {
1012       SetErrorCode("The given edges belong to several faces of the given shape");
1013       return NULL;
1014     } else {
1015       TopTools_IndexedMapOfShape anIndices;
1016       TopExp::MapShapes(aBlockOrComp, anIndices);
1017       Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
1018       anArray->SetValue(1, anIndices.FindIndex(aShape));
1019       aResult = GetEngine()->AddSubShape(theShape, anArray);
1020     }
1021   }
1022   catch (Standard_Failure) {
1023     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1024     SetErrorCode(aFail->GetMessageString());
1025     return NULL;
1026   }
1027
1028   Handle(GEOM_Function) aFunction = aResult->GetLastFunction();
1029
1030   //Make a Python command
1031   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetFaceByEdges("
1032     << theShape << ", " << theEdge1 << ", " << theEdge2 << ")";
1033
1034   SetErrorCode(OK);
1035   return aResult;
1036 }
1037
1038 //=============================================================================
1039 /*!
1040  *  GetOppositeFace
1041  */
1042 //=============================================================================
1043 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::GetOppositeFace
1044                                                 (Handle(GEOM_Object) theShape,
1045                                                  Handle(GEOM_Object) theFace)
1046 {
1047   SetErrorCode(KO);
1048
1049   //New object
1050   Handle(GEOM_Object) aResult;
1051
1052   // Arguments
1053   if (theShape.IsNull() || theFace.IsNull()) return NULL;
1054
1055   TopoDS_Shape aBlockOrComp = theShape->GetValue();
1056   if (aBlockOrComp.IsNull()) {
1057     SetErrorCode("Block is null");
1058     return NULL;
1059   }
1060   if (aBlockOrComp.ShapeType() != TopAbs_SOLID) {
1061     SetErrorCode("Shape is not a block");
1062     return NULL;
1063   }
1064
1065   TopoDS_Shape anArg = theFace->GetValue();
1066   if (anArg.IsNull()) {
1067     SetErrorCode("Null shape is given as argument");
1068     return NULL;
1069   }
1070   if (anArg.ShapeType() != TopAbs_FACE) {
1071     SetErrorCode("Element for face identification is not a face");
1072     return NULL;
1073   }
1074
1075   //Compute the Face value
1076   try {
1077     TopoDS_Shape aShape;
1078
1079     GEOMImpl_Block6Explorer aBlockTool;
1080     aBlockTool.InitByBlockAndFace(aBlockOrComp, anArg);
1081     aShape = aBlockTool.GetFace(2);
1082
1083     TopTools_IndexedMapOfShape anIndices;
1084     TopExp::MapShapes(aBlockOrComp, anIndices);
1085     Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
1086     anArray->SetValue(1, anIndices.FindIndex(aShape));
1087     aResult = GetEngine()->AddSubShape(theShape, anArray);
1088   }
1089   catch (Standard_Failure) {
1090     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1091     SetErrorCode(aFail->GetMessageString());
1092     return NULL;
1093   }
1094
1095   Handle(GEOM_Function) aFunction = aResult->GetLastFunction();
1096
1097   //Make a Python command
1098   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetOppositeFace("
1099                                << theShape << ", " << theFace << ")";
1100
1101   SetErrorCode(OK);
1102   return aResult;
1103 }
1104
1105 //=============================================================================
1106 /*!
1107  *  GetFaceNearPoint
1108  */
1109 //=============================================================================
1110 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::GetFaceNearPoint
1111                                                 (Handle(GEOM_Object) theShape,
1112                                                  Handle(GEOM_Object) thePoint)
1113 {
1114   SetErrorCode(KO);
1115
1116   //New object
1117   Handle(GEOM_Object) aResult;
1118
1119   // Arguments
1120   if (theShape.IsNull() || thePoint.IsNull()) return NULL;
1121
1122   TopoDS_Shape aBlockOrComp = theShape->GetValue();
1123   if (aBlockOrComp.IsNull()) {
1124     SetErrorCode("Block or compound is null");
1125     return NULL;
1126   }
1127   if (aBlockOrComp.ShapeType() != TopAbs_SOLID &&
1128       aBlockOrComp.ShapeType() != TopAbs_COMPOUND &&
1129       aBlockOrComp.ShapeType() != TopAbs_COMPSOLID) {
1130     SetErrorCode("Shape is neither a block, nor a compound of blocks");
1131     return NULL;
1132   }
1133
1134   TopoDS_Shape anArg = thePoint->GetValue();
1135   if (anArg.IsNull()) {
1136     SetErrorCode("Null shape is given as argument");
1137     return NULL;
1138   }
1139   if (anArg.ShapeType() != TopAbs_VERTEX) {
1140     SetErrorCode("Element for face identification is not a vertex");
1141     return NULL;
1142   }
1143
1144   //Compute the Face value
1145   try {
1146     TopoDS_Shape aShape;
1147
1148     TopoDS_Vertex aVert = TopoDS::Vertex(anArg);
1149     gp_Pnt aPnt = BRep_Tool::Pnt(aVert);
1150     Standard_Real PX, PY, PZ;
1151     aPnt.Coord(PX, PY, PZ);
1152
1153     // 1. Classify the point relatively each face
1154     Standard_Integer nearest = 2, nbFound = 0;
1155     TopTools_DataMapOfShapeInteger mapShapeDist;
1156     TopExp_Explorer exp (aBlockOrComp, TopAbs_FACE);
1157     for (; exp.More(); exp.Next()) {
1158       TopoDS_Shape aFace = exp.Current();
1159
1160       if (!mapShapeDist.IsBound(aFace)) {
1161         Standard_Integer aDistance = 2;
1162
1163         // 1.a. Classify relatively Surface
1164         Handle(Geom_Surface) aSurf = BRep_Tool::Surface(TopoDS::Face(aFace));
1165         Handle(ShapeAnalysis_Surface) aSurfAna = new ShapeAnalysis_Surface (aSurf);
1166         gp_Pnt2d p2dOnSurf = aSurfAna->ValueOfUV(aPnt, Precision::Confusion());
1167         gp_Pnt p3dOnSurf = aSurfAna->Value(p2dOnSurf);
1168         Standard_Real aDist = p3dOnSurf.Distance(aPnt);
1169         if (aDist > Precision::Confusion()) {
1170           // OUT of Surface
1171           aDistance = 1;
1172         } else {
1173           // 1.b. Classify relatively the face itself
1174           BRepClass_FaceClassifier FC (TopoDS::Face(aFace), p2dOnSurf, Precision::Confusion());
1175           if (FC.State() == TopAbs_IN) {
1176             aDistance = -1;
1177           } else if (FC.State() == TopAbs_ON) {
1178             aDistance = 0;
1179           } else { // OUT
1180             aDistance = 1;
1181           }
1182         }
1183
1184         if (aDistance < nearest) {
1185           nearest = aDistance;
1186           aShape = aFace;
1187           nbFound = 1;
1188
1189           // A first found face, containing the point inside, will be returned.
1190           // It is the solution, if there are no
1191           // coincident or intersecting faces in the compound.
1192           if (nearest == -1) break;
1193
1194         } else if (aDistance == nearest) {
1195           nbFound++;
1196         } else {
1197         }
1198
1199         mapShapeDist.Bind(aFace, aDistance);
1200       } // if (!mapShapeDist.IsBound(aFace))
1201     }
1202
1203     // 2. Define face, containing the point or having minimum distance to it
1204     if (nbFound > 1) {
1205       if (nearest == 0) {
1206         // The point is on boundary of some faces and there are
1207         // no faces, having the point inside
1208         SetErrorCode("Multiple faces near the given point are found");
1209         return NULL;
1210
1211       } else if (nearest == 1) {
1212         // The point is outside some faces and there are
1213         // no faces, having the point inside or on boundary.
1214         // We will get a nearest face
1215         Standard_Real bigReal = RealLast();
1216         Standard_Real minDist = bigReal;
1217         TopTools_DataMapIteratorOfDataMapOfShapeInteger mapShapeDistIter (mapShapeDist);
1218         for (; mapShapeDistIter.More(); mapShapeDistIter.Next()) {
1219           if (mapShapeDistIter.Value() == 1) {
1220             TopoDS_Shape aFace = mapShapeDistIter.Key();
1221             Standard_Real aDist = bigReal;
1222
1223             // 2.a. Fast check of distance - if point projection on surface is on face
1224             Handle(Geom_Surface) aSurf = BRep_Tool::Surface(TopoDS::Face(aFace));
1225             Handle(ShapeAnalysis_Surface) aSurfAna = new ShapeAnalysis_Surface (aSurf);
1226             gp_Pnt2d p2dOnSurf = aSurfAna->ValueOfUV(aPnt, Precision::Confusion());
1227             gp_Pnt p3dOnSurf = aSurfAna->Value(p2dOnSurf);
1228             aDist = p3dOnSurf.Distance(aPnt);
1229
1230             BRepClass_FaceClassifier FC (TopoDS::Face(aFace), p2dOnSurf, Precision::Confusion());
1231             if (FC.State() == TopAbs_OUT) {
1232               if (aDist < minDist) {
1233                 // 2.b. Slow check - if point projection on surface is outside of face
1234                 BRepExtrema_DistShapeShape aDistTool (aVert, aFace);
1235                 if (!aDistTool.IsDone()) {
1236                   SetErrorCode("Can not find a distance from the given point to one of faces");
1237                   return NULL;
1238                 }
1239                 aDist = aDistTool.Value();
1240               } else {
1241                 aDist = bigReal;
1242               }
1243             }
1244
1245             if (aDist < minDist) {
1246               minDist = aDist;
1247               aShape = aFace;
1248             }
1249           }
1250         }
1251       } else { // nearest == -1
1252 //        // The point is inside some faces.
1253 //        // We will get a face with nearest center
1254 //        Standard_Real minDist = RealLast();
1255 //        TopTools_DataMapIteratorOfDataMapOfShapeInteger mapShapeDistIter (mapShapeDist);
1256 //        for (; mapShapeDistIter.More(); mapShapeDistIter.Next()) {
1257 //          if (mapShapeDistIter.Value() == -1) {
1258 //            TopoDS_Shape aFace = mapShapeDistIter.Key();
1259 //            GProp_GProps aSystem;
1260 //            BRepGProp::SurfaceProperties(aFace, aSystem);
1261 //            gp_Pnt aCenterMass = aSystem.CentreOfMass();
1262 //
1263 //            Standard_Real aDist = aCenterMass.Distance(aPnt);
1264 //            if (aDist < minDist) {
1265 //              minDist = aDist;
1266 //              aShape = aFace;
1267 //            }
1268 //          }
1269 //        }
1270       }
1271     } // if (nbFound > 1)
1272
1273     if (nbFound == 0) {
1274       SetErrorCode("There are no faces near the given point");
1275       return NULL;
1276     } else {
1277       TopTools_IndexedMapOfShape anIndices;
1278       TopExp::MapShapes(aBlockOrComp, anIndices);
1279       Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
1280       anArray->SetValue(1, anIndices.FindIndex(aShape));
1281       aResult = GetEngine()->AddSubShape(theShape, anArray);
1282     }
1283   }
1284   catch (Standard_Failure) {
1285     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1286     SetErrorCode(aFail->GetMessageString());
1287     return NULL;
1288   }
1289
1290   Handle(GEOM_Function) aFunction = aResult->GetLastFunction();
1291
1292   //Make a Python command
1293   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetFaceNearPoint("
1294                                << theShape << ", " << thePoint << ")";
1295
1296   SetErrorCode(OK);
1297   return aResult;
1298 }
1299
1300 //=============================================================================
1301 /*!
1302  *  GetFaceByNormale
1303  */
1304 //=============================================================================
1305 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::GetFaceByNormale
1306                                                 (Handle(GEOM_Object) theShape,
1307                                                  Handle(GEOM_Object) theVector)
1308 {
1309   SetErrorCode(KO);
1310
1311   //New object
1312   Handle(GEOM_Object) aResult;
1313
1314   // Arguments
1315   if (theShape.IsNull() || theVector.IsNull()) return NULL;
1316
1317   TopoDS_Shape aBlockOrComp = theShape->GetValue();
1318   if (aBlockOrComp.IsNull()) {
1319     SetErrorCode("Block or compound is null");
1320     return NULL;
1321   }
1322   if (aBlockOrComp.ShapeType() != TopAbs_SOLID &&
1323       aBlockOrComp.ShapeType() != TopAbs_COMPOUND &&
1324       aBlockOrComp.ShapeType() != TopAbs_COMPSOLID) {
1325     SetErrorCode("Shape is neither a block, nor a compound of blocks");
1326     return NULL;
1327   }
1328
1329   TopoDS_Shape anArg = theVector->GetValue();
1330   if (anArg.IsNull()) {
1331     SetErrorCode("Null shape is given as argument");
1332     return NULL;
1333   }
1334   if (anArg.ShapeType() != TopAbs_EDGE) {
1335     SetErrorCode("Element for normale identification is not an edge");
1336     return NULL;
1337   }
1338
1339   //Compute the Face value
1340   try {
1341     TopoDS_Shape aShape;
1342
1343     TopoDS_Edge anEdge = TopoDS::Edge(anArg);
1344     TopoDS_Vertex V1, V2;
1345     TopExp::Vertices(anEdge, V1, V2, Standard_True);
1346     gp_Pnt P1 = BRep_Tool::Pnt(V1);
1347     gp_Pnt P2 = BRep_Tool::Pnt(V2);
1348     gp_Vec aVec (P1, P2);
1349     if (aVec.Magnitude() < Precision::Confusion()) {
1350       SetErrorCode("Vector with null magnitude is given");
1351       return NULL;
1352     }
1353
1354     Standard_Real minAngle = RealLast();
1355     TopTools_MapOfShape mapShape;
1356     TopExp_Explorer exp (aBlockOrComp, TopAbs_FACE);
1357     for (; exp.More(); exp.Next()) {
1358       if (mapShape.Add(exp.Current())) {
1359         TopoDS_Face aFace = TopoDS::Face(exp.Current());
1360         BRepAdaptor_Surface SF (aFace);
1361
1362         Standard_Real u, v, x;
1363
1364         // find a point on the surface to get normal direction in
1365         u = SF.FirstUParameter();
1366         x = SF.LastUParameter();
1367         if (Precision::IsInfinite(u)) {
1368           u =  (Precision::IsInfinite(x)) ? 0. : x;
1369         } else if (!Precision::IsInfinite(x)) {
1370           u = (u+x) / 2.;
1371         }
1372
1373         v = SF.FirstVParameter();
1374         x = SF.LastVParameter();
1375         if (Precision::IsInfinite(v)) {
1376           v =  (Precision::IsInfinite(x)) ? 0. : x;
1377         } else if (!Precision::IsInfinite(x)) {
1378           v = (v+x) / 2.;
1379         }
1380
1381         // compute the normal direction
1382         gp_Vec Vec1,Vec2;
1383         SF.D1(u,v,P1,Vec1,Vec2);
1384         gp_Vec V = Vec1.Crossed(Vec2);
1385         x = V.Magnitude();
1386         if (V.Magnitude() < Precision::Confusion()) {
1387           SetErrorCode("Normal vector of a face has null magnitude");
1388           return NULL;
1389         }
1390
1391         // consider the face orientation
1392         if (aFace.Orientation() == TopAbs_REVERSED ||
1393             aFace.Orientation() == TopAbs_INTERNAL) {
1394           V = - V;
1395         }
1396
1397         // compute the angle and compare with the minimal one
1398         Standard_Real anAngle = aVec.Angle(V);
1399         if (anAngle < minAngle) {
1400           minAngle = anAngle;
1401           aShape = aFace;
1402         }
1403       }
1404     }
1405
1406     if (aShape.IsNull()) {
1407       SetErrorCode("Failed to find a face by the given normale");
1408       return NULL;
1409     } else {
1410       TopTools_IndexedMapOfShape anIndices;
1411       TopExp::MapShapes(aBlockOrComp, anIndices);
1412       Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
1413       anArray->SetValue(1, anIndices.FindIndex(aShape));
1414       aResult = GetEngine()->AddSubShape(theShape, anArray);
1415     }
1416   }
1417   catch (Standard_Failure) {
1418     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1419     SetErrorCode(aFail->GetMessageString());
1420     return NULL;
1421   }
1422
1423   Handle(GEOM_Function) aFunction = aResult->GetLastFunction();
1424
1425   //Make a Python command
1426   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetFaceByNormale("
1427     << theShape << ", " << theVector << ")";
1428
1429   SetErrorCode(OK);
1430   return aResult;
1431 }
1432
1433 //=============================================================================
1434 /*!
1435  *  IsCompoundOfBlocks
1436  */
1437 //=============================================================================
1438 Standard_Boolean GEOMImpl_IBlocksOperations::IsCompoundOfBlocks
1439                                                 (Handle(GEOM_Object)    theCompound,
1440                                                  const Standard_Integer theMinNbFaces,
1441                                                  const Standard_Integer theMaxNbFaces,
1442                                                  Standard_Integer&      theNbBlocks)
1443 {
1444   SetErrorCode(KO);
1445   Standard_Boolean isCompOfBlocks = Standard_False;
1446   theNbBlocks = 0;
1447
1448   if (theCompound.IsNull()) return isCompOfBlocks;
1449   TopoDS_Shape aBlockOrComp = theCompound->GetValue();
1450
1451   //Check
1452   isCompOfBlocks = Standard_True;
1453   try {
1454     TopTools_MapOfShape mapShape;
1455     TopExp_Explorer exp (aBlockOrComp, TopAbs_SOLID);
1456     for (; exp.More(); exp.Next()) {
1457       if (mapShape.Add(exp.Current())) {
1458         TopoDS_Shape aSolid = exp.Current();
1459
1460         TopTools_MapOfShape mapFaces;
1461         TopExp_Explorer expF (aSolid, TopAbs_FACE);
1462         Standard_Integer nbFaces = 0;
1463         for (; expF.More(); expF.Next()) {
1464           if (mapFaces.Add(expF.Current())) {
1465             nbFaces++;
1466             if (nbFaces > theMaxNbFaces) {
1467               isCompOfBlocks = Standard_False;
1468               break;
1469             }
1470           }
1471         }
1472         if (nbFaces < theMinNbFaces || theMaxNbFaces < nbFaces) {
1473           isCompOfBlocks = Standard_False;
1474         } else {
1475           theNbBlocks++;
1476         }
1477       }
1478     }
1479   }
1480   catch (Standard_Failure) {
1481     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
1482     SetErrorCode(aFail->GetMessageString());
1483     return isCompOfBlocks;
1484   }
1485
1486   SetErrorCode(OK);
1487   return isCompOfBlocks;
1488 }
1489
1490 //=============================================================================
1491 /*!
1492  *  Set of functions, used by CheckCompoundOfBlocks() method
1493  */
1494 //=============================================================================
1495 void GEOMImpl_IBlocksOperations::AddBlocksFrom (const TopoDS_Shape&   theShape,
1496                                                 TopTools_ListOfShape& BLO,
1497                                                 TopTools_ListOfShape& NOT,
1498                                                 TopTools_ListOfShape& EXT)
1499 {
1500   TopAbs_ShapeEnum aType = theShape.ShapeType();
1501   switch (aType) {
1502   case TopAbs_COMPOUND:
1503   case TopAbs_COMPSOLID:
1504     {
1505       TopoDS_Iterator It (theShape);
1506       for (; It.More(); It.Next()) {
1507         AddBlocksFrom(It.Value(), BLO, NOT, EXT);
1508       }
1509     }
1510     break;
1511   case TopAbs_SOLID:
1512     {
1513       // Check, if there are seam or degenerated edges
1514       BlockFix_CheckTool aTool;
1515       aTool.SetShape(theShape);
1516       aTool.Perform();
1517       if (aTool.NbPossibleBlocks() > 0) {
1518         EXT.Append(theShape);
1519       } else {
1520         // Count faces and edges in each face to recognize blocks
1521         TopTools_MapOfShape mapFaces;
1522         Standard_Integer nbFaces = 0;
1523         Standard_Boolean hasNonQuadr = Standard_False;
1524         TopExp_Explorer expF (theShape, TopAbs_FACE);
1525
1526         for (; expF.More(); expF.Next()) {
1527           if (mapFaces.Add(expF.Current())) {
1528             nbFaces++;
1529             if (nbFaces > 6) break;
1530
1531             // get wire
1532             TopoDS_Shape aF = expF.Current();
1533             TopExp_Explorer wires (aF, TopAbs_WIRE);
1534             if (!wires.More()) {
1535               // no wire in the face
1536               hasNonQuadr = Standard_True;
1537               break;
1538             }
1539             TopoDS_Shape aWire = wires.Current();
1540             wires.Next();
1541             if (wires.More()) {
1542               // multiple wires in the face
1543               hasNonQuadr = Standard_True;
1544               break;
1545             }
1546
1547             // Check number of edges in the face
1548             Standard_Integer nbEdges = 0;
1549             TopTools_MapOfShape mapEdges;
1550             TopExp_Explorer expW (aWire, TopAbs_EDGE);
1551             for (; expW.More(); expW.Next()) {
1552               if (mapEdges.Add(expW.Current())) {
1553                 nbEdges++;
1554                 if (nbEdges > 4) break;
1555               }
1556             }
1557             if (nbEdges != 4) {
1558               hasNonQuadr = Standard_True;
1559             }
1560           }
1561         }
1562
1563         if (nbFaces == 6 && !hasNonQuadr) {
1564           BLO.Append(theShape);
1565         } else {
1566           NOT.Append(theShape);
1567         }
1568       }
1569     }
1570     break;
1571   default:
1572     NOT.Append(theShape);
1573   }
1574 }
1575
1576 void AddBlocksFromOld (const TopoDS_Shape&   theShape,
1577                        TopTools_ListOfShape& BLO,
1578                        TopTools_ListOfShape& NOT,
1579                        TopTools_ListOfShape& DEG,
1580                        TopTools_ListOfShape& SEA)
1581 {
1582   TopAbs_ShapeEnum aType = theShape.ShapeType();
1583   switch (aType) {
1584   case TopAbs_COMPOUND:
1585   case TopAbs_COMPSOLID:
1586     {
1587       TopoDS_Iterator It (theShape);
1588       for (; It.More(); It.Next()) {
1589         AddBlocksFromOld(It.Value(), BLO, NOT, DEG, SEA);
1590       }
1591     }
1592     break;
1593   case TopAbs_SOLID:
1594     {
1595       TopTools_MapOfShape mapFaces;
1596       TopExp_Explorer expF (theShape, TopAbs_FACE);
1597       Standard_Integer nbFaces = 0;
1598       Standard_Boolean hasNonQuadr = Standard_False;
1599       Standard_Boolean hasDegenerated = Standard_False;
1600       Standard_Boolean hasSeam = Standard_False;
1601       for (; expF.More(); expF.Next()) {
1602         if (mapFaces.Add(expF.Current())) {
1603           nbFaces++;
1604           if (nbFaces > 6) break;
1605
1606           // Check number of edges in the face
1607           Standard_Integer nbEdges = 0;
1608           TopTools_MapOfShape mapEdges;
1609
1610           // get wire
1611           TopoDS_Shape aF = expF.Current();
1612           TopExp_Explorer wires (aF, TopAbs_WIRE);
1613           if (!wires.More()) {
1614             // no wire in the face
1615             hasNonQuadr = Standard_True;
1616             break;
1617           }
1618           TopoDS_Shape aWire = wires.Current();
1619           wires.Next();
1620           if (wires.More()) {
1621             // multiple wires in the face
1622             hasNonQuadr = Standard_True;
1623             break;
1624           }
1625
1626           // iterate on wire
1627           BRepTools_WireExplorer aWE (TopoDS::Wire(aWire), TopoDS::Face(aF));
1628           for (; aWE.More(); aWE.Next(), nbEdges++) {
1629             if (BRep_Tool::Degenerated(aWE.Current())) {
1630               // degenerated edge found
1631               hasDegenerated = Standard_True;
1632 //              break;
1633             }
1634             if (mapEdges.Contains(aWE.Current())) {
1635               // seam edge found
1636               hasSeam = Standard_True;
1637 //              break;
1638             }
1639             mapEdges.Add(aWE.Current());
1640           }
1641           if (nbEdges != 4) {
1642             hasNonQuadr = Standard_True;
1643           }
1644         }
1645       }
1646       if (nbFaces == 6) {
1647         if (hasDegenerated || hasSeam) {
1648           if (hasDegenerated) {
1649             DEG.Append(theShape);
1650           }
1651           if (hasSeam) {
1652             SEA.Append(theShape);
1653           }
1654         } else if (hasNonQuadr) {
1655           NOT.Append(theShape);
1656         } else {
1657           BLO.Append(theShape);
1658         }
1659       } else {
1660         NOT.Append(theShape);
1661       }
1662     }
1663     break;
1664   default:
1665     NOT.Append(theShape);
1666   }
1667 }
1668
1669 #define REL_NOT_CONNECTED 0
1670 #define REL_OK            1
1671 #define REL_NOT_GLUED     2
1672 #define REL_COLLISION_VV  3
1673 #define REL_COLLISION_FF  4
1674 #define REL_COLLISION_EE  5
1675 #define REL_UNKNOWN       6
1676
1677 Standard_Integer BlocksRelation (const TopoDS_Shape& theBlock1,
1678                                  const TopoDS_Shape& theBlock2)
1679 {
1680   // Compare bounding boxes before calling BRepExtrema_DistShapeShape
1681   Standard_Real Xmin1, Ymin1, Zmin1, Xmax1, Ymax1, Zmax1;
1682   Standard_Real Xmin2, Ymin2, Zmin2, Xmax2, Ymax2, Zmax2;
1683   Bnd_Box B1, B2;
1684   BRepBndLib::Add(theBlock1, B1);
1685   BRepBndLib::Add(theBlock2, B2);
1686   B1.Get(Xmin1, Ymin1, Zmin1, Xmax1, Ymax1, Zmax1);
1687   B2.Get(Xmin2, Ymin2, Zmin2, Xmax2, Ymax2, Zmax2);
1688   if (Xmax2 < Xmin1 || Xmax1 < Xmin2 ||
1689       Ymax2 < Ymin1 || Ymax1 < Ymin2 ||
1690       Zmax2 < Zmin1 || Zmax1 < Zmin2) {
1691     return REL_NOT_CONNECTED;
1692   }
1693
1694   BRepExtrema_DistShapeShape dst (theBlock1, theBlock2);
1695   if (!dst.IsDone()) {
1696     return REL_UNKNOWN;
1697   }
1698
1699   if (dst.Value() > Precision::Confusion()) {
1700     return REL_NOT_CONNECTED;
1701   }
1702
1703   if (dst.InnerSolution()) {
1704     return REL_COLLISION_VV;
1705   }
1706
1707   Standard_Integer nbSol = dst.NbSolution();
1708   Standard_Integer relation = REL_OK;
1709   Standard_Integer nbVerts = 0;
1710   Standard_Integer nbEdges = 0;
1711   Standard_Integer sol = 1;
1712   for (; sol <= nbSol; sol++) {
1713     BRepExtrema_SupportType supp1 = dst.SupportTypeShape1(sol);
1714     BRepExtrema_SupportType supp2 = dst.SupportTypeShape2(sol);
1715     if (supp1 == BRepExtrema_IsVertex && supp2 == BRepExtrema_IsVertex) {
1716       nbVerts++;
1717     } else if (supp1 == BRepExtrema_IsInFace || supp2 == BRepExtrema_IsInFace) {
1718       return REL_COLLISION_FF;
1719     } else if (supp1 == BRepExtrema_IsOnEdge && supp2 == BRepExtrema_IsOnEdge) {
1720       nbEdges++;
1721     } else if ((supp1 == BRepExtrema_IsOnEdge && supp2 == BRepExtrema_IsVertex) ||
1722                (supp2 == BRepExtrema_IsOnEdge && supp1 == BRepExtrema_IsVertex)) {
1723       relation = REL_COLLISION_EE;
1724     } else {
1725     }
1726   }
1727
1728   if (relation != REL_OK) {
1729     return relation;
1730   }
1731
1732   TColStd_Array1OfInteger vertSol (1, nbVerts);
1733   TopTools_Array1OfShape V1 (1, nbVerts);
1734   TopTools_Array1OfShape V2 (1, nbVerts);
1735   Standard_Integer ivs = 0;
1736   for (sol = 1; sol <= nbSol; sol++) {
1737     if (dst.SupportTypeShape1(sol) == BRepExtrema_IsVertex &&
1738         dst.SupportTypeShape2(sol) == BRepExtrema_IsVertex) {
1739       TopoDS_Vertex Vcur = TopoDS::Vertex(dst.SupportOnShape1(sol));
1740       // Check, that this vertex is far enough from other solution vertices.
1741       Standard_Integer ii = 1;
1742       for (; ii <= ivs; ii++) {
1743         if (BRepTools::Compare(TopoDS::Vertex(V1(ii)), Vcur)) {
1744           continue;
1745         }
1746       }
1747       ivs++;
1748       vertSol(ivs) = sol;
1749       V1(ivs) = Vcur;
1750       V2(ivs) = dst.SupportOnShape2(sol);
1751     }
1752   }
1753
1754   // As we deal only with quadrangles,
1755   // 2, 3 or 4 vertex solutions can be found.
1756   if (ivs <= 1) {
1757     if (nbEdges > 0) {
1758       return REL_COLLISION_FF;
1759     }
1760     return REL_NOT_CONNECTED;
1761   }
1762   if (ivs > 4) {
1763     return REL_UNKNOWN;
1764   }
1765
1766   // Check sharing of coincident entities.
1767   if (ivs == 2 || ivs == 3) {
1768     // Map vertices and edges of the blocks
1769     TopTools_IndexedDataMapOfShapeListOfShape MVE1, MVE2;
1770     GEOMImpl_Block6Explorer::MapShapesAndAncestors
1771       (theBlock1, TopAbs_VERTEX, TopAbs_EDGE, MVE1);
1772     GEOMImpl_Block6Explorer::MapShapesAndAncestors
1773       (theBlock2, TopAbs_VERTEX, TopAbs_EDGE, MVE2);
1774
1775     if (ivs == 2) {
1776       // Find common edge
1777       TopoDS_Shape anEdge1, anEdge2;
1778       GEOMImpl_Block6Explorer::FindEdge(anEdge1, V1(1), V1(2), MVE1);
1779       if (anEdge1.IsNull()) return REL_UNKNOWN;
1780
1781       GEOMImpl_Block6Explorer::FindEdge(anEdge2, V2(1), V2(2), MVE2);
1782       if (anEdge2.IsNull()) return REL_UNKNOWN;
1783
1784       if (!anEdge1.IsSame(anEdge2)) return REL_NOT_GLUED;
1785
1786     } else { // ivs == 3
1787       // Find common edges
1788       Standard_Integer e1_v1 = 1;
1789       Standard_Integer e1_v2 = 2;
1790       Standard_Integer e2_v1 = 3;
1791       Standard_Integer e2_v2 = 1;
1792
1793       TopoDS_Shape anEdge11, anEdge12;
1794       GEOMImpl_Block6Explorer::FindEdge(anEdge11, V1(e1_v1), V1(e1_v2), MVE1);
1795       if (anEdge11.IsNull()) {
1796         e1_v2 = 3;
1797         e2_v1 = 2;
1798         GEOMImpl_Block6Explorer::FindEdge(anEdge11, V1(e1_v1), V1(e1_v2), MVE1);
1799         if (anEdge11.IsNull()) return REL_UNKNOWN;
1800       }
1801       GEOMImpl_Block6Explorer::FindEdge(anEdge12, V1(e2_v1), V1(e2_v2), MVE1);
1802       if (anEdge12.IsNull()) {
1803         e2_v2 = 5 - e2_v1;
1804         GEOMImpl_Block6Explorer::FindEdge(anEdge12, V1(e2_v1), V1(e2_v2), MVE1);
1805         if (anEdge12.IsNull()) return REL_UNKNOWN;
1806       }
1807
1808       TopoDS_Shape anEdge21, anEdge22;
1809       GEOMImpl_Block6Explorer::FindEdge(anEdge21, V2(e1_v1), V2(e1_v2), MVE2);
1810       if (anEdge21.IsNull()) return REL_UNKNOWN;
1811       GEOMImpl_Block6Explorer::FindEdge(anEdge22, V2(e2_v1), V2(e2_v2), MVE2);
1812       if (anEdge22.IsNull()) return REL_UNKNOWN;
1813
1814       // Check of edges coincidence (with some precision) have to be done here
1815       // if (!anEdge11.IsEqual(anEdge21)) return REL_UNKNOWN;
1816       // if (!anEdge12.IsEqual(anEdge22)) return REL_UNKNOWN;
1817
1818       // Check of edges sharing
1819       if (!anEdge11.IsSame(anEdge21)) return REL_NOT_GLUED;
1820       if (!anEdge12.IsSame(anEdge22)) return REL_NOT_GLUED;
1821     }
1822   }
1823
1824   if (ivs == 4) {
1825     // Map vertices and faces of the blocks
1826     TopTools_IndexedDataMapOfShapeListOfShape MVF1, MVF2;
1827     GEOMImpl_Block6Explorer::MapShapesAndAncestors
1828       (theBlock1, TopAbs_VERTEX, TopAbs_FACE, MVF1);
1829     GEOMImpl_Block6Explorer::MapShapesAndAncestors
1830       (theBlock2, TopAbs_VERTEX, TopAbs_FACE, MVF2);
1831
1832     TopoDS_Shape aFace1, aFace2;
1833     GEOMImpl_Block6Explorer::FindFace(aFace1, V1(1), V1(2), V1(3), V1(4), MVF1);
1834     if (aFace1.IsNull()) return REL_UNKNOWN;
1835     GEOMImpl_Block6Explorer::FindFace(aFace2, V2(1), V2(2), V2(3), V2(4), MVF2);
1836     if (aFace2.IsNull()) return REL_UNKNOWN;
1837
1838     // Check of faces coincidence (with some precision) have to be done here
1839     // if (!aFace1.IsEqual(aFace2)) return REL_UNKNOWN;
1840
1841     // Check of faces sharing
1842     if (!aFace1.IsSame(aFace2)) return REL_NOT_GLUED;
1843   }
1844
1845   return REL_OK;
1846 }
1847
1848 void FindConnected (const Standard_Integer         theBlockIndex,
1849                     const TColStd_Array2OfInteger& theRelations,
1850                     TColStd_MapOfInteger&          theProcessedMap,
1851                     TColStd_MapOfInteger&          theConnectedMap)
1852 {
1853   theConnectedMap.Add(theBlockIndex);
1854   theProcessedMap.Add(theBlockIndex);
1855
1856   Standard_Integer nbBlocks = theRelations.ColLength();
1857   Standard_Integer col = 1;
1858   for (; col <= nbBlocks; col++) {
1859     if (theRelations(theBlockIndex, col) == REL_OK ||
1860         theRelations(theBlockIndex, col) == REL_NOT_GLUED) {
1861       if (!theProcessedMap.Contains(col)) {
1862         FindConnected(col, theRelations, theProcessedMap, theConnectedMap);
1863       }
1864     }
1865   }
1866 }
1867
1868 Standard_Boolean HasAnyConnection (const Standard_Integer         theBlockIndex,
1869                                    const TColStd_MapOfInteger&    theWith,
1870                                    const TColStd_Array2OfInteger& theRelations,
1871                                    TColStd_MapOfInteger&          theProcessedMap)
1872 {
1873   theProcessedMap.Add(theBlockIndex);
1874
1875   Standard_Integer nbBlocks = theRelations.ColLength();
1876   Standard_Integer col = 1;
1877   for (; col <= nbBlocks; col++) {
1878     if (theRelations(theBlockIndex, col) != REL_NOT_CONNECTED) {
1879       if (!theProcessedMap.Contains(col)) {
1880         if (theWith.Contains(col))
1881           return Standard_True;
1882         if (HasAnyConnection(col, theWith, theRelations, theProcessedMap))
1883           return Standard_True;
1884       }
1885     }
1886   }
1887
1888   return Standard_False;
1889 }
1890
1891 //=============================================================================
1892 /*!
1893  *  CheckCompoundOfBlocksOld
1894  */
1895 //=============================================================================
1896 Standard_Boolean GEOMImpl_IBlocksOperations::CheckCompoundOfBlocksOld
1897                                                 (Handle(GEOM_Object) theCompound,
1898                                                  list<BCError>&      theErrors)
1899 {
1900   SetErrorCode(KO);
1901
1902   if (theCompound.IsNull()) return Standard_False;
1903   TopoDS_Shape aBlockOrComp = theCompound->GetValue();
1904
1905   Standard_Boolean isCompOfBlocks = Standard_True;
1906
1907   // Map sub-shapes and their indices
1908   TopTools_IndexedMapOfShape anIndices;
1909   TopExp::MapShapes(aBlockOrComp, anIndices);
1910
1911   // 1. Report non-blocks
1912   TopTools_ListOfShape NOT; // Not blocks
1913   TopTools_ListOfShape DEG; // Hexahedral solids, having degenerated edges
1914   TopTools_ListOfShape SEA; // Hexahedral solids, having seam edges
1915   TopTools_ListOfShape BLO; // All blocks from the given compound
1916   AddBlocksFromOld(aBlockOrComp, BLO, NOT, DEG, SEA);
1917
1918   if (NOT.Extent() > 0) {
1919     isCompOfBlocks = Standard_False;
1920     BCError anErr;
1921     anErr.error = NOT_BLOCK;
1922     TopTools_ListIteratorOfListOfShape it (NOT);
1923     for (; it.More(); it.Next()) {
1924       anErr.incriminated.push_back(anIndices.FindIndex(it.Value()));
1925     }
1926     theErrors.push_back(anErr);
1927   }
1928
1929   if (DEG.Extent() > 0 || SEA.Extent() > 0) {
1930     isCompOfBlocks = Standard_False;
1931     BCError anErr;
1932     anErr.error = EXTRA_EDGE;
1933
1934     TopTools_ListIteratorOfListOfShape itDEG (DEG);
1935     for (; itDEG.More(); itDEG.Next()) {
1936       anErr.incriminated.push_back(anIndices.FindIndex(itDEG.Value()));
1937     }
1938
1939     TopTools_ListIteratorOfListOfShape itSEA (SEA);
1940     for (; itSEA.More(); itSEA.Next()) {
1941       anErr.incriminated.push_back(anIndices.FindIndex(itSEA.Value()));
1942     }
1943
1944     theErrors.push_back(anErr);
1945   }
1946
1947   Standard_Integer nbBlocks = BLO.Extent();
1948   if (nbBlocks == 0) {
1949     isCompOfBlocks = Standard_False;
1950     SetErrorCode(OK);
1951     return isCompOfBlocks;
1952   }
1953   if (nbBlocks == 1) {
1954     SetErrorCode(OK);
1955     return isCompOfBlocks;
1956   }
1957
1958   // Convert list of blocks into array for easy and fast access
1959   Standard_Integer ibl = 1;
1960   TopTools_Array1OfShape aBlocks (1, nbBlocks);
1961   TopTools_ListIteratorOfListOfShape BLOit (BLO);
1962   for (; BLOit.More(); BLOit.Next(), ibl++) {
1963     aBlocks.SetValue(ibl, BLOit.Value());
1964   }
1965
1966   // 2. Find relations between all blocks,
1967   //    report connection errors (NOT_GLUED and INVALID_CONNECTION)
1968   TColStd_Array2OfInteger aRelations (1, nbBlocks, 1, nbBlocks);
1969   aRelations.Init(REL_NOT_CONNECTED);
1970
1971   Standard_Integer row = 1;
1972   for (row = 1; row <= nbBlocks; row++) {
1973     TopoDS_Shape aBlock = aBlocks.Value(row);
1974
1975     Standard_Integer col = row + 1;
1976     for (; col <= nbBlocks; col++) {
1977       Standard_Integer aRel = BlocksRelation(aBlock, aBlocks.Value(col));
1978       if (aRel != REL_NOT_CONNECTED) {
1979         aRelations.SetValue(row, col, aRel);
1980         aRelations.SetValue(col, row, aRel);
1981         if (aRel == REL_NOT_GLUED) {
1982           // report connection error
1983           isCompOfBlocks = Standard_False;
1984           BCError anErr;
1985           anErr.error = NOT_GLUED;
1986           anErr.incriminated.push_back(anIndices.FindIndex(aBlocks.Value(row)));
1987           anErr.incriminated.push_back(anIndices.FindIndex(aBlocks.Value(col)));
1988           theErrors.push_back(anErr);
1989         } else if (aRel == REL_COLLISION_VV ||
1990                    aRel == REL_COLLISION_FF ||
1991                    aRel == REL_COLLISION_EE ||
1992                    aRel == REL_UNKNOWN) {
1993           // report connection error
1994           isCompOfBlocks = Standard_False;
1995           BCError anErr;
1996           anErr.error = INVALID_CONNECTION;
1997           anErr.incriminated.push_back(anIndices.FindIndex(aBlocks.Value(row)));
1998           anErr.incriminated.push_back(anIndices.FindIndex(aBlocks.Value(col)));
1999           theErrors.push_back(anErr);
2000         } else {
2001         }
2002       }
2003     }
2004   }
2005
2006   // 3. Find largest set of connected (good connection or not glued) blocks
2007   TColStd_MapOfInteger aProcessedMap;
2008   TColStd_MapOfInteger aLargestSet;
2009   TColStd_MapOfInteger aCurrentSet;
2010   for (ibl = 1; ibl <= nbBlocks; ibl++) {
2011     if (!aProcessedMap.Contains(ibl)) {
2012       aCurrentSet.Clear();
2013       FindConnected(ibl, aRelations, aProcessedMap, aCurrentSet);
2014       if (aCurrentSet.Extent() > aLargestSet.Extent()) {
2015         aLargestSet = aCurrentSet;
2016       }
2017     }
2018   }
2019
2020   // 4. Report all blocks, isolated from <aLargestSet>
2021   BCError anErr;
2022   anErr.error = NOT_CONNECTED;
2023   Standard_Boolean hasIsolated = Standard_False;
2024   for (ibl = 1; ibl <= nbBlocks; ibl++) {
2025     if (!aLargestSet.Contains(ibl)) {
2026       aProcessedMap.Clear();
2027       if (!HasAnyConnection(ibl, aLargestSet, aRelations, aProcessedMap)) {
2028         // report connection absence
2029         hasIsolated = Standard_True;
2030         anErr.incriminated.push_back(anIndices.FindIndex(aBlocks.Value(ibl)));
2031       }
2032     }
2033   }
2034   if (hasIsolated) {
2035     isCompOfBlocks = Standard_False;
2036     theErrors.push_back(anErr);
2037   }
2038
2039   SetErrorCode(OK);
2040   return isCompOfBlocks;
2041 }
2042
2043 //=============================================================================
2044 /*!
2045  *  PrintBCErrors
2046  */
2047 //=============================================================================
2048 TCollection_AsciiString GEOMImpl_IBlocksOperations::PrintBCErrors
2049                                                 (Handle(GEOM_Object)  theCompound,
2050                                                  const list<BCError>& theErrors)
2051 {
2052   TCollection_AsciiString aDescr;
2053
2054   list<BCError>::const_iterator errIt = theErrors.begin();
2055   int i = 0;
2056   for (; errIt != theErrors.end(); i++, errIt++) {
2057     BCError errStruct = *errIt;
2058
2059     switch (errStruct.error) {
2060     case NOT_BLOCK:
2061       aDescr += "\n\tNot a Blocks: ";
2062       break;
2063     case EXTRA_EDGE:
2064       aDescr += "\n\tHexahedral solids with degenerated and/or seam edges: ";
2065       break;
2066     case INVALID_CONNECTION:
2067       aDescr += "\n\tInvalid connection between two blocks: ";
2068       break;
2069     case NOT_CONNECTED:
2070       aDescr += "\n\tBlocks, not connected with main body: ";
2071       break;
2072     case NOT_GLUED:
2073       aDescr += "\n\tNot glued blocks: ";
2074       break;
2075     default:
2076       break;
2077     }
2078
2079     list<int> sshList = errStruct.incriminated;
2080     list<int>::iterator sshIt = sshList.begin();
2081     int jj = 0;
2082     for (; sshIt != sshList.end(); jj++, sshIt++) {
2083       if (jj > 0)
2084         aDescr += ", ";
2085       aDescr += TCollection_AsciiString(*sshIt);
2086     }
2087   }
2088
2089   return aDescr;
2090 }
2091
2092 //=============================================================================
2093 /*!
2094  *  CheckCompoundOfBlocks
2095  */
2096 //=============================================================================
2097 Standard_Boolean GEOMImpl_IBlocksOperations::CheckCompoundOfBlocks
2098                                               (Handle(GEOM_Object) theCompound,
2099                                                list<BCError>&      theErrors)
2100 {
2101   SetErrorCode(KO);
2102
2103   if (theCompound.IsNull()) return Standard_False;
2104   TopoDS_Shape aBlockOrComp = theCompound->GetValue();
2105
2106   Standard_Boolean isCompOfBlocks = Standard_True;
2107
2108   // Map sub-shapes and their indices
2109   TopTools_IndexedMapOfShape anIndices;
2110   TopExp::MapShapes(aBlockOrComp, anIndices);
2111
2112   // 1. Separate blocks from non-blocks
2113   TopTools_ListOfShape NOT; // Not blocks
2114   TopTools_ListOfShape EXT; // Hexahedral solids, having degenerated and/or seam edges
2115   TopTools_ListOfShape BLO; // All blocks from the given compound
2116   AddBlocksFrom(aBlockOrComp, BLO, NOT, EXT);
2117
2118   // Report non-blocks
2119   if (NOT.Extent() > 0) {
2120     isCompOfBlocks = Standard_False;
2121     BCError anErr;
2122     anErr.error = NOT_BLOCK;
2123     TopTools_ListIteratorOfListOfShape it (NOT);
2124     for (; it.More(); it.Next()) {
2125       anErr.incriminated.push_back(anIndices.FindIndex(it.Value()));
2126     }
2127     theErrors.push_back(anErr);
2128   }
2129
2130   // Report solids, having degenerated and/or seam edges
2131   if (EXT.Extent() > 0) {
2132     isCompOfBlocks = Standard_False;
2133     BCError anErr;
2134     anErr.error = EXTRA_EDGE;
2135     TopTools_ListIteratorOfListOfShape it (EXT);
2136     for (; it.More(); it.Next()) {
2137       anErr.incriminated.push_back(anIndices.FindIndex(it.Value()));
2138     }
2139     theErrors.push_back(anErr);
2140   }
2141
2142   Standard_Integer nbBlocks = BLO.Extent();
2143   if (nbBlocks == 0) {
2144     isCompOfBlocks = Standard_False;
2145     SetErrorCode(OK);
2146     return isCompOfBlocks;
2147   }
2148   if (nbBlocks == 1) {
2149     SetErrorCode(OK);
2150     return isCompOfBlocks;
2151   }
2152
2153   // Prepare data for 2. and 3.
2154   TColStd_Array2OfInteger aRelations (1, nbBlocks, 1, nbBlocks);
2155   aRelations.Init(REL_NOT_CONNECTED);
2156
2157   TopTools_IndexedMapOfShape mapBlocks;
2158
2159   BRep_Builder BB;
2160   TopoDS_Compound aComp;
2161   BB.MakeCompound(aComp);
2162
2163   TopTools_ListIteratorOfListOfShape BLOit (BLO);
2164   for (; BLOit.More(); BLOit.Next()) {
2165     mapBlocks.Add(BLOit.Value());
2166     BB.Add(aComp, BLOit.Value());
2167   }
2168
2169   // 2. Find glued blocks (having shared faces)
2170   TopTools_IndexedDataMapOfShapeListOfShape mapFaceBlocks;
2171   GEOMImpl_Block6Explorer::MapShapesAndAncestors
2172     (aComp, TopAbs_FACE, TopAbs_SOLID, mapFaceBlocks);
2173
2174   Standard_Integer prevInd = 0, curInd = 0;
2175   Standard_Integer ind = 1, nbFaces = mapFaceBlocks.Extent();
2176   for (; ind <= nbFaces; ind++) {
2177     const TopTools_ListOfShape& aGluedBlocks = mapFaceBlocks.FindFromIndex(ind);
2178     if (aGluedBlocks.Extent() > 1) { // Shared face found
2179       TopTools_ListIteratorOfListOfShape aGluedBlocksIt (aGluedBlocks);
2180       TopoDS_Shape prevBlock, curBlock;
2181       for (; aGluedBlocksIt.More(); aGluedBlocksIt.Next()) {
2182         curBlock = aGluedBlocksIt.Value();
2183         if (!prevBlock.IsNull()) {
2184           prevInd = mapBlocks.FindIndex(prevBlock);
2185           curInd  = mapBlocks.FindIndex(curBlock);
2186           aRelations.SetValue(prevInd, curInd, REL_OK);
2187           aRelations.SetValue(curInd, prevInd, REL_OK);
2188         }
2189         prevBlock = curBlock;
2190       }
2191     }
2192   }
2193
2194   // 3. Find not glued blocks
2195   GEOMAlgo_GlueAnalyser aGD;
2196
2197   aGD.SetShape(aComp);
2198   aGD.SetTolerance(Precision::Confusion());
2199   aGD.SetCheckGeometry(Standard_True);
2200   aGD.Perform();
2201
2202   Standard_Integer iErr, iWrn;
2203   iErr = aGD.ErrorStatus();
2204   if (iErr) {
2205     SetErrorCode("Error in GEOMAlgo_GlueAnalyser");
2206     return isCompOfBlocks;
2207   }
2208   iWrn = aGD.WarningStatus();
2209   if (iWrn) {
2210     MESSAGE("Warning in GEOMAlgo_GlueAnalyser");
2211   }
2212
2213   // Report not glued blocks
2214   if (aGD.HasSolidsToGlue()) {
2215     isCompOfBlocks = Standard_False;
2216     Standard_Integer aSx1Ind, aSx2Ind;
2217
2218     const GEOMAlgo_ListOfCoupleOfShapes& aLCS = aGD.SolidsToGlue();
2219     GEOMAlgo_ListIteratorOfListOfCoupleOfShapes aItCS (aLCS);
2220     for (; aItCS.More(); aItCS.Next()) {
2221       const GEOMAlgo_CoupleOfShapes& aCS = aItCS.Value();
2222       const TopoDS_Shape& aSx1 = aCS.Shape1();
2223       const TopoDS_Shape& aSx2 = aCS.Shape2();
2224
2225       aSx1Ind = mapBlocks.FindIndex(aSx1);
2226       aSx2Ind = mapBlocks.FindIndex(aSx2);
2227       aRelations.SetValue(aSx1Ind, aSx2Ind, NOT_GLUED);
2228       aRelations.SetValue(aSx2Ind, aSx1Ind, NOT_GLUED);
2229
2230       BCError anErr;
2231       anErr.error = NOT_GLUED;
2232       anErr.incriminated.push_back(anIndices.FindIndex(aSx1));
2233       anErr.incriminated.push_back(anIndices.FindIndex(aSx2));
2234       theErrors.push_back(anErr);
2235     }
2236   }
2237
2238   // 4. Find largest set of connected (good connection or not glued) blocks
2239   Standard_Integer ibl = 1;
2240   TColStd_MapOfInteger aProcessedMap;
2241   TColStd_MapOfInteger aLargestSet;
2242   TColStd_MapOfInteger aCurrentSet;
2243   for (ibl = 1; ibl <= nbBlocks; ibl++) {
2244     if (!aProcessedMap.Contains(ibl)) {
2245       aCurrentSet.Clear();
2246       FindConnected(ibl, aRelations, aProcessedMap, aCurrentSet);
2247       if (aCurrentSet.Extent() > aLargestSet.Extent()) {
2248         aLargestSet = aCurrentSet;
2249       }
2250     }
2251   }
2252
2253   // 5. Report all blocks, isolated from <aLargestSet>
2254   BCError anErr;
2255   anErr.error = NOT_CONNECTED;
2256   Standard_Boolean hasIsolated = Standard_False;
2257   for (ibl = 1; ibl <= nbBlocks; ibl++) {
2258     if (!aLargestSet.Contains(ibl)) {
2259       aProcessedMap.Clear();
2260       if (!HasAnyConnection(ibl, aLargestSet, aRelations, aProcessedMap)) {
2261         // report connection absence
2262         hasIsolated = Standard_True;
2263         anErr.incriminated.push_back(anIndices.FindIndex(mapBlocks.FindKey(ibl)));
2264       }
2265     }
2266   }
2267   if (hasIsolated) {
2268     isCompOfBlocks = Standard_False;
2269     theErrors.push_back(anErr);
2270   }
2271
2272   SetErrorCode(OK);
2273   return isCompOfBlocks;
2274 }
2275
2276 //=============================================================================
2277 /*!
2278  *  RemoveExtraEdges
2279  */
2280 //=============================================================================
2281 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::RemoveExtraEdges
2282                                              (Handle(GEOM_Object) theObject)
2283 {
2284   SetErrorCode(KO);
2285
2286   if (theObject.IsNull()) return NULL;
2287
2288   Handle(GEOM_Function) aLastFunction = theObject->GetLastFunction();
2289   if (aLastFunction.IsNull()) return NULL; //There is no function which creates an object to be fixed
2290
2291   //Add a new Copy object
2292   Handle(GEOM_Object) aCopy = GetEngine()->AddObject(GetDocID(), GEOM_COPY);
2293
2294   //Add a function
2295   Handle(GEOM_Function) aFunction =
2296     aCopy->AddFunction(GEOMImpl_BlockDriver::GetID(), BLOCK_REMOVE_EXTRA);
2297
2298   //Check if the function is set correctly
2299   if (aFunction->GetDriverGUID() != GEOMImpl_BlockDriver::GetID()) return NULL;
2300
2301   GEOMImpl_IBlockTrsf aTI (aFunction);
2302   aTI.SetOriginal(aLastFunction);
2303
2304   //Compute the fixed shape
2305   try {
2306     if (!GetSolver()->ComputeFunction(aFunction)) {
2307       SetErrorCode("Block driver failed to remove extra edges of the given shape");
2308       return NULL;
2309     }
2310   }
2311   catch (Standard_Failure) {
2312     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2313     SetErrorCode(aFail->GetMessageString());
2314     return NULL;
2315   }
2316
2317   //Make a Python command
2318   GEOM::TPythonDump(aFunction) << aCopy
2319     << " = geompy.RemoveExtraEdges(" << theObject << ")";
2320
2321   SetErrorCode(OK);
2322   return aCopy;
2323 }
2324
2325 //=============================================================================
2326 /*!
2327  *  CheckAndImprove
2328  */
2329 //=============================================================================
2330 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::CheckAndImprove
2331                                              (Handle(GEOM_Object) theObject)
2332 {
2333   SetErrorCode(KO);
2334
2335   if (theObject.IsNull()) return NULL;
2336
2337   Handle(GEOM_Function) aLastFunction = theObject->GetLastFunction();
2338   if (aLastFunction.IsNull()) return NULL; //There is no function which creates an object to be fixed
2339
2340   //Add a new Copy object
2341   Handle(GEOM_Object) aCopy = GetEngine()->AddObject(GetDocID(), GEOM_COPY);
2342
2343   //Add a function
2344   Handle(GEOM_Function) aFunction =
2345     aCopy->AddFunction(GEOMImpl_BlockDriver::GetID(), BLOCK_COMPOUND_IMPROVE);
2346
2347   //Check if the function is set correctly
2348   if (aFunction->GetDriverGUID() != GEOMImpl_BlockDriver::GetID()) return NULL;
2349
2350   GEOMImpl_IBlockTrsf aTI (aFunction);
2351   aTI.SetOriginal(aLastFunction);
2352
2353   //Compute the fixed shape
2354   try {
2355     if (!GetSolver()->ComputeFunction(aFunction)) {
2356       SetErrorCode("Block driver failed to improve the given blocks compound");
2357       return NULL;
2358     }
2359   }
2360   catch (Standard_Failure) {
2361     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2362     SetErrorCode(aFail->GetMessageString());
2363     return NULL;
2364   }
2365
2366   //Make a Python command
2367   GEOM::TPythonDump(aFunction) << aCopy
2368     << " = geompy.CheckAndImprove(" << theObject << ")";
2369
2370   SetErrorCode(OK);
2371   return aCopy;
2372 }
2373
2374 //=============================================================================
2375 /*!
2376  *  ExplodeCompoundOfBlocks
2377  */
2378 //=============================================================================
2379 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IBlocksOperations::ExplodeCompoundOfBlocks
2380                                                 (Handle(GEOM_Object)    theCompound,
2381                                                  const Standard_Integer theMinNbFaces,
2382                                                  const Standard_Integer theMaxNbFaces)
2383 {
2384   SetErrorCode(KO);
2385
2386   if (theCompound.IsNull()) return NULL;
2387   TopoDS_Shape aBlockOrComp = theCompound->GetValue();
2388   if (aBlockOrComp.IsNull()) return NULL;
2389
2390   Handle(TColStd_HSequenceOfTransient) aBlocks = new TColStd_HSequenceOfTransient;
2391   Handle(GEOM_Object) anObj;
2392   Handle(GEOM_Function) aFunction;
2393
2394   TopTools_MapOfShape mapShape;
2395   TCollection_AsciiString anAsciiList, anEntry;
2396
2397   // Map shapes
2398   TopTools_IndexedMapOfShape anIndices;
2399   TopExp::MapShapes(aBlockOrComp, anIndices);
2400   Handle(TColStd_HArray1OfInteger) anArray;
2401
2402   // Explode
2403   try {
2404     TopExp_Explorer exp (aBlockOrComp, TopAbs_SOLID);
2405     for (; exp.More(); exp.Next()) {
2406       if (mapShape.Add(exp.Current())) {
2407         TopoDS_Shape aSolid = exp.Current();
2408
2409         TopTools_MapOfShape mapFaces;
2410         TopExp_Explorer expF (aSolid, TopAbs_FACE);
2411         Standard_Integer nbFaces = 0;
2412         for (; expF.More(); expF.Next()) {
2413           if (mapFaces.Add(expF.Current())) {
2414             nbFaces++;
2415           }
2416         }
2417
2418         if (theMinNbFaces <= nbFaces && nbFaces <= theMaxNbFaces) {
2419           anArray = new TColStd_HArray1OfInteger(1,1);
2420           anArray->SetValue(1, anIndices.FindIndex(aSolid));
2421           anObj = GetEngine()->AddSubShape(theCompound, anArray);
2422           aBlocks->Append(anObj);
2423
2424           //Make a Python command
2425           TDF_Tool::Entry(anObj->GetEntry(), anEntry);
2426           anAsciiList += anEntry + ", ";
2427         }
2428       }
2429     }
2430   }
2431   catch (Standard_Failure) {
2432     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2433     SetErrorCode(aFail->GetMessageString());
2434     return aBlocks;
2435   }
2436
2437   if (aBlocks->IsEmpty()) {
2438     SetErrorCode("There are no specified blocks in the given shape");
2439     return aBlocks;
2440   }
2441
2442   anAsciiList.Trunc(anAsciiList.Length() - 2);
2443
2444   //The explode doesn't change object so no new function is required.
2445   aFunction = theCompound->GetLastFunction();
2446   TCollection_AsciiString anOldDescr = aFunction->GetDescription();
2447
2448   //Make a Python command
2449   GEOM::TPythonDump(aFunction) << anOldDescr.ToCString() << "\n\t["
2450     << anAsciiList.ToCString() << "] = geompy.MakeBlockExplode("
2451       << theCompound << ", " << theMinNbFaces << ", " << theMaxNbFaces << ")";
2452
2453   SetErrorCode(OK);
2454   return aBlocks;
2455 }
2456
2457 //=============================================================================
2458 /*!
2459  *  GetBlockNearPoint
2460  */
2461 //=============================================================================
2462 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::GetBlockNearPoint
2463                                                 (Handle(GEOM_Object) theCompound,
2464                                                  Handle(GEOM_Object) thePoint)
2465 {
2466   SetErrorCode(KO);
2467
2468   //New object
2469   Handle(GEOM_Object) aResult;
2470
2471   // Arguments
2472   if (theCompound.IsNull() || thePoint.IsNull()) return NULL;
2473
2474   TopoDS_Shape aBlockOrComp = theCompound->GetValue();
2475   if (aBlockOrComp.IsNull()) {
2476     SetErrorCode("Compound is null");
2477     return NULL;
2478   }
2479   if (aBlockOrComp.ShapeType() != TopAbs_COMPOUND &&
2480       aBlockOrComp.ShapeType() != TopAbs_COMPSOLID) {
2481     SetErrorCode("Shape to find block in is not a compound");
2482     return NULL;
2483   }
2484
2485   TopoDS_Shape anArg = thePoint->GetValue();
2486   if (anArg.IsNull()) {
2487     SetErrorCode("Point is null");
2488     return NULL;
2489   }
2490   if (anArg.ShapeType() != TopAbs_VERTEX) {
2491     SetErrorCode("Shape for block identification is not a vertex");
2492     return NULL;
2493   }
2494
2495   //Compute the Block value
2496   try {
2497     TopoDS_Shape aShape;
2498
2499     TopoDS_Vertex aVert = TopoDS::Vertex(anArg);
2500     gp_Pnt aPnt = BRep_Tool::Pnt(aVert);
2501     Standard_Real PX, PY, PZ;
2502     aPnt.Coord(PX, PY, PZ);
2503
2504     // 1. Classify the point relatively each block
2505     Standard_Integer nearest = 2, nbFound = 0;
2506     TopTools_DataMapOfShapeInteger mapShapeDist;
2507     TopExp_Explorer exp (aBlockOrComp, TopAbs_SOLID);
2508     for (; exp.More(); exp.Next()) {
2509       TopoDS_Shape aSolid = exp.Current();
2510
2511       if (!mapShapeDist.IsBound(aSolid)) {
2512         Standard_Integer aDistance = 2;
2513
2514         // 1.a. Classify relatively Bounding box
2515         Standard_Real Xmin, Ymin, Zmin, Xmax, Ymax, Zmax;
2516         Bnd_Box BB;
2517         BRepBndLib::Add(aSolid, BB);
2518         BB.Get(Xmin, Ymin, Zmin, Xmax, Ymax, Zmax);
2519         if (PX < Xmin || Xmax < PX ||
2520             PY < Ymin || Ymax < PY ||
2521             PZ < Zmin || Zmax < PZ) {
2522           // OUT of bounding box
2523           aDistance = 1;
2524         } else {
2525           // 1.b. Classify relatively the solid itself
2526           BRepClass3d_SolidClassifier SC (aSolid, aPnt, Precision::Confusion());
2527           if (SC.State() == TopAbs_IN) {
2528             aDistance = -1;
2529           } else if (SC.State() == TopAbs_ON) {
2530             aDistance = 0;
2531           } else { // OUT
2532             aDistance = 1;
2533           }
2534         }
2535
2536         if (aDistance < nearest) {
2537           nearest = aDistance;
2538           aShape = aSolid;
2539           nbFound = 1;
2540
2541           // A first found block, containing the point inside, will be returned.
2542           // It is the solution, if there are no intersecting blocks in the compound.
2543           if (nearest == -1) break;
2544
2545         } else if (aDistance == nearest) {
2546           nbFound++;
2547         } else {
2548         }
2549
2550         mapShapeDist.Bind(aSolid, aDistance);
2551       } // if (!mapShapeDist.IsBound(aSolid))
2552     }
2553
2554     // 2. Define block, containing the point or having minimum distance to it
2555     if (nbFound > 1) {
2556       if (nearest == 0) {
2557         // The point is on boundary of some blocks and there are
2558         // no blocks, having the point inside their volume
2559         SetErrorCode("Multiple blocks near the given point are found");
2560         return NULL;
2561
2562       } else if (nearest == 1) {
2563         // The point is outside some blocks and there are
2564         // no blocks, having the point inside or on boundary.
2565         // We will get a nearest block
2566         Standard_Real minDist = RealLast();
2567         TopTools_DataMapIteratorOfDataMapOfShapeInteger mapShapeDistIter (mapShapeDist);
2568         for (; mapShapeDistIter.More(); mapShapeDistIter.Next()) {
2569           if (mapShapeDistIter.Value() == 1) {
2570             TopoDS_Shape aSolid = mapShapeDistIter.Key();
2571             BRepExtrema_DistShapeShape aDistTool (aVert, aSolid);
2572             if (!aDistTool.IsDone()) {
2573               SetErrorCode("Can not find a distance from the given point to one of blocks");
2574               return NULL;
2575             }
2576             Standard_Real aDist = aDistTool.Value();
2577             if (aDist < minDist) {
2578               minDist = aDist;
2579               aShape = aSolid;
2580             }
2581           }
2582         }
2583       } else { // nearest == -1
2584 //        // The point is inside some blocks.
2585 //        // We will get a block with nearest center
2586 //        Standard_Real minDist = RealLast();
2587 //        TopTools_DataMapIteratorOfDataMapOfShapeInteger mapShapeDistIter (mapShapeDist);
2588 //        for (; mapShapeDistIter.More(); mapShapeDistIter.Next()) {
2589 //          if (mapShapeDistIter.Value() == -1) {
2590 //            TopoDS_Shape aSolid = mapShapeDistIter.Key();
2591 //            GProp_GProps aSystem;
2592 //            BRepGProp::VolumeProperties(aSolid, aSystem);
2593 //            gp_Pnt aCenterMass = aSystem.CentreOfMass();
2594 //
2595 //            Standard_Real aDist = aCenterMass.Distance(aPnt);
2596 //            if (aDist < minDist) {
2597 //              minDist = aDist;
2598 //              aShape = aSolid;
2599 //            }
2600 //          }
2601 //        }
2602       }
2603     } // if (nbFound > 1)
2604
2605     if (nbFound == 0) {
2606       SetErrorCode("There are no blocks near the given point");
2607       return NULL;
2608     } else {
2609       TopTools_IndexedMapOfShape anIndices;
2610       TopExp::MapShapes(aBlockOrComp, anIndices);
2611       Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
2612       anArray->SetValue(1, anIndices.FindIndex(aShape));
2613       aResult = GetEngine()->AddSubShape(theCompound, anArray);
2614     }
2615   }
2616   catch (Standard_Failure) {
2617     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2618     SetErrorCode(aFail->GetMessageString());
2619     return NULL;
2620   }
2621
2622   Handle(GEOM_Function) aFunction = aResult->GetLastFunction();
2623
2624   //Make a Python command
2625   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetBlockNearPoint("
2626                                << theCompound << ", " << thePoint << ")";
2627
2628   SetErrorCode(OK);
2629   return aResult;
2630 }
2631
2632 //=============================================================================
2633 /*!
2634  *  GetBlockByParts
2635  */
2636 //=============================================================================
2637 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::GetBlockByParts
2638                       (Handle(GEOM_Object)                         theCompound,
2639                        const Handle(TColStd_HSequenceOfTransient)& theParts)
2640 {
2641   SetErrorCode(KO);
2642
2643   Handle(GEOM_Object) aResult;
2644
2645   if (theCompound.IsNull() || theParts.IsNull()) return NULL;
2646   TopoDS_Shape aBlockOrComp = theCompound->GetValue();
2647   if (aBlockOrComp.IsNull()) return NULL;
2648
2649   //Get the parts
2650   Standard_Integer argi, aLen = theParts->Length();
2651   TopTools_Array1OfShape anArgs (1, aLen);
2652   TCollection_AsciiString anEntry, aPartsDescr;
2653   for (argi = 1; argi <= aLen; argi++) {
2654     Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast(theParts->Value(argi));
2655     Handle(GEOM_Function) aRef = anObj->GetLastFunction();
2656     if (aRef.IsNull()) return NULL;
2657
2658     TopoDS_Shape anArg = aRef->GetValue();
2659     if (anArg.IsNull()) {
2660       SetErrorCode("Null shape is given as argument");
2661       return NULL;
2662     }
2663     anArgs(argi) = anArg;
2664
2665     // For Python command
2666     TDF_Tool::Entry(anObj->GetEntry(), anEntry);
2667     if (argi > 1) aPartsDescr += ", ";
2668     aPartsDescr += anEntry;
2669   }
2670
2671   //Compute the Block value
2672   try {
2673     // 1. Explode compound on solids
2674     TopTools_MapOfShape mapShape;
2675     Standard_Integer nbSolids = 0;
2676     TopExp_Explorer exp (aBlockOrComp, TopAbs_SOLID);
2677     for (; exp.More(); exp.Next()) {
2678       if (mapShape.Add(exp.Current())) {
2679         nbSolids++;
2680       }
2681     }
2682
2683     mapShape.Clear();
2684     Standard_Integer ind = 1;
2685     TopTools_Array1OfShape aSolids (1, nbSolids);
2686     TColStd_Array1OfInteger aNbParts (1, nbSolids);
2687     for (exp.Init(aBlockOrComp, TopAbs_SOLID); exp.More(); exp.Next(), ind++) {
2688       if (mapShape.Add(exp.Current())) {
2689         TopoDS_Shape aSolid = exp.Current();
2690         aSolids(ind) = aSolid;
2691         aNbParts(ind) = 0;
2692
2693         // 2. Define quantity of parts, contained in each solid
2694         TopTools_IndexedMapOfShape aSubShapes;
2695         TopExp::MapShapes(aSolid, aSubShapes);
2696         for (argi = 1; argi <= aLen; argi++) {
2697           if (aSubShapes.Contains(anArgs(argi))) {
2698             aNbParts(ind)++;
2699           }
2700         }
2701       }
2702     }
2703
2704     // 3. Define solid, containing maximum quantity of parts
2705     Standard_Integer maxNb = 0, nbFound = 0;
2706     TopoDS_Shape aShape;
2707     for (ind = 1; ind <= nbSolids; ind++) {
2708       if (aNbParts(ind) > maxNb) {
2709         maxNb = aNbParts(ind);
2710         aShape = aSolids(ind);
2711         nbFound = 1;
2712       } else if (aNbParts(ind) == maxNb) {
2713         nbFound++;
2714       } else {
2715       }
2716     }
2717     if (nbFound > 1) {
2718       SetErrorCode("Multiple blocks, containing maximum quantity of the given parts, are found");
2719       return NULL;
2720     } else if (nbFound == 0) {
2721       SetErrorCode("There are no blocks, containing the given parts");
2722       return NULL;
2723     } else {
2724       TopTools_IndexedMapOfShape anIndices;
2725       TopExp::MapShapes(aBlockOrComp, anIndices);
2726       Handle(TColStd_HArray1OfInteger) anArray = new TColStd_HArray1OfInteger(1,1);
2727       anArray->SetValue(1, anIndices.FindIndex(aShape));
2728       aResult = GetEngine()->AddSubShape(theCompound, anArray);
2729     }
2730   } catch (Standard_Failure) {
2731     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2732     SetErrorCode(aFail->GetMessageString());
2733     return NULL;
2734   }
2735
2736   Handle(GEOM_Function) aFunction = aResult->GetLastFunction();
2737
2738   //Make a Python command
2739   GEOM::TPythonDump(aFunction) << aResult << " = geompy.GetBlockByParts("
2740     << theCompound << ", [" << aPartsDescr.ToCString() << "])";
2741
2742   SetErrorCode(OK);
2743   return aResult;
2744 }
2745
2746 //=============================================================================
2747 /*!
2748  *  GetBlocksByParts
2749  */
2750 //=============================================================================
2751 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IBlocksOperations::GetBlocksByParts
2752                       (Handle(GEOM_Object)                         theCompound,
2753                        const Handle(TColStd_HSequenceOfTransient)& theParts)
2754 {
2755   SetErrorCode(KO);
2756
2757   if (theCompound.IsNull() || theParts.IsNull()) return NULL;
2758   TopoDS_Shape aBlockOrComp = theCompound->GetValue();
2759   if (aBlockOrComp.IsNull()) return NULL;
2760
2761   Handle(TColStd_HSequenceOfTransient) aBlocks = new TColStd_HSequenceOfTransient;
2762   Handle(GEOM_Object) anObj;
2763   Handle(GEOM_Function) aFunction;
2764
2765   //Get the parts
2766   Standard_Integer argi, aLen = theParts->Length();
2767   TopTools_Array1OfShape anArgs (1, aLen);
2768   TCollection_AsciiString anEntry, aPartsDescr, anAsciiList;
2769
2770   for (argi = 1; argi <= aLen; argi++) {
2771     Handle(GEOM_Object) anObj = Handle(GEOM_Object)::DownCast(theParts->Value(argi));
2772     Handle(GEOM_Function) aRef = anObj->GetLastFunction();
2773     if (aRef.IsNull()) return NULL;
2774
2775     TopoDS_Shape anArg = aRef->GetValue();
2776     if (anArg.IsNull()) {
2777       SetErrorCode("Null shape is given as argument");
2778       return NULL;
2779     }
2780     anArgs(argi) = anArg;
2781
2782     // For Python command
2783     TDF_Tool::Entry(anObj->GetEntry(), anEntry);
2784     aPartsDescr += anEntry + ", ";
2785   }
2786
2787   //Get the Blocks
2788   try {
2789     TopTools_MapOfShape mapShape;
2790     Standard_Integer nbSolids = 0;
2791     TopExp_Explorer exp (aBlockOrComp, TopAbs_SOLID);
2792     for (; exp.More(); exp.Next()) {
2793       if (mapShape.Add(exp.Current())) {
2794         nbSolids++;
2795       }
2796     }
2797
2798     mapShape.Clear();
2799     Standard_Integer ind = 1;
2800     TopTools_Array1OfShape aSolids (1, nbSolids);
2801     TColStd_Array1OfInteger aNbParts (1, nbSolids);
2802     for (exp.Init(aBlockOrComp, TopAbs_SOLID); exp.More(); exp.Next(), ind++) {
2803       if (mapShape.Add(exp.Current())) {
2804         TopoDS_Shape aSolid = exp.Current();
2805         aSolids(ind) = aSolid;
2806         aNbParts(ind) = 0;
2807
2808         // 2. Define quantity of parts, contained in each solid
2809         TopTools_IndexedMapOfShape aSubShapes;
2810         TopExp::MapShapes(aSolid, aSubShapes);
2811         for (argi = 1; argi <= aLen; argi++) {
2812           if (aSubShapes.Contains(anArgs(argi))) {
2813             aNbParts(ind)++;
2814           }
2815         }
2816       }
2817     }
2818
2819     // 3. Define solid, containing maximum quantity of parts
2820     Standard_Integer maxNb = 0, nbFound = 0;
2821     for (ind = 1; ind <= nbSolids; ind++) {
2822       if (aNbParts(ind) > maxNb) {
2823         maxNb = aNbParts(ind);
2824         nbFound = 1;
2825       } else if (aNbParts(ind) == maxNb) {
2826         nbFound++;
2827       } else {
2828       }
2829     }
2830     if (nbFound == 0) {
2831       SetErrorCode("There are no blocks, containing the given parts");
2832       return NULL;
2833     }
2834
2835     // Map shapes
2836     TopTools_IndexedMapOfShape anIndices;
2837     TopExp::MapShapes(aBlockOrComp, anIndices);
2838     Handle(TColStd_HArray1OfInteger) anArray;
2839
2840     for (ind = 1; ind <= nbSolids; ind++) {
2841       if (aNbParts(ind) == maxNb) {
2842         anArray = new TColStd_HArray1OfInteger(1,1);
2843         anArray->SetValue(1, anIndices.FindIndex(aSolids(ind)));
2844         anObj = GetEngine()->AddSubShape(theCompound, anArray);
2845         aBlocks->Append(anObj);
2846
2847         // For Python command
2848         TDF_Tool::Entry(anObj->GetEntry(), anEntry);
2849         anAsciiList += anEntry + ", ";
2850         if (aFunction.IsNull())
2851           aFunction = anObj->GetLastFunction();
2852       }
2853     }
2854   }
2855   catch (Standard_Failure) {
2856     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2857     SetErrorCode(aFail->GetMessageString());
2858     return NULL;
2859   }
2860
2861   //Make a Python command
2862   aPartsDescr.Trunc(aPartsDescr.Length() - 2);
2863   anAsciiList.Trunc(anAsciiList.Length() - 2);
2864
2865   GEOM::TPythonDump(aFunction) << "[" << anAsciiList.ToCString()
2866     << "] = geompy.GetBlocksByParts(" << theCompound
2867       << ", [" << aPartsDescr.ToCString() << "])";
2868
2869   SetErrorCode(OK);
2870   return aBlocks;
2871 }
2872
2873 //=============================================================================
2874 /*!
2875  *  MakeMultiTransformation1D
2876  */
2877 //=============================================================================
2878 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::MakeMultiTransformation1D
2879                                              (Handle(GEOM_Object)    theObject,
2880                                               const Standard_Integer theDirFace1,
2881                                               const Standard_Integer theDirFace2,
2882                                               const Standard_Integer theNbTimes)
2883 {
2884   SetErrorCode(KO);
2885
2886   if (theObject.IsNull()) return NULL;
2887
2888   Handle(GEOM_Function) aLastFunction = theObject->GetLastFunction();
2889   if (aLastFunction.IsNull()) return NULL; //There is no function which creates an object to be moved
2890
2891   //Add a new Copy object
2892   Handle(GEOM_Object) aCopy = GetEngine()->AddObject(GetDocID(), GEOM_COPY);
2893
2894   //Add a translate function
2895   Handle(GEOM_Function) aFunction =
2896     aCopy->AddFunction(GEOMImpl_BlockDriver::GetID(), BLOCK_MULTI_TRANSFORM_1D);
2897
2898   //Check if the function is set correctly
2899   if (aFunction->GetDriverGUID() != GEOMImpl_BlockDriver::GetID()) return NULL;
2900
2901   GEOMImpl_IBlockTrsf aTI (aFunction);
2902   aTI.SetOriginal(aLastFunction);
2903   aTI.SetFace1U(theDirFace1);
2904   aTI.SetFace2U(theDirFace2);
2905   aTI.SetNbIterU(theNbTimes);
2906
2907   //Compute the transformation
2908   try {
2909     if (!GetSolver()->ComputeFunction(aFunction)) {
2910       SetErrorCode("Block driver failed to make multi-transformation");
2911       return NULL;
2912     }
2913   }
2914   catch (Standard_Failure) {
2915     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2916     SetErrorCode(aFail->GetMessageString());
2917     return NULL;
2918   }
2919
2920   //Make a Python command
2921   GEOM::TPythonDump(aFunction) << aCopy << " = geompy.MakeMultiTransformation1D("
2922     << theObject << ", " << theDirFace1 << ", " << theDirFace2 << ", " << theNbTimes << ")";
2923
2924   SetErrorCode(OK);
2925   return aCopy;
2926 }
2927
2928 //=============================================================================
2929 /*!
2930  *  MakeMultiTransformation2D
2931  */
2932 //=============================================================================
2933 Handle(GEOM_Object) GEOMImpl_IBlocksOperations::MakeMultiTransformation2D
2934                                              (Handle(GEOM_Object)    theObject,
2935                                               const Standard_Integer theDirFace1U,
2936                                               const Standard_Integer theDirFace2U,
2937                                               const Standard_Integer theNbTimesU,
2938                                               const Standard_Integer theDirFace1V,
2939                                               const Standard_Integer theDirFace2V,
2940                                               const Standard_Integer theNbTimesV)
2941 {
2942   SetErrorCode(KO);
2943
2944   if (theObject.IsNull()) return NULL;
2945
2946   Handle(GEOM_Function) aLastFunction = theObject->GetLastFunction();
2947   if (aLastFunction.IsNull()) return NULL; //There is no function which creates an object to be moved
2948
2949   //Add a new Copy object
2950   Handle(GEOM_Object) aCopy = GetEngine()->AddObject(GetDocID(), GEOM_COPY);
2951
2952   //Add a translate function
2953   Handle(GEOM_Function) aFunction =
2954     aCopy->AddFunction(GEOMImpl_BlockDriver::GetID(), BLOCK_MULTI_TRANSFORM_2D);
2955
2956   //Check if the function is set correctly
2957   if (aFunction->GetDriverGUID() != GEOMImpl_BlockDriver::GetID()) return NULL;
2958
2959   GEOMImpl_IBlockTrsf aTI (aFunction);
2960   aTI.SetOriginal(aLastFunction);
2961   aTI.SetFace1U(theDirFace1U);
2962   aTI.SetFace2U(theDirFace2U);
2963   aTI.SetNbIterU(theNbTimesU);
2964   aTI.SetFace1V(theDirFace1V);
2965   aTI.SetFace2V(theDirFace2V);
2966   aTI.SetNbIterV(theNbTimesV);
2967
2968   //Compute the transformation
2969   try {
2970     if (!GetSolver()->ComputeFunction(aFunction)) {
2971       SetErrorCode("Block driver failed to make multi-transformation");
2972       return NULL;
2973     }
2974   }
2975   catch (Standard_Failure) {
2976     Handle(Standard_Failure) aFail = Standard_Failure::Caught();
2977     SetErrorCode(aFail->GetMessageString());
2978     return NULL;
2979   }
2980
2981   //Make a Python command
2982   GEOM::TPythonDump(aFunction) << aCopy << " = geompy.MakeMultiTransformation2D("
2983     << theObject << ", " << theDirFace1U << ", " << theDirFace2U << ", " << theNbTimesU
2984       << ", " << theDirFace1V << ", " << theDirFace2V << ", " << theNbTimesV << ")";
2985
2986   SetErrorCode(OK);
2987   return aCopy;
2988 }
2989
2990 //=============================================================================
2991 /*!
2992  *  Propagate
2993  */
2994 //=============================================================================
2995 Handle(TColStd_HSequenceOfTransient) GEOMImpl_IBlocksOperations::Propagate
2996                                                  (Handle(GEOM_Object) theShape)
2997 {
2998   SetErrorCode(KO);
2999
3000   if (theShape.IsNull()) return NULL;
3001
3002   TopoDS_Shape aShape = theShape->GetValue();
3003   if (aShape.IsNull()) return NULL;
3004
3005   TopTools_IndexedMapOfShape anIndices;
3006   TopExp::MapShapes(aShape, anIndices);
3007
3008   TopTools_IndexedDataMapOfShapeListOfShape MEW;
3009   GEOMImpl_Block6Explorer::MapShapesAndAncestors
3010     (aShape, TopAbs_EDGE, TopAbs_WIRE, MEW);
3011   Standard_Integer ie, nbEdges = MEW.Extent();
3012
3013   // Result
3014   Handle(TColStd_HSequenceOfTransient) aSeq = new TColStd_HSequenceOfTransient;
3015
3016   TopTools_MapOfShape mapAcceptedEdges;
3017   TCollection_AsciiString aListRes, anEntry;
3018
3019   for (ie = 1; ie <= nbEdges; ie++) {
3020     TopoDS_Shape curE = MEW.FindKey(ie);
3021
3022     if (mapAcceptedEdges.Contains(curE)) continue;
3023
3024     // Build the chain
3025     TopTools_ListOfShape currentChain;
3026     TopTools_ListOfShape listPrevEdges;
3027
3028     currentChain.Append(curE);
3029     listPrevEdges.Append(curE);
3030     mapAcceptedEdges.Add(curE);
3031
3032     // Collect all edges pass by pass
3033     while (listPrevEdges.Extent() > 0) {
3034       // List of edges, added to chain on this cycle pass
3035       TopTools_ListOfShape listCurEdges;
3036
3037       // Find the next portion of edges
3038       TopTools_ListIteratorOfListOfShape itE (listPrevEdges);
3039       for (; itE.More(); itE.Next()) {
3040         TopoDS_Shape anE = itE.Value();
3041
3042         // Iterate on faces, having edge <anE>
3043         TopTools_ListIteratorOfListOfShape itW (MEW.FindFromKey(anE));
3044         for (; itW.More(); itW.Next()) {
3045           TopoDS_Shape aW = itW.Value();
3046           TopoDS_Shape anOppE;
3047
3048           BRepTools_WireExplorer aWE (TopoDS::Wire(aW));
3049           Standard_Integer nb = 1, found = 0;
3050           TopTools_Array1OfShape anEdges (1,4);
3051           for (; aWE.More(); aWE.Next(), nb++) {
3052             if (nb > 4) {
3053               found = 0;
3054               break;
3055             }
3056             anEdges(nb) = aWE.Current();
3057             if (anEdges(nb).IsSame(anE)) found = nb;
3058           }
3059
3060           if (nb == 5 && found > 0) {
3061             // Quadrangle face found, get an opposite edge
3062             Standard_Integer opp = found + 2;
3063             if (opp > 4) opp -= 4;
3064             anOppE = anEdges(opp);
3065
3066             if (!mapAcceptedEdges.Contains(anOppE)) {
3067               // Add found edge to the chain
3068               currentChain.Append(anOppE);
3069               listCurEdges.Append(anOppE);
3070               mapAcceptedEdges.Add(anOppE);
3071             }
3072           } // if (nb == 5 && found > 0)
3073         } // for (; itF.More(); itF.Next())
3074       } // for (; itE.More(); itE.Next())
3075
3076       listPrevEdges = listCurEdges;
3077     } // while (listPrevEdges.Extent() > 0)
3078
3079     // Store the chain in the document
3080     Handle(TColStd_HArray1OfInteger) anArray =
3081       new TColStd_HArray1OfInteger (1, currentChain.Extent());
3082
3083     // Fill array of sub-shape indices
3084     TopTools_ListIteratorOfListOfShape itSub (currentChain);
3085     for (int index = 1; itSub.More(); itSub.Next(), ++index) {
3086       int id = anIndices.FindIndex(itSub.Value());
3087       anArray->SetValue(index, id);
3088     }
3089
3090     // Add a new group object
3091     Handle(GEOM_Object) aChain = GetEngine()->AddSubShape(theShape, anArray);
3092
3093     // Set a GROUP type
3094     aChain->SetType(GEOM_GROUP);
3095
3096     // Set a sub shape type
3097     TDF_Label aFreeLabel = aChain->GetFreeLabel();
3098     TDataStd_Integer::Set(aFreeLabel, (Standard_Integer)TopAbs_EDGE);
3099
3100     // Add the chain to the result
3101     aSeq->Append(aChain);
3102
3103     //Make a Python command
3104     TDF_Tool::Entry(aChain->GetEntry(), anEntry);
3105     aListRes += anEntry + ", ";
3106   }
3107
3108   if (aSeq->IsEmpty()) {
3109     SetErrorCode("There are no quadrangle faces in the shape");
3110     return aSeq;
3111   }
3112
3113   aListRes.Trunc(aListRes.Length() - 2);
3114
3115   // The Propagation doesn't change object so no new function is required.
3116   Handle(GEOM_Function) aFunction = theShape->GetLastFunction();
3117   TCollection_AsciiString anOldDescr = aFunction->GetDescription();
3118
3119   // Make a Python command
3120   GEOM::TPythonDump(aFunction) << anOldDescr.ToCString() << "\n\t["
3121     << aListRes.ToCString() << "] = geompy.Propagate(" << theShape << ")";
3122
3123   SetErrorCode(OK);
3124   return aSeq;
3125 }