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