1 2 #ifdef OCCGEOMETRY 3 4 #include <mystdlib.h> 5 #include <occgeom.hpp> 6 #include <core/register_archive.hpp> 7 #include <cstdio> 8 #include "ShapeAnalysis_ShapeTolerance.hxx" 9 #include "ShapeAnalysis_ShapeContents.hxx" 10 #include "ShapeAnalysis_CheckSmallFace.hxx" 11 #include "ShapeAnalysis_DataMapOfShapeListOfReal.hxx" 12 #include "ShapeAnalysis_Surface.hxx" 13 14 #include "BRepCheck_Analyzer.hxx" 15 #include "BRepLib.hxx" 16 #include "ShapeBuild_ReShape.hxx" 17 #include "ShapeFix.hxx" 18 #include "ShapeFix_FixSmallFace.hxx" 19 #include "Partition_Spliter.hxx" 20 #include "BRepAlgoAPI_Fuse.hxx" 21 #include "Interface_InterfaceModel.hxx" 22 23 #include "XSControl_WorkSession.hxx" 24 #include "XSControl_TransferReader.hxx" 25 #include "StepRepr_RepresentationItem.hxx" 26 #include "StepBasic_ProductDefinitionRelationship.hxx" 27 #include "Transfer_TransientProcess.hxx" 28 #include "TransferBRep.hxx" 29 #ifndef _Standard_Version_HeaderFile 30 #include <Standard_Version.hxx> 31 #endif 32 33 #if OCC_VERSION_HEX < 0x070000 34 // pass 35 #elif OCC_VERSION_HEX < 0x070200 36 #include "StlTransfer.hxx" 37 #include "TopoDS_Iterator.hxx" 38 #else 39 #include "TopoDS_Iterator.hxx" 40 #endif 41 42 namespace netgen 43 { 44 45 // std::map<Handle(TopoDS_TShape), string> OCCGeometry::global_shape_names; 46 // std::map<Handle(TopoDS_TShape), Vec<3>> OCCGeometry::global_shape_cols; 47 std::map<Handle(TopoDS_TShape), ShapeProperties> OCCGeometry::global_shape_properties; 48 std::map<Handle(TopoDS_TShape), std::vector<OCCIdentification>> OCCGeometry::identifications; 49 OCCGeometry(const TopoDS_Shape & _shape,int aoccdim)50 OCCGeometry::OCCGeometry(const TopoDS_Shape& _shape, int aoccdim) 51 { 52 shape = _shape; 53 occdim = aoccdim; 54 changed = true; 55 BuildFMap(); 56 CalcBoundingBox(); 57 58 snames.SetSize(somap.Size()); 59 for(auto i : Range(snames)) 60 snames[i] = "domain_" + ToString(i+1); 61 for (TopExp_Explorer e(shape, TopAbs_SOLID); e.More(); e.Next()) 62 { 63 TopoDS_Solid solid = TopoDS::Solid(e.Current()); 64 if (auto name = global_shape_properties[solid.TShape()].name) 65 snames[somap.FindIndex(solid)-1] = *name; 66 } 67 68 fnames.SetSize(fmap.Size()); 69 for(auto i : Range(fnames)) 70 fnames[i] = "bc_" + ToString(i+1); 71 for(TopExp_Explorer e(shape, TopAbs_FACE); e.More(); e.Next()) 72 { 73 TopoDS_Face face = TopoDS::Face(e.Current()); 74 if (auto name = global_shape_properties[face.TShape()].name) 75 fnames[fmap.FindIndex(face)-1] = *name; 76 } 77 enames.SetSize(emap.Size()); 78 enames = ""; 79 for (TopExp_Explorer e(shape, TopAbs_EDGE); e.More(); e.Next()) 80 { 81 TopoDS_Edge edge = TopoDS::Edge(e.Current()); 82 if (auto name = global_shape_properties[edge.TShape()].name) 83 enames[emap.FindIndex(edge)-1] = *name; 84 } 85 } 86 STEP_GetEntityName(const TopoDS_Shape & theShape,STEPCAFControl_Reader * aReader)87 string STEP_GetEntityName(const TopoDS_Shape & theShape, STEPCAFControl_Reader * aReader) 88 { 89 const Handle(XSControl_WorkSession)& theSession = aReader->Reader().WS(); 90 const Handle(XSControl_TransferReader)& aTransferReader = 91 theSession->TransferReader(); 92 93 Handle(Standard_Transient) anEntity = 94 aTransferReader->EntityFromShapeResult(theShape, 1); 95 96 if (anEntity.IsNull()) // as just mapped 97 anEntity = aTransferReader->EntityFromShapeResult (theShape,-1); 98 99 if (anEntity.IsNull()) // as anything 100 anEntity = aTransferReader->EntityFromShapeResult (theShape,4); 101 102 if (anEntity.IsNull()) 103 { 104 cout<<"Warning: cannot get entity from shape" <<endl; 105 return "none"; 106 } 107 108 auto aReprItem = Handle(StepRepr_RepresentationItem)::DownCast(anEntity); 109 if(!aReprItem.IsNull()) 110 return aReprItem->Name()->ToCString();; 111 112 auto bReprItem = Handle(StepBasic_ProductDefinitionRelationship)::DownCast(anEntity); 113 if (!bReprItem.IsNull()) 114 return bReprItem->Description()->ToCString(); 115 116 cout<<"Warning: unknown entity type " << anEntity->DynamicType() << endl; 117 return "none"; 118 } 119 Analyse(Mesh & mesh,const MeshingParameters & mparam) const120 void OCCGeometry :: Analyse(Mesh& mesh, 121 const MeshingParameters& mparam) const 122 { 123 OCCSetLocalMeshSize(*this, mesh, mparam, occparam); 124 } 125 FindEdges(Mesh & mesh,const MeshingParameters & mparam) const126 void OCCGeometry :: FindEdges(Mesh& mesh, 127 const MeshingParameters& mparam) const 128 { 129 OCCFindEdges(*this, mesh, mparam); 130 } 131 MeshSurface(Mesh & mesh,const MeshingParameters & mparam) const132 void OCCGeometry :: MeshSurface(Mesh& mesh, 133 const MeshingParameters& mparam) const 134 { 135 OCCMeshSurface(*this, mesh, mparam); 136 } 137 FinalizeMesh(Mesh & mesh) const138 void OCCGeometry :: FinalizeMesh(Mesh& mesh) const 139 { 140 for (int i = 0; i < mesh.GetNDomains(); i++) 141 if (snames.Size()) 142 mesh.SetMaterial (i+1, snames[i]); 143 } 144 PrintNrShapes()145 void OCCGeometry :: PrintNrShapes () 146 { 147 TopExp_Explorer e; 148 int count = 0; 149 for (e.Init(shape, TopAbs_COMPSOLID); e.More(); e.Next()) count++; 150 cout << "CompSolids: " << count << endl; 151 152 cout << "Solids : " << somap.Extent() << endl; 153 cout << "Shells : " << shmap.Extent() << endl; 154 cout << "Faces : " << fmap.Extent() << endl; 155 cout << "Edges : " << emap.Extent() << endl; 156 cout << "Vertices : " << vmap.Extent() << endl; 157 } 158 159 160 161 PrintContents(OCCGeometry * geom)162 void PrintContents (OCCGeometry * geom) 163 { 164 ShapeAnalysis_ShapeContents cont; 165 cont.Clear(); 166 cont.Perform(geom->shape); 167 168 (*testout) << "OCC CONTENTS" << endl; 169 (*testout) << "============" << endl; 170 (*testout) << "SOLIDS : " << cont.NbSolids() << endl; 171 (*testout) << "SHELLS : " << cont.NbShells() << endl; 172 (*testout) << "FACES : " << cont.NbFaces() << endl; 173 (*testout) << "WIRES : " << cont.NbWires() << endl; 174 (*testout) << "EDGES : " << cont.NbEdges() << endl; 175 (*testout) << "VERTICES : " << cont.NbVertices() << endl; 176 177 TopExp_Explorer e; 178 int count = 0; 179 for (e.Init(geom->shape, TopAbs_COMPOUND); e.More(); e.Next()) 180 count++; 181 (*testout) << "Compounds: " << count << endl; 182 183 count = 0; 184 for (e.Init(geom->shape, TopAbs_COMPSOLID); e.More(); e.Next()) 185 count++; 186 (*testout) << "CompSolids: " << count << endl; 187 188 (*testout) << endl; 189 190 cout << "Highest entry in topology hierarchy: " << endl; 191 if (count) 192 cout << count << " composite solid(s)" << endl; 193 else 194 if (geom->somap.Extent()) 195 cout << geom->somap.Extent() << " solid(s)" << endl; 196 else 197 if (geom->shmap.Extent()) 198 cout << geom->shmap.Extent() << " shells(s)" << endl; 199 else 200 if (geom->fmap.Extent()) 201 cout << geom->fmap.Extent() << " face(s)" << endl; 202 else 203 if (geom->wmap.Extent()) 204 cout << geom->wmap.Extent() << " wire(s)" << endl; 205 else 206 if (geom->emap.Extent()) 207 cout << geom->emap.Extent() << " edge(s)" << endl; 208 else 209 if (geom->vmap.Extent()) 210 cout << geom->vmap.Extent() << " vertices(s)" << endl; 211 else 212 cout << "no entities" << endl; 213 214 } 215 GlueGeometry()216 void OCCGeometry :: GlueGeometry() 217 { 218 PrintMessage(1, "OCC Glue Geometry"); 219 /* 220 // 221 BRep_Builder builder; 222 TopoDS_Shape my_fuse; 223 int cnt = 0; 224 for (TopExp_Explorer exp_solid(shape, TopAbs_SOLID); exp_solid.More(); exp_solid.Next()) 225 { 226 cout << "cnt = " << cnt << endl; 227 if (cnt == 0) 228 my_fuse = exp_solid.Current(); 229 else 230 // my_fuse = BRepAlgoAPI_Fuse (my_fuse, exp_solid.Current()); 231 my_fuse = QANewModTopOpe_Glue::QANewModTopOpe_Glue(my_fuse, exp_solid.Current()); 232 cnt++; 233 } 234 cout << "remove" << endl; 235 // for (int i = 1; i <= somap.Size(); i++) 236 // builder.Remove (shape, somap(i)); 237 cout << "now add" << endl; 238 // builder.Add (shape, my_fuse); 239 shape = my_fuse; 240 cout << "build fmap" << endl; 241 BuildFMap(); 242 */ 243 244 245 // from 246 // https://www.opencascade.com/doc/occt-7.4.0/overview/html/occt_user_guides__boolean_operations.html 247 BOPAlgo_Builder aBuilder; 248 249 // Setting arguments 250 TopTools_ListOfShape aLSObjects; 251 for (TopExp_Explorer exp_solid(shape, TopAbs_SOLID); exp_solid.More(); exp_solid.Next()) 252 aLSObjects.Append (exp_solid.Current()); 253 aBuilder.SetArguments(aLSObjects); 254 255 // Setting options for GF 256 // Set parallel processing mode (default is false) 257 // Standard_Boolean bRunParallel = Standard_True; 258 // aBuilder.SetRunParallel(bRunParallel); 259 260 // Set Fuzzy value (default is Precision::Confusion()) 261 // Standard_Real aFuzzyValue = 1.e-5; 262 // aBuilder.SetFuzzyValue(aFuzzyValue); 263 264 // Set safe processing mode (default is false) 265 // Standard_Boolean bSafeMode = Standard_True; 266 // aBuilder.SetNonDestructive(bSafeMode); 267 268 // Set Gluing mode for coinciding arguments (default is off) 269 // BOPAlgo_GlueEnum aGlue = BOPAlgo_GlueShift; 270 // aBuilder.SetGlue(aGlue); 271 272 // Disabling/Enabling the check for inverted solids (default is true) 273 // Standard Boolean bCheckInverted = Standard_False; 274 // aBuilder.SetCheckInverted(bCheckInverted); 275 276 // Set OBB usage (default is false) 277 // Standard_Boolean bUseOBB = Standard_True; 278 // aBuilder.SetUseOBB(buseobb); 279 280 // Perform the operation 281 aBuilder.Perform(); 282 // Check for the errors 283 #if OCC_VERSION_HEX >= 0x070200 284 if (aBuilder.HasErrors()) 285 { 286 cout << "builder has errors" << endl; 287 return; 288 } 289 // Check for the warnings 290 if (aBuilder.HasWarnings()) 291 { 292 // treatment of the warnings 293 ; 294 } 295 #endif 296 297 #ifdef OCC_HAVE_HISTORY 298 Handle(BRepTools_History) history = aBuilder.History (); 299 300 for (TopExp_Explorer e(shape, TopAbs_SOLID); e.More(); e.Next()) 301 { 302 if (auto name = OCCGeometry::global_shape_properties[e.Current().TShape()].name) 303 for (auto mods : history->Modified(e.Current())) 304 OCCGeometry::global_shape_properties[mods.TShape()].name = *name; 305 } 306 #endif // OCC_HAVE_HISTORY 307 308 // result of the operation 309 shape = aBuilder.Shape(); 310 BuildFMap(); 311 } 312 HealGeometry()313 void OCCGeometry :: HealGeometry () 314 { 315 int nrc = 0, nrcs = 0, 316 nrso = somap.Extent(), 317 nrsh = shmap.Extent(), 318 nrf = fmap.Extent(), 319 nrw = wmap.Extent(), 320 nre = emap.Extent(), 321 nrv = vmap.Extent(); 322 323 TopExp_Explorer exp0; 324 TopExp_Explorer exp1; 325 326 327 for (exp0.Init(shape, TopAbs_COMPOUND); exp0.More(); exp0.Next()) nrc++; 328 for (exp0.Init(shape, TopAbs_COMPSOLID); exp0.More(); exp0.Next()) nrcs++; 329 330 double surfacecont = 0; 331 332 { 333 Handle(ShapeBuild_ReShape) rebuild = new ShapeBuild_ReShape; 334 rebuild->Apply(shape); 335 for (exp1.Init (shape, TopAbs_EDGE); exp1.More(); exp1.Next()) 336 { 337 TopoDS_Edge edge = TopoDS::Edge(exp1.Current()); 338 if ( BRep_Tool::Degenerated(edge) ) 339 rebuild->Remove(edge); 340 } 341 shape = rebuild->Apply(shape); 342 } 343 344 BuildFMap(); 345 346 347 for (exp0.Init (shape, TopAbs_FACE); exp0.More(); exp0.Next()) 348 { 349 TopoDS_Face face = TopoDS::Face(exp0.Current()); 350 351 GProp_GProps system; 352 BRepGProp::SurfaceProperties(face, system); 353 surfacecont += system.Mass(); 354 } 355 356 357 cout << "Starting geometry healing procedure (tolerance: " << tolerance << ")" << endl 358 << "-----------------------------------" << endl; 359 360 { 361 cout << endl << "- repairing faces" << endl; 362 363 Handle(ShapeFix_Face) sff; 364 Handle(ShapeBuild_ReShape) rebuild = new ShapeBuild_ReShape; 365 rebuild->Apply(shape); 366 367 368 for (exp0.Init (shape, TopAbs_FACE); exp0.More(); exp0.Next()) 369 { 370 // Variable to hold the colour (if there exists one) of 371 // the current face being processed 372 Quantity_Color face_colour; 373 374 TopoDS_Face face = TopoDS::Face (exp0.Current()); 375 376 if(face_colours.IsNull() 377 || (!(face_colours->GetColor(face,XCAFDoc_ColorSurf,face_colour)))) 378 { 379 // Set the default face colour to green (Netgen Standard) 380 // if no colour has been defined for the face 381 face_colour = Quantity_Color(0.0,1.0,0.0,Quantity_TOC_RGB); 382 } 383 384 sff = new ShapeFix_Face (face); 385 sff->FixAddNaturalBoundMode() = Standard_True; 386 sff->FixSmallAreaWireMode() = Standard_True; 387 sff->Perform(); 388 389 if(sff->Status(ShapeExtend_DONE1) || 390 sff->Status(ShapeExtend_DONE2) || 391 sff->Status(ShapeExtend_DONE3) || 392 sff->Status(ShapeExtend_DONE4) || 393 sff->Status(ShapeExtend_DONE5)) 394 { 395 cout << "repaired face " << fmap.FindIndex(face) << " "; 396 if(sff->Status(ShapeExtend_DONE1)) 397 cout << "(some wires are fixed)" <<endl; 398 else if(sff->Status(ShapeExtend_DONE2)) 399 cout << "(orientation of wires fixed)" <<endl; 400 else if(sff->Status(ShapeExtend_DONE3)) 401 cout << "(missing seam added)" <<endl; 402 else if(sff->Status(ShapeExtend_DONE4)) 403 cout << "(small area wire removed)" <<endl; 404 else if(sff->Status(ShapeExtend_DONE5)) 405 cout << "(natural bounds added)" <<endl; 406 TopoDS_Face newface = sff->Face(); 407 408 rebuild->Replace(face, newface); 409 } 410 411 // Set the original colour of the face to the newly created 412 // face (after the healing process) 413 face = TopoDS::Face (exp0.Current()); 414 face_colours->SetColor(face,face_colour,XCAFDoc_ColorSurf); 415 } 416 shape = rebuild->Apply(shape); 417 } 418 419 420 { 421 Handle(ShapeBuild_ReShape) rebuild = new ShapeBuild_ReShape; 422 rebuild->Apply(shape); 423 for (exp1.Init (shape, TopAbs_EDGE); exp1.More(); exp1.Next()) 424 { 425 TopoDS_Edge edge = TopoDS::Edge(exp1.Current()); 426 if ( BRep_Tool::Degenerated(edge) ) 427 rebuild->Remove(edge); 428 } 429 shape = rebuild->Apply(shape); 430 } 431 432 433 if (fixsmalledges) 434 { 435 cout << endl << "- fixing small edges" << endl; 436 437 Handle(ShapeFix_Wire) sfw; 438 Handle(ShapeBuild_ReShape) rebuild = new ShapeBuild_ReShape; 439 rebuild->Apply(shape); 440 441 442 for (exp0.Init (shape, TopAbs_FACE); exp0.More(); exp0.Next()) 443 { 444 TopoDS_Face face = TopoDS::Face(exp0.Current()); 445 446 for (exp1.Init (face, TopAbs_WIRE); exp1.More(); exp1.Next()) 447 { 448 TopoDS_Wire oldwire = TopoDS::Wire(exp1.Current()); 449 sfw = new ShapeFix_Wire (oldwire, face ,tolerance); 450 sfw->ModifyTopologyMode() = Standard_True; 451 452 sfw->ClosedWireMode() = Standard_True; 453 454 bool replace = false; 455 456 replace = sfw->FixReorder() || replace; 457 458 replace = sfw->FixConnected() || replace; 459 460 461 462 if (sfw->FixSmall (Standard_False, tolerance) && ! (sfw->StatusSmall(ShapeExtend_FAIL1) || 463 sfw->StatusSmall(ShapeExtend_FAIL2) || 464 sfw->StatusSmall(ShapeExtend_FAIL3))) 465 { 466 cout << "Fixed small edge in wire " << wmap.FindIndex (oldwire) << endl; 467 replace = true; 468 469 } 470 else if (sfw->StatusSmall(ShapeExtend_FAIL1)) 471 cerr << "Failed to fix small edge in wire " << wmap.FindIndex (oldwire) 472 << ", edge cannot be checked (no 3d curve and no pcurve)" << endl; 473 else if (sfw->StatusSmall(ShapeExtend_FAIL2)) 474 cerr << "Failed to fix small edge in wire " << wmap.FindIndex (oldwire) 475 << ", edge is null-length and has different vertives at begin and end, and lockvtx is True or ModifiyTopologyMode is False" << endl; 476 else if (sfw->StatusSmall(ShapeExtend_FAIL3)) 477 cerr << "Failed to fix small edge in wire " << wmap.FindIndex (oldwire) 478 << ", CheckConnected has failed" << endl; 479 480 replace = sfw->FixEdgeCurves() || replace; 481 482 replace = sfw->FixDegenerated() || replace; 483 484 replace = sfw->FixSelfIntersection() || replace; 485 486 replace = sfw->FixLacking(Standard_True) || replace; 487 488 if(replace) 489 { 490 TopoDS_Wire newwire = sfw->Wire(); 491 rebuild->Replace(oldwire, newwire); 492 } 493 494 //delete sfw; sfw = NULL; 495 496 } 497 } 498 499 shape = rebuild->Apply(shape); 500 501 502 503 { 504 BuildFMap(); 505 Handle(ShapeBuild_ReShape) rebuild = new ShapeBuild_ReShape; 506 rebuild->Apply(shape); 507 508 for (exp1.Init (shape, TopAbs_EDGE); exp1.More(); exp1.Next()) 509 { 510 TopoDS_Edge edge = TopoDS::Edge(exp1.Current()); 511 if (vmap.FindIndex(TopExp::FirstVertex (edge)) == 512 vmap.FindIndex(TopExp::LastVertex (edge))) 513 { 514 GProp_GProps system; 515 BRepGProp::LinearProperties(edge, system); 516 if (system.Mass() < tolerance) 517 { 518 cout << "removing degenerated edge " << emap.FindIndex(edge) 519 << " from vertex " << vmap.FindIndex(TopExp::FirstVertex (edge)) 520 << " to vertex " << vmap.FindIndex(TopExp::LastVertex (edge)) << endl; 521 rebuild->Remove(edge); 522 } 523 } 524 } 525 shape = rebuild->Apply(shape); 526 527 //delete rebuild; rebuild = NULL; 528 } 529 530 531 532 { 533 Handle(ShapeBuild_ReShape) rebuild = new ShapeBuild_ReShape; 534 rebuild->Apply(shape); 535 for (exp1.Init (shape, TopAbs_EDGE); exp1.More(); exp1.Next()) 536 { 537 TopoDS_Edge edge = TopoDS::Edge(exp1.Current()); 538 if ( BRep_Tool::Degenerated(edge) ) 539 rebuild->Remove(edge); 540 } 541 shape = rebuild->Apply(shape); 542 } 543 544 545 546 547 Handle(ShapeFix_Wireframe) sfwf = new ShapeFix_Wireframe; 548 sfwf->SetPrecision(tolerance); 549 sfwf->Load (shape); 550 sfwf->ModeDropSmallEdges() = Standard_True; 551 552 sfwf->SetPrecision(boundingbox.Diam()); 553 554 if (sfwf->FixWireGaps()) 555 { 556 cout << endl << "- fixing wire gaps" << endl; 557 if (sfwf->StatusWireGaps(ShapeExtend_OK)) cout << "no gaps found" << endl; 558 if (sfwf->StatusWireGaps(ShapeExtend_DONE1)) cout << "some 2D gaps fixed" << endl; 559 if (sfwf->StatusWireGaps(ShapeExtend_DONE2)) cout << "some 3D gaps fixed" << endl; 560 if (sfwf->StatusWireGaps(ShapeExtend_FAIL1)) cout << "failed to fix some 2D gaps" << endl; 561 if (sfwf->StatusWireGaps(ShapeExtend_FAIL2)) cout << "failed to fix some 3D gaps" << endl; 562 } 563 564 sfwf->SetPrecision(tolerance); 565 566 567 { 568 for (exp1.Init (shape, TopAbs_EDGE); exp1.More(); exp1.Next()) 569 { 570 TopoDS_Edge edge = TopoDS::Edge(exp1.Current()); 571 if ( BRep_Tool::Degenerated(edge) ) 572 cout << "degenerated edge at position 4" << endl; 573 } 574 } 575 576 577 578 if (sfwf->FixSmallEdges()) 579 { 580 cout << endl << "- fixing wire frames" << endl; 581 if (sfwf->StatusSmallEdges(ShapeExtend_OK)) cout << "no small edges found" << endl; 582 if (sfwf->StatusSmallEdges(ShapeExtend_DONE1)) cout << "some small edges fixed" << endl; 583 if (sfwf->StatusSmallEdges(ShapeExtend_FAIL1)) cout << "failed to fix some small edges" << endl; 584 } 585 586 587 588 shape = sfwf->Shape(); 589 590 //delete sfwf; sfwf = NULL; 591 //delete rebuild; rebuild = NULL; 592 593 } 594 595 596 597 598 599 { 600 for (exp1.Init (shape, TopAbs_EDGE); exp1.More(); exp1.Next()) 601 { 602 TopoDS_Edge edge = TopoDS::Edge(exp1.Current()); 603 if ( BRep_Tool::Degenerated(edge) ) 604 cout << "degenerated edge at position 5" << endl; 605 } 606 } 607 608 609 610 611 if (fixspotstripfaces) 612 { 613 614 cout << endl << "- fixing spot and strip faces" << endl; 615 Handle(ShapeFix_FixSmallFace) sffsm = new ShapeFix_FixSmallFace(); 616 sffsm -> Init (shape); 617 sffsm -> SetPrecision (tolerance); 618 sffsm -> Perform(); 619 620 shape = sffsm -> FixShape(); 621 //delete sffsm; sffsm = NULL; 622 } 623 624 625 { 626 for (exp1.Init (shape, TopAbs_EDGE); exp1.More(); exp1.Next()) 627 { 628 TopoDS_Edge edge = TopoDS::Edge(exp1.Current()); 629 if ( BRep_Tool::Degenerated(edge) ) 630 cout << "degenerated edge at position 6" << endl; 631 } 632 } 633 634 635 636 if (sewfaces) 637 { 638 cout << endl << "- sewing faces" << endl; 639 640 BRepOffsetAPI_Sewing sewedObj(tolerance); 641 642 for (exp0.Init (shape, TopAbs_FACE); exp0.More(); exp0.Next()) 643 { 644 TopoDS_Face face = TopoDS::Face (exp0.Current()); 645 sewedObj.Add (face); 646 } 647 648 sewedObj.Perform(); 649 650 if (!sewedObj.SewedShape().IsNull()) 651 shape = sewedObj.SewedShape(); 652 else 653 cout << " not possible"; 654 } 655 656 657 658 { 659 Handle(ShapeBuild_ReShape) rebuild = new ShapeBuild_ReShape; 660 rebuild->Apply(shape); 661 for (exp1.Init (shape, TopAbs_EDGE); exp1.More(); exp1.Next()) 662 { 663 TopoDS_Edge edge = TopoDS::Edge(exp1.Current()); 664 if ( BRep_Tool::Degenerated(edge) ) 665 rebuild->Remove(edge); 666 } 667 shape = rebuild->Apply(shape); 668 } 669 670 671 if (makesolids) 672 { 673 cout << endl << "- making solids" << endl; 674 675 BRepBuilderAPI_MakeSolid ms; 676 int count = 0; 677 for (exp0.Init(shape, TopAbs_SHELL); exp0.More(); exp0.Next()) 678 { 679 count++; 680 ms.Add (TopoDS::Shell(exp0.Current())); 681 } 682 683 if (!count) 684 { 685 cout << " not possible (no shells)" << endl; 686 } 687 else 688 { 689 BRepCheck_Analyzer ba(ms); 690 if (ba.IsValid ()) 691 { 692 Handle(ShapeFix_Shape) sfs = new ShapeFix_Shape; 693 sfs->Init (ms); 694 sfs->SetPrecision(tolerance); 695 sfs->SetMaxTolerance(tolerance); 696 sfs->Perform(); 697 shape = sfs->Shape(); 698 699 for (exp0.Init(shape, TopAbs_SOLID); exp0.More(); exp0.Next()) 700 { 701 TopoDS_Solid solid = TopoDS::Solid(exp0.Current()); 702 TopoDS_Solid newsolid = solid; 703 BRepLib::OrientClosedSolid (newsolid); 704 Handle(ShapeBuild_ReShape) rebuild = new ShapeBuild_ReShape; 705 // rebuild->Apply(shape); 706 rebuild->Replace(solid, newsolid); 707 TopoDS_Shape newshape = rebuild->Apply(shape, TopAbs_COMPSOLID);//, 1); 708 // TopoDS_Shape newshape = rebuild->Apply(shape); 709 shape = newshape; 710 } 711 712 //delete sfs; sfs = NULL; 713 } 714 else 715 cout << " not possible" << endl; 716 } 717 } 718 719 720 721 if (splitpartitions) 722 { 723 cout << "- running SALOME partition splitter" << endl; 724 725 TopExp_Explorer e2; 726 Partition_Spliter ps; 727 int count = 0; 728 729 for (e2.Init (shape, TopAbs_SOLID); 730 e2.More(); e2.Next()) 731 { 732 count++; 733 ps.AddShape (e2.Current()); 734 } 735 736 ps.Compute(); 737 shape = ps.Shape(); 738 739 cout << " before: " << count << " solids" << endl; 740 741 count = 0; 742 for (e2.Init (shape, TopAbs_SOLID); 743 e2.More(); e2.Next()) count++; 744 745 cout << " after : " << count << " solids" << endl; 746 } 747 748 BuildFMap(); 749 750 751 752 { 753 for (exp1.Init (shape, TopAbs_EDGE); exp1.More(); exp1.Next()) 754 { 755 TopoDS_Edge edge = TopoDS::Edge(exp1.Current()); 756 if ( BRep_Tool::Degenerated(edge) ) 757 cout << "degenerated edge at position 8" << endl; 758 } 759 } 760 761 762 double newsurfacecont = 0; 763 764 765 for (exp0.Init (shape, TopAbs_FACE); exp0.More(); exp0.Next()) 766 { 767 TopoDS_Face face = TopoDS::Face(exp0.Current()); 768 GProp_GProps system; 769 BRepGProp::SurfaceProperties(face, system); 770 newsurfacecont += system.Mass(); 771 } 772 773 774 int nnrc = 0, nnrcs = 0, 775 nnrso = somap.Extent(), 776 nnrsh = shmap.Extent(), 777 nnrf = fmap.Extent(), 778 nnrw = wmap.Extent(), 779 nnre = emap.Extent(), 780 nnrv = vmap.Extent(); 781 782 for (exp0.Init(shape, TopAbs_COMPOUND); exp0.More(); exp0.Next()) nnrc++; 783 for (exp0.Init(shape, TopAbs_COMPSOLID); exp0.More(); exp0.Next()) nnrcs++; 784 785 cout << "-----------------------------------" << endl; 786 cout << "Compounds : " << nnrc << " (" << nrc << ")" << endl; 787 cout << "Composite solids: " << nnrcs << " (" << nrcs << ")" << endl; 788 cout << "Solids : " << nnrso << " (" << nrso << ")" << endl; 789 cout << "Shells : " << nnrsh << " (" << nrsh << ")" << endl; 790 cout << "Wires : " << nnrw << " (" << nrw << ")" << endl; 791 cout << "Faces : " << nnrf << " (" << nrf << ")" << endl; 792 cout << "Edges : " << nnre << " (" << nre << ")" << endl; 793 cout << "Vertices : " << nnrv << " (" << nrv << ")" << endl; 794 cout << endl; 795 cout << "Totol surface area : " << newsurfacecont << " (" << surfacecont << ")" << endl; 796 cout << endl; 797 } 798 799 800 801 BuildFMap()802 void OCCGeometry :: BuildFMap() 803 { 804 somap.Clear(); 805 shmap.Clear(); 806 fmap.Clear(); 807 wmap.Clear(); 808 emap.Clear(); 809 vmap.Clear(); 810 811 TopExp_Explorer exp0, exp1, exp2, exp3, exp4, exp5; 812 813 for (exp0.Init(shape, TopAbs_COMPOUND); 814 exp0.More(); exp0.Next()) 815 { 816 TopoDS_Compound compound = TopoDS::Compound (exp0.Current()); 817 (*testout) << "compound" << endl; 818 int i = 0; 819 for (exp1.Init(compound, TopAbs_SHELL); 820 exp1.More(); exp1.Next()) 821 { 822 (*testout) << "shell " << ++i << endl; 823 } 824 } 825 826 for (exp0.Init(shape, TopAbs_SOLID); 827 exp0.More(); exp0.Next()) 828 { 829 TopoDS_Solid solid = TopoDS::Solid (exp0.Current()); 830 831 if (somap.FindIndex(solid) < 1) 832 { 833 somap.Add (solid); 834 835 for (exp1.Init(solid, TopAbs_SHELL); 836 exp1.More(); exp1.Next()) 837 { 838 TopoDS_Shell shell = TopoDS::Shell (exp1.Current()); 839 if (shmap.FindIndex(shell) < 1) 840 { 841 shmap.Add (shell); 842 843 for (exp2.Init(shell, TopAbs_FACE); 844 exp2.More(); exp2.Next()) 845 { 846 TopoDS_Face face = TopoDS::Face(exp2.Current()); 847 if (fmap.FindIndex(face) < 1) 848 { 849 fmap.Add (face); 850 (*testout) << "face " << fmap.FindIndex(face) << " "; 851 (*testout) << ((face.Orientation() == TopAbs_REVERSED) ? "-" : "+") << ", "; 852 (*testout) << ((exp2.Current().Orientation() == TopAbs_REVERSED) ? "-" : "+") << endl; 853 for (exp3.Init(exp2.Current(), TopAbs_WIRE); 854 exp3.More(); exp3.Next()) 855 { 856 TopoDS_Wire wire = TopoDS::Wire (exp3.Current()); 857 if (wmap.FindIndex(wire) < 1) 858 { 859 wmap.Add (wire); 860 861 for (exp4.Init(exp3.Current(), TopAbs_EDGE); 862 exp4.More(); exp4.Next()) 863 { 864 TopoDS_Edge edge = TopoDS::Edge(exp4.Current()); 865 if (emap.FindIndex(edge) < 1) 866 { 867 emap.Add (edge); 868 for (exp5.Init(exp4.Current(), TopAbs_VERTEX); 869 exp5.More(); exp5.Next()) 870 { 871 TopoDS_Vertex vertex = TopoDS::Vertex(exp5.Current()); 872 if (vmap.FindIndex(vertex) < 1) 873 vmap.Add (vertex); 874 } 875 } 876 } 877 } 878 } 879 } 880 } 881 } 882 } 883 } 884 } 885 886 // Free Shells 887 for (exp1.Init(shape, TopAbs_SHELL, TopAbs_SOLID); exp1.More(); exp1.Next()) 888 { 889 TopoDS_Shell shell = TopoDS::Shell(exp1.Current()); 890 if (shmap.FindIndex(shell) < 1) 891 { 892 shmap.Add (shell); 893 894 (*testout) << "shell " << shmap.FindIndex(shell) << " "; 895 (*testout) << ((shell.Orientation() == TopAbs_REVERSED) ? "-" : "+") << ", "; 896 (*testout) << ((exp1.Current().Orientation() == TopAbs_REVERSED) ? "-" : "+") << endl; 897 898 for (exp2.Init(shell, TopAbs_FACE); exp2.More(); exp2.Next()) 899 { 900 TopoDS_Face face = TopoDS::Face(exp2.Current()); 901 if (fmap.FindIndex(face) < 1) 902 { 903 fmap.Add (face); 904 905 for (exp3.Init(face, TopAbs_WIRE); exp3.More(); exp3.Next()) 906 { 907 TopoDS_Wire wire = TopoDS::Wire (exp3.Current()); 908 if (wmap.FindIndex(wire) < 1) 909 { 910 wmap.Add (wire); 911 912 for (exp4.Init(wire, TopAbs_EDGE); exp4.More(); exp4.Next()) 913 { 914 TopoDS_Edge edge = TopoDS::Edge(exp4.Current()); 915 if (emap.FindIndex(edge) < 1) 916 { 917 emap.Add (edge); 918 for (exp5.Init(edge, TopAbs_VERTEX); exp5.More(); exp5.Next()) 919 { 920 TopoDS_Vertex vertex = TopoDS::Vertex(exp5.Current()); 921 if (vmap.FindIndex(vertex) < 1) 922 vmap.Add (vertex); 923 } 924 } 925 } 926 } 927 } 928 } 929 } 930 } 931 } 932 933 934 // Free Faces 935 936 for (exp2.Init(shape, TopAbs_FACE, TopAbs_SHELL); exp2.More(); exp2.Next()) 937 { 938 TopoDS_Face face = TopoDS::Face(exp2.Current()); 939 if (fmap.FindIndex(face) < 1) 940 { 941 fmap.Add (face); 942 943 for (exp3.Init(exp2.Current(), TopAbs_WIRE); exp3.More(); exp3.Next()) 944 { 945 TopoDS_Wire wire = TopoDS::Wire (exp3.Current()); 946 if (wmap.FindIndex(wire) < 1) 947 { 948 wmap.Add (wire); 949 950 for (exp4.Init(exp3.Current(), TopAbs_EDGE); exp4.More(); exp4.Next()) 951 { 952 TopoDS_Edge edge = TopoDS::Edge(exp4.Current()); 953 if (emap.FindIndex(edge) < 1) 954 { 955 emap.Add (edge); 956 for (exp5.Init(exp4.Current(), TopAbs_VERTEX); exp5.More(); exp5.Next()) 957 { 958 TopoDS_Vertex vertex = TopoDS::Vertex(exp5.Current()); 959 if (vmap.FindIndex(vertex) < 1) 960 vmap.Add (vertex); 961 } 962 } 963 } 964 } 965 } 966 } 967 } 968 969 970 // Free Wires 971 972 for (exp3.Init(shape, TopAbs_WIRE, TopAbs_FACE); exp3.More(); exp3.Next()) 973 { 974 TopoDS_Wire wire = TopoDS::Wire (exp3.Current()); 975 if (wmap.FindIndex(wire) < 1) 976 { 977 wmap.Add (wire); 978 979 for (exp4.Init(exp3.Current(), TopAbs_EDGE); exp4.More(); exp4.Next()) 980 { 981 TopoDS_Edge edge = TopoDS::Edge(exp4.Current()); 982 if (emap.FindIndex(edge) < 1) 983 { 984 emap.Add (edge); 985 for (exp5.Init(exp4.Current(), TopAbs_VERTEX); exp5.More(); exp5.Next()) 986 { 987 TopoDS_Vertex vertex = TopoDS::Vertex(exp5.Current()); 988 if (vmap.FindIndex(vertex) < 1) 989 vmap.Add (vertex); 990 } 991 } 992 } 993 } 994 } 995 996 997 // Free Edges 998 999 for (exp4.Init(shape, TopAbs_EDGE, TopAbs_WIRE); exp4.More(); exp4.Next()) 1000 { 1001 TopoDS_Edge edge = TopoDS::Edge(exp4.Current()); 1002 if (emap.FindIndex(edge) < 1) 1003 { 1004 emap.Add (edge); 1005 for (exp5.Init(exp4.Current(), TopAbs_VERTEX); exp5.More(); exp5.Next()) 1006 { 1007 TopoDS_Vertex vertex = TopoDS::Vertex(exp5.Current()); 1008 if (vmap.FindIndex(vertex) < 1) 1009 vmap.Add (vertex); 1010 } 1011 } 1012 } 1013 1014 1015 // Free Vertices 1016 1017 for (exp5.Init(shape, TopAbs_VERTEX, TopAbs_EDGE); exp5.More(); exp5.Next()) 1018 { 1019 TopoDS_Vertex vertex = TopoDS::Vertex(exp5.Current()); 1020 if (vmap.FindIndex(vertex) < 1) 1021 vmap.Add (vertex); 1022 } 1023 1024 1025 1026 1027 facemeshstatus.DeleteAll(); 1028 facemeshstatus.SetSize (fmap.Extent()); 1029 facemeshstatus = 0; 1030 1031 // Philippose - 15/01/2009 1032 face_maxh.DeleteAll(); 1033 face_maxh.SetSize (fmap.Extent()); 1034 face_maxh = 1e99; // mparam.maxh; 1035 1036 // Philippose - 15/01/2010 1037 face_maxh_modified.DeleteAll(); 1038 face_maxh_modified.SetSize(fmap.Extent()); 1039 face_maxh_modified = 0; 1040 1041 1042 // Philippose - 17/01/2009 1043 face_sel_status.DeleteAll(); 1044 face_sel_status.SetSize (fmap.Extent()); 1045 face_sel_status = 0; 1046 1047 fvispar.SetSize (fmap.Extent()); 1048 evispar.SetSize (emap.Extent()); 1049 vvispar.SetSize (vmap.Extent()); 1050 1051 fsingular.SetSize (fmap.Extent()); 1052 esingular.SetSize (emap.Extent()); 1053 vsingular.SetSize (vmap.Extent()); 1054 1055 fsingular = esingular = vsingular = false; 1056 } 1057 1058 1059 SewFaces()1060 void OCCGeometry :: SewFaces () 1061 { 1062 (*testout) << "Trying to sew faces ..." << endl; 1063 cout << "Trying to sew faces ..." << flush; 1064 1065 BRepOffsetAPI_Sewing sewedObj(1); 1066 1067 for (int i = 1; i <= fmap.Extent(); i++) 1068 { 1069 TopoDS_Face face = TopoDS::Face (fmap(i)); 1070 sewedObj.Add (face); 1071 } 1072 1073 sewedObj.Perform(); 1074 1075 if (!sewedObj.SewedShape().IsNull()) 1076 { 1077 shape = sewedObj.SewedShape(); 1078 cout << " done" << endl; 1079 } 1080 else 1081 cout << " not possible"; 1082 } 1083 1084 1085 1086 1087 MakeSolid()1088 void OCCGeometry :: MakeSolid () 1089 { 1090 TopExp_Explorer exp0; 1091 1092 (*testout) << "Trying to build solids ..." << endl; 1093 cout << "Trying to build solids ..." << flush; 1094 1095 BRepBuilderAPI_MakeSolid ms; 1096 int count = 0; 1097 for (exp0.Init(shape, TopAbs_SHELL); exp0.More(); exp0.Next()) 1098 { 1099 count++; 1100 ms.Add (TopoDS::Shell(exp0.Current())); 1101 } 1102 1103 if (!count) 1104 { 1105 cout << " not possible (no shells)" << endl; 1106 return; 1107 } 1108 1109 BRepCheck_Analyzer ba(ms); 1110 if (ba.IsValid ()) 1111 { 1112 Handle(ShapeFix_Shape) sfs = new ShapeFix_Shape; 1113 sfs->Init (ms); 1114 1115 sfs->SetPrecision(1e-5); 1116 sfs->SetMaxTolerance(1e-5); 1117 1118 sfs->Perform(); 1119 1120 shape = sfs->Shape(); 1121 1122 for (exp0.Init(shape, TopAbs_SOLID); exp0.More(); exp0.Next()) 1123 { 1124 TopoDS_Solid solid = TopoDS::Solid(exp0.Current()); 1125 TopoDS_Solid newsolid = solid; 1126 BRepLib::OrientClosedSolid (newsolid); 1127 Handle(ShapeBuild_ReShape) rebuild = new ShapeBuild_ReShape; 1128 rebuild->Replace(solid, newsolid); 1129 1130 TopoDS_Shape newshape = rebuild->Apply(shape, TopAbs_SHAPE, 1); 1131 shape = newshape; 1132 } 1133 1134 cout << " done" << endl; 1135 } 1136 else 1137 cout << " not possible" << endl; 1138 } 1139 1140 1141 1142 BuildVisualizationMesh(double deflection)1143 void OCCGeometry :: BuildVisualizationMesh (double deflection) 1144 { 1145 cout << "Preparing visualization (deflection = " << deflection << ") ... " << flush; 1146 1147 BRepTools::Clean (shape); 1148 // BRepMesh_IncrementalMesh:: 1149 BRepMesh_IncrementalMesh (shape, deflection, true); 1150 cout << "done" << endl; 1151 } 1152 1153 1154 1155 CalcBoundingBox()1156 void OCCGeometry :: CalcBoundingBox () 1157 { 1158 Bnd_Box bb; 1159 #if OCC_VERSION_HEX < 0x070000 1160 BRepBndLib::Add (shape, bb); 1161 #else 1162 BRepBndLib::Add ((const TopoDS_Shape) shape, bb,(Standard_Boolean)true); 1163 #endif 1164 1165 double x1,y1,z1,x2,y2,z2; 1166 bb.Get (x1,y1,z1,x2,y2,z2); 1167 Point<3> p1 = Point<3> (x1,y1,z1); 1168 Point<3> p2 = Point<3> (x2,y2,z2); 1169 1170 (*testout) << "Bounding Box = [" << p1 << " - " << p2 << "]" << endl; 1171 boundingbox = Box<3> (p1,p2); 1172 SetCenter(); 1173 } 1174 ProjectPoint(int surfi,Point<3> & p) const1175 PointGeomInfo OCCGeometry :: ProjectPoint(int surfi, Point<3> & p) const 1176 { 1177 static int cnt = 0; 1178 if (++cnt % 1000 == 0) cout << "Project cnt = " << cnt << endl; 1179 1180 gp_Pnt pnt(p(0), p(1), p(2)); 1181 1182 double u,v; 1183 Handle( Geom_Surface ) thesurf = BRep_Tool::Surface(TopoDS::Face(fmap(surfi))); 1184 Handle( ShapeAnalysis_Surface ) su = new ShapeAnalysis_Surface( thesurf ); 1185 gp_Pnt2d suval = su->ValueOfUV ( pnt, BRep_Tool::Tolerance( TopoDS::Face(fmap(surfi)) ) ); 1186 suval.Coord( u, v); 1187 pnt = thesurf->Value( u, v ); 1188 1189 PointGeomInfo gi; 1190 gi.trignum = surfi; 1191 gi.u = u; 1192 gi.v = v; 1193 p = Point<3> (pnt.X(), pnt.Y(), pnt.Z()); 1194 return gi; 1195 } 1196 ProjectPointGI(int surfind,Point<3> & p,PointGeomInfo & gi) const1197 bool OCCGeometry :: ProjectPointGI(int surfind, Point<3>& p, PointGeomInfo& gi) const 1198 { 1199 double u = gi.u; 1200 double v = gi.v; 1201 1202 Point<3> hp = p; 1203 if (FastProject (surfind, hp, u, v)) 1204 { 1205 p = hp; 1206 return 1; 1207 } 1208 ProjectPoint (surfind, p); 1209 return CalcPointGeomInfo (surfind, gi, p); 1210 } 1211 ProjectPointEdge(int surfind,INDEX surfind2,Point<3> & p,EdgePointGeomInfo * gi) const1212 void OCCGeometry :: ProjectPointEdge(int surfind, INDEX surfind2, 1213 Point<3> & p, EdgePointGeomInfo* gi) const 1214 { 1215 TopExp_Explorer exp0, exp1; 1216 bool done = false; 1217 Handle(Geom_Curve) c; 1218 1219 for (exp0.Init(fmap(surfind), TopAbs_EDGE); !done && exp0.More(); exp0.Next()) 1220 for (exp1.Init(fmap(surfind2), TopAbs_EDGE); !done && exp1.More(); exp1.Next()) 1221 { 1222 if (TopoDS::Edge(exp0.Current()).IsSame(TopoDS::Edge(exp1.Current()))) 1223 { 1224 done = true; 1225 double s0, s1; 1226 c = BRep_Tool::Curve(TopoDS::Edge(exp0.Current()), s0, s1); 1227 } 1228 } 1229 1230 gp_Pnt pnt(p(0), p(1), p(2)); 1231 GeomAPI_ProjectPointOnCurve proj(pnt, c); 1232 pnt = proj.NearestPoint(); 1233 p(0) = pnt.X(); 1234 p(1) = pnt.Y(); 1235 p(2) = pnt.Z(); 1236 1237 } 1238 FastProject(int surfi,Point<3> & ap,double & u,double & v) const1239 bool OCCGeometry :: FastProject (int surfi, Point<3> & ap, double& u, double& v) const 1240 { 1241 gp_Pnt p(ap(0), ap(1), ap(2)); 1242 1243 Handle(Geom_Surface) surface = BRep_Tool::Surface(TopoDS::Face(fmap(surfi))); 1244 1245 gp_Pnt x = surface->Value (u,v); 1246 1247 if (p.SquareDistance(x) <= sqr(PROJECTION_TOLERANCE)) return true; 1248 1249 gp_Vec du, dv; 1250 1251 surface->D1(u,v,x,du,dv); 1252 1253 int count = 0; 1254 1255 gp_Pnt xold; 1256 gp_Vec n; 1257 double det, lambda, mu; 1258 1259 do { 1260 count++; 1261 1262 n = du^dv; 1263 1264 det = Det3 (n.X(), du.X(), dv.X(), 1265 n.Y(), du.Y(), dv.Y(), 1266 n.Z(), du.Z(), dv.Z()); 1267 1268 if (det < 1e-15) return false; 1269 1270 lambda = Det3 (n.X(), p.X()-x.X(), dv.X(), 1271 n.Y(), p.Y()-x.Y(), dv.Y(), 1272 n.Z(), p.Z()-x.Z(), dv.Z())/det; 1273 1274 mu = Det3 (n.X(), du.X(), p.X()-x.X(), 1275 n.Y(), du.Y(), p.Y()-x.Y(), 1276 n.Z(), du.Z(), p.Z()-x.Z())/det; 1277 1278 u += lambda; 1279 v += mu; 1280 1281 xold = x; 1282 surface->D1(u,v,x,du,dv); 1283 1284 } while (xold.SquareDistance(x) > sqr(PROJECTION_TOLERANCE) && count < 50); 1285 1286 // (*testout) << "FastProject count: " << count << endl; 1287 1288 if (count == 50) return false; 1289 1290 ap = Point<3> (x.X(), x.Y(), x.Z()); 1291 1292 return true; 1293 } 1294 GetNormal(int surfind,const Point<3> & p,const PointGeomInfo * geominfo) const1295 Vec<3> OCCGeometry :: GetNormal(int surfind, const Point<3> & p, const PointGeomInfo* geominfo) const 1296 { 1297 if(geominfo) 1298 { 1299 gp_Pnt pnt; 1300 gp_Vec du, dv; 1301 1302 Handle(Geom_Surface) occface; 1303 occface = BRep_Tool::Surface(TopoDS::Face(fmap(surfind))); 1304 1305 occface->D1(geominfo->u,geominfo->v,pnt,du,dv); 1306 1307 auto n = Cross (Vec<3>(du.X(), du.Y(), du.Z()), 1308 Vec<3>(dv.X(), dv.Y(), dv.Z())); 1309 n.Normalize(); 1310 1311 if (fmap(surfind).Orientation() == TopAbs_REVERSED) n *= -1; 1312 return n; 1313 } 1314 Standard_Real u,v; 1315 1316 gp_Pnt pnt(p(0), p(1), p(2)); 1317 1318 Handle(Geom_Surface) occface; 1319 occface = BRep_Tool::Surface(TopoDS::Face(fmap(surfind))); 1320 1321 /* 1322 GeomAPI_ProjectPointOnSurf proj(pnt, occface); 1323 1324 if (proj.NbPoints() < 1) 1325 { 1326 cout << "ERROR: OCCSurface :: GetNormalVector: GeomAPI_ProjectPointOnSurf failed!" 1327 << endl; 1328 cout << p << endl; 1329 return; 1330 } 1331 1332 proj.LowerDistanceParameters (u, v); 1333 */ 1334 1335 Handle( ShapeAnalysis_Surface ) su = new ShapeAnalysis_Surface( occface ); 1336 gp_Pnt2d suval = su->ValueOfUV ( pnt, BRep_Tool::Tolerance( TopoDS::Face(fmap(surfind)) ) ); 1337 suval.Coord( u, v); 1338 pnt = occface->Value( u, v ); 1339 1340 gp_Vec du, dv; 1341 occface->D1(u,v,pnt,du,dv); 1342 1343 /* 1344 if (!occface->IsCNu (1) || !occface->IsCNv (1)) 1345 (*testout) << "SurfOpt: Differentiation FAIL" << endl; 1346 */ 1347 1348 auto n = Cross (Vec3d(du.X(), du.Y(), du.Z()), 1349 Vec3d(dv.X(), dv.Y(), dv.Z())); 1350 n.Normalize(); 1351 1352 if (fmap(surfind).Orientation() == TopAbs_REVERSED) n *= -1; 1353 return n; 1354 } 1355 CalcPointGeomInfo(int surfind,PointGeomInfo & gi,const Point<3> & p) const1356 bool OCCGeometry :: CalcPointGeomInfo(int surfind, PointGeomInfo& gi, const Point<3> & p) const 1357 { 1358 Standard_Real u,v; 1359 1360 gp_Pnt pnt(p(0), p(1), p(2)); 1361 1362 Handle(Geom_Surface) occface; 1363 occface = BRep_Tool::Surface(TopoDS::Face(fmap(surfind))); 1364 1365 /* 1366 GeomAPI_ProjectPointOnSurf proj(pnt, occface); 1367 1368 if (proj.NbPoints() < 1) 1369 { 1370 cout << "ERROR: OCCSurface :: GetNormalVector: GeomAPI_ProjectPointOnSurf failed!" 1371 << endl; 1372 cout << p << endl; 1373 return 0; 1374 } 1375 1376 proj.LowerDistanceParameters (u, v); 1377 */ 1378 1379 Handle( ShapeAnalysis_Surface ) su = new ShapeAnalysis_Surface( occface ); 1380 gp_Pnt2d suval = su->ValueOfUV ( pnt, BRep_Tool::Tolerance( TopoDS::Face(fmap(surfind)) ) ); 1381 suval.Coord( u, v); 1382 //pnt = occface->Value( u, v ); 1383 1384 1385 gi.u = u; 1386 gi.v = v; 1387 return true; 1388 } 1389 PointBetween(const Point<3> & p1,const Point<3> & p2,double secpoint,int surfi,const PointGeomInfo & gi1,const PointGeomInfo & gi2,Point<3> & newp,PointGeomInfo & newgi) const1390 void OCCGeometry :: PointBetween(const Point<3> & p1, const Point<3> & p2, double secpoint, 1391 int surfi, 1392 const PointGeomInfo & gi1, 1393 const PointGeomInfo & gi2, 1394 Point<3> & newp, PointGeomInfo & newgi) const 1395 { 1396 Point<3> hnewp; 1397 hnewp = p1+secpoint*(p2-p1); 1398 1399 if (surfi > 0) 1400 { 1401 double u = gi1.u+secpoint*(gi2.u-gi1.u); 1402 double v = gi1.v+secpoint*(gi2.v-gi1.v); 1403 1404 auto savept = hnewp; 1405 if (!FastProject(surfi, hnewp, u, v) || Dist(hnewp, savept) > Dist(p1,p2)) 1406 { 1407 // cout << "Fast projection to surface fails! Using OCC projection" << endl; 1408 hnewp = savept; 1409 ProjectPoint(surfi, hnewp); 1410 } 1411 newgi.trignum = 1; 1412 newgi.u = u; 1413 newgi.v = v; 1414 } 1415 newp = hnewp; 1416 } 1417 1418 PointBetweenEdge(const Point<3> & p1,const Point<3> & p2,double secpoint,int surfi1,int surfi2,const EdgePointGeomInfo & ap1,const EdgePointGeomInfo & ap2,Point<3> & newp,EdgePointGeomInfo & newgi) const1419 void OCCGeometry :: PointBetweenEdge(const Point<3> & p1, 1420 const Point<3> & p2, double secpoint, 1421 int surfi1, int surfi2, 1422 const EdgePointGeomInfo & ap1, 1423 const EdgePointGeomInfo & ap2, 1424 Point<3> & newp, EdgePointGeomInfo & newgi) const 1425 { 1426 double s0, s1; 1427 1428 newp = p1+secpoint*(p2-p1); 1429 if(ap1.edgenr > emap.Size() || ap1.edgenr == 0) 1430 return; 1431 gp_Pnt pnt(newp(0), newp(1), newp(2)); 1432 GeomAPI_ProjectPointOnCurve proj(pnt, BRep_Tool::Curve(TopoDS::Edge(emap(ap1.edgenr)), s0, s1)); 1433 pnt = proj.NearestPoint(); 1434 newp = Point<3> (pnt.X(), pnt.Y(), pnt.Z()); 1435 newgi = ap1; 1436 }; 1437 1438 1439 // void OCCGeometry :: WriteOCC_STL(char * filename) 1440 // { 1441 // cout << "writing stl..."; cout.flush(); 1442 // StlAPI_Writer writer; 1443 // writer.RelativeMode() = Standard_False; 1444 // 1445 // writer.SetDeflection(0.02); 1446 // writer.Write(shape,filename); 1447 // 1448 // cout << "done" << endl; 1449 // } 1450 1451 LoadOCCInto(OCCGeometry * occgeo,const char * filename)1452 void LoadOCCInto(OCCGeometry* occgeo, const char* filename) 1453 { 1454 static Timer timer_all("LoadOCC"); RegionTimer rtall(timer_all); 1455 static Timer timer_readfile("LoadOCC-ReadFile"); 1456 static Timer timer_transfer("LoadOCC-Transfer"); 1457 static Timer timer_getnames("LoadOCC-get names"); 1458 1459 // Initiate a dummy XCAF Application to handle the STEP XCAF Document 1460 static Handle(XCAFApp_Application) dummy_app = XCAFApp_Application::GetApplication(); 1461 1462 // Create an XCAF Document to contain the STEP file itself 1463 Handle(TDocStd_Document) step_doc; 1464 1465 // Check if a STEP File is already open under this handle, if so, close it to prevent 1466 // Segmentation Faults when trying to create a new document 1467 if(dummy_app->NbDocuments() > 0) 1468 { 1469 dummy_app->GetDocument(1,step_doc); 1470 dummy_app->Close(step_doc); 1471 } 1472 dummy_app->NewDocument ("STEP-XCAF",step_doc); 1473 1474 timer_readfile.Start(); 1475 STEPCAFControl_Reader reader; 1476 1477 // Enable transfer of colours 1478 reader.SetColorMode(Standard_True); 1479 reader.SetNameMode(Standard_True); 1480 Standard_Integer stat = reader.ReadFile((char*)filename); 1481 timer_readfile.Stop(); 1482 1483 timer_transfer.Start(); 1484 if(stat != IFSelect_RetDone) 1485 { 1486 throw NgException("Couldn't load OCC geometry"); 1487 } 1488 1489 reader.Transfer(step_doc); 1490 timer_transfer.Stop(); 1491 1492 // Read in the shape(s) and the colours present in the STEP File 1493 Handle(XCAFDoc_ShapeTool) step_shape_contents = XCAFDoc_DocumentTool::ShapeTool(step_doc->Main()); 1494 Handle(XCAFDoc_ColorTool) step_colour_contents = XCAFDoc_DocumentTool::ColorTool(step_doc->Main()); 1495 1496 TDF_LabelSequence step_shapes; 1497 step_shape_contents->GetShapes(step_shapes); 1498 1499 // List out the available colours in the STEP File as Colour Names 1500 TDF_LabelSequence all_colours; 1501 step_colour_contents->GetColors(all_colours); 1502 PrintMessage(1,"Number of colours in STEP File: ",all_colours.Length()); 1503 for(int i = 1; i <= all_colours.Length(); i++) 1504 { 1505 Quantity_Color col; 1506 stringstream col_rgb; 1507 step_colour_contents->GetColor(all_colours.Value(i),col); 1508 col_rgb << " : (" << col.Red() << "," << col.Green() << "," << col.Blue() << ")"; 1509 PrintMessage(1, "Colour [", i, "] = ",col.StringName(col.Name()),col_rgb.str()); 1510 } 1511 1512 1513 // For the STEP File Reader in OCC, the 1st Shape contains the entire 1514 // compound geometry as one shape 1515 occgeo->shape = step_shape_contents->GetShape(step_shapes.Value(1)); 1516 occgeo->face_colours = step_colour_contents; 1517 occgeo->changed = 1; 1518 occgeo->BuildFMap(); 1519 1520 occgeo->CalcBoundingBox(); 1521 PrintContents (occgeo); 1522 string name; 1523 TopExp_Explorer exp0,exp1; 1524 1525 1526 1527 std::map<Handle(TopoDS_TShape), string> shape_names; 1528 { 1529 static Timer t("file shape_names"); RegionTimer r(t); 1530 // code inspired from 1531 // https://www.opencascade.com/content/reading-step-entity-id-slow 1532 const Handle(XSControl_WorkSession) workSession = reader.Reader().WS(); 1533 const Handle(Interface_InterfaceModel) model = workSession->Model(); 1534 const Handle(XSControl_TransferReader) transferReader = workSession->TransferReader(); 1535 Handle(Transfer_TransientProcess) transProc = transferReader->TransientProcess(); 1536 1537 Standard_Integer nb = model->NbEntities(); 1538 for (Standard_Integer i = 1; i < nb; i++) 1539 { 1540 Handle(Standard_Transient) entity = model->Value(i); 1541 1542 // if (!entity->DynamicType()->SubType("StepShape_OpenShell")) continue; 1543 1544 Handle(StepRepr_RepresentationItem) SRRI = 1545 Handle(StepRepr_RepresentationItem)::DownCast(entity); 1546 1547 if (SRRI.IsNull()) { 1548 // cout << "no StepRepr_RepresentationItem found in " << entity->DynamicType()->Name(); 1549 continue; 1550 } 1551 Handle(TCollection_HAsciiString) hName = SRRI->Name(); 1552 string shapeName = hName->ToCString(); 1553 1554 // cout << "STEP " << i << " " << entity->DynamicType()->Name() << ", shapename = " << shapeName; 1555 Handle(Transfer_Binder) binder; 1556 if (!transProc->IsBound(SRRI)) { 1557 // cout << "found unbound entity " << shapeName; 1558 continue; 1559 } 1560 binder = transProc->Find(SRRI); 1561 TopoDS_Shape shape = TransferBRep::ShapeResult(binder); 1562 // if (!shape.IsNull()) 1563 shape_names[shape.TShape()] = shapeName; 1564 /* 1565 if (!shape.IsNull()) 1566 cout << " shapetype = " << shape.ShapeType() << endl; 1567 else 1568 cout << "is-Null" << endl; 1569 */ 1570 } 1571 // for (auto pair : shape_names) 1572 // cout << "name = " << pair.second << endl; 1573 } 1574 1575 for (auto [s,n] : shape_names) 1576 OCCGeometry::global_shape_properties[s].name = n; 1577 1578 1579 timer_getnames.Start(); 1580 for (exp0.Init(occgeo->shape, TopAbs_SOLID); exp0.More(); exp0.Next()) 1581 { 1582 TopoDS_Solid solid = TopoDS::Solid(exp0.Current()); 1583 // name = STEP_GetEntityName(solid,&reader); 1584 // cout << "solidname = " << name << ", mapname = " << shape_names[solid.TShape()] << endl; 1585 name = shape_names[solid.TShape()]; 1586 if (name == "") 1587 name = string("domain_") + ToString(occgeo->snames.Size()); 1588 occgeo->snames.Append(name); 1589 } 1590 1591 for (exp0.Init(occgeo->shape, TopAbs_FACE); exp0.More(); exp0.Next()) 1592 { 1593 TopoDS_Face face = TopoDS::Face(exp0.Current()); 1594 // name = STEP_GetEntityName(face,&reader); 1595 // cout << "getname = " << name << ", mapname = " << shape_names[face.TShape()] << endl; 1596 name = shape_names[face.TShape()]; 1597 if (name == "") 1598 name = string("bc_") + ToString(occgeo->fnames.Size()); 1599 occgeo->fnames.Append(name); 1600 for (exp1.Init(face, TopAbs_EDGE); exp1.More(); exp1.Next()) 1601 { 1602 TopoDS_Edge edge = TopoDS::Edge(exp1.Current()); 1603 // name = STEP_GetEntityName(edge,&reader); 1604 // cout << "getname = " << name << ", mapname = " << shape_names[edge.TShape()] << endl; 1605 name = shape_names[edge.TShape()]; 1606 occgeo->enames.Append(name); 1607 } 1608 } 1609 timer_getnames.Stop(); 1610 } 1611 1612 // Philippose - 23/02/2009 1613 /* Special IGES File load function including the ability 1614 to extract individual surface colours via the extended 1615 OpenCascade XDE and XCAF Feature set. 1616 */ LoadOCC_IGES(const char * filename)1617 OCCGeometry *LoadOCC_IGES(const char *filename) 1618 { 1619 OCCGeometry *occgeo; 1620 occgeo = new OCCGeometry; 1621 // Initiate a dummy XCAF Application to handle the IGES XCAF Document 1622 static Handle(XCAFApp_Application) dummy_app = XCAFApp_Application::GetApplication(); 1623 1624 // Create an XCAF Document to contain the IGES file itself 1625 Handle(TDocStd_Document) iges_doc; 1626 1627 // Check if a IGES File is already open under this handle, if so, close it to prevent 1628 // Segmentation Faults when trying to create a new document 1629 if(dummy_app->NbDocuments() > 0) 1630 { 1631 dummy_app->GetDocument(1,iges_doc); 1632 dummy_app->Close(iges_doc); 1633 } 1634 dummy_app->NewDocument ("IGES-XCAF",iges_doc); 1635 1636 IGESCAFControl_Reader reader; 1637 1638 Standard_Integer stat = reader.ReadFile((char*)filename); 1639 1640 if(stat != IFSelect_RetDone) 1641 { 1642 throw NgException("Couldn't load occ"); 1643 } 1644 1645 // Enable transfer of colours 1646 reader.SetColorMode(Standard_True); 1647 1648 reader.Transfer(iges_doc); 1649 1650 // Read in the shape(s) and the colours present in the IGES File 1651 Handle(XCAFDoc_ShapeTool) iges_shape_contents = XCAFDoc_DocumentTool::ShapeTool(iges_doc->Main()); 1652 Handle(XCAFDoc_ColorTool) iges_colour_contents = XCAFDoc_DocumentTool::ColorTool(iges_doc->Main()); 1653 1654 TDF_LabelSequence iges_shapes; 1655 iges_shape_contents->GetShapes(iges_shapes); 1656 1657 // List out the available colours in the IGES File as Colour Names 1658 TDF_LabelSequence all_colours; 1659 iges_colour_contents->GetColors(all_colours); 1660 PrintMessage(1,"Number of colours in IGES File: ",all_colours.Length()); 1661 for(int i = 1; i <= all_colours.Length(); i++) 1662 { 1663 Quantity_Color col; 1664 stringstream col_rgb; 1665 iges_colour_contents->GetColor(all_colours.Value(i),col); 1666 col_rgb << " : (" << col.Red() << "," << col.Green() << "," << col.Blue() << ")"; 1667 PrintMessage(1, "Colour [", i, "] = ",col.StringName(col.Name()),col_rgb.str()); 1668 } 1669 1670 1671 // For the IGES Reader, all the shapes can be exported as one compound shape 1672 // using the "OneShape" member 1673 occgeo->shape = reader.OneShape(); 1674 occgeo->face_colours = iges_colour_contents; 1675 occgeo->changed = 1; 1676 occgeo->BuildFMap(); 1677 1678 occgeo->CalcBoundingBox(); 1679 PrintContents (occgeo); 1680 return occgeo; 1681 } 1682 1683 1684 1685 1686 1687 // Philippose - 29/01/2009 1688 /* Special STEP File load function including the ability 1689 to extract individual surface colours via the extended 1690 OpenCascade XDE and XCAF Feature set. 1691 */ LoadOCC_STEP(const char * filename)1692 OCCGeometry * LoadOCC_STEP (const char * filename) 1693 { 1694 OCCGeometry * occgeo; 1695 occgeo = new OCCGeometry; 1696 1697 LoadOCCInto(occgeo, filename); 1698 return occgeo; 1699 } 1700 1701 1702 1703 LoadOCC_BREP(const char * filename)1704 OCCGeometry *LoadOCC_BREP (const char *filename) 1705 { 1706 OCCGeometry * occgeo; 1707 occgeo = new OCCGeometry; 1708 1709 BRep_Builder aBuilder; 1710 Standard_Boolean result = BRepTools::Read(occgeo->shape, const_cast<char*> (filename),aBuilder); 1711 1712 if(!result) 1713 { 1714 delete occgeo; 1715 return NULL; 1716 } 1717 1718 // Philippose - 23/02/2009 1719 // Fixed a bug in the OpenCascade XDE Colour handling when 1720 // opening BREP Files, since BREP Files have no colour data. 1721 // Hence, the face_colours Handle needs to be created as a NULL handle. 1722 occgeo->face_colours = Handle(XCAFDoc_ColorTool)(); 1723 occgeo->face_colours.Nullify(); 1724 occgeo->changed = 1; 1725 occgeo->BuildFMap(); 1726 1727 occgeo->CalcBoundingBox(); 1728 PrintContents (occgeo); 1729 1730 return occgeo; 1731 } 1732 1733 Save(string sfilename) const1734 void OCCGeometry :: Save (string sfilename) const 1735 { 1736 const char * filename = sfilename.c_str(); 1737 if (strlen(filename) < 4) 1738 throw NgException ("illegal filename"); 1739 1740 if (strcmp (&filename[strlen(filename)-3], "igs") == 0) 1741 { 1742 IGESControl_Writer writer("millimeters", 1); 1743 writer.AddShape (shape); 1744 writer.Write (filename); 1745 } 1746 else if (strcmp (&filename[strlen(filename)-3], "stp") == 0) 1747 { 1748 STEPControl_Writer writer; 1749 writer.Transfer (shape, STEPControl_AsIs); 1750 writer.Write (filename); 1751 } 1752 else if (strcmp (&filename[strlen(filename)-3], "stl") == 0) 1753 { 1754 StlAPI_Writer writer; 1755 writer.ASCIIMode() = Standard_True; 1756 writer.Write (shape, filename); 1757 } 1758 else if (strcmp (&filename[strlen(filename)-4], "stlb") == 0) 1759 { 1760 StlAPI_Writer writer; 1761 writer.ASCIIMode() = Standard_False; 1762 writer.Write (shape, filename); 1763 } 1764 } 1765 DoArchive(Archive & ar)1766 void OCCGeometry :: DoArchive(Archive& ar) 1767 { 1768 if(ar.Output()) 1769 { 1770 std::stringstream ss; 1771 STEPControl_Writer writer; 1772 writer.Transfer(shape, STEPControl_AsIs); 1773 auto filename = ".tmpfile_out.step"; 1774 writer.Write(filename); 1775 std::ifstream is(filename); 1776 ss << is.rdbuf(); 1777 ar << ss.str(); 1778 std::remove(filename); 1779 } 1780 else 1781 { 1782 std::string str; 1783 ar & str; 1784 1785 auto filename = ".tmpfile.step"; 1786 auto tmpfile = std::fopen(filename, "w"); 1787 std::fputs(str.c_str(), tmpfile); 1788 std::fclose(tmpfile); 1789 LoadOCCInto(this, filename); 1790 std::remove(filename); 1791 } 1792 } 1793 1794 const char * shapesname[] = 1795 {" ", "CompSolids", "Solids", "Shells", 1796 1797 "Faces", "Wires", "Edges", "Vertices"}; 1798 1799 const char * shapename[] = 1800 {" ", "CompSolid", "Solid", "Shell", 1801 "Face", "Wire", "Edge", "Vertex"}; 1802 1803 const char * orientationstring[] = 1804 {"+", "-"}; 1805 1806 1807 1808 RecursiveTopologyTree(const TopoDS_Shape & sh,stringstream & str,TopAbs_ShapeEnum l,bool isfree,const char * lname)1809 void OCCGeometry :: RecursiveTopologyTree (const TopoDS_Shape & sh, 1810 stringstream & str, 1811 TopAbs_ShapeEnum l, 1812 bool isfree, 1813 const char * lname) 1814 { 1815 if (l > TopAbs_VERTEX) return; 1816 1817 TopExp_Explorer e; 1818 int count = 0; 1819 int count2 = 0; 1820 1821 if (isfree) 1822 e.Init(sh, l, TopAbs_ShapeEnum(l-1)); 1823 else 1824 e.Init(sh, l); 1825 1826 for (; e.More(); e.Next()) 1827 { 1828 count++; 1829 1830 stringstream lname2; 1831 lname2 << lname << "/" << shapename[l] << count; 1832 str << lname2.str() << " "; 1833 1834 switch (e.Current().ShapeType()) 1835 { 1836 case TopAbs_SOLID: 1837 count2 = somap.FindIndex(TopoDS::Solid(e.Current())); break; 1838 case TopAbs_SHELL: 1839 count2 = shmap.FindIndex(TopoDS::Shell(e.Current())); break; 1840 case TopAbs_FACE: 1841 count2 = fmap.FindIndex(TopoDS::Face(e.Current())); break; 1842 case TopAbs_WIRE: 1843 count2 = wmap.FindIndex(TopoDS::Wire(e.Current())); break; 1844 case TopAbs_EDGE: 1845 count2 = emap.FindIndex(TopoDS::Edge(e.Current())); break; 1846 case TopAbs_VERTEX: 1847 count2 = vmap.FindIndex(TopoDS::Vertex(e.Current())); break; 1848 default: 1849 cout << "RecursiveTopologyTree: Case " << e.Current().ShapeType() << " not handeled" << endl; 1850 } 1851 1852 int nrsubshapes = 0; 1853 1854 if (l <= TopAbs_WIRE) 1855 { 1856 TopExp_Explorer e2; 1857 for (e2.Init (e.Current(), TopAbs_ShapeEnum (l+1)); 1858 e2.More(); e2.Next()) 1859 nrsubshapes++; 1860 } 1861 1862 str << "{" << shapename[l] << " " << count2; 1863 1864 if (l <= TopAbs_EDGE) 1865 { 1866 str << " (" << orientationstring[e.Current().Orientation()]; 1867 if (nrsubshapes != 0) str << ", " << nrsubshapes; 1868 str << ") } "; 1869 } 1870 else 1871 str << " } "; 1872 1873 RecursiveTopologyTree (e.Current(), str, TopAbs_ShapeEnum (l+1), 1874 false, (char*)lname2.str().c_str()); 1875 1876 } 1877 } 1878 1879 1880 1881 GetTopologyTree(stringstream & str)1882 void OCCGeometry :: GetTopologyTree (stringstream & str) 1883 { 1884 cout << "Building topology tree ... " << flush; 1885 RecursiveTopologyTree (shape, str, TopAbs_COMPSOLID, false, "CompSolids"); 1886 RecursiveTopologyTree (shape, str, TopAbs_SOLID, true, "FreeSolids"); 1887 RecursiveTopologyTree (shape, str, TopAbs_SHELL, true, "FreeShells"); 1888 RecursiveTopologyTree (shape, str, TopAbs_FACE, true, "FreeFaces"); 1889 RecursiveTopologyTree (shape, str, TopAbs_WIRE, true, "FreeWires"); 1890 RecursiveTopologyTree (shape, str, TopAbs_EDGE, true, "FreeEdges"); 1891 RecursiveTopologyTree (shape, str, TopAbs_VERTEX, true, "FreeVertices"); 1892 str << flush; 1893 // cout << "done" << endl; 1894 } 1895 1896 1897 1898 CheckIrregularEntities(stringstream & str)1899 void OCCGeometry :: CheckIrregularEntities(stringstream & str) 1900 { 1901 ShapeAnalysis_CheckSmallFace csm; 1902 1903 csm.SetTolerance (1e-6); 1904 1905 TopTools_DataMapOfShapeListOfShape mapEdges; 1906 ShapeAnalysis_DataMapOfShapeListOfReal mapParam; 1907 TopoDS_Compound theAllVert; 1908 1909 int spotfaces = 0; 1910 int stripsupportfaces = 0; 1911 int singlestripfaces = 0; 1912 int stripfaces = 0; 1913 int facessplitbyvertices = 0; 1914 int stretchedpinfaces = 0; 1915 int smoothpinfaces = 0; 1916 int twistedfaces = 0; 1917 // int edgessamebutnotidentified = 0; 1918 1919 cout << "checking faces ... " << flush; 1920 1921 int i; 1922 for (i = 1; i <= fmap.Extent(); i++) 1923 { 1924 TopoDS_Face face = TopoDS::Face (fmap(i)); 1925 TopoDS_Edge e1, e2; 1926 1927 if (csm.CheckSpotFace (face)) 1928 { 1929 if (!spotfaces++) 1930 str << "SpotFace {Spot face} "; 1931 1932 (*testout) << "Face " << i << " is a spot face" << endl; 1933 str << "SpotFace/Face" << i << " "; 1934 str << "{Face " << i << " } "; 1935 } 1936 1937 if (csm.IsStripSupport (face)) 1938 { 1939 if (!stripsupportfaces++) 1940 str << "StripSupportFace {Strip support face} "; 1941 1942 (*testout) << "Face " << i << " has strip support" << endl; 1943 str << "StripSupportFace/Face" << i << " "; 1944 str << "{Face " << i << " } "; 1945 } 1946 1947 if (csm.CheckSingleStrip(face, e1, e2)) 1948 { 1949 if (!singlestripfaces++) 1950 str << "SingleStripFace {Single strip face} "; 1951 1952 (*testout) << "Face " << i << " is a single strip (edge " << emap.FindIndex(e1) 1953 << " and edge " << emap.FindIndex(e2) << " are identical)" << endl; 1954 str << "SingleStripFace/Face" << i << " "; 1955 str << "{Face " << i << " (edge " << emap.FindIndex(e1) 1956 << " and edge " << emap.FindIndex(e2) << " are identical)} "; 1957 } 1958 1959 if (csm.CheckStripFace(face, e1, e2)) 1960 { 1961 if (!stripfaces++) 1962 str << "StripFace {Strip face} "; 1963 1964 (*testout) << "Face " << i << " is a strip (edge " << emap.FindIndex(e1) 1965 << " and edge " << emap.FindIndex(e2) 1966 << " are identical)" << endl; 1967 str << "StripFace/Face" << i << " "; 1968 str << "{Face " << i << " (edge " << emap.FindIndex(e1) 1969 << " and edge " << emap.FindIndex(e2) << " are identical)} "; 1970 } 1971 1972 if (int count = csm.CheckSplittingVertices(face, mapEdges, mapParam, theAllVert)) 1973 { 1974 if (!facessplitbyvertices++) 1975 str << "FaceSplitByVertices {Face split by vertices} "; 1976 1977 (*testout) << "Face " << i << " is split by " << count 1978 << " vertex/vertices " << endl; 1979 str << "FaceSplitByVertices/Face" << i << " "; 1980 str << "{Face " << i << " (split by " << count << "vertex/vertices)} "; 1981 } 1982 1983 int whatrow, sens; 1984 if (int type = csm.CheckPin (face, whatrow, sens)) 1985 { 1986 if (type == 1) 1987 { 1988 if (!smoothpinfaces++) 1989 str << "SmoothPinFace {Smooth pin face} "; 1990 1991 (*testout) << "Face " << i << " is a smooth pin" << endl; 1992 str << "SmoothPinFace/Face" << i << " "; 1993 str << "{Face " << i << " } "; 1994 } 1995 else 1996 { 1997 if (!stretchedpinfaces++) 1998 str << "StretchedPinFace {Stretched pin face} "; 1999 2000 (*testout) << "Face " << i << " is a stretched pin" << endl; 2001 str << "StretchedPinFace/Face" << i << " "; 2002 str << "{Face " << i << " } "; 2003 } 2004 } 2005 2006 double paramu, paramv; 2007 if (csm.CheckTwisted (face, paramu, paramv)) 2008 { 2009 if (!twistedfaces++) 2010 str << "TwistedFace {Twisted face} "; 2011 2012 (*testout) << "Face " << i << " is twisted" << endl; 2013 str << "TwistedFace/Face" << i << " "; 2014 str << "{Face " << i << " } "; 2015 } 2016 } 2017 2018 cout << "done" << endl; 2019 cout << "checking edges ... " << flush; 2020 2021 // double dmax; 2022 // int cnt = 0; 2023 NgArray <double> edgeLengths; 2024 NgArray <int> order; 2025 edgeLengths.SetSize (emap.Extent()); 2026 order.SetSize (emap.Extent()); 2027 2028 for (i = 1; i <= emap.Extent(); i++) 2029 { 2030 TopoDS_Edge edge1 = TopoDS::Edge (emap(i)); 2031 GProp_GProps system; 2032 BRepGProp::LinearProperties(edge1, system); 2033 edgeLengths[i-1] = system.Mass(); 2034 } 2035 2036 Sort (edgeLengths, order); 2037 2038 str << "ShortestEdges {Shortest edges} "; 2039 for (i = 1; i <= min(20, emap.Extent()); i++) 2040 { 2041 str << "ShortestEdges/Edge" << i; 2042 str << " {Edge " << order[i-1] << " (L=" << edgeLengths[order[i-1]-1] << ")} "; 2043 } 2044 2045 str << flush; 2046 2047 cout << "done" << endl; 2048 } 2049 2050 2051 2052 GetUnmeshedFaceInfo(stringstream & str)2053 void OCCGeometry :: GetUnmeshedFaceInfo (stringstream & str) 2054 { 2055 for (int i = 1; i <= fmap.Extent(); i++) 2056 { 2057 if (facemeshstatus[i-1] == -1) 2058 str << "Face" << i << " {Face " << i << " } "; 2059 } 2060 str << flush; 2061 } 2062 2063 2064 2065 GetNotDrawableFaces(stringstream & str)2066 void OCCGeometry :: GetNotDrawableFaces (stringstream & str) 2067 { 2068 for (int i = 1; i <= fmap.Extent(); i++) 2069 { 2070 if (!fvispar[i-1].IsDrawable()) 2071 str << "Face" << i << " {Face " << i << " } "; 2072 } 2073 str << flush; 2074 } 2075 2076 2077 2078 ErrorInSurfaceMeshing()2079 bool OCCGeometry :: ErrorInSurfaceMeshing () 2080 { 2081 for (int i = 1; i <= fmap.Extent(); i++) 2082 if (facemeshstatus[i-1] == -1) 2083 return true; 2084 2085 return false; 2086 } 2087 Print(ostream & ost) const2088 void OCCParameters :: Print(ostream & ost) const 2089 { 2090 ost << "OCC Parameters:" << endl 2091 << "minimum edge length: " << resthminedgelenenable 2092 << ", min len = " << resthminedgelen << endl; 2093 } 2094 2095 DLL_HEADER extern OCCParameters occparam; 2096 OCCParameters occparam; 2097 2098 // int OCCGeometry :: GenerateMesh (shared_ptr<Mesh> & mesh, MeshingParameters & mparam) 2099 // { 2100 // return OCCGenerateMesh (*this, mesh, mparam, occparam); 2101 // } 2102 static RegisterClassForArchive<OCCGeometry, NetgenGeometry> regnggeo; 2103 } 2104 2105 2106 #endif 2107