From b780824be6a37305c3df96ee6d614a852211a885 Mon Sep 17 00:00:00 2001 From: untereiner Date: Wed, 13 Jun 2012 15:40:14 +0200 Subject: [PATCH] - adding multilevel extrusion to importMeshExtrud - volumetric decimation --- Apps/deprecated/tutoriel.cpp | 2 +- include/Algo/Decimation/decimation.hpp | 2 +- include/Algo/Decimation/edgeSelector.hpp | 10 +- include/Algo/DecimationVolumes/approximator.h | 88 ++++- include/Algo/DecimationVolumes/decimator.hpp | 60 +++- include/Algo/DecimationVolumes/edgeSelector.h | 166 +++------ .../Algo/DecimationVolumes/edgeSelector.hpp | 337 ++++++++++-------- .../DecimationVolumes/geometryApproximator.h | 22 +- .../geometryApproximator.hpp | 93 +++-- include/Algo/DecimationVolumes/operator.h | 19 +- include/Algo/DecimationVolumes/operator.hpp | 121 +++---- include/Algo/DecimationVolumes/selector.h | 11 +- include/Algo/Import/import.h | 5 +- include/Algo/Import/importMesh.hpp | 176 +++++++-- .../Algo/Modelisation/tetrahedralization.h | 5 + .../Algo/Modelisation/tetrahedralization.hpp | 79 ++++ include/Algo/Selection/raySelector.h | 18 + include/Geometry/intersection.h | 14 + include/Geometry/intersection.hpp | 82 +++++ src/Topology/map/embeddedMap3.cpp | 73 +++- 20 files changed, 914 insertions(+), 469 deletions(-) diff --git a/Apps/deprecated/tutoriel.cpp b/Apps/deprecated/tutoriel.cpp index 31d48ed55..d8d8664be 100644 --- a/Apps/deprecated/tutoriel.cpp +++ b/Apps/deprecated/tutoriel.cpp @@ -621,7 +621,7 @@ void MyQT::cb_keyPress(int keycode) // Sélectionne des faces case 'f': d_faces.clear(); - Algo::Selection:: facesRaySelection(myMap, position, SelectorTrue(), rayA, AB, d_faces); + Algo::Selection::facesRaySelection(myMap, position, SelectorTrue(), rayA, AB, d_faces); if (!d_faces.empty()) { diff --git a/include/Algo/Decimation/decimation.hpp b/include/Algo/Decimation/decimation.hpp index 91b350475..6cbc490a0 100644 --- a/include/Algo/Decimation/decimation.hpp +++ b/include/Algo/Decimation/decimation.hpp @@ -126,7 +126,7 @@ void decimate( selector = new EdgeSelector_Curvature(map, position, approximators, selected) ; break ; case S_MinDetail : - selector = new EdgeSelector_Random(map, position, approximators, selected) ; + selector = new EdgeSelector_MinDetail(map, position, approximators, selected) ; break ; case S_hLightfield : selector = new HalfEdgeSelector_Lightfield(map, position, approximators, selected) ; diff --git a/include/Algo/Decimation/edgeSelector.hpp b/include/Algo/Decimation/edgeSelector.hpp index 13800cb0e..91ce97e8f 100644 --- a/include/Algo/Decimation/edgeSelector.hpp +++ b/include/Algo/Decimation/edgeSelector.hpp @@ -946,14 +946,10 @@ bool EdgeSelector_MinDetail::init() edges.clear() ; - CellMarker eMark(m) ; - for(Dart d = m.begin(); d != m.end(); m.next(d)) + TraversorE travE(m); + for(Dart dit = travE.begin() ; dit != travE.end() ; dit = travE.next()) { - if(!eMark.isMarked(d)) - { - initEdgeInfo(d) ; // init the edges with their optimal position - eMark.mark(d) ; // and insert them in the multimap according to their error - } + initEdgeInfo(dit); } cur = edges.begin() ; // init the current edge to the first one diff --git a/include/Algo/DecimationVolumes/approximator.h b/include/Algo/DecimationVolumes/approximator.h index 202228290..92596bb6f 100644 --- a/include/Algo/DecimationVolumes/approximator.h +++ b/include/Algo/DecimationVolumes/approximator.h @@ -26,6 +26,7 @@ #define __APPROXIMATOR_VOLUMES_H__ #include "Algo/DecimationVolumes/operator.h" +#include "Algo/DecimationVolumes/predictor.h" namespace CGoGN { @@ -38,7 +39,7 @@ namespace DecimationVolumes enum ApproximatorType { - A_Centroid //barycenter of the n-cells + A_QEM }; template @@ -60,9 +61,12 @@ public: virtual const std::string& getApproximatedAttributeName() const = 0 ; virtual ApproximatorType getType() const = 0 ; virtual bool init() = 0 ; + virtual const PredictorGen* getPredictor() const = 0 ; virtual void approximate(Operator* op) = 0 ; virtual void saveApprox(Operator* op) = 0 ; virtual void affectApprox(Operator* op) = 0 ; + + virtual void approximate(Dart d) = 0; } ; @@ -74,39 +78,91 @@ public: typedef typename PFP::REAL REAL; protected: + Predictor* m_predictor ; - //TODO ajouter un predictor - - VertexAttribute& m_attrV; // vertex attribute to be approximated - - //TODO Attribute to store approximation result - //TODO attribute to store detail information for reconstruction + //Vertex attribute to be approximated + VertexAttribute& m_attrV; + //Attribute to store approximation result + EdgeAttribute m_approx; + // attribute to store detail information for reconstruction + EdgeAttribute m_detail ; T m_app; public: - Approximator(MAP& m, VertexAttribute& a): - ApproximatorGen(m), m_attrV(a) - {} + Approximator(MAP& m, VertexAttribute& a, Predictor* predictor) : + ApproximatorGen(m), m_predictor(predictor), m_attrV(a) + { + std::stringstream aname ; + aname << "approx_" << m_attrV.name() ; + m_approx = this->m_map.template addAttribute(aname.str()) ; + + if(m_predictor) // if a predictor is associated to the approximator + { // create an attribute to store the detail needed for reconstruction + std::stringstream dname ; + dname << "detail_" << m_attrV.name() ; + m_detail = this->m_map.template addAttribute(dname.str()) ; + } + } - //virtual ~Approximator(); + virtual ~Approximator() + { + this->m_map.template removeAttribute(m_approx) ; + + if(m_predictor) + this->m_map.template removeAttribute(m_detail) ; + } const std::string& getApproximatedAttributeName() const { return m_attrV.name() ; } - void saveApprox(Operator* op) + const T& getApprox(Dart d) const { - Dart d = op->getEdge(); - //m_app = m_approx[d] ; + return m_approx[d] ; + } + + const Predictor* getPredictor() const + { + return m_predictor ; } - void affectApprox(Operator* op) + const T& getDetail(Dart d) const + { + assert(m_predictor || !"Trying to get detail on a non-predictive scheme") ; + return m_detail[d] ; + } + + void setDetail(Dart d, T& val) + { + assert(m_predictor || !"Trying to set detail on a non-predictive scheme") ; + m_detail[d] = val ; + } + + virtual void saveApprox(Operator* op) { Dart d = op->getEdge(); - m_attrV[d] = m_app ; + m_app = m_approx[d]; } + virtual void affectApprox(Operator* op) + { + Dart d = op->getEdge(); + m_attrV[d] = m_app; + } + + + // void saveApprox(Dart d) + // { + // m_app = m_approx[d] ; + // } + // + // void affectApprox(Dart d) + // { + // m_attrV[d] = m_app ; + // } + + }; diff --git a/include/Algo/DecimationVolumes/decimator.hpp b/include/Algo/DecimationVolumes/decimator.hpp index 5e861c388..ecbb622d5 100644 --- a/include/Algo/DecimationVolumes/decimator.hpp +++ b/include/Algo/DecimationVolumes/decimator.hpp @@ -38,13 +38,13 @@ void decimate( ) { std::vector*> approximators ; - Selector* selector = NULL ; + EdgeSelector* selector = NULL ; //choose the Approximator switch(a) { - case A_Centroid : - approximators.push_back(new Approximator_Centroid(map, position)) ; + case A_QEM : + approximators.push_back(new Approximator_QEM(map, position)) ; break ; default : CGoGNout << "not yet implemented" << CGoGNendl; @@ -54,11 +54,8 @@ void decimate( //choose the Selector switch(s) { - case S_MapOrder : - selector = new EdgeSelector_MapOrder(map, position, approximators, selected); - break ; - case S_Random : - //selector = new EdgeSelector_Random(map, position, approximators, selected); + case S_QEM : + selector = new EdgeSelector_QEM(map, position, approximators, selected) ; break ; default: CGoGNout << "not yet implemented" << CGoGNendl; @@ -77,33 +74,38 @@ void decimate( CGoGNout << " decimate (" << nbVertices << " vertices).." << /* flush */ CGoGNendl ; bool finished = false ; + Dart d; + while(!finished) { //Next Operator to perform Operator *op; - if(!selector->nextOperator(&op)) + //if(!selector->nextOperator(&op)) + if((op = selector->nextOperator()) == NULL) break; // compute approximated attributes for(typename std::vector*>::iterator it = approximators.begin(); it != approximators.end(); ++it) { (*it)->approximate(op) ; + (*it)->saveApprox(op) ; } //Update the selector before performing operation - selector->updateBeforeOperation(&op); + if(!selector->updateBeforeOperation(op)) + break; //Perform the topological operation and //compute the number of resulting cells nbVertices -= op->perform(map, position); - //for(typename std::vector*>::iterator it = approximators.begin(); it != approximators.end(); ++it) - // (*it)->affectApprox(op); // affect data to the resulting vertex + for(typename std::vector*>::iterator it = approximators.begin(); it != approximators.end(); ++it) + (*it)->affectApprox(op); // affect data to the resulting vertex //Update the embedded position and //search the next operation to perform - selector->updateAfterOperation(&op); + selector->updateAfterOperation(op); if(nbVertices <= nbWantedVertices) finished = true ; @@ -126,3 +128,35 @@ void decimate( } //namespace Algo } //namespace CGoGN + +// if(!selector->nextEdge(d)) +// break ; +// +// +// Dart d2 = map.phi2(map.phi_1(d)) ; +// Dart dd2 = map.phi2(map.phi_1(map.phi2(d))) ; +// +// std::cout << "bin a contracter : " << d << std::endl; +// std::cout << "voisin d2 : " << d2 << std::endl; +// std::cout << "voisin dd2 : " << dd2 << std::endl; +// +// --nbVertices ; +// +// for(typename std::vector*>::iterator it = approximators.begin(); it != approximators.end(); ++it) +// { +// (*it)->approximate(d) ; // compute approximated attributes +// (*it)->saveApprox(d) ; +// } +// +// selector->updateBeforeCollapse(d) ; // update selector +// +// map.collapseEdge(d) ; // collapse edge +// +// for(typename std::vector*>::iterator it = approximators.begin(); it != approximators.end(); ++it) +// (*it)->affectApprox(d2); // affect data to the resulting vertex +// +// selector->updateAfterCollapse(d2, dd2) ;// update selector +// +// if(nbVertices <= nbWantedVertices) +// finished = true ; + diff --git a/include/Algo/DecimationVolumes/edgeSelector.h b/include/Algo/DecimationVolumes/edgeSelector.h index da35a100f..12a7138cc 100644 --- a/include/Algo/DecimationVolumes/edgeSelector.h +++ b/include/Algo/DecimationVolumes/edgeSelector.h @@ -4,6 +4,8 @@ #include "Container/fakeAttribute.h" #include "Algo/DecimationVolumes/selector.h" #include "Algo/DecimationVolumes/operator.h" +#include "Utils/qem.h" +#include "Topology/generic/traversorCell.h" namespace CGoGN { @@ -35,140 +37,62 @@ public: }; -/******************************************************************************** - * Map Order Edge Selector * - ********************************************************************************/ template -class EdgeSelector_MapOrder : public EdgeSelector +class EdgeSelector_QEM : public EdgeSelector { public: - typedef typename PFP::MAP MAP ; - typedef typename PFP::VEC3 VEC3 ; + typedef typename PFP::MAP MAP; + typedef typename PFP::VEC3 VEC3; typedef typename PFP::REAL REAL; private: - Dart cur; + typedef struct + { + typename std::multimap::iterator it; + bool valid; + static std::string CGoGNnameOfType() { return "QEMedgeInfo" ; } + } QEMedgeInfo; -public: - EdgeSelector_MapOrder(MAP& m, VertexAttribute& pos, std::vector*>& approx, const FunctorSelect& select) : - EdgeSelector(m, pos, approx, select) - {} + typedef NoMathIOAttribute EdgeInfo; - ~EdgeSelector_MapOrder() - {} + EdgeAttribute edgeInfo; + VertexAttribute > quadric; + Quadric tmpQ; - SelectorType getType() { return S_MapOrder; } - Dart nextCell() { return cur; } + std::multimap edges; + typename std::multimap::iterator cur; - bool init(); - bool nextOperator(Operator** op); - void updateBeforeOperation(Operator** op) - { } - void updateAfterOperation(Operator** op); - void finish() - { } + Approximator* m_positionApproximator; + void initEdgeInfo(Dart d); + void updateEdgeInfo(Dart d, bool recompute); + void computeEdgeInfo(Dart d, EdgeInfo& einfo); + +public: + EdgeSelector_QEM(MAP& m, VertexAttribute& pos, + std::vector* >& approx, const FunctorSelect& select) : + EdgeSelector(m, pos, approx, select) + { + edgeInfo = m.template addAttribute("edgeInfo"); + quadric = m.template addAttribute, VERTEX>("QEMquadric"); + } + + ~EdgeSelector_QEM() + { + this->m_map.removeAttribute(quadric); + this->m_map.removeAttribute(edgeInfo); + } + + SelectorType getType() { return S_QEM; } + bool init(); + Operator* nextOperator(); + bool updateBeforeOperation(Operator* op); + void updateAfterOperation(Operator* op); + void finish() { } }; -/******************************************************************************** - * Random Selector * - ********************************************************************************/ -//template -//class EdgeSelector_Random : public EdgeSelector -//{ -//public: -// typedef typename PFP::MAP MAP ; -// typedef typename PFP::VEC3 VEC3 ; -// typedef typename PFP::REAL REAL ; -// -//private: -// std::vector darts ; -// unsigned int cur ; -// bool allSkipped ; -// -//public: -// EdgeSelector_Random(MAP& m, VertexAttribute& pos, std::vector*>& approx, const FunctorSelect& select) : -// EdgeSelector(m, pos, approx, select) -// {} -// -// ~EdgeSelector_Random() -// {} -// -// SelectorType getType() { return S_Random; } -//// Dart nextCell() { return darts[cur]; } -// -// bool init(); -// bool nextOperator(Operator** op); -// void updateBeforeOperation(Operator* op) -// { } -// void updateAfterOperation(Operator* op); -// void finish() -// { } -// -//} ; - -///******************************************************************************** -// * Length Edge Selector * -// ********************************************************************************/ -//template -//class EdgeSelector_Length : public EdgeSelector -//{ -//public: -// typedef typename PFP::MAP MAP; -// typedef typename PFP::VEC3 VEC3; -// typedef typename PFP::REAL REAL; -// -//protected: -// /* -// * Update the Edge Informations -// */ -// void updateEdgeInfo(Dart d, typename PFP::VEC3 e); -// /* -// * Init the Edge informations -// */ -// void initEdgeInfo(Dart d); -// /* -// * Erase The Edge informations -// */ -// void eraseEdgeInfo(Dart d); -// -//private: -// //New embedding type for storing informations -// typedef struct { -// typename std::multimap::iterator it ; -// bool valid ; -// } LengthEdgeInfo ; -// -// typedef NoMathIOAttribute EdgeInfo ; -// -// EdgeAttribute edgeInfo ; -// -// std::multimap edges ; -// typename std::multimap::iterator cur ; -// -// -//public: -// EdgeSelector_Length(MAP& m, VertexAttribute& pos) : -// EdgeSelector(m, pos) -// { -// edgeInfo = m.template addAttribute(EDGE, "edgeInfo") ; -// } -// -// ~EdgeSelector_Length() -// { -// this->m_map.removeAttribute(edgeInfo) ; -// } -// -// SelectorType getType() { return S_EdgeLength; } -// Dart nexCell() { return (*cur).second; } -// -// void init(Algo::DecimationVolumique::Approximator* approx); -// bool nextOperator(); -// void updateBeforeOperation(Algo::DecimationVolumique::Operator* op); -// void updateAfterOperation(Algo::DecimationVolumique::Operator* op); -// void finish() -// { } -//}; + + } //end namespace DecimationVolumes diff --git a/include/Algo/DecimationVolumes/edgeSelector.hpp b/include/Algo/DecimationVolumes/edgeSelector.hpp index 604a6fb26..6025e76cb 100644 --- a/include/Algo/DecimationVolumes/edgeSelector.hpp +++ b/include/Algo/DecimationVolumes/edgeSelector.hpp @@ -35,179 +35,222 @@ namespace DecimationVolumes { /************************************************************************************ - * Map Order Edge Selector * + * QUADRIC ERROR METRIC * ************************************************************************************/ template -bool EdgeSelector_MapOrder::init() +bool EdgeSelector_QEM::init() { - cur = this->m_map.begin() ; + MAP& m = this->m_map ; + + //searching the Geometry Approximator + bool ok = false ; + for(typename std::vector*>::iterator it = this->m_approximators.begin(); + it != this->m_approximators.end() && !ok; + ++it) + { + if((*it)->getApproximatedAttributeName() == "position") + { + m_positionApproximator = reinterpret_cast* >(*it) ; + ok = true ; + } + } + + if(!ok) + return false ; + + // init the Quadric error metric for cells + edges.clear(); + TraversorV travV(m); + for(Dart dit = travV.begin() ; dit != travV.end() ; dit = travV.next()) + { + //create one quadric per vertex + Quadric q; + quadric[dit] = q; + } + + // init the quadric for each boundary triangle + TraversorF travF(m); + for(Dart dit = travF.begin() ; dit != travF.end() ; dit = travF.next()) + { + Dart d1 = m.phi1(dit) ; // for each triangle, + Dart d_1 = m.phi_1(dit) ; // initialize the quadric of the triangle + Quadric q(this->m_position[dit], this->m_position[d1], this->m_position[d_1]) ; + quadric[dit] += q ; // and add the contribution of + quadric[d1] += q ; // this quadric to the ones + quadric[d_1] += q ; // of the 3 incident vertices + + } + + // init the edges with their optimal position + // and insert them in the multimap according to their error + TraversorE travE(m); + for(Dart dit = travE.begin() ; dit != travE.end() ; dit = travE.next()) + { + initEdgeInfo(dit); + } - //TODO Choix de nextOperator - //this->nextOp = Algo::DecimationVolumique::O_CEdge; + // init the current edge to the first one + cur = edges.begin(); return true; } template -bool EdgeSelector_MapOrder::nextOperator(Operator** op) +Operator* EdgeSelector_QEM::nextOperator() { - if(cur == this->m_map.end()) - return false ; + MAP& m = this->m_map ; + Operator* op = NULL; + + if(cur != edges.end() && !edges.empty()) + { + op = new CollapseEdgeOperator((*cur).second); - *op = new CollapseEdgeOperator(cur); + op->setFirstIncidentEdge(m.phi2(m.phi_1((*cur).second))); + op->setSecondIncidentEdge(m.phi2(m.phi_1(m.phi2((*cur).second)))); + } - return true ; + return op ; } template -void EdgeSelector_MapOrder::updateAfterOperation(Operator** op) +bool EdgeSelector_QEM::updateBeforeOperation(Operator* op) { MAP& m = this->m_map ; + Dart d = op->getEdge(); - //Mise à jour de nextOperator + EdgeInfo& edgeE = edgeInfo[d] ; + if(edgeE.valid) + edges.erase(edgeE.it) ; - cur = m.begin() ; + edgeE = edgeInfo[m.phi1(d)] ; + if(edgeE.valid) // remove all + edges.erase(edgeE.it) ; - while(!(*op)->canPerform(m,cur,this->m_position)) + edgeE = edgeInfo[m.phi_1(d)] ; // the concerned edges + if(edgeE.valid) + edges.erase(edgeE.it) ; + // from the multimap + Dart dd = m.phi2(d) ; + if(dd != d) { - //CGoGNout << "update " << cur << CGoGNendl; - m.next(cur) ; - if(cur == m.end()) - break ; + edgeE = edgeInfo[m.phi1(dd)] ; + if(edgeE.valid) + edges.erase(edgeE.it) ; + + edgeE = edgeInfo[m.phi_1(dd)] ; + if(edgeE.valid) + edges.erase(edgeE.it) ; } + + tmpQ.zero() ; // compute quadric for the new + tmpQ += quadric[d] ; // vertex as the sum of those + tmpQ += quadric[dd] ; // of the contracted vertices + + return true; } -/************************************************************************************ - * Random Edge Selector * - ************************************************************************************/ +template +void EdgeSelector_QEM::updateAfterOperation(Operator* op) +{ + MAP& m = this->m_map ; + Dart d2 = op->getFirstIncidentEdge(); + Dart dd2 = op->getSecondIncidentEdge(); + + quadric[d2] = tmpQ ; + + Dart vit = d2 ; + do + { + updateEdgeInfo(m.phi1(vit), false) ; // must recompute some edge infos in the + if(vit == d2 || vit == dd2) // neighborhood of the collapsed edge + { + initEdgeInfo(vit) ; // various optimizations are applied here by + // treating differently : + Dart vit2 = m.phi1(m.phi2(m.phi1(vit))) ; // - edges for which the criteria must be recomputed + Dart stop = m.phi2(vit) ; // - edges that must be re-embedded + do // - edges for which only the collapsibility must be re-tested + { + updateEdgeInfo(vit2, false) ; + updateEdgeInfo(m.phi1(vit2), false) ; + vit2 = m.phi1(m.phi2(vit2)) ; + } while(vit2 != stop) ; + } + else + updateEdgeInfo(vit, true) ; + + vit = m.phi2(m.phi_1(vit)) ; + } while(vit != d2) ; + + cur = edges.begin() ; // set the current edge to the first one + +} + +template +void EdgeSelector_QEM::initEdgeInfo(Dart d) +{ + MAP& m = this->m_map ; + EdgeInfo einfo ; + if(m.edgeCanCollapse(d)) + computeEdgeInfo(d, einfo) ; + else + einfo.valid = false ; + edgeInfo[d] = einfo ; +} + +template +void EdgeSelector_QEM::updateEdgeInfo(Dart d, bool recompute) +{ + MAP& m = this->m_map ; + EdgeInfo& einfo = edgeInfo[d] ; + + if(recompute) + { + if(einfo.valid) + edges.erase(einfo.it) ; // remove the edge from the multimap + if(m.edgeCanCollapse(d)) + computeEdgeInfo(d, einfo) ; + else + einfo.valid = false ; + } + else + { + if(m.edgeCanCollapse(d)) + { // if the edge can be collapsed now + if(!einfo.valid) // but it was not before + computeEdgeInfo(d, einfo) ; + } + else + { // if the edge cannot be collapsed now + if(einfo.valid) // and it was before + { + edges.erase(einfo.it) ; + einfo.valid = false ; + } + } + } +} + +template +void EdgeSelector_QEM::computeEdgeInfo(Dart d, EdgeInfo& einfo) +{ + MAP& m = this->m_map ; + Dart dd = m.phi2(d) ; + + Quadric quad ; + quad += quadric[d] ; // compute the sum of the + quad += quadric[dd] ; // two vertices quadrics + + m_positionApproximator->approximate(d) ; + + REAL err = std::max(REAL(0),REAL(quad(m_positionApproximator->getApprox(d)))) ; + + einfo.it = edges.insert(std::make_pair(err, d)) ; + einfo.valid = true ; +} + + -//template -//bool EdgeSelector_Random::init() -//{ -// MAP& m = this->m_map ; -// -// darts.reserve(m.getNbDarts()) ; -// darts.clear() ; -// -// for(Dart d = m.begin(); d != m.end(); m.next(d)) -// darts.push_back(d) ; -// -// srand(time(NULL)) ; -// int remains = darts.size() ; -// for(unsigned int i = 0; i < darts.size()-1; ++i) // generate the random permutation -// { -// int r = (rand() % remains) + i ; -// // swap ith and rth elements -// Dart tmp = darts[i] ; -// darts[i] = darts[r] ; -// darts[r] = tmp ; -// --remains ; -// } -// cur = 0 ; -// allSkipped = true ; -// -// //TODO Choix de nextOperator -// //this->nextOp = O_CEdge; -// -// return true; -//} -// -//template -//bool EdgeSelector_Random::nextOperator(Operator** op) -//{ -// if(cur == darts.size() && allSkipped) -// return false ; -// -// return true ; -//} -// -// -//template -//void EdgeSelector_Random::updateAfterOperation(Operator* op) -//{ -// MAP& m = this->m_map ; -// allSkipped = false ; -// -// //if(this->nextOp == O_CEdge) -// //{ -//// for(typename std::vector::iterator it = darts.begin() ; it != darts.end() ; ++it) -//// CGoGNout << *it << CGoGNendl; -// -// do -// { -// ++cur ; -// if(cur == darts.size()) -// { -// //cur = 0 ; -// allSkipped = true ; -// } -// } while(!allSkipped && !op->canPerform(m,darts[cur],this->m_position)) ; -// //} -//} - -///************************************************************************************ -// * Edge Length Selector * -// ************************************************************************************/ -//template -//void EdgeSelector_Length::init() -//{ -// MAP& m = this->m_map ; -// -// edges.clear() ; -// -// CellMarker eMark(m, EDGE) ; -// for(Dart d = m.begin(); d != m.end(); m.next(d)) -// { -// if(!eMark.isMarked(d)) -// { -// initEdgeInfo(d) ; -// eMark.mark(d) ; -// } -// } -// -// cur = edges.begin() ; // init the current edge to the first one -// -// //TODO Choix de nextOperator -// this->nextOp = Algo::DecimationVolumique::O_CEdge; -//} -// -//template -//bool EdgeSelector_Length::nextOperator(Operator* op) -//{ -// if(cur == edges.end() || edges.empty()) -// return false ; -// -// return true ; -//} -// -//template -//void EdgeSelector_Length::updateBeforeOperation(Operator* op) -//{ -// MAP& m = this->m_map; -// Dart d = op->getEdge(); -// -// -// //Remove all the concerned darts from the multimap -//} -// -//template -//void EdgeSelector_Length::updateAfterOperation(Operator* op) -//{ -// -//} -// -//template -//void EdgeSelector_Length::initEdgeInfo(Dart d) -//{ -// MAP& m = this->m_map ; -// EdgeInfo einfo ; -// if(m.edgeCanCollapse(d)) -// computeEdgeInfo(d, einfo) ; -// else -// einfo.valid = false ; -// edgeInfo[d] = einfo ; -//} } //end namespace DecimationVolumique } //end namespace Algo diff --git a/include/Algo/DecimationVolumes/geometryApproximator.h b/include/Algo/DecimationVolumes/geometryApproximator.h index c49854da2..48dd04cea 100644 --- a/include/Algo/DecimationVolumes/geometryApproximator.h +++ b/include/Algo/DecimationVolumes/geometryApproximator.h @@ -37,24 +37,26 @@ namespace DecimationVolumes { template -class Approximator_Centroid : public Approximator +class Approximator_QEM : public Approximator { public: typedef typename PFP::MAP MAP ; typedef typename PFP::VEC3 VEC3 ; - typedef typename PFP::REAL REAL; + typedef typename PFP::REAL REAL ; - Approximator_Centroid(MAP& m, VertexAttribute& pos): - Approximator(m, pos) - {} +protected: + VertexAttribute > m_quadric ; - ~Approximator_Centroid() +public: + Approximator_QEM(MAP& m, VertexAttribute& pos, Predictor* pred = NULL) : + Approximator(m, pos, pred) {} - - ApproximatorType getType() const { return A_Centroid; } - bool init(); + ~Approximator_QEM() + {} + ApproximatorType getType() const { return A_QEM ; } + bool init() ; void approximate(Operator *op); - + void approximate(Dart d); } ; } //namespace DecimationVolumes diff --git a/include/Algo/DecimationVolumes/geometryApproximator.hpp b/include/Algo/DecimationVolumes/geometryApproximator.hpp index dcce94e21..c1bb93ec9 100644 --- a/include/Algo/DecimationVolumes/geometryApproximator.hpp +++ b/include/Algo/DecimationVolumes/geometryApproximator.hpp @@ -35,52 +35,85 @@ namespace DecimationVolumes { /************************************************************************************ - * Centroid * + * QUADRIC ERROR METRIC * ************************************************************************************/ + template -bool Approximator_Centroid::init() +bool Approximator_QEM::init() { - return true; + m_quadric = this->m_map.template getAttribute, VERTEX>("QEMquadric") ; + + if(this->m_predictor) + { + return false ; + } + return true ; } template -void Approximator_Centroid::approximate(Operator *op) +void Approximator_QEM::approximate(Operator *op) { Dart d = op->getEdge(); - VEC3 a = this->m_attrV[d]; + approximate(d); +} - //CGoGNout << "d=" << d << " a=" << a << CGoGNendl; +template +void Approximator_QEM::approximate(Dart d) +{ + MAP& m = this->m_map ; // get some darts - Dart dd = this->m_map.phi2(d) ; + Dart dd = m.phi2(d) ; - switch(op->getType()) + Quadric q1, q2 ; + if(!m_quadric.isValid()) // if the selector is not QEM, compute local error quadrics { - case O_CVolume : { - a = Algo::Geometry::volumeCentroid(this->m_map, d, this->m_attrV); - break; - } - case O_CFace : { - a = Algo::Geometry::faceCentroid(this->m_map, d, this->m_attrV); - break; - } - case O_CEdge : { - // get the contracted edge vertices positions - VEC3 v2 = this->m_attrV[dd] ; - - // Compute the approximated position - a = (a + v2) / REAL(2) ; - break; - } - default : { - CGoGNerr << "Approximate : reducer does not exist" << CGoGNendl; - break; - } + // compute the error quadric associated to v1 + Dart it = d ; + do + { + Quadric q(this->m_attrV[it], this->m_attrV[m.phi1(it)], this->m_attrV[m.phi_1(it)]) ; + q1 += q ; + it = m.phi2_1(it) ; + } while(it != d) ; + + // compute the error quadric associated to v2 + it = dd ; + do + { + Quadric q(this->m_attrV[it], this->m_attrV[m.phi1(it)], this->m_attrV[m.phi_1(it)]) ; + q2 += q ; + it = m.phi2_1(it) ; + } while(it != dd) ; + } + else // if the selector is QEM, use the error quadrics computed by the selector + { + q1 = m_quadric[d] ; + q2 = m_quadric[dd] ; } - this->m_app = a; -} + Quadric quad ; + quad += q1 ; // compute the sum of the + quad += q2 ; // two vertices quadrics + VEC3 res ; + bool opt = quad.findOptimizedPos(res) ; // try to compute an optimized position for the contraction of this edge + if(!opt) + { + VEC3 p1 = this->m_attrV[d] ; // let the new vertex lie + VEC3 p2 = this->m_attrV[dd] ; // on either one of the two endpoints + VEC3 p12 = (p1 + p2) / 2.0f ; // or the middle of the edge + REAL e1 = quad(p1) ; + REAL e2 = quad(p2) ; + REAL e12 = quad(p12) ; + REAL minerr = std::min(std::min(e1, e2), e12) ; // consider only the one for + if(minerr == e12) this->m_approx[d] = p12 ; // which the error is minimal + else if(minerr == e1) this->m_approx[d] = p1 ; + else this->m_approx[d] = p2 ; + } + else + this->m_approx[d] = res ; +} } //end namespace DecimationVolumes diff --git a/include/Algo/DecimationVolumes/operator.h b/include/Algo/DecimationVolumes/operator.h index 100240ea2..982287231 100644 --- a/include/Algo/DecimationVolumes/operator.h +++ b/include/Algo/DecimationVolumes/operator.h @@ -61,6 +61,9 @@ protected: // Dart m_edge; + Dart d2; + Dart dd2; + /** * need a pointer to the current approximator if the current selector needs * the future result of a collapse to estimate its cost @@ -82,6 +85,14 @@ public: virtual unsigned int perform(MAP& m, VertexAttribute& position) = 0; virtual bool canPerform(MAP &m ,Dart d, VertexAttribute& position) = 0; + + + + void setFirstIncidentEdge(Dart d) { d2 = d; } + void setSecondIncidentEdge(Dart d) { dd2 = d; } + + Dart getFirstIncidentEdge() { return d2; } + Dart getSecondIncidentEdge() { return dd2; } }; template @@ -113,13 +124,12 @@ public: typedef typename PFP::VEC3 VEC3 ; typedef typename PFP::REAL REAL; -protected: - - public: CollapseEdgeOperator(Dart d) ://, ApproximatorGen* approx) : CollapseOperator(d)//, approx) - {} + { + + } ~CollapseEdgeOperator() { } @@ -127,6 +137,7 @@ public: OperatorType getType() { return O_CEdge; } unsigned int perform(MAP &m, VertexAttribute& position); bool canPerform(MAP &m ,Dart d, VertexAttribute& position); + }; diff --git a/include/Algo/DecimationVolumes/operator.hpp b/include/Algo/DecimationVolumes/operator.hpp index 9ffebef97..1c7e26d46 100644 --- a/include/Algo/DecimationVolumes/operator.hpp +++ b/include/Algo/DecimationVolumes/operator.hpp @@ -39,10 +39,7 @@ unsigned int CollapseEdgeOperator::perform(MAP& m, VertexAttributem_edge] = this->m_approximator->getApproximation(); m.collapseEdge(this->m_edge); - ++nbCell; return nbCell; @@ -53,77 +50,53 @@ bool CollapseEdgeOperator::canPerform(MAP &m ,Dart d, VertexAttributem_approximator->getApproximation(); - typename PFP::VEC3 p2, p3, p4; - - typename PFP::VEC3::DATA_TYPE v1; - typename PFP::VEC3::DATA_TYPE v2; - - p2 = position[m.phi1(a)]; - p3 = position[m.phi_1(a)]; - p4 = position[m.phi_1(m.phi2(a))]; - - v1 = Geom::tetraSignedVolume(p1, p2, p3, p4); - if (v1 < 0) - canCollapse = false; - - //CGoGNout << " v1 = " << v1; - - //tetrahedre du bas - p2 = position[m.phi1(b)]; - p3 = position[m.phi_1(b)]; - p4 = position[m.phi_1(m.phi2(b))]; - - v2 = Geom::tetraSignedVolume(p1, p2, p3, p4); - - if (v2 < 0) - canCollapse = false; - - //CGoGNout << " v2 = " << v2; - } - } - - //CGoGNout << CGoGNendl; - - e = m.alpha2(e); - }while ( e != d && canCollapse); - - //CGoGNout << CGoGNendl << CGoGNendl; - - //CGoGNout << "is collapsable ? " << canCollapse << CGoGNendl; - +// Dart e = d; +// do +// { +// //isBoundaryVolume +// if(m.isBoundaryVolume(e)) +// { +// canCollapse = false; +// } +// //l'un de ces 2 voisins existe +// else if(m.phi3(m.phi2(m.phi1(e))) != m.phi2(m.phi1(e)) && m.phi3(m.phi2(m.phi_1(e))) != m.phi2(m.phi_1(e))) +// { +// //l'un de ces 2 voisins est au bord +// if(m.isBoundaryVolume(m.phi3(m.phi2(m.phi1(e)))) || m.isBoundaryVolume(m.phi3(m.phi2(m.phi_1(e))))) +// { +// canCollapse = false; +// +// } +// } +// else +// { +// //Edge Criteria Valide +// if(m.edgeDegree(m.phi1(m.phi2(m.phi_1(e)))) < 3) +// canCollapse = false; +// elseframe +// { +// //Is inverted +// Dart a = m.phi3(m.phi2(m.phi1(e))); +// Dart b = m.phi1(m.phi3(m.phi2(m.phi_1(e)))); +// +// typename PFP::VEC3::DATA_TYPE v1; +// typename PFP::VEC3::DATA_TYPE v2; +// +// v1 = Algo::Geometry::tetrahedronSignedVolume(m,a,position); +// +// if (v1 < 0) +// canCollapse = false; +// +// v2 = Algo::Geometry::tetrahedronSignedVolume(m,b,position); +// if (v2 < 0) +// canCollapse = false; +// +// //CGoGNout << " v2 = " << v2; +// } +// } +// +// e = m.alpha2(e); +// }while ( e != d && canCollapse); return canCollapse; } diff --git a/include/Algo/DecimationVolumes/selector.h b/include/Algo/DecimationVolumes/selector.h index e4f8ed32d..e28dc19be 100644 --- a/include/Algo/DecimationVolumes/selector.h +++ b/include/Algo/DecimationVolumes/selector.h @@ -15,9 +15,7 @@ namespace DecimationVolumes enum SelectorType { - S_MapOrder, - S_Random, - S_EdgeLength + S_QEM } ; template class ApproximatorGen ; @@ -85,24 +83,23 @@ public: /** * */ - virtual bool nextOperator(Operator** op) = 0 ; + virtual Operator* nextOperator() = 0 ; /** * */ - virtual void updateBeforeOperation(Operator** op) = 0 ; + virtual bool updateBeforeOperation(Operator* op) = 0 ; /** * */ - virtual void updateAfterOperation(Operator** op) = 0 ; + virtual void updateAfterOperation(Operator* op) = 0 ; /** * */ virtual void finish() = 0 ; - // /** // * // */ diff --git a/include/Algo/Import/import.h b/include/Algo/Import/import.h index fcf9a62a3..f9580faff 100644 --- a/include/Algo/Import/import.h +++ b/include/Algo/Import/import.h @@ -71,7 +71,10 @@ bool importMeshV(typename PFP::MAP& map, const std::string& filename, std::vecto * @return a boolean indicating if import was successful */ template -bool importMeshToExtrude(typename PFP::MAP& map, const std::string& filename, std::vector& attrNames); +bool importMeshToExtrude(typename PFP::MAP& map, const std::string& filename, std::vector& attrNames, float scale = 5.0f, unsigned int nbStage = 1); + +template +bool importMeshSAsV(typename PFP::MAP& map, const std::string& filename, std::vector& attrNames); /* * import a MOKA file diff --git a/include/Algo/Import/importMesh.hpp b/include/Algo/Import/importMesh.hpp index 12f92cf87..e11b2fbd0 100644 --- a/include/Algo/Import/importMesh.hpp +++ b/include/Algo/Import/importMesh.hpp @@ -237,7 +237,7 @@ bool importMeshSToV(typename PFP::MAP& map, MeshTablesSurface& mts, float d } template -bool importMeshSurfToVol(typename PFP::MAP& map, MeshTablesSurface& mts, float dist) +bool importMeshSurfToVol(typename PFP::MAP& map, MeshTablesSurface& mts, float scale, unsigned int nbStage) { VertexAutoAttribute< NoMathIONameAttribute< std::vector > > vecDartsPerVertex(map); unsigned nbf = mts.getNbFaces(); @@ -248,8 +248,10 @@ bool importMeshSurfToVol(typename PFP::MAP& map, MeshTablesSurface& mts, fl DartMarkerNoUnmark m(map) ; + unsigned int nbVertices = mts.getNbVertices(); + VertexAttribute position = map.template getAttribute("position"); - std::vector backEdgesBuffer(mts.getNbVertices(), EMBNULL); + std::vector backEdgesBuffer(nbVertices*nbStage, EMBNULL); // for each face of table -> create a prism for(unsigned int i = 0; i < nbf; ++i) @@ -275,37 +277,55 @@ bool importMeshSurfToVol(typename PFP::MAP& map, MeshTablesSurface& mts, fl nbe = edgesBuffer.size(); if (nbe > 2) { - Dart d = Algo::Modelisation::createPrism(map, nbe); + Dart dprev = NIL; - //Embed the base faces - for (unsigned int j = 0; j < nbe; ++j) + for(unsigned int k = 0 ; k < nbStage ; ++k) { - unsigned int em = edgesBuffer[j]; // get embedding + Dart d = Algo::Modelisation::createPrism(map, nbe); - if(backEdgesBuffer[em] == EMBNULL) + //Embed the base faces + for (unsigned int j = 0; j < nbe; ++j) { - unsigned int emn = map.template newCell(); - map.template copyCell(emn, em); - backEdgesBuffer[em] = emn; - position[emn] += typename PFP::VEC3(0,0,dist); + unsigned int em = edgesBuffer[j]; // get embedding + Dart d2 = map.phi1(map.phi1(map.phi2(d))); + + if(k==0) + { + FunctorSetEmb fsetemb(map, em); + map.template foreach_dart_of_orbit(d, fsetemb); + vecDartsPerVertex[em].push_back(d); // store incident darts for fast adjacency reconstruction + m.mark(d) ; // mark on the fly to unmark on second loop + } + else + { + unsigned int emn = backEdgesBuffer[((k-1)*nbVertices) + em]; + FunctorSetEmb fsetemb(map, emn); + map.template foreach_dart_of_orbit(d, fsetemb); + vecDartsPerVertex[emn].push_back(d); // store incident darts for fast adjacency reconstruction + m.mark(d) ; // mark on the fly to unmark on second loop + } + + if(backEdgesBuffer[(k*nbVertices) + em] == EMBNULL) + { + unsigned int emn = map.template newCell(); + map.template copyCell(emn, em); + backEdgesBuffer[(k*nbVertices) + em] = emn; + position[emn] += typename PFP::VEC3(0,0, (k+1) * scale); + } + + unsigned int em2 = backEdgesBuffer[(k*nbVertices) + em]; + FunctorSetEmb fsetemb(map, em2); + map.template foreach_dart_of_orbit(d2, fsetemb); + + d = map.phi_1(d); } - FunctorSetEmb fsetemb(map, em); - //foreach_dart_of_orbit_in_parent(&map, VERTEX, d, fsetemb) ; - map.template foreach_dart_of_orbit(d, fsetemb); - //Embed the other base face - Dart d2 = map.phi1(map.phi1(map.phi2(d))); - unsigned int em2 = backEdgesBuffer[em]; - FunctorSetEmb fsetemb2(map, em2); - //foreach_dart_of_orbit_in_parent(&map, VERTEX, d2, fsetemb2) ; - map.template foreach_dart_of_orbit(d2, fsetemb2); + if(dprev != NIL) + map.sewVolumes(d, map.phi2(map.phi1(map.phi1(map.phi2(dprev)))), false); - m.mark(d) ; // mark on the fly to unmark on second loop - vecDartsPerVertex[em].push_back(d); // store incident darts for fast adjacency reconstruction - d = map.phi_1(d); + dprev = d; } - } } @@ -344,6 +364,98 @@ bool importMeshSurfToVol(typename PFP::MAP& map, MeshTablesSurface& mts, fl return true ; } +template +bool importMeshSAsV(typename PFP::MAP& map, MeshTablesSurface& mts) +{ + VertexAutoAttribute< NoMathIONameAttribute< std::vector > > vecDartsPerVertex(map, "incidents"); + + unsigned nbf = mts.getNbFaces(); + int index = 0; + // buffer for tempo faces (used to remove degenerated edges) + std::vector edgesBuffer; + edgesBuffer.reserve(16); + + DartMarkerNoUnmark m(map) ; + + // for each face of table + for(unsigned int i = 0; i < nbf; ++i) + { + // store face in buffer, removing degenerated edges + unsigned int nbe = mts.getNbEdgesFace(i); + edgesBuffer.clear(); + unsigned int prec = EMBNULL; + for (unsigned int j = 0; j < nbe; ++j) + { + unsigned int em = mts.getEmbIdx(index++); + if (em != prec) + { + prec = em; + edgesBuffer.push_back(em); + } + } + // check first/last vertices + if (edgesBuffer.front() == edgesBuffer.back()) + edgesBuffer.pop_back(); + + // create only non degenerated faces + nbe = edgesBuffer.size(); + if (nbe > 2) + { + Dart d = map.newFace(nbe, false); + for (unsigned int j = 0; j < nbe; ++j) + { + unsigned int em = edgesBuffer[j]; // get embedding + + FunctorSetEmb fsetemb(map, em); +// foreach_dart_of_orbit_in_parent(&map, VERTEX, d, fsetemb) ; + map.template foreach_dart_of_orbit(d, fsetemb); + + m.mark(d) ; // mark on the fly to unmark on second loop + vecDartsPerVertex[em].push_back(d); // store incident darts for fast adjacency reconstruction + d = map.phi1(d); + } + } + } + + // reconstruct neighbourhood + unsigned int nbBoundaryEdges = 0; + for (Dart d = map.begin(); d != map.end(); map.next(d)) + { + if (m.isMarked(d)) + { + // darts incident to end vertex of edge + std::vector& vec = vecDartsPerVertex[map.phi1(d)]; + + unsigned int embd = map.template getEmbedding(d); + Dart good_dart = NIL; + for (typename std::vector::iterator it = vec.begin(); it != vec.end() && good_dart == NIL; ++it) + { + if (map.template getEmbedding(map.phi1(*it)) == embd) + good_dart = *it; + } + + if (good_dart != NIL) + { + map.sewFaces(d, good_dart, false); + m.unmarkOrbit(d); + } + else + { + m.unmark(d); + ++nbBoundaryEdges; + } + } + } + + unsigned int nbH = map.closeMap(); + CGoGNout << "Map closed (" << map.template getNbOrbits() << " boundary faces / " << nbH << " holes)" << CGoGNendl; + std::cout << "nb darts : " << map.getNbDarts() << std::endl ; + // ensure bijection between topo and embedding + //map.template bijectiveOrbitEmbedding(); + + return true ; +} + template bool importMesh(typename PFP::MAP& map, MeshTablesVolume& mtv) { @@ -416,15 +528,25 @@ bool importMeshV(typename PFP::MAP& map, const std::string& filename, std::vecto } template -bool importMeshToExtrude(typename PFP::MAP& map, const std::string& filename, std::vector& attrNames) +bool importMeshToExtrude(typename PFP::MAP& map, const std::string& filename, std::vector& attrNames, float scale, unsigned int nbStage) +{ + MeshTablesSurface mts(map); + + if(!mts.importMesh(filename, attrNames)) + return false; + + return importMeshSurfToVol(map, mts, scale, nbStage); +} + +template +bool importMeshSAsV(typename PFP::MAP& map, const std::string& filename, std::vector& attrNames) { - float dist = 1.0f; MeshTablesSurface mts(map); if(!mts.importMesh(filename, attrNames)) return false; - return importMeshSurfToVol(map, mts, dist); + return importMeshSAsV(map, mts); } } // namespace Import diff --git a/include/Algo/Modelisation/tetrahedralization.h b/include/Algo/Modelisation/tetrahedralization.h index 27d1b2b31..11ac8f5b2 100644 --- a/include/Algo/Modelisation/tetrahedralization.h +++ b/include/Algo/Modelisation/tetrahedralization.h @@ -50,6 +50,9 @@ void hexahedronToTetrahedron(typename PFP::MAP& map, Dart d); template void hexahedronsToTetrahedrons(typename PFP::MAP& map); +template +void tetrahedrizeVolume(typename PFP::MAP& map, VertexAttribute& position); + /************************************************************************************************ * Collapse / Split Operators ************************************************************************************************/ @@ -77,6 +80,8 @@ bool isTetrahedron(typename PFP::MAP& the_map, Dart d); template bool isTetrahedralization(typename PFP::MAP& map, const FunctorSelect& selected = allDarts); + + /************************************************************************************************ * Swap Functions * ************************************************************************************************/ diff --git a/include/Algo/Modelisation/tetrahedralization.hpp b/include/Algo/Modelisation/tetrahedralization.hpp index d0fda8d7d..ab6242441 100644 --- a/include/Algo/Modelisation/tetrahedralization.hpp +++ b/include/Algo/Modelisation/tetrahedralization.hpp @@ -91,6 +91,85 @@ void hexahedronsToTetrahedrons(typename PFP::MAP& map) } } +template +void tetrahedrizeVolume(typename PFP::MAP& map, VertexAttribute& position) +{ + //mark bad edges + DartMarker me(map); + + int i = 0; + + TraversorE travE(map); + + for(Dart dit = travE.begin() ; dit != travE.end() ; dit = travE.next()) + { + //check if this edge is an "ear-edge" + if(!me.isMarked(dit)) + { + //search three positions + typename PFP::VEC3 tris1[3]; + tris1[0] = position[dit]; + tris1[1] = position[map.phi_1(dit)]; + tris1[2] = position[map.phi_1(map.phi2(dit))]; + + //search if the triangle formed by these three points intersect the rest of the mesh (intersection triangle/triangle) + TraversorF travF(map); + for(Dart ditF = travF.begin() ; ditF != travF.end() ; ditF = travF.next()) + { + //get vertices position + typename PFP::VEC3 tris2[3]; + tris2[0] = position[ditF]; + tris2[1] = position[map.phi1(ditF)]; + tris2[2] = position[map.phi_1(ditF)]; + + bool intersection = false; + + for (unsigned int i = 0; i < 3 && !intersection; ++i) + { + typename PFP::VEC3 inter; + intersection = Geom::intersectionSegmentTriangle(tris1[i], tris1[(i+1)%3], tris2[0], tris2[1], tris2[2], inter); + } + + if(!intersection) + { + for (unsigned int i = 0; i < 3 && !intersection; ++i) + { + typename PFP::VEC3 inter; + intersection = Geom::intersectionSegmentTriangle(tris2[i], tris2[(i+1)%3], tris1[0], tris1[1], tris1[2], inter); + } + } + + std::cout << "intersection ? " << (intersection ? "true" : "false") << std::endl; + + if(intersection) + { + me.markOrbit(dit); + } + else //cut a tetrahedron + { + Dart dring = map.phi_1(dit); + std::vector vPath; + + vPath.push_back(map.phi_1(dring)); + vPath.push_back(map.phi1(map.phi2(dring))); + vPath.push_back(map.phi_1(map.phi2(dring))); + vPath.push_back(map.phi1(dring)); + + map.splitVolume(vPath); + + //map.splitFace(map.phi2(map.phi1(dit)), map.phi2(map.phi1(map.phi2(dit)))); + } + + ++i; + + if(i == 16) + return; + } + } + } + +} + /************************************************************************************************ * Collapse / Split Operators diff --git a/include/Algo/Selection/raySelector.h b/include/Algo/Selection/raySelector.h index d6570339b..2192321e2 100644 --- a/include/Algo/Selection/raySelector.h +++ b/include/Algo/Selection/raySelector.h @@ -264,6 +264,24 @@ void dartsRaySelection(typename PFP::MAP& map, const VertexAttribute +void facesPlanSelection(typename PFP::MAP& map, const VertexAttribute& position, + const typename Geom::Plane3D& plan, std::vector& vecDarts, + const FunctorSelect& good = allDarts) +{ + TraversorF travF(map); + + for(Dart dit = travF.begin() ; dit != travF.end() ; dit = travF.next() ) + { + if(Geom::intersectionTrianglePlan(position[dit], position[map.phi1(dit)], position[map.phi_1(dit)],plan.d(), plan.normal()) == Geom::FACE_INTERSECTION) + { + vecDarts.push_back(dit); + } + } + + std::cout << "nb faces = " << vecDarts.size() << std::endl; +} + } //namespace Selection } //namespace Algo diff --git a/include/Geometry/intersection.h b/include/Geometry/intersection.h index f3d524823..48851e7bc 100644 --- a/include/Geometry/intersection.h +++ b/include/Geometry/intersection.h @@ -163,6 +163,20 @@ Intersection intersectionPlaneRay(const PLANE3D& pl,const VEC3& p1,const VEC3& d template Intersection intersection2DSegmentSegment(const VEC3& PA, const VEC3& PB, const VEC3& QA, const VEC3& QB, VEC3& Inter) ; +template +Intersection intersectionSegmentPlan(const VEC3& PA, const VEC3& PB, const VEC3& PlaneP, const VEC3& NormP); //, VEC3& Inter) ; + +template +Intersection intersectionTrianglePlan(const VEC3& Ta, const VEC3& Tb, const VEC3& Tc, const VEC3& PlaneP, const VEC3& NormP);//, VEC3& Inter) ; + +template +Intersection intersectionSegmentHalfPlan(const VEC3& PA, const VEC3& PB, + const VEC3& P, const VEC3& DirP, const VEC3& OrientP);//, VEC3& Inter) + +template +Intersection intersectionTriangleHalfPlan(const VEC3& Ta, const VEC3& Tb, const VEC3& Tc, + const VEC3& P, const VEC3& DirP, const VEC3& OrientP); //, VEC3& Inter) ; + } // namespace Geom } // namespace CGoGN diff --git a/include/Geometry/intersection.hpp b/include/Geometry/intersection.hpp index d1afa409e..7d5619624 100644 --- a/include/Geometry/intersection.hpp +++ b/include/Geometry/intersection.hpp @@ -475,6 +475,88 @@ Intersection intersection2DSegmentSegment(const VEC3& PA, const VEC3& PB, const return EDGE_INTERSECTION; } +template +Intersection intersectionSegmentPlan(const VEC3& PA, const VEC3& PB, const VEC3& PlaneP, const VEC3& NormP)//, VEC3& Inter) +{ + typename VEC3::DATA_TYPE panp = NormP * PA; + typename VEC3::DATA_TYPE pbnp = NormP * PB; + + if(panp == 0 || pbnp == 0) + return VERTEX_INTERSECTION; + else if((panp < 0 && pbnp > 0) || (panp > 0 && pbnp < 0)) + return EDGE_INTERSECTION; + else + return NO_INTERSECTION; + +} + + +template +Intersection intersectionTrianglePlan(const VEC3& Ta, const VEC3& Tb, const VEC3& Tc, const VEC3& PlaneP, const VEC3& NormP) //, VEC3& Inter) ; +{ + if((intersectionSegmentPlan(Ta,Tb,PlaneP, NormP) == EDGE_INTERSECTION) + || (intersectionSegmentPlan(Ta,Tc,PlaneP, NormP) == EDGE_INTERSECTION) + || (intersectionSegmentPlan(Tb,Tc,PlaneP, NormP) == EDGE_INTERSECTION)) + { + return FACE_INTERSECTION; + } + else if((intersectionSegmentPlan(Ta,Tb,PlaneP, NormP) == VERTEX_INTERSECTION) + || (intersectionSegmentPlan(Ta,Tc,PlaneP, NormP) == VERTEX_INTERSECTION) + || (intersectionSegmentPlan(Tb,Tc,PlaneP, NormP) == VERTEX_INTERSECTION)) + { + return VERTEX_INTERSECTION; + } + else + { + return NO_INTERSECTION; + } +} + +template +Intersection intersectionSegmentHalfPlan(const VEC3& PA, const VEC3& PB, + const VEC3& P, const VEC3& DirP, const VEC3& OrientP)//, VEC3& Inter) +{ + VEC3 NormP = (DirP-P) ^ (OrientP-P) ; + NormP.normalize() ; + + //intersection SegmentPlan + Intersection inter = intersectionSegmentPlan(PA,PB,P,NormP); + if(inter == EDGE_INTERSECTION) + { + //and one of the two points must be in the right side of the line + return intersectionSegmentPlan(PA,PB, P, OrientP); + } + else + { + return inter; + } + + + +} + +template +Intersection intersectionTriangleHalfPlan(const VEC3& Ta, const VEC3& Tb, const VEC3& Tc, + const VEC3& P, const VEC3& DirP, const VEC3& OrientP) //, VEC3& Inter) +{ + if((intersectionSegmentHalfPlan(Ta,Tb,P, DirP, OrientP) == EDGE_INTERSECTION) + || (intersectionSegmentHalfPlan(Ta,Tc,P, DirP, OrientP) == EDGE_INTERSECTION) + || (intersectionSegmentHalfPlan(Tb,Tc,P, DirP, OrientP) == EDGE_INTERSECTION)) + { + return FACE_INTERSECTION; + } + else if((intersectionSegmentHalfPlan(Ta,Tb,P, DirP, OrientP) == VERTEX_INTERSECTION) + || (intersectionSegmentHalfPlan(Ta,Tc,P, DirP, OrientP) == VERTEX_INTERSECTION) + || (intersectionSegmentHalfPlan(Tb,Tb,P, DirP, OrientP) == VERTEX_INTERSECTION)) + { + return FACE_INTERSECTION; + } + else + { + return NO_INTERSECTION; + } +} + } } diff --git a/src/Topology/map/embeddedMap3.cpp b/src/Topology/map/embeddedMap3.cpp index 66f930966..98c42c792 100644 --- a/src/Topology/map/embeddedMap3.cpp +++ b/src/Topology/map/embeddedMap3.cpp @@ -23,6 +23,8 @@ *******************************************************************************/ #include "Topology/map/embeddedMap3.h" +#include +#include namespace CGoGN { @@ -153,21 +155,71 @@ Dart EmbeddedMap3::deleteEdge(Dart d) bool EmbeddedMap3::edgeCanCollapse(Dart d) { - //Criteres sur le bord - if(isBoundaryEdge(d)) +// //Criteres sur le bord +// if(isBoundaryEdge(d)) +// { +// //fusion de deux bords +// +// //deconnection du bord +// } +// +// return true; + + +// if(Map2::isBoundaryVertex(d) || Map2::isBoundaryVertex(phi1(d))) +// return false ; + + unsigned int val_v1 = Map2::vertexDegree(d) ; + unsigned int val_v2 = Map2::vertexDegree(phi1(d)) ; + + if(val_v1 + val_v2 < 8 || val_v1 + val_v2 > 14) + return false ; + + if(Map2::faceDegree(d) == 3) { - //fusion de deux bords + if(Map2::vertexDegree(phi_1(d)) < 4) + return false ; + } - //deconnection du bord + Dart dd = phi2(d) ; + if(Map2::faceDegree(dd) == 3) + { + if(Map2::vertexDegree(phi_1(dd)) < 4) + return false ; } - return false; + // Check vertex sharing condition + std::vector vu1 ; + vu1.reserve(256) ; + Dart vit1 = phi2(phi_1(phi2(phi_1(d)))) ; + Dart end = phi1(dd) ; + do + { + unsigned int ve = getEmbedding(phi2(vit1)) ; + vu1.push_back(ve) ; + vit1 = phi2(phi_1(vit1)) ; + } while(vit1 != end) ; + end = phi1(d) ; + Dart vit2 = phi2(phi_1(phi2(phi_1(dd)))) ; + do + { + unsigned int ve = getEmbedding(phi2(vit2)) ; + std::vector::iterator it = std::find(vu1.begin(), vu1.end(), ve) ; + if(it != vu1.end()) + return false ; + vit2 = phi2(phi_1(vit2)) ; + } while(vit2 != end) ; + + return true ; } Dart EmbeddedMap3::collapseEdge(Dart d, bool delDegenerateVolumes) { unsigned int vEmb = getEmbedding(d) ; + Dart d2 = phi2(phi_1(d)) ; + Dart dd2 = phi2(phi_1(phi2(d))) ; + Dart resV = Map3::collapseEdge(d, delDegenerateVolumes); if(resV != NIL) @@ -176,16 +228,17 @@ Dart EmbeddedMap3::collapseEdge(Dart d, bool delDegenerateVolumes) { embedOrbit(resV, vEmb); } + + if(isOrbitEmbedded()) + { + embedOrbit(d2, getEmbedding(d2)); + embedOrbit(dd2, getEmbedding(dd2)); + } } return resV; } -//bool EmbeddedMap3::collapseDegeneratedFace(Dart d) -//{ -// return Map3::collapseDegeneratedFace(d); -//} - void EmbeddedMap3::splitFace(Dart d, Dart e) { Dart dd = phi1(phi3(d)); -- GitLab