00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024 #include "PreCompiled.h"
00025
00026 #ifndef _PreComp_
00027 # include <algorithm>
00028 #endif
00029
00030 #include "Algorithm.h"
00031 #include "Approximation.h"
00032 #include "Elements.h"
00033 #include "Iterator.h"
00034 #include "Grid.h"
00035 #include "Triangulation.h"
00036
00037 #include <Base/Console.h>
00038 #include <Base/Sequencer.h>
00039
00040 using namespace MeshCore;
00041 using Base::BoundBox3f;
00042 using Base::BoundBox2D;
00043 using Base::Polygon2D;
00044
00045
00046 bool MeshAlgorithm::IsVertexVisible (const Base::Vector3f &rcVertex, const Base::Vector3f &rcView, const MeshFacetGrid &rclGrid ) const
00047 {
00048 Base::Vector3f cDirection = rcVertex-rcView;
00049 float fDistance = cDirection.Length();
00050 Base::Vector3f cIntsct; unsigned long uInd;
00051
00052
00053 if ( NearestFacetOnRay( rcView, cDirection, rclGrid, cIntsct, uInd) )
00054 {
00055
00056 float fLen = Base::Distance( rcView, cIntsct );
00057 if ( fLen < fDistance )
00058 {
00059
00060 if ( Base::Distance(rcVertex, cIntsct) > 0.001f )
00061 {
00062
00063 return false;
00064 }
00065 }
00066 }
00067
00068 return true;
00069 }
00070
00071 bool MeshAlgorithm::NearestFacetOnRay (const Base::Vector3f &rclPt, const Base::Vector3f &rclDir, Base::Vector3f &rclRes,
00072 unsigned long &rulFacet) const
00073 {
00074 Base::Vector3f clProj, clRes;
00075 bool bSol = false;
00076 unsigned long ulInd = 0;
00077
00078
00079 MeshFacetIterator clFIter(_rclMesh);
00080 for (clFIter.Init(); clFIter.More(); clFIter.Next()) {
00081 if (clFIter->Foraminate( rclPt, rclDir, clRes ) == true) {
00082 if (bSol == false) {
00083 bSol = true;
00084 clProj = clRes;
00085 ulInd = clFIter.Position();
00086 }
00087 else {
00088 if ((clRes - rclPt).Length() < (clProj - rclPt).Length()) {
00089 clProj = clRes;
00090 ulInd = clFIter.Position();
00091 }
00092 }
00093 }
00094 }
00095
00096 if (bSol) {
00097 rclRes = clProj;
00098 rulFacet = ulInd;
00099 }
00100
00101 return bSol;
00102 }
00103
00104 bool MeshAlgorithm::NearestFacetOnRay (const Base::Vector3f &rclPt, const Base::Vector3f &rclDir, const MeshFacetGrid &rclGrid,
00105 Base::Vector3f &rclRes, unsigned long &rulFacet) const
00106 {
00107 std::vector<unsigned long> aulFacets;
00108 MeshGridIterator clGridIter(rclGrid);
00109
00110 if (clGridIter.InitOnRay(rclPt, rclDir, aulFacets) == true) {
00111 if (RayNearestField(rclPt, rclDir, aulFacets, rclRes, rulFacet) == false) {
00112 aulFacets.clear();
00113 while (clGridIter.NextOnRay(aulFacets) == true) {
00114 if (RayNearestField(rclPt, rclDir, aulFacets, rclRes, rulFacet) == true)
00115 return true;
00116 }
00117 }
00118 else
00119 return true;
00120 }
00121
00122 return false;
00123 }
00124
00125 bool MeshAlgorithm::NearestFacetOnRay (const Base::Vector3f &rclPt, const Base::Vector3f &rclDir, float fMaxSearchArea,
00126 const MeshFacetGrid &rclGrid, Base::Vector3f &rclRes, unsigned long &rulFacet) const
00127 {
00128 std::vector<unsigned long> aulFacets;
00129 MeshGridIterator clGridIter(rclGrid);
00130
00131 if (clGridIter.InitOnRay(rclPt, rclDir, fMaxSearchArea, aulFacets) == true) {
00132 if (RayNearestField(rclPt, rclDir, aulFacets, rclRes, rulFacet, 1.75f) == false) {
00133 aulFacets.clear();
00134 while (clGridIter.NextOnRay(aulFacets) == true) {
00135 if (RayNearestField(rclPt, rclDir, aulFacets, rclRes, rulFacet, 1.75f) == true)
00136 return true;
00137 }
00138 }
00139 else
00140 return true;
00141 }
00142
00143 return false;
00144 }
00145
00146 bool MeshAlgorithm::NearestFacetOnRay (const Base::Vector3f &rclPt, const Base::Vector3f &rclDir, const std::vector<unsigned long> &raulFacets,
00147 Base::Vector3f &rclRes, unsigned long &rulFacet) const
00148 {
00149 Base::Vector3f clProj, clRes;
00150 bool bSol = false;
00151 unsigned long ulInd = 0;
00152
00153 for (std::vector<unsigned long>::const_iterator pI = raulFacets.begin(); pI != raulFacets.end(); pI++) {
00154 MeshGeomFacet rclSFacet = _rclMesh.GetFacet(*pI);
00155 if (rclSFacet.Foraminate(rclPt, rclDir, clRes) == true) {
00156 if (bSol == false) {
00157 bSol = true;
00158 clProj = clRes;
00159 ulInd = *pI;
00160 }
00161 else {
00162 if ((clRes - rclPt).Length() < (clProj - rclPt).Length()) {
00163 clProj = clRes;
00164 ulInd = *pI;
00165 }
00166 }
00167 }
00168 }
00169
00170 if (bSol) {
00171 rclRes = clProj;
00172 rulFacet = ulInd;
00173 }
00174
00175 return bSol;
00176 }
00177
00178 bool MeshAlgorithm::RayNearestField (const Base::Vector3f &rclPt, const Base::Vector3f &rclDir, const std::vector<unsigned long> &raulFacets,
00179 Base::Vector3f &rclRes, unsigned long &rulFacet, float fMaxAngle) const
00180 {
00181 Base::Vector3f clProj, clRes;
00182 bool bSol = false;
00183 unsigned long ulInd = 0;
00184
00185 for (std::vector<unsigned long>::const_iterator pF = raulFacets.begin(); pF != raulFacets.end(); pF++) {
00186 if (_rclMesh.GetFacet(*pF).Foraminate(rclPt, rclDir, clRes) == true) {
00187 if (bSol == false) {
00188 bSol = true;
00189 clProj = clRes;
00190 ulInd = *pF;
00191 }
00192 else {
00193 if ((clRes - rclPt).Length() < (clProj - rclPt).Length()) {
00194 clProj = clRes;
00195 ulInd = *pF;
00196 }
00197 }
00198 }
00199 }
00200
00201 if (bSol) {
00202 rclRes = clProj;
00203 rulFacet = ulInd;
00204 }
00205
00206 return bSol;
00207 }
00208
00209 bool MeshAlgorithm::FirstFacetToVertex(const Base::Vector3f &rPt, float fMaxDistance, const MeshFacetGrid &rGrid, unsigned long &uIndex) const
00210 {
00211 const float fEps = 0.001f;
00212
00213 bool found = false;
00214 std::vector<unsigned long> facets;
00215
00216
00217 rGrid.GetElements(rPt, facets);
00218
00219
00220 for (std::vector<unsigned long>::iterator it = facets.begin(); it != facets.end(); ++it) {
00221 MeshGeomFacet cFacet = this->_rclMesh.GetFacet(*it);
00222 if (cFacet.IsPointOfFace(rPt, fMaxDistance)) {
00223 found = true;
00224 uIndex = *it;
00225 break;
00226 }
00227 else {
00228
00229 Base::Vector3f res = rPt;
00230 float fDist;
00231 unsigned short uSide;
00232 cFacet.ProjectPointToPlane(res);
00233 cFacet.NearestEdgeToPoint(res, fDist, uSide);
00234 if (fDist < fEps) {
00235 found = true;
00236 uIndex = *it;
00237 break;
00238 }
00239 }
00240 }
00241
00242 return found;
00243 }
00244
00245 float MeshAlgorithm::GetAverageEdgeLength() const
00246 {
00247 float fLen = 0.0f;
00248 MeshFacetIterator cF(_rclMesh);
00249 for (cF.Init(); cF.More(); cF.Next()) {
00250 for (int i=0; i<3; i++)
00251 fLen += Base::Distance(cF->_aclPoints[i], cF->_aclPoints[(i+1)%3]);
00252 }
00253
00254 fLen = fLen / (3.0f * _rclMesh.CountFacets() );
00255 return fLen;
00256 }
00257
00258 void MeshAlgorithm::GetMeshBorders (std::list<std::vector<Base::Vector3f> > &rclBorders) const
00259 {
00260 std::vector<unsigned long> aulAllFacets(_rclMesh.CountFacets());
00261 unsigned long k = 0;
00262 for (std::vector<unsigned long>::iterator pI = aulAllFacets.begin(); pI != aulAllFacets.end(); pI++)
00263 *pI = k++;
00264
00265 GetFacetBorders(aulAllFacets, rclBorders);
00266 }
00267
00268 void MeshAlgorithm::GetMeshBorders (std::list<std::vector<unsigned long> > &rclBorders) const
00269 {
00270 std::vector<unsigned long> aulAllFacets(_rclMesh.CountFacets());
00271 unsigned long k = 0;
00272 for (std::vector<unsigned long>::iterator pI = aulAllFacets.begin(); pI != aulAllFacets.end(); pI++)
00273 *pI = k++;
00274
00275 GetFacetBorders(aulAllFacets, rclBorders, true);
00276 }
00277
00278 void MeshAlgorithm::GetFacetBorders (const std::vector<unsigned long> &raulInd, std::list<std::vector<Base::Vector3f> > &rclBorders) const
00279 {
00280 #if 1
00281 const MeshPointArray &rclPAry = _rclMesh._aclPointArray;
00282 std::list<std::vector<unsigned long> > aulBorders;
00283
00284 GetFacetBorders (raulInd, aulBorders, true);
00285 for ( std::list<std::vector<unsigned long> >::iterator it = aulBorders.begin(); it != aulBorders.end(); ++it )
00286 {
00287 std::vector<Base::Vector3f> boundary;
00288 boundary.reserve( it->size() );
00289
00290 for ( std::vector<unsigned long>::iterator jt = it->begin(); jt != it->end(); ++jt )
00291 boundary.push_back(rclPAry[*jt]);
00292
00293 rclBorders.push_back( boundary );
00294 }
00295
00296 #else
00297 const MeshFacetArray &rclFAry = _rclMesh._aclFacetArray;
00298
00299
00300 ResetFacetFlag(MeshFacet::VISIT);
00301 for (std::vector<unsigned long>::const_iterator pIter = raulInd.begin(); pIter != raulInd.end(); pIter++)
00302 rclFAry[*pIter].SetFlag(MeshFacet::VISIT);
00303
00304 std::list<std::pair<unsigned long, unsigned long> > aclEdges;
00305
00306 for (std::vector<unsigned long>::const_iterator pIter2 = raulInd.begin(); pIter2 != raulInd.end(); pIter2++)
00307 {
00308 const MeshFacet &rclFacet = rclFAry[*pIter2];
00309 for (int i = 0; i < 3; i++)
00310 {
00311 unsigned long ulNB = rclFacet._aulNeighbours[i];
00312 if (ulNB != ULONG_MAX)
00313 {
00314 if (rclFAry[ulNB].IsFlag(MeshFacet::VISIT) == true)
00315 continue;
00316 }
00317
00318 aclEdges.push_back(rclFacet.GetEdge(i));
00319 }
00320 }
00321
00322 if (aclEdges.size() == 0)
00323 return;
00324
00325
00326 const MeshPointArray &rclPAry = _rclMesh._aclPointArray;
00327 unsigned long ulFirst, ulLast;
00328 std::list<Base::Vector3f> clBorder;
00329 ulFirst = aclEdges.begin()->first;
00330 ulLast = aclEdges.begin()->second;
00331
00332 aclEdges.erase(aclEdges.begin());
00333 clBorder.push_back(rclPAry[ulFirst]);
00334 clBorder.push_back(rclPAry[ulLast]);
00335
00336 while (aclEdges.size() > 0)
00337 {
00338
00339 std::list<std::pair<unsigned long, unsigned long> >::iterator pEI;
00340 for (pEI = aclEdges.begin(); pEI != aclEdges.end(); pEI++)
00341 {
00342 if (pEI->first == ulLast)
00343 {
00344 ulLast = pEI->second;
00345 clBorder.push_back(rclPAry[ulLast]);
00346 aclEdges.erase(pEI);
00347 break;
00348 }
00349 else if (pEI->second == ulLast)
00350 {
00351 ulLast = pEI->first;
00352 clBorder.push_back(rclPAry[ulLast]);
00353 aclEdges.erase(pEI);
00354 break;
00355 }
00356 else if (pEI->first == ulFirst)
00357 {
00358 ulFirst = pEI->second;
00359 clBorder.push_front(rclPAry[ulFirst]);
00360 aclEdges.erase(pEI);
00361 break;
00362 }
00363 else if (pEI->second == ulFirst)
00364 {
00365 ulFirst = pEI->first;
00366 clBorder.push_front(rclPAry[ulFirst]);
00367 aclEdges.erase(pEI);
00368 break;
00369 }
00370 }
00371 if ((pEI == aclEdges.end()) || (ulLast == ulFirst))
00372 {
00373 rclBorders.push_back(std::vector<Base::Vector3f>(clBorder.begin(), clBorder.end()));
00374 clBorder.clear();
00375
00376 if (aclEdges.size() > 0)
00377 {
00378 ulFirst = aclEdges.begin()->first;
00379 ulLast = aclEdges.begin()->second;
00380 aclEdges.erase(aclEdges.begin());
00381 clBorder.push_back(rclPAry[ulFirst]);
00382 clBorder.push_back(rclPAry[ulLast]);
00383 }
00384 }
00385 }
00386 #endif
00387 }
00388
00389 void MeshAlgorithm::GetFacetBorders (const std::vector<unsigned long> &raulInd,
00390 std::list<std::vector<unsigned long> > &rclBorders,
00391 bool ignoreOrientation) const
00392 {
00393 const MeshFacetArray &rclFAry = _rclMesh._aclFacetArray;
00394
00395
00396 ResetFacetFlag(MeshFacet::VISIT);
00397 for (std::vector<unsigned long>::const_iterator it = raulInd.begin(); it != raulInd.end(); ++it)
00398 rclFAry[*it].SetFlag(MeshFacet::VISIT);
00399
00400
00401 std::list<std::pair<unsigned long, unsigned long> > aclEdges;
00402 for (std::vector<unsigned long>::const_iterator it = raulInd.begin(); it != raulInd.end(); ++it) {
00403 const MeshFacet &rclFacet = rclFAry[*it];
00404 for (int i = 0; i < 3; i++) {
00405 unsigned long ulNB = rclFacet._aulNeighbours[i];
00406 if (ulNB != ULONG_MAX) {
00407 if (rclFAry[ulNB].IsFlag(MeshFacet::VISIT) == true)
00408 continue;
00409 }
00410
00411 aclEdges.push_back(rclFacet.GetEdge(i));
00412 }
00413 }
00414
00415 if (aclEdges.size() == 0)
00416 return;
00417
00418
00419 unsigned long ulFirst, ulLast;
00420 std::list<unsigned long> clBorder;
00421 ulFirst = aclEdges.begin()->first;
00422 ulLast = aclEdges.begin()->second;
00423
00424 aclEdges.erase(aclEdges.begin());
00425 clBorder.push_back(ulFirst);
00426 clBorder.push_back(ulLast);
00427
00428 while (aclEdges.size() > 0) {
00429
00430 std::list<std::pair<unsigned long, unsigned long> >::iterator pEI;
00431 for (pEI = aclEdges.begin(); pEI != aclEdges.end(); ++pEI) {
00432 if (pEI->first == ulLast) {
00433 ulLast = pEI->second;
00434 clBorder.push_back(ulLast);
00435 aclEdges.erase(pEI);
00436 pEI = aclEdges.begin();
00437 break;
00438 }
00439 else if (pEI->second == ulFirst) {
00440 ulFirst = pEI->first;
00441 clBorder.push_front(ulFirst);
00442 aclEdges.erase(pEI);
00443 pEI = aclEdges.begin();
00444 break;
00445 }
00446
00447
00448
00449 else if (pEI->second == ulLast && ignoreOrientation) {
00450 ulLast = pEI->first;
00451 clBorder.push_back(ulLast);
00452 aclEdges.erase(pEI);
00453 pEI = aclEdges.begin();
00454 break;
00455 }
00456 else if (pEI->first == ulFirst && ignoreOrientation) {
00457 ulFirst = pEI->second;
00458 clBorder.push_front(ulFirst);
00459 aclEdges.erase(pEI);
00460 pEI = aclEdges.begin();
00461 break;
00462 }
00463 }
00464
00465
00466
00467 if ((pEI == aclEdges.end()) || aclEdges.empty() || (ulLast == ulFirst)) {
00468
00469 rclBorders.push_back(std::vector<unsigned long>(clBorder.begin(), clBorder.end()));
00470 clBorder.clear();
00471
00472 if (aclEdges.size() > 0) {
00473
00474 ulFirst = aclEdges.begin()->first;
00475 ulLast = aclEdges.begin()->second;
00476 aclEdges.erase(aclEdges.begin());
00477 clBorder.push_back(ulFirst);
00478 clBorder.push_back(ulLast);
00479 }
00480 }
00481 }
00482 }
00483
00484 void MeshAlgorithm::GetMeshBorder(unsigned long uFacet, std::list<unsigned long>& rBorder) const
00485 {
00486 const MeshFacetArray &rFAry = _rclMesh._aclFacetArray;
00487 std::list<std::pair<unsigned long, unsigned long> > openEdges;
00488 if (uFacet >= rFAry.size())
00489 return;
00490
00491 MeshFacetArray::_TConstIterator face = rFAry.begin() + uFacet;
00492 for (int i = 0; i < 3; i++)
00493 {
00494 if (face->_aulNeighbours[i] == ULONG_MAX)
00495 openEdges.push_back(face->GetEdge(i));
00496 }
00497
00498 if (openEdges.empty())
00499 return;
00500
00501 for (MeshFacetArray::_TConstIterator it = rFAry.begin(); it != rFAry.end(); ++it)
00502 {
00503 if (it == face)
00504 continue;
00505 for (int i = 0; i < 3; i++)
00506 {
00507 if (it->_aulNeighbours[i] == ULONG_MAX)
00508 openEdges.push_back(it->GetEdge(i));
00509 }
00510 }
00511
00512
00513 unsigned long ulFirst = openEdges.begin()->first;
00514 unsigned long ulLast = openEdges.begin()->second;
00515
00516 openEdges.erase(openEdges.begin());
00517 rBorder.push_back(ulFirst);
00518 rBorder.push_back(ulLast);
00519
00520 while (ulLast != ulFirst)
00521 {
00522
00523 std::list<std::pair<unsigned long, unsigned long> >::iterator pEI;
00524 for (pEI = openEdges.begin(); pEI != openEdges.end(); pEI++)
00525 {
00526 if (pEI->first == ulLast)
00527 {
00528 ulLast = pEI->second;
00529 rBorder.push_back(ulLast);
00530 openEdges.erase(pEI);
00531 break;
00532 }
00533 else if (pEI->second == ulFirst)
00534 {
00535 ulFirst = pEI->first;
00536 rBorder.push_front(ulFirst);
00537 openEdges.erase(pEI);
00538 break;
00539 }
00540 }
00541
00542
00543 if (pEI == openEdges.end())
00544 break;
00545 }
00546 }
00547
00548 void MeshAlgorithm::SplitBoundaryLoops( std::list<std::vector<unsigned long> >& aBorders )
00549 {
00550
00551 std::map<unsigned long, int> openPointDegree;
00552 for (MeshFacetArray::_TConstIterator jt = _rclMesh._aclFacetArray.begin();
00553 jt != _rclMesh._aclFacetArray.end(); ++jt) {
00554 for (int i=0; i<3; i++) {
00555 if (jt->_aulNeighbours[i] == ULONG_MAX) {
00556 openPointDegree[jt->_aulPoints[i]]++;
00557 openPointDegree[jt->_aulPoints[(i+1)%3]]++;
00558 }
00559 }
00560 }
00561
00562
00563 std::list<std::vector<unsigned long> > aSplitBorders;
00564 for (std::list<std::vector<unsigned long> >::iterator it = aBorders.begin();
00565 it != aBorders.end(); ++it) {
00566 bool split=false;
00567 for (std::vector<unsigned long>::iterator jt = it->begin(); jt != it->end(); ++jt) {
00568
00569 if (openPointDegree[*jt] > 2) {
00570 split = true;
00571 break;
00572 }
00573 }
00574
00575 if (!split)
00576 aSplitBorders.push_back( *it );
00577 else
00578 SplitBoundaryLoops( *it, aSplitBorders );
00579 }
00580
00581 aBorders = aSplitBorders;
00582 }
00583
00584 void MeshAlgorithm::SplitBoundaryLoops(const std::vector<unsigned long>& rBound,
00585 std::list<std::vector<unsigned long> >& aBorders)
00586 {
00587 std::map<unsigned long, int> aPtDegree;
00588 std::vector<unsigned long> cBound;
00589 for (std::vector<unsigned long>::const_iterator it = rBound.begin(); it != rBound.end(); ++it) {
00590 int deg = (aPtDegree[*it]++);
00591 if (deg > 0) {
00592 for (std::vector<unsigned long>::iterator jt = cBound.begin(); jt != cBound.end(); ++jt) {
00593 if (*jt == *it) {
00594 std::vector<unsigned long> cBoundLoop;
00595 cBoundLoop.insert(cBoundLoop.end(), jt, cBound.end());
00596 cBoundLoop.push_back(*it);
00597 cBound.erase(jt, cBound.end());
00598 aBorders.push_back(cBoundLoop);
00599 (aPtDegree[*it]--);
00600 break;
00601 }
00602 }
00603 }
00604
00605 cBound.push_back(*it);
00606 }
00607 }
00608
00609 bool MeshAlgorithm::FillupHole(const std::vector<unsigned long>& boundary,
00610 AbstractPolygonTriangulator& cTria,
00611 MeshFacetArray& rFaces, MeshPointArray& rPoints,
00612 int level, const MeshRefPointToFacets* pP2FStructure) const
00613 {
00614 if (boundary.front() == boundary.back()) {
00615
00616 if (boundary.size() < 4)
00617 return false;
00618 }
00619 else if (boundary.size() < 3) {
00620 return false;
00621 }
00622
00623
00624 MeshGeomFacet rTriangle;
00625 MeshFacet rFace;
00626 unsigned long refPoint0 = *(boundary.begin());
00627 unsigned long refPoint1 = *(boundary.begin()+1);
00628 if (pP2FStructure) {
00629 const std::set<unsigned long>& ring = (*pP2FStructure)[refPoint0];
00630 bool found = false;
00631 for (std::set<unsigned long>::const_iterator it = ring.begin(); it != ring.end() && !found; ++it) {
00632 for (int i=0; i<3; i++) {
00633 if (pP2FStructure->getFacet(*it)->_aulPoints[i] == refPoint1) {
00634 rFace = *pP2FStructure->getFacet(*it);
00635 rTriangle = _rclMesh.GetFacet(rFace);
00636 found = true;
00637 break;
00638 }
00639 }
00640 }
00641 } else {
00642 bool ready = false;
00643 for (MeshFacetArray::_TConstIterator it = _rclMesh._aclFacetArray.begin(); it != _rclMesh._aclFacetArray.end(); ++it) {
00644 for (int i=0; i<3; i++) {
00645 if ((it->_aulPoints[i] == refPoint0) && (it->_aulPoints[(i+1)%3] == refPoint1)) {
00646 rFace = *it;
00647 rTriangle = _rclMesh.GetFacet(*it);
00648 ready = true;
00649 break;
00650 }
00651 }
00652
00653 if (ready)
00654 break;
00655 }
00656 }
00657
00658
00659 std::vector<Base::Vector3f> polygon;
00660 for (std::vector<unsigned long>::const_iterator jt = boundary.begin(); jt != boundary.end(); ++jt) {
00661 polygon.push_back(_rclMesh._aclPointArray[*jt]);
00662 rPoints.push_back(_rclMesh._aclPointArray[*jt]);
00663 }
00664
00665
00666 if (boundary.front() == boundary.back()) {
00667 polygon.pop_back();
00668 rPoints.pop_back();
00669 }
00670
00671
00672
00673
00674 cTria.SetPolygon(polygon);
00675
00676 std::vector<Base::Vector3f> surf_pts = cTria.GetPolygon();
00677 if (pP2FStructure && level > 0) {
00678 std::set<unsigned long> index = pP2FStructure->NeighbourPoints(boundary, level);
00679 for (std::set<unsigned long>::iterator it = index.begin(); it != index.end(); ++it) {
00680 Base::Vector3f pt(_rclMesh._aclPointArray[*it]);
00681 surf_pts.push_back(pt);
00682 }
00683 }
00684
00685 if (cTria.TriangulatePolygon()) {
00686
00687
00688 cTria.ProjectOntoSurface(surf_pts);
00689
00690 rFaces.insert(rFaces.end(), cTria.GetFacets().begin(), cTria.GetFacets().end());
00691 std::vector<Base::Vector3f> newVertices = cTria.AddedPoints();
00692 for (std::vector<Base::Vector3f>::iterator pt = newVertices.begin(); pt != newVertices.end(); ++pt) {
00693 rPoints.push_back((*pt));
00694 }
00695
00696
00697
00698
00699 std::vector<MeshFacet> faces = cTria.GetFacets();
00700
00701
00702
00703 if (faces.size()==1){
00704 MeshFacet first;
00705 first._aulPoints[0] = boundary[faces.front()._aulPoints[0]];
00706 first._aulPoints[1] = boundary[faces.front()._aulPoints[1]];
00707 first._aulPoints[2] = boundary[faces.front()._aulPoints[2]];
00708 if ( first.IsEqual(rFace) ) {
00709 rFaces.clear();
00710 rPoints.clear();
00711 cTria.Discard();
00712 return false;
00713 }
00714 }
00715
00716
00717 MeshFacet facet;
00718 unsigned short ref_side = rFace.Side(refPoint0, refPoint1);
00719 unsigned short tri_side = USHRT_MAX;
00720 if (ref_side < USHRT_MAX) {
00721 for (std::vector<MeshFacet>::iterator it = faces.begin(); it != faces.end(); ++it) {
00722 tri_side = it->Side(0, 1);
00723 if (tri_side < USHRT_MAX) {
00724 facet = *it;
00725 break;
00726 }
00727 }
00728 }
00729
00730
00731 if (ref_side == USHRT_MAX || tri_side == USHRT_MAX) {
00732 Base::Console().Log("MeshAlgorithm::FillupHole: Expected open edge for facet <%d, %d, %d>\n",
00733 rFace._aulPoints[0], rFace._aulPoints[1], rFace._aulPoints[2]);
00734 rFaces.clear();
00735 rPoints.clear();
00736 cTria.Discard();
00737 return false;
00738 }
00739
00740 MeshGeomFacet triangle;
00741 triangle._aclPoints[0] = rPoints[facet._aulPoints[0]];
00742 triangle._aclPoints[1] = rPoints[facet._aulPoints[1]];
00743 triangle._aclPoints[2] = rPoints[facet._aulPoints[2]];
00744
00745
00746
00747 Base::Vector3f planeNormal = rTriangle.GetNormal() % (rTriangle._aclPoints[(ref_side+1)%3]-rTriangle._aclPoints[ref_side]);
00748 Base::Vector3f planeBase = rTriangle._aclPoints[ref_side];
00749 Base::Vector3f ref_point = rTriangle._aclPoints[(ref_side+2)%3];
00750 Base::Vector3f tri_point = triangle._aclPoints[(tri_side+2)%3];
00751
00752 float ref_dist = (ref_point - planeBase) * planeNormal;
00753 float tri_dist = (tri_point - planeBase) * planeNormal;
00754 if (ref_dist * tri_dist > 0.0f) {
00755 rFaces.clear();
00756 rPoints.clear();
00757 cTria.Discard();
00758 return false;
00759 }
00760
00761
00762 if ( triangle.GetNormal() * rTriangle.GetNormal() <= 0.0f ) {
00763 for (MeshFacetArray::_TIterator it = rFaces.begin(); it != rFaces.end(); ++it)
00764 it->FlipNormal();
00765 }
00766
00767 return true;
00768 }
00769
00770 return false;
00771 }
00772
00773 void MeshAlgorithm::SetFacetsProperty(const std::vector<unsigned long> &raulInds, const std::vector<unsigned long> &raulProps) const
00774 {
00775 if (raulInds.size() != raulProps.size()) return;
00776
00777 std::vector<unsigned long>::const_iterator iP = raulProps.begin();
00778 for (std::vector<unsigned long>::const_iterator i = raulInds.begin(); i != raulInds.end(); i++, iP++)
00779 _rclMesh._aclFacetArray[*i].SetProperty(*iP);
00780 }
00781
00782 void MeshAlgorithm::SetFacetsFlag (const std::vector<unsigned long> &raulInds, MeshFacet::TFlagType tF) const
00783 {
00784 for (std::vector<unsigned long>::const_iterator i = raulInds.begin(); i != raulInds.end(); i++)
00785 _rclMesh._aclFacetArray[*i].SetFlag(tF);
00786 }
00787
00788 void MeshAlgorithm::SetPointsFlag (const std::vector<unsigned long> &raulInds, MeshPoint::TFlagType tF) const
00789 {
00790 for (std::vector<unsigned long>::const_iterator i = raulInds.begin(); i != raulInds.end(); i++)
00791 _rclMesh._aclPointArray[*i].SetFlag(tF);
00792 }
00793
00794 void MeshAlgorithm::GetFacetsFlag (std::vector<unsigned long> &raulInds, MeshFacet::TFlagType tF) const
00795 {
00796 raulInds.reserve(raulInds.size() + CountFacetFlag(tF));
00797 MeshFacetArray::_TConstIterator beg = _rclMesh._aclFacetArray.begin();
00798 MeshFacetArray::_TConstIterator end = _rclMesh._aclFacetArray.end();
00799 for (MeshFacetArray::_TConstIterator it = beg; it != end; ++it) {
00800 if (it->IsFlag(tF))
00801 raulInds.push_back(it-beg);
00802 }
00803 }
00804
00805 void MeshAlgorithm::GetPointsFlag (std::vector<unsigned long> &raulInds, MeshPoint::TFlagType tF) const
00806 {
00807 raulInds.reserve(raulInds.size() + CountPointFlag(tF));
00808 MeshPointArray::_TConstIterator beg = _rclMesh._aclPointArray.begin();
00809 MeshPointArray::_TConstIterator end = _rclMesh._aclPointArray.end();
00810 for (MeshPointArray::_TConstIterator it = beg; it != end; ++it) {
00811 if (it->IsFlag(tF))
00812 raulInds.push_back(it-beg);
00813 }
00814 }
00815
00816 void MeshAlgorithm::ResetFacetsFlag (const std::vector<unsigned long> &raulInds, MeshFacet::TFlagType tF) const
00817 {
00818 for (std::vector<unsigned long>::const_iterator i = raulInds.begin(); i != raulInds.end(); i++)
00819 _rclMesh._aclFacetArray[*i].ResetFlag(tF);
00820 }
00821
00822 void MeshAlgorithm::ResetPointsFlag (const std::vector<unsigned long> &raulInds, MeshPoint::TFlagType tF) const
00823 {
00824 for (std::vector<unsigned long>::const_iterator i = raulInds.begin(); i != raulInds.end(); i++)
00825 _rclMesh._aclPointArray[*i].ResetFlag(tF);
00826 }
00827
00828 void MeshAlgorithm::SetFacetFlag (MeshFacet::TFlagType tF) const
00829 {
00830 _rclMesh._aclFacetArray.SetFlag(tF);
00831 }
00832
00833 void MeshAlgorithm::SetPointFlag (MeshPoint::TFlagType tF) const
00834 {
00835 _rclMesh._aclPointArray.SetFlag(tF);
00836 }
00837
00838 void MeshAlgorithm::ResetFacetFlag (MeshFacet::TFlagType tF) const
00839 {
00840 _rclMesh._aclFacetArray.ResetFlag(tF);
00841 }
00842
00843 void MeshAlgorithm::ResetPointFlag (MeshPoint::TFlagType tF) const
00844 {
00845 _rclMesh._aclPointArray.ResetFlag(tF);
00846 }
00847
00848 unsigned long MeshAlgorithm::CountFacetFlag (MeshFacet::TFlagType tF) const
00849 {
00850 return std::count_if(_rclMesh._aclFacetArray.begin(), _rclMesh._aclFacetArray.end(),
00851 std::bind2nd(MeshIsFlag<MeshFacet>(), tF));
00852 }
00853
00854 unsigned long MeshAlgorithm::CountPointFlag (MeshPoint::TFlagType tF) const
00855 {
00856 return std::count_if(_rclMesh._aclPointArray.begin(), _rclMesh._aclPointArray.end(),
00857 std::bind2nd(MeshIsFlag<MeshPoint>(), tF));
00858 }
00859
00860 void MeshAlgorithm::GetFacetsFromToolMesh( const MeshKernel& rToolMesh, const Base::Vector3f& rcDir, std::vector<unsigned long> &raclCutted ) const
00861 {
00862 MeshFacetIterator cFIt(_rclMesh);
00863 MeshFacetIterator cTIt(rToolMesh);
00864
00865 BoundBox3f cBB = rToolMesh.GetBoundBox();
00866
00867 Base::SequencerLauncher seq("Check facets...", _rclMesh.CountFacets());
00868
00869
00870 Base::Vector3f tmp;
00871 for (cFIt.Init(); cFIt.More(); cFIt.Next()) {
00872
00873 for (int i=0; i<3; i++) {
00874
00875 if (cBB.IsInBox( cFIt->_aclPoints[i])) {
00876
00877 int ct=0;
00878 for (cTIt.Init(); cTIt.More(); cTIt.Next()) {
00879 if (cTIt->IsPointOfFace( cFIt->_aclPoints[i], MeshPoint::epsilon())) {
00880 ct=1;
00881 break;
00882 }
00883 else if (cTIt->Foraminate( cFIt->_aclPoints[i], rcDir, tmp)) {
00884
00885 if ((tmp - cFIt->_aclPoints[i]) * rcDir > 0)
00886 ct++;
00887 }
00888 }
00889
00890
00891 if (ct % 2 == 1) {
00892 raclCutted.push_back( cFIt.Position() );
00893 break;
00894 }
00895 }
00896 }
00897
00898 seq.next();
00899 }
00900 }
00901
00902 void MeshAlgorithm::GetFacetsFromToolMesh(const MeshKernel& rToolMesh, const Base::Vector3f& rcDir,
00903 const MeshFacetGrid& rGrid, std::vector<unsigned long> &raclCutted) const
00904 {
00905
00906 MeshGridIterator clGridIter(rGrid);
00907 BoundBox3f cBB = rToolMesh.GetBoundBox();
00908 Base::Vector3f tmp;
00909
00910 MeshFacetIterator cFIt(_rclMesh);
00911 MeshFacetIterator cTIt(rToolMesh);
00912 MeshAlgorithm cToolAlg(rToolMesh);
00913
00914
00915
00916
00917
00918
00919
00920
00921 std::vector<unsigned long> aulInds;
00922 for (clGridIter.Init(); clGridIter.More(); clGridIter.Next()) {
00923 int ret = cToolAlg.Surround(clGridIter.GetBoundBox(), rcDir);
00924
00925
00926 if (ret == 1) {
00927
00928 clGridIter.GetElements(raclCutted);
00929 }
00930
00931 else if (ret == 0) {
00932
00933 clGridIter.GetElements(aulInds);
00934 }
00935
00936
00937 else if (ret == -1) {
00938
00939 clGridIter.GetElements(aulInds);
00940 }
00941 }
00942
00943
00944 std::sort(aulInds.begin(), aulInds.end());
00945 aulInds.erase(std::unique(aulInds.begin(), aulInds.end()), aulInds.end());
00946 std::sort(raclCutted.begin(), raclCutted.end());
00947 raclCutted.erase(std::unique(raclCutted.begin(), raclCutted.end()), raclCutted.end());
00948
00949 Base::SequencerLauncher seq("Check facets...", aulInds.size());
00950
00951
00952 for (std::vector<unsigned long>::iterator it = aulInds.begin(); it != aulInds.end(); ++it) {
00953 cFIt.Set(*it);
00954
00955
00956 for (int i=0; i<3; i++) {
00957
00958 if (cBB.IsInBox(cFIt->_aclPoints[i])) {
00959
00960 int ct=0;
00961 for (cTIt.Init(); cTIt.More(); cTIt.Next()) {
00962 if (cTIt->IsPointOfFace(cFIt->_aclPoints[i], MeshPoint::epsilon())) {
00963 ct=1;
00964 break;
00965 }
00966 else if (cTIt->Foraminate(cFIt->_aclPoints[i], rcDir, tmp)) {
00967
00968 if ((tmp - cFIt->_aclPoints[i]) * rcDir > 0)
00969 ct++;
00970 }
00971 }
00972
00973
00974 if (ct % 2 == 1) {
00975 raclCutted.push_back(cFIt.Position());
00976 break;
00977 }
00978 }
00979 }
00980
00981 seq.next();
00982 }
00983
00984
00985 std::sort(raclCutted.begin(), raclCutted.end());
00986 raclCutted.erase(std::unique(raclCutted.begin(), raclCutted.end()), raclCutted.end());
00987 }
00988
00989 int MeshAlgorithm::Surround(const Base::BoundBox3f& rBox, const Base::Vector3f& rcDir)
00990 {
00991 Base::Vector3f pt1, pt2, tmp;
00992 const BoundBox3f& cBB = _rclMesh.GetBoundBox();
00993
00994
00995 if (cBB && rBox) {
00996
00997 Base::Vector3f cCorner[8] = {
00998 Base::Vector3f(rBox.MinX,rBox.MinY,rBox.MinZ), Base::Vector3f(rBox.MaxX,rBox.MinY,rBox.MinZ),
00999 Base::Vector3f(rBox.MaxX,rBox.MaxY,rBox.MinZ), Base::Vector3f(rBox.MinX,rBox.MaxY,rBox.MinZ),
01000 Base::Vector3f(rBox.MinX,rBox.MinY,rBox.MaxZ), Base::Vector3f(rBox.MaxX,rBox.MinY,rBox.MaxZ),
01001 Base::Vector3f(rBox.MaxX,rBox.MaxY,rBox.MaxZ), Base::Vector3f(rBox.MinX,rBox.MaxY,rBox.MaxZ)};
01002
01003 MeshFacetIterator cTIt(_rclMesh);
01004
01005
01006 int triangles[36] = {
01007 0,1,2,0,2,3,
01008 0,1,5,0,5,4,
01009 0,4,7,0,7,3,
01010 6,7,4,6,4,5,
01011 6,2,3,6,3,7,
01012 6,1,2,6,5,1
01013 };
01014
01015 std::vector<MeshGeomFacet> cFacet(12);
01016 int id=0;
01017 for (int ii=0; ii<12; ii++) {
01018 cFacet[ii]._aclPoints[0]=cCorner[triangles[id++]];
01019 cFacet[ii]._aclPoints[1]=cCorner[triangles[id++]];
01020 cFacet[ii]._aclPoints[2]=cCorner[triangles[id++]];
01021 }
01022
01023
01024 for (std::vector<MeshGeomFacet>::iterator it = cFacet.begin(); it != cFacet.end(); ++it) {
01025 for (cTIt.Init(); cTIt.More(); cTIt.Next()) {
01026 int ret = cTIt->IntersectWithFacet(*it,pt1, pt2);
01027
01028
01029 if (ret != 0)
01030 return 0;
01031 }
01032 }
01033
01034
01035
01036
01037 int ct=0;
01038 for (cTIt.Init(); cTIt.More(); cTIt.Next()) {
01039 if (cTIt->IsPointOfFace(cCorner[0], MeshPoint::epsilon())) {
01040 ct=1;
01041 break;
01042 }
01043 else if (cTIt->Foraminate(cCorner[0], rcDir, tmp)) {
01044
01045 if ((tmp - cCorner[0]) * rcDir > 0)
01046 ct++;
01047 }
01048 }
01049
01050
01051 return (ct % 2 == 1) ? 1 : -1;
01052 }
01053
01054
01055 return -1;
01056 }
01057
01058 void MeshAlgorithm::CheckFacets(const MeshFacetGrid& rclGrid, const Base::ViewProjMethod* pclProj, const Base::Polygon2D& rclPoly,
01059 bool bInner, std::vector<unsigned long> &raulFacets) const
01060 {
01061 std::vector<unsigned long>::iterator it;
01062 MeshFacetIterator clIter(_rclMesh, 0);
01063 Base::Vector3f clPt2d;
01064 Base::Vector3f clGravityOfFacet;
01065 bool bNoPointInside;
01066
01067
01068 if (bInner)
01069 {
01070 BoundBox3f clBBox3d;
01071 BoundBox2D clViewBBox, clPolyBBox;
01072 std::vector<unsigned long> aulAllElements;
01073
01074
01075 clPolyBBox = rclPoly.CalcBoundBox();
01076
01077 MeshGridIterator clGridIter(rclGrid);
01078
01079 for (clGridIter.Init(); clGridIter.More(); clGridIter.Next())
01080 {
01081 clBBox3d = clGridIter.GetBoundBox();
01082 clViewBBox = clBBox3d.ProjectBox(pclProj);
01083 if (clViewBBox || clPolyBBox)
01084 {
01085
01086 clGridIter.GetElements(aulAllElements);
01087 }
01088 }
01089
01090 std::sort(aulAllElements.begin(), aulAllElements.end());
01091 aulAllElements.erase(std::unique(aulAllElements.begin(), aulAllElements.end()), aulAllElements.end());
01092
01093 Base::SequencerLauncher seq( "Check facets", aulAllElements.size() );
01094
01095 for (it = aulAllElements.begin(); it != aulAllElements.end(); it++)
01096 {
01097 bNoPointInside = true;
01098 clGravityOfFacet.Set(0.0f, 0.0f, 0.0f);
01099 MeshGeomFacet rclFacet = _rclMesh.GetFacet(*it);
01100 for (int j=0; j<3; j++)
01101 {
01102 clPt2d = pclProj->operator()(rclFacet._aclPoints[j]);
01103 clGravityOfFacet += clPt2d;
01104 if (rclPoly.Contains(Base::Vector2D(clPt2d.x, clPt2d.y)) == bInner)
01105 {
01106 raulFacets.push_back(*it);
01107 bNoPointInside = false;
01108 break;
01109 }
01110 }
01111
01112
01113 if (bNoPointInside == true)
01114 {
01115 clGravityOfFacet *= 1.0f/3.0f;
01116
01117 if (rclPoly.Contains(Base::Vector2D(clGravityOfFacet.x, clGravityOfFacet.y)) == bInner)
01118 raulFacets.push_back(*it);
01119 }
01120
01121 seq.next();
01122 }
01123 }
01124
01125 else
01126 {
01127 Base::SequencerLauncher seq("Check facets", _rclMesh.CountFacets());
01128 for (clIter.Init(); clIter.More(); clIter.Next())
01129 {
01130 for (int j=0; j<3; j++)
01131 {
01132 clPt2d = pclProj->operator()(clIter->_aclPoints[j]);
01133 if (rclPoly.Contains(Base::Vector2D(clPt2d.x, clPt2d.y)) == bInner)
01134 {
01135 raulFacets.push_back(clIter.Position());
01136 break;
01137 }
01138 }
01139 seq.next();
01140 }
01141 }
01142 }
01143
01144 void MeshAlgorithm::CheckFacets(const Base::ViewProjMethod* pclProj, const Base::Polygon2D& rclPoly,
01145 bool bInner, std::vector<unsigned long> &raulFacets) const
01146 {
01147 const MeshPointArray& p = _rclMesh.GetPoints();
01148 const MeshFacetArray& f = _rclMesh.GetFacets();
01149 Base::SequencerLauncher seq("Check facets", f.size());
01150 Base::Vector3f pt2d;
01151 unsigned long index=0;
01152 for (MeshFacetArray::_TConstIterator it = f.begin(); it != f.end(); ++it,++index) {
01153 for (int i = 0; i < 3; i++) {
01154 pt2d = (*pclProj)(p[it->_aulPoints[i]]);
01155 if (rclPoly.Contains(Base::Vector2D(pt2d.x, pt2d.y)) == bInner) {
01156 raulFacets.push_back(index);
01157 break;
01158 }
01159 }
01160 seq.next();
01161 }
01162 }
01163
01164 float MeshAlgorithm::Surface (void) const
01165 {
01166 float fTotal = 0.0f;
01167 MeshFacetIterator clFIter(_rclMesh);
01168
01169 for (clFIter.Init(); clFIter.More(); clFIter.Next())
01170 fTotal += clFIter->Area();
01171
01172 return fTotal;
01173 }
01174
01175 void MeshAlgorithm::SubSampleByDist (float fDist, std::vector<Base::Vector3f> &rclPoints) const
01176 {
01177 rclPoints.clear();
01178 MeshFacetIterator clFIter(_rclMesh);
01179 for (clFIter.Init(); clFIter.More(); clFIter.Next()) {
01180 size_t k = rclPoints.size();
01181 clFIter->SubSample(fDist, rclPoints);
01182 if (rclPoints.size() == k)
01183 rclPoints.push_back(clFIter->GetGravityPoint());
01184 }
01185 }
01186
01187 void MeshAlgorithm::SubSampleAllPoints (std::vector<Base::Vector3f> &rclPoints) const
01188 {
01189 rclPoints.clear();
01190
01191
01192
01193 MeshPointIterator clPIter(_rclMesh);
01194 for (clPIter.Init(); clPIter.More(); clPIter.Next()) {
01195 rclPoints.push_back(*clPIter);
01196 }
01197 }
01198
01199 void MeshAlgorithm::SubSampleByCount (unsigned long ulCtPoints, std::vector<Base::Vector3f> &rclPoints) const
01200 {
01201 float fDist = float(sqrt(Surface() / float(ulCtPoints)));
01202 SubSampleByDist(fDist, rclPoints);
01203 }
01204
01205 void MeshAlgorithm::SearchFacetsFromPolyline (const std::vector<Base::Vector3f> &rclPolyline, float fRadius,
01206 const MeshFacetGrid& rclGrid, std::vector<unsigned long> &rclResultFacetsIndices) const
01207 {
01208 rclResultFacetsIndices.clear();
01209 if ( rclPolyline.size() < 3 )
01210 return;
01211
01212 std::set<unsigned long> aclFacets;
01213 for (std::vector<Base::Vector3f>::const_iterator pV = rclPolyline.begin(); pV < (rclPolyline.end() - 1); pV++)
01214 {
01215 const Base::Vector3f &rclP0 = *pV, &rclP1 = *(pV + 1);
01216
01217
01218 BoundBox3f clSegmBB(rclP0.x, rclP0.y, rclP0.z, rclP0.x, rclP0.y, rclP0.z);
01219 clSegmBB &= rclP1;
01220 clSegmBB.Enlarge(fRadius);
01221
01222 std::vector<unsigned long> aclBBFacets;
01223 unsigned long k = rclGrid.Inside(clSegmBB, aclBBFacets, false);
01224 for (unsigned long i = 0; i < k; i++)
01225 {
01226 if (_rclMesh.GetFacet(aclBBFacets[i]).DistanceToLineSegment(rclP0, rclP1) < fRadius)
01227 aclFacets.insert(aclBBFacets[i]);
01228 }
01229 }
01230
01231 rclResultFacetsIndices.insert(rclResultFacetsIndices.begin(), aclFacets.begin(), aclFacets.end());
01232 }
01233
01234 void MeshAlgorithm::CutBorderFacets (std::vector<unsigned long> &raclFacetIndices, unsigned short usLevel) const
01235 {
01236 std::vector<unsigned long> aclToDelete;
01237
01238 CheckBorderFacets(raclFacetIndices, aclToDelete, usLevel);
01239
01240
01241 std::vector<unsigned long> aclResult;
01242 std::set<unsigned long> aclTmp(aclToDelete.begin(), aclToDelete.end());
01243
01244 for (std::vector<unsigned long>::iterator pI = raclFacetIndices.begin(); pI != raclFacetIndices.end(); pI++)
01245 {
01246 if (aclTmp.find(*pI) == aclTmp.end())
01247 aclResult.push_back(*pI);
01248 }
01249
01250 raclFacetIndices = aclResult;
01251 }
01252
01253 unsigned long MeshAlgorithm::CountBorderEdges() const
01254 {
01255 unsigned long cnt=0;
01256 const MeshFacetArray &rclFAry = _rclMesh._aclFacetArray;
01257 MeshFacetArray::_TConstIterator end = rclFAry.end();
01258 for (MeshFacetArray::_TConstIterator it = rclFAry.begin(); it != end; ++it) {
01259 for (int i=0; i<3; i++) {
01260 if (it->_aulNeighbours[i] == ULONG_MAX)
01261 cnt++;
01262 }
01263 }
01264
01265 return cnt;
01266 }
01267
01268 void MeshAlgorithm::CheckBorderFacets (const std::vector<unsigned long> &raclFacetIndices, std::vector<unsigned long> &raclResultIndices, unsigned short usLevel) const
01269 {
01270 ResetFacetFlag(MeshFacet::TMP0);
01271 SetFacetsFlag(raclFacetIndices, MeshFacet::TMP0);
01272
01273 const MeshFacetArray &rclFAry = _rclMesh._aclFacetArray;
01274
01275 for (unsigned short usL = 0; usL < usLevel; usL++)
01276 {
01277 for (std::vector<unsigned long>::const_iterator pF = raclFacetIndices.begin(); pF != raclFacetIndices.end(); pF++)
01278 {
01279 for (int i = 0; i < 3; i++)
01280 {
01281 unsigned long ulNB = rclFAry[*pF]._aulNeighbours[i];
01282 if (ulNB == ULONG_MAX)
01283 {
01284 raclResultIndices.push_back(*pF);
01285 rclFAry[*pF].ResetFlag(MeshFacet::TMP0);
01286 continue;
01287 }
01288 if (rclFAry[ulNB].IsFlag(MeshFacet::TMP0) == false)
01289 {
01290 raclResultIndices.push_back(*pF);
01291 rclFAry[*pF].ResetFlag(MeshFacet::TMP0);
01292 continue;
01293 }
01294 }
01295 }
01296 }
01297 }
01298
01299 void MeshAlgorithm::GetBorderPoints (const std::vector<unsigned long> &raclFacetIndices, std::set<unsigned long> &raclResultPointsIndices) const
01300 {
01301 ResetFacetFlag(MeshFacet::TMP0);
01302 SetFacetsFlag(raclFacetIndices, MeshFacet::TMP0);
01303
01304 const MeshFacetArray &rclFAry = _rclMesh._aclFacetArray;
01305
01306 for (std::vector<unsigned long>::const_iterator pF = raclFacetIndices.begin(); pF != raclFacetIndices.end(); pF++)
01307 {
01308 for (int i = 0; i < 3; i++)
01309 {
01310 const MeshFacet &rclFacet = rclFAry[*pF];
01311 unsigned long ulNB = rclFacet._aulNeighbours[i];
01312 if (ulNB == ULONG_MAX)
01313 {
01314 raclResultPointsIndices.insert(rclFacet._aulPoints[i]);
01315 raclResultPointsIndices.insert(rclFacet._aulPoints[(i+1)%3]);
01316 continue;
01317 }
01318 if (rclFAry[ulNB].IsFlag(MeshFacet::TMP0) == false)
01319 {
01320 raclResultPointsIndices.insert(rclFacet._aulPoints[i]);
01321 raclResultPointsIndices.insert(rclFacet._aulPoints[(i+1)%3]);
01322 continue;
01323 }
01324 }
01325 }
01326 }
01327
01328 bool MeshAlgorithm::NearestPointFromPoint (const Base::Vector3f &rclPt, unsigned long &rclResFacetIndex, Base::Vector3f &rclResPoint) const
01329 {
01330 if (_rclMesh.CountFacets() == 0)
01331 return false;
01332
01333
01334 float fMinDist = FLOAT_MAX;
01335 unsigned long ulInd = ULONG_MAX;
01336 MeshFacetIterator pF(_rclMesh);
01337 for (pF.Init(); pF.More(); pF.Next())
01338 {
01339 float fDist = pF->DistanceToPoint(rclPt);
01340 if (fDist < fMinDist)
01341 {
01342 fMinDist = fDist;
01343 ulInd = pF.Position();
01344 }
01345 }
01346
01347 MeshGeomFacet rclSFacet = _rclMesh.GetFacet(ulInd);
01348 rclSFacet.DistanceToPoint(rclPt, rclResPoint);
01349 rclResFacetIndex = ulInd;
01350
01351 return true;
01352 }
01353
01354 bool MeshAlgorithm::NearestPointFromPoint (const Base::Vector3f &rclPt, const MeshFacetGrid& rclGrid, unsigned long &rclResFacetIndex, Base::Vector3f &rclResPoint) const
01355 {
01356 unsigned long ulInd = rclGrid.SearchNearestFromPoint(rclPt);
01357
01358 if (ulInd == ULONG_MAX)
01359 {
01360 return false;
01361 }
01362
01363 MeshGeomFacet rclSFacet = _rclMesh.GetFacet(ulInd);
01364 rclSFacet.DistanceToPoint(rclPt, rclResPoint);
01365 rclResFacetIndex = ulInd;
01366
01367 return true;
01368 }
01369
01370 bool MeshAlgorithm::NearestPointFromPoint (const Base::Vector3f &rclPt, const MeshFacetGrid& rclGrid, float fMaxSearchArea,
01371 unsigned long &rclResFacetIndex, Base::Vector3f &rclResPoint) const
01372 {
01373 unsigned long ulInd = rclGrid.SearchNearestFromPoint(rclPt, fMaxSearchArea);
01374
01375 if (ulInd == ULONG_MAX)
01376 return false;
01377
01378 MeshGeomFacet rclSFacet = _rclMesh.GetFacet(ulInd);
01379 rclSFacet.DistanceToPoint(rclPt, rclResPoint);
01380 rclResFacetIndex = ulInd;
01381
01382 return true;
01383 }
01384
01385 bool MeshAlgorithm::CutWithPlane (const Base::Vector3f &clBase, const Base::Vector3f &clNormal, const MeshFacetGrid &rclGrid,
01386 std::list<std::vector<Base::Vector3f> > &rclResult, float fMinEps, bool bConnectPolygons) const
01387 {
01388 std::vector<unsigned long> aulFacets;
01389
01390
01391 MeshGridIterator clGridIter(rclGrid);
01392 for (clGridIter.Init(); clGridIter.More(); clGridIter.Next())
01393 {
01394
01395 if (clGridIter.GetBoundBox().IsCutPlane(clBase, clNormal) == true)
01396 clGridIter.GetElements(aulFacets);
01397 }
01398
01399
01400 std::sort(aulFacets.begin(), aulFacets.end());
01401 aulFacets.erase(std::unique(aulFacets.begin(), aulFacets.end()), aulFacets.end());
01402
01403
01404 std::list<std::pair<Base::Vector3f, Base::Vector3f> > clTempPoly;
01405
01406 for (std::vector<unsigned long>::iterator pF = aulFacets.begin(); pF != aulFacets.end(); pF++)
01407 {
01408 Base::Vector3f clE1, clE2;
01409 const MeshGeomFacet clF(_rclMesh.GetFacet(*pF));
01410
01411
01412 if (clF.IntersectWithPlane(clBase, clNormal, clE1, clE2) == true)
01413 clTempPoly.push_back(std::pair<Base::Vector3f, Base::Vector3f>(clE1, clE2));
01414 }
01415
01416 if(bConnectPolygons)
01417 {
01418
01419 std::list<std::pair<Base::Vector3f, Base::Vector3f> > rclResultLines(clTempPoly.begin(),clTempPoly.end());
01420 std::list<std::vector<Base::Vector3f> > tempList;
01421 ConnectLines(clTempPoly, tempList, fMinEps);
01422 ConnectPolygons(tempList, clTempPoly);
01423
01424 for(std::list<std::pair<Base::Vector3f, Base::Vector3f> >::iterator iter = clTempPoly.begin(); iter != clTempPoly.end(); iter++)
01425 {
01426 rclResultLines.push_front(*iter);
01427 }
01428
01429 return ConnectLines(rclResultLines, rclResult, fMinEps);
01430 }
01431
01432 return ConnectLines(clTempPoly, rclResult, fMinEps);
01433 }
01434
01435 bool MeshAlgorithm::ConnectLines (std::list<std::pair<Base::Vector3f, Base::Vector3f> > &rclLines,
01436 std::list<std::vector<Base::Vector3f> > &rclPolylines, float fMinEps) const
01437 {
01438 typedef std::list<std::pair<Base::Vector3f, Base::Vector3f> >::iterator TCIter;
01439
01440
01441
01442 fMinEps = fMinEps * fMinEps;
01443
01444
01445 std::list<TCIter> _clToDelete;
01446 float fToDelDist = fMinEps / 10.0f;
01447 for (TCIter pF = rclLines.begin(); pF != rclLines.end(); pF++)
01448 {
01449 if (Base::DistanceP2(pF->first, pF->second) < fToDelDist)
01450 _clToDelete.push_back(pF);
01451 }
01452 for (std::list<TCIter>::iterator pI = _clToDelete.begin(); pI != _clToDelete.end(); pI++)
01453 rclLines.erase(*pI);
01454
01455
01456 while (rclLines.size() > 0)
01457 {
01458 TCIter pF;
01459
01460
01461 std::list<Base::Vector3f> clPoly;
01462
01463
01464 Base::Vector3f clFront = rclLines.begin()->first;
01465 Base::Vector3f clEnd = rclLines.begin()->second;
01466 clPoly.push_back(clFront);
01467 clPoly.push_back(clEnd);
01468 rclLines.erase(rclLines.begin());
01469
01470
01471 TCIter pFront, pEnd;
01472 do
01473 {
01474 float fFrontMin = fMinEps, fEndMin = fMinEps;
01475 bool bFrontFirst=false, bEndFirst=false;
01476
01477 pFront = rclLines.end();
01478 pEnd = rclLines.end();
01479 for (pF = rclLines.begin(); pF != rclLines.end(); pF++)
01480 {
01481 if (Base::DistanceP2(clFront, pF->first) < fFrontMin)
01482 {
01483 fFrontMin = Base::DistanceP2(clFront, pF->first);
01484 pFront = pF;
01485 bFrontFirst = true;
01486 }
01487 else if (Base::DistanceP2(clEnd, pF->first) < fEndMin)
01488 {
01489 fEndMin = Base::DistanceP2(clEnd, pF->first);
01490 pEnd = pF;
01491 bEndFirst = true;
01492 }
01493 else if (Base::DistanceP2(clFront, pF->second) < fFrontMin)
01494 {
01495 fFrontMin = Base::DistanceP2(clFront, pF->second);
01496 pFront = pF;
01497 bFrontFirst = false;
01498 }
01499 else if (Base::DistanceP2(clEnd, pF->second) < fEndMin)
01500 {
01501 fEndMin = Base::DistanceP2(clEnd, pF->second);
01502 pEnd = pF;
01503 bEndFirst = false;
01504 }
01505 }
01506
01507 if (pFront != rclLines.end())
01508 {
01509 if (bFrontFirst == true)
01510 {
01511 clPoly.push_front(pFront->second);
01512 clFront = pFront->second;
01513 }
01514 else
01515 {
01516 clPoly.push_front(pFront->first);
01517 clFront = pFront->first;
01518 }
01519 rclLines.erase(pFront);
01520 }
01521
01522 if (pEnd != rclLines.end())
01523 {
01524 if (bEndFirst == true)
01525 {
01526 clPoly.push_back(pEnd->second);
01527 clEnd = pEnd->second;
01528 }
01529 else
01530 {
01531 clPoly.push_back(pEnd->first);
01532 clEnd = pEnd->first;
01533 }
01534 rclLines.erase(pEnd);
01535 }
01536 }
01537 while ((pFront != rclLines.end()) || (pEnd != rclLines.end()));
01538
01539 rclPolylines.push_back(std::vector<Base::Vector3f>(clPoly.begin(), clPoly.end()));
01540 }
01541
01542
01543 typedef std::list<std::vector<Base::Vector3f> >::iterator TPIter;
01544 std::list<TPIter> _clPolyToDelete;
01545 for (TPIter pJ = rclPolylines.begin(); pJ != rclPolylines.end(); pJ++)
01546 {
01547 if (pJ->size() == 2)
01548 {
01549 if (Base::DistanceP2(*pJ->begin(), *(pJ->begin() + 1)) <= fMinEps)
01550 _clPolyToDelete.push_back(pJ);
01551 }
01552 }
01553 for (std::list<TPIter>::iterator pK = _clPolyToDelete.begin(); pK != _clPolyToDelete.end(); pK++)
01554 rclPolylines.erase(*pK);
01555
01556 return true;
01557 }
01558
01559 bool MeshAlgorithm::ConnectPolygons(std::list<std::vector<Base::Vector3f> > &clPolyList,
01560 std::list<std::pair<Base::Vector3f, Base::Vector3f> > &rclLines) const
01561 {
01562
01563 for(std::list< std::vector<Base::Vector3f> >::iterator OutIter = clPolyList.begin(); OutIter != clPolyList.end() ; OutIter++)
01564 {
01565 std::pair<Base::Vector3f,Base::Vector3f> currentSort;
01566 float fDist = Base::Distance(OutIter->front(),OutIter->back());
01567 currentSort.first = OutIter->front();
01568 currentSort.second = OutIter->back();
01569
01570 for(std::list< std::vector<Base::Vector3f> >::iterator InnerIter = clPolyList.begin(); InnerIter != clPolyList.end(); InnerIter++)
01571 {
01572 if(OutIter == InnerIter) continue;
01573
01574 if(Base::Distance(OutIter->front(),InnerIter->front()) < fDist)
01575 {
01576 currentSort.second = InnerIter->front();
01577 fDist = Base::Distance(OutIter->front(),InnerIter->front());
01578 }
01579 if(Base::Distance(OutIter->front(),InnerIter->back()) < fDist)
01580 {
01581 currentSort.second = InnerIter->back();
01582 fDist = Base::Distance(OutIter->front(),InnerIter->back());
01583 }
01584 }
01585
01586 rclLines.push_front(currentSort);
01587
01588 }
01589
01590 return true;
01591 }
01592
01593 void MeshAlgorithm::GetFacetsFromPlane (const MeshFacetGrid &rclGrid, const Base::Vector3f& clNormal, float d, const Base::Vector3f &rclLeft,
01594 const Base::Vector3f &rclRight, std::vector<unsigned long> &rclRes) const
01595 {
01596 std::vector<unsigned long> aulFacets;
01597
01598 Base::Vector3f clBase = d * clNormal;
01599
01600 Base::Vector3f clPtNormal(rclLeft - rclRight);
01601 clPtNormal.Normalize();
01602
01603
01604 MeshGridIterator clGridIter(rclGrid);
01605 for (clGridIter.Init(); clGridIter.More(); clGridIter.Next())
01606 {
01607
01608 if (clGridIter.GetBoundBox().IsCutPlane(clBase, clNormal) == true)
01609 clGridIter.GetElements(aulFacets);
01610 }
01611
01612
01613 for (std::vector<unsigned long>::iterator pI = aulFacets.begin(); pI != aulFacets.end(); pI++)
01614 {
01615 MeshGeomFacet clSFacet = _rclMesh.GetFacet(*pI);
01616 if (clSFacet.IntersectWithPlane(clBase, clNormal) == true)
01617 {
01618 bool bInner = false;
01619 for (int i = 0; (i < 3) && (bInner == false); i++)
01620 {
01621 Base::Vector3f clPt = clSFacet._aclPoints[i];
01622 if ((clPt.DistanceToPlane(rclLeft, clPtNormal) <= 0.0f) && (clPt.DistanceToPlane(rclRight, clPtNormal) >= 0.0f))
01623 bInner = true;
01624 }
01625
01626 if (bInner == true)
01627 rclRes.push_back(*pI);
01628 }
01629 }
01630 }
01631
01632 void MeshAlgorithm::PointsFromFacetsIndices (const std::vector<unsigned long> &rvecIndices, std::vector<Base::Vector3f> &rvecPoints) const
01633 {
01634 const MeshFacetArray &rclFAry = _rclMesh._aclFacetArray;
01635 const MeshPointArray &rclPAry = _rclMesh._aclPointArray;
01636
01637 std::set<unsigned long> setPoints;
01638
01639 for (std::vector<unsigned long>::const_iterator itI = rvecIndices.begin(); itI != rvecIndices.end(); itI++)
01640 {
01641 for (int i = 0; i < 3; i++)
01642 setPoints.insert(rclFAry[*itI]._aulPoints[i]);
01643 }
01644
01645 rvecPoints.clear();
01646 for (std::set<unsigned long>::iterator itP = setPoints.begin(); itP != setPoints.end(); itP++)
01647 rvecPoints.push_back(rclPAry[*itP]);
01648 }
01649
01650 bool MeshAlgorithm::Distance (const Base::Vector3f &rclPt, unsigned long ulFacetIdx, float fMaxDistance, float &rfDistance) const
01651 {
01652 const MeshFacetArray &rclFAry = _rclMesh._aclFacetArray;
01653 const MeshPointArray &rclPAry = _rclMesh._aclPointArray;
01654 const unsigned long *pulIdx = rclFAry[ulFacetIdx]._aulPoints;
01655
01656 BoundBox3f clBB;
01657 clBB &= rclPAry[*(pulIdx++)];
01658 clBB &= rclPAry[*(pulIdx++)];
01659 clBB &= rclPAry[*pulIdx];
01660 clBB.Enlarge(fMaxDistance);
01661
01662 if (clBB.IsInBox(rclPt) == false)
01663 return false;
01664
01665 rfDistance = _rclMesh.GetFacet(ulFacetIdx).DistanceToPoint(rclPt);
01666
01667 return rfDistance < fMaxDistance;
01668 }
01669
01670 float MeshAlgorithm::CalculateMinimumGridLength(float fLength, const Base::BoundBox3f& rBBox, unsigned long maxElements) const
01671 {
01672
01673 float fMaxGridElements=(float)maxElements;
01674
01675
01676 float fMinGridLen = (float)pow((rBBox.LengthX()*rBBox.LengthY()*rBBox.LengthZ()/fMaxGridElements), 0.3333f);
01677 return std::max<float>(fMinGridLen, fLength);
01678 }
01679
01680
01681
01682 void MeshRefPointToFacets::Rebuild (void)
01683 {
01684 _map.clear();
01685
01686 const MeshPointArray& rPoints = _rclMesh.GetPoints();
01687 const MeshFacetArray& rFacets = _rclMesh.GetFacets();
01688 _map.resize(rPoints.size());
01689
01690 MeshFacetArray::_TConstIterator pFBegin = rFacets.begin();
01691 for (MeshFacetArray::_TConstIterator pFIter = rFacets.begin(); pFIter != rFacets.end(); ++pFIter) {
01692 _map[pFIter->_aulPoints[0]].insert(pFIter - pFBegin);
01693 _map[pFIter->_aulPoints[1]].insert(pFIter - pFBegin);
01694 _map[pFIter->_aulPoints[2]].insert(pFIter - pFBegin);
01695 }
01696 }
01697
01698 Base::Vector3f MeshRefPointToFacets::GetNormal(unsigned long pos) const
01699 {
01700 const std::set<unsigned long>& n = _map[pos];
01701 Base::Vector3f normal;
01702 MeshGeomFacet f;
01703 for (std::set<unsigned long>::const_iterator it = n.begin(); it != n.end(); ++it) {
01704 f = _rclMesh.GetFacet(*it);
01705 normal += f.Area() * f.GetNormal();
01706 }
01707
01708 normal.Normalize();
01709 return normal;
01710 }
01711
01712 std::set<unsigned long> MeshRefPointToFacets::NeighbourPoints(const std::vector<unsigned long>& pt, int level) const
01713 {
01714 std::set<unsigned long> cp,nb,lp;
01715 cp.insert(pt.begin(), pt.end());
01716 lp.insert(pt.begin(), pt.end());
01717 MeshFacetArray::_TConstIterator f_it = _rclMesh.GetFacets().begin();
01718 for (int i=0; i < level; i++) {
01719 std::set<unsigned long> cur;
01720 for (std::set<unsigned long>::iterator it = lp.begin(); it != lp.end(); ++it) {
01721 const std::set<unsigned long>& ft = (*this)[*it];
01722 for (std::set<unsigned long>::const_iterator jt = ft.begin(); jt != ft.end(); ++jt) {
01723 for (int j = 0; j < 3; j++) {
01724 unsigned long index = f_it[*jt]._aulPoints[j];
01725 if (cp.find(index) == cp.end() && nb.find(index) == nb.end()) {
01726 nb.insert(index);
01727 cur.insert(index);
01728 }
01729 }
01730 }
01731 }
01732
01733 lp = cur;
01734 if (lp.empty())
01735 break;
01736 }
01737 return nb;
01738 }
01739
01740
01741
01743 void MeshRefPointToFacets::Neighbours (unsigned long ulFacetInd, float fMpxDist, std::vector<MeshFacetArray::_TConstIterator> &rclNb)
01744 {
01745 rclNb.clear();
01746 Base::Vector3f clCenter = _rclMesh.GetFacet(ulFacetInd).GetGravityPoint();
01747
01748 const MeshFacetArray& rFacets = _rclMesh.GetFacets();
01749 SearchNeighbours(rFacets.begin() + ulFacetInd, clCenter, fMpxDist * fMpxDist, rclNb);
01750
01751 for (std::vector<MeshFacetArray::_TConstIterator>::iterator i = rclNb.begin(); i != rclNb.end(); i++)
01752 (*i)->ResetFlag(MeshFacet::VISIT);
01753 }
01754
01756 void MeshRefPointToFacets::SearchNeighbours(MeshFacetArray::_TConstIterator f_it, const Base::Vector3f &rclCenter, float fMpxDist, std::vector<MeshFacetArray::_TConstIterator> &rclNb)
01757 {
01758 if (f_it->IsFlag(MeshFacet::VISIT) == true)
01759 return;
01760
01761 if (Base::DistanceP2(rclCenter, _rclMesh.GetFacet(*f_it).GetGravityPoint()) > fMpxDist)
01762 return;
01763
01764 rclNb.push_back(f_it);
01765 f_it->SetFlag(MeshFacet::VISIT);
01766
01767 MeshFacetArray::_TConstIterator f_beg = _rclMesh.GetFacets().begin();
01768
01769 for (int i = 0; i < 3; i++) {
01770 const std::set<unsigned long> &f = (*this)[f_it->_aulPoints[i]];
01771
01772 for (std::set<unsigned long>::const_iterator j = f.begin(); j != f.end(); ++j) {
01773 SearchNeighbours(f_beg+*j, rclCenter, fMpxDist, rclNb);
01774 }
01775 }
01776 }
01777
01778 MeshFacetArray::_TConstIterator
01779 MeshRefPointToFacets::getFacet (unsigned long index) const
01780 {
01781 return _rclMesh.GetFacets().begin() + index;
01782 }
01783
01784 const std::set<unsigned long>&
01785 MeshRefPointToFacets::operator[] (unsigned long pos) const
01786 {
01787 return _map[pos];
01788 }
01789
01790
01791
01792 void MeshRefFacetToFacets::Rebuild (void)
01793 {
01794 _map.clear();
01795
01796 const MeshFacetArray& rFacets = _rclMesh.GetFacets();
01797 _map.resize(rFacets.size());
01798
01799 MeshRefPointToFacets vertexFace(_rclMesh);
01800 MeshFacetArray::_TConstIterator pFBegin = rFacets.begin();
01801 for (MeshFacetArray::_TConstIterator pFIter = pFBegin; pFIter != rFacets.end(); ++pFIter) {
01802 for (int i = 0; i < 3; i++) {
01803 const std::set<unsigned long>& faces = vertexFace[pFIter->_aulPoints[i]];
01804 for (std::set<unsigned long>::const_iterator it = faces.begin(); it != faces.end(); ++it)
01805 _map[pFIter - pFBegin].insert(*it);
01806 }
01807 }
01808 }
01809
01810 const std::set<unsigned long>&
01811 MeshRefFacetToFacets::operator[] (unsigned long pos) const
01812 {
01813 return _map[pos];
01814 }
01815
01816
01817
01818 void MeshRefPointToPoints::Rebuild (void)
01819 {
01820 _map.clear();
01821
01822 const MeshPointArray& rPoints = _rclMesh.GetPoints();
01823 _map.resize(rPoints.size());
01824
01825 const MeshFacetArray& rFacets = _rclMesh.GetFacets();
01826 for (MeshFacetArray::_TConstIterator pFIter = rFacets.begin(); pFIter != rFacets.end(); ++pFIter) {
01827 unsigned long ulP0 = pFIter->_aulPoints[0];
01828 unsigned long ulP1 = pFIter->_aulPoints[1];
01829 unsigned long ulP2 = pFIter->_aulPoints[2];
01830
01831 _map[ulP0].insert(ulP1);
01832 _map[ulP0].insert(ulP2);
01833 _map[ulP1].insert(ulP0);
01834 _map[ulP1].insert(ulP2);
01835 _map[ulP2].insert(ulP0);
01836 _map[ulP2].insert(ulP1);
01837 }
01838 }
01839
01840 Base::Vector3f MeshRefPointToPoints::GetNormal(unsigned long pos) const
01841 {
01842 const MeshPointArray& rPoints = _rclMesh.GetPoints();
01843 MeshCore::PlaneFit pf;
01844 pf.AddPoint(rPoints[pos]);
01845 MeshCore::MeshPoint center = rPoints[pos];
01846 const std::set<unsigned long>& cv = _map[pos];
01847 for (std::set<unsigned long>::const_iterator cv_it = cv.begin(); cv_it !=cv.end(); ++cv_it) {
01848 pf.AddPoint(rPoints[*cv_it]);
01849 center += rPoints[*cv_it];
01850 }
01851
01852 pf.Fit();
01853
01854 Base::Vector3f normal = pf.GetNormal();
01855 normal.Normalize();
01856 return normal;
01857 }
01858
01859 float MeshRefPointToPoints::GetAverageEdgeLength(unsigned long index) const
01860 {
01861 const MeshPointArray& rPoints = _rclMesh.GetPoints();
01862 float len=0.0f;
01863 const std::set<unsigned long>& n = (*this)[index];
01864 const Base::Vector3f& p = rPoints[index];
01865 for (std::set<unsigned long>::const_iterator it = n.begin(); it != n.end(); ++it) {
01866 len += Base::Distance(p, rPoints[*it]);
01867 }
01868 return (len/n.size());
01869 }
01870
01871 const std::set<unsigned long>&
01872 MeshRefPointToPoints::operator[] (unsigned long pos) const
01873 {
01874 return _map[pos];
01875 }
01876
01877
01878
01879 void MeshRefEdgeToFacets::Rebuild (void)
01880 {
01881 _map.clear();
01882
01883 const MeshFacetArray& rFacets = _rclMesh.GetFacets();
01884 unsigned long index = 0;
01885 for (MeshFacetArray::_TConstIterator it = rFacets.begin(); it != rFacets.end(); ++it, ++index) {
01886 for (int i=0; i<3; i++) {
01887 MeshEdge e;
01888 e.first = it->_aulPoints[i];
01889 e.second = it->_aulPoints[(i+1)%3];
01890 std::map<MeshEdge, MeshFacetPair, EdgeOrder>::iterator jt =
01891 _map.find(e);
01892 if (jt == _map.end()) {
01893 _map[e].first = index;
01894 _map[e].second = ULONG_MAX;
01895 }
01896 else {
01897 _map[e].second = index;
01898 }
01899 }
01900 }
01901 }
01902
01903 const std::pair<unsigned long, unsigned long>&
01904 MeshRefEdgeToFacets::operator[] (const MeshEdge& edge) const
01905 {
01906 return _map.find(edge)->second;
01907 }
01908
01909
01910
01911 void MeshRefNormalToPoints::Rebuild (void)
01912 {
01913 _norm.clear();
01914
01915 const MeshPointArray& rPoints = _rclMesh.GetPoints();
01916 _norm.resize(rPoints.size());
01917
01918 const MeshFacetArray& rFacets = _rclMesh.GetFacets();
01919 for (MeshFacetArray::_TConstIterator pF = rFacets.begin(); pF != rFacets.end(); ++pF) {
01920 const MeshPoint &p0 = rPoints[pF->_aulPoints[0]];
01921 const MeshPoint &p1 = rPoints[pF->_aulPoints[1]];
01922 const MeshPoint &p2 = rPoints[pF->_aulPoints[2]];
01923 float l2p01 = Base::DistanceP2(p0,p1);
01924 float l2p12 = Base::DistanceP2(p1,p2);
01925 float l2p20 = Base::DistanceP2(p2,p0);
01926
01927 Base::Vector3f facenormal = _rclMesh.GetFacet(*pF).GetNormal();
01928 _norm[pF->_aulPoints[0]] += facenormal * (1.0f / (l2p01 * l2p20));
01929 _norm[pF->_aulPoints[1]] += facenormal * (1.0f / (l2p12 * l2p01));
01930 _norm[pF->_aulPoints[2]] += facenormal * (1.0f / (l2p20 * l2p12));
01931 }
01932 for (std::vector<Base::Vector3f>::iterator it = _norm.begin(); it != _norm.end(); ++it)
01933 it->Normalize();
01934 }
01935
01936 const Base::Vector3f&
01937 MeshRefNormalToPoints::operator[] (unsigned long pos) const
01938 {
01939 return _norm[pos];
01940 }