Commit 87f29bb0 authored by untereiner's avatar untereiner

compute dual app, some changes in Map3MR

parent eb68fc91
...@@ -98,13 +98,11 @@ public: ...@@ -98,13 +98,11 @@ public:
bool m_drawTopo ; bool m_drawTopo ;
VertexAttribute<VEC3> position ; VertexAttribute<VEC3> position ;
VertexAttribute<VEC3> normal ;
Algo::Render::GL2::MapRender* m_render ; Algo::Render::GL2::MapRender* m_render ;
Algo::Render::GL2::TopoRender* m_topoRender ; Algo::Render::GL2::TopoRender* m_topoRender ;
Utils::VBO* m_positionVBO ; Utils::VBO* m_positionVBO ;
Utils::VBO* m_normalVBO ;
Utils::ShaderPhong* m_phongShader ; Utils::ShaderPhong* m_phongShader ;
Utils::ShaderFlat* m_flatShader ; Utils::ShaderFlat* m_flatShader ;
......
...@@ -57,15 +57,24 @@ int main(int argc, char **argv) ...@@ -57,15 +57,24 @@ int main(int argc, char **argv)
std::string filename(argv[1]); std::string filename(argv[1]);
size_t pos = filename.rfind("."); // position of "." in filename
std::string extension = filename.substr(pos);
// declaration of the map // declaration of the map
PFP::MAP myMap; PFP::MAP myMap;
std::vector<std::string> attrNames ; std::vector<std::string> attrNames ;
Algo::Volume::Import::importMesh<PFP>(myMap, argv[1], attrNames);
if (extension == std::string(".map"))
myMap.loadMapBin(filename);
else
Algo::Volume::Import::importMesh<PFP>(myMap, filename, attrNames);
// get a handler to the 3D vector attribute created by the import // get a handler to the 3D vector attribute created by the import
VertexAttribute<PFP::VEC3> position = myMap.getAttribute<PFP::VEC3, VERTEX>(attrNames[0]); VertexAttribute<PFP::VEC3> position = myMap.getAttribute<PFP::VEC3, VERTEX>(attrNames[0]);
// Les faces vont devenir des aretes -> echange de FACE ET EDGE
VolumeAttribute<PFP::VEC3> positionV = myMap.getAttribute<PFP::VEC3, VOLUME>("position") ; VolumeAttribute<PFP::VEC3> positionV = myMap.getAttribute<PFP::VEC3, VOLUME>("position") ;
if(!positionV.isValid()) if(!positionV.isValid())
positionV = myMap.addAttribute<PFP::VEC3, VOLUME>("position") ; positionV = myMap.addAttribute<PFP::VEC3, VOLUME>("position") ;
...@@ -81,21 +90,73 @@ int main(int argc, char **argv) ...@@ -81,21 +90,73 @@ int main(int argc, char **argv)
break; break;
} }
} }
if(dsave != NIL)
positionV[dsave] = PFP::VEC3(0.0);
//Algo::Modelisation::computeDual3<PFP>(myMap,allDarts) ; Dart dcenter = myMap.explodBorderTopo(dsave);
myMap.computeDual();
DartMarker mf(myMap);
for(Dart dit = myMap.begin() ; dit != myMap.end() ; myMap.next(dit))
{
if(myMap.isBoundaryMarked3(dit) && !mf.isMarked(dit))
{
mf.markOrbit<FACE>(dit);
positionV[dit] = Algo::Surface::Geometry::faceCentroid<PFP>(myMap,dit,position);
}
}
for(Dart dit = myMap.begin() ; dit != myMap.end() ; myMap.next(dit))
{
if(myMap.isBoundaryMarked3(dit))
{
myMap.fillHole(dit);
}
}
//
//Compute Dual Test -- begin
//
DartAttribute<Dart> old_phi1 = myMap.getAttribute<Dart, DART>("phi1") ;
DartAttribute<Dart> old_phi_1 = myMap.getAttribute<Dart, DART>("phi_1") ;
DartAttribute<Dart> new_phi1 = myMap.addAttribute<Dart, DART>("new_phi1") ;
DartAttribute<Dart> new_phi_1 = myMap.addAttribute<Dart, DART>("new_phi_1") ;
DartAttribute<Dart> old_phi2 = myMap.getAttribute<Dart, DART>("phi2") ;
DartAttribute<Dart> new_phi2 = myMap.addAttribute<Dart, DART>("new_phi2") ;
for(Dart d = myMap.begin(); d != myMap.end(); myMap.next(d))
{
Dart dd = myMap.phi2(myMap.phi3(d)) ;
new_phi1[d] = dd ;
new_phi_1[dd] = d ;
Dart ddd = myMap.phi1(myMap.phi3(d));
new_phi2[d] = ddd;
new_phi2[ddd] = d;
}
myMap.swapAttributes<Dart>(old_phi1, new_phi1) ;
myMap.swapAttributes<Dart>(old_phi_1, new_phi_1) ;
myMap.swapAttributes<Dart>(old_phi2, new_phi2) ;
myMap.removeAttribute(new_phi1) ;
myMap.removeAttribute(new_phi_1) ;
myMap.removeAttribute(new_phi2) ;
myMap.swapEmbeddingContainers(VERTEX, VOLUME) ;
//
//ComputeDualTest -- end
//
myMap.createHole(dcenter);
VolumeAttribute<PFP::VEC3> del;
del = position;
position = positionV ; position = positionV ;
//turn_to<VertexAttribute<PFP::VEC3>, FaceAttribute<PFP::VEC3> >(&positionF); myMap.removeAttribute(del);
//turn_to<FaceAttribute<PFP::VEC3>, VertexAttribute<PFP::VEC3> >(position);
//const std::type_info &t1 = typeid(&positionF); myMap.dumpAttributesAndMarkers();
//std::cout << "type name : " << t1.name() << std::endl;
//Algo::Export::exportOFF<PFP>(myMap, position, "result.off");
myMap.saveMapBin("result.map"); myMap.saveMapBin("result.map");
std::cout << "Exported" << std::endl; std::cout << "Exported" << std::endl;
......
...@@ -46,7 +46,108 @@ namespace Primal ...@@ -46,7 +46,108 @@ namespace Primal
namespace Masks namespace Masks
{ {
//
//template <typename PFP>
//class MJ96VertexVertexFunctor : public FunctorType
//{
//protected:
// typename PFP::MAP& m_map ;
// VertexAttribute<typename PFP::VEC3>& m_position ;
//
//public:
// MJ96VertexVertexFunctor(typename PFP::MAP& m, VertexAttribute<typename PFP::VEC3>& p) : m_map(m), m_position(p)
// {}
//
// bool operator() (Dart d)
// {
// if(m_map.isBoundaryVertex(d))
// {
// std::cout << "boundary" << std::endl;
//
// m_map.decCurrentLevel() ;
//
// Dart db = m_map.findBoundaryFaceOfVertex(d);
//
// typename PFP::VEC3 np1(0) ;
// typename PFP::VEC3 np2(0) ;
// unsigned int degree1 = 0 ;
// unsigned int degree2 = 0 ;
// Dart it = db ;
// do
// {
// ++degree1 ;
// Dart dd = m_map.phi1(it) ;
// np1 += m_position[dd] ;
// Dart end = m_map.phi_1(it) ;
// dd = m_map.phi1(dd) ;
// do
// {
// ++degree2 ;
// np2 += m_position[dd] ;
// dd = m_map.phi1(dd) ;
// } while(dd != end) ;
// it = m_map.phi2(m_map.phi_1(it)) ;
// } while(it != db) ;
//
// float beta = 3.0 / (2.0 * degree1) ;
// float gamma = 1.0 / (4.0 * degree2) ;
// np1 *= beta / degree1 ;
// np2 *= gamma / degree2 ;
//
// typename PFP::VEC3 vp = m_position[db] ;
// vp *= 1.0 - beta - gamma ;
//
// m_map.incCurrentLevel() ;
//
// m_position[d] = np1 + np2 + vp ;
//
// }
// else
// {
// m_map.decCurrentLevel() ;
// typename PFP::VEC3 p = m_position[d] ;
// m_map.incCurrentLevel() ;
//
// m_position[d] = p ;
//// m_map.decCurrentLevel() ;
////
//// typename PFP::VEC3 P = m_position[d];
////
//// //vertex points
//// typename PFP::VEC3 Cavg = typename PFP::VEC3(0);
//// unsigned int degree = m_map.template degree<typename PFP::MAP, VERTEX, VOLUME>(d);
//// Cavg = Algo::Surface::Geometry::volumeCentroid<PFP>(m_map, d, m_position);
//// Cavg /= degree;
////
//// //P /= degree;
////
//// typename PFP::VEC3 Aavg = typename PFP::VEC3(0);
//// degree = m_map.template degree<typename PFP::MAP, VERTEX, FACE>(d);
//// Aavg = Algo::Surface::Geometry::faceCentroid<PFP>(m_map, d, m_position);
//// Aavg /= degree;
////
//// typename PFP::VEC3 Mavg = typename PFP::VEC3(0);
//// degree = m_map.template degree<typename PFP::MAP, VERTEX, EDGE>(d);
//// Dart d2 = m_map.phi2(d);
//// Mavg = (m_position[d] + m_position[d2]) * typename PFP::REAL(0.5);
//// Mavg /= degree;
////
//// typename PFP::VEC3 vp = Cavg + Aavg * 3 + Mavg * 3 + P;
//// vp /= 8;
////
//// m_map.incCurrentLevel() ;
////
//// m_position[d] = P; // += vp; //P;
// }
//
// return false;
// }
//};
/* MJ96 basic functions : polyhedral meshes
*********************************************************************************/
template <typename PFP> template <typename PFP>
class MJ96VertexVertexFunctor : public FunctorType class MJ96VertexVertexFunctor : public FunctorType
{ {
...@@ -62,8 +163,6 @@ public: ...@@ -62,8 +163,6 @@ public:
{ {
if(m_map.isBoundaryVertex(d)) if(m_map.isBoundaryVertex(d))
{ {
std::cout << "boundary" << std::endl;
m_map.decCurrentLevel() ; m_map.decCurrentLevel() ;
Dart db = m_map.findBoundaryFaceOfVertex(d); Dart db = m_map.findBoundaryFaceOfVertex(d);
...@@ -110,21 +209,34 @@ public: ...@@ -110,21 +209,34 @@ public:
//vertex points //vertex points
typename PFP::VEC3 Cavg = typename PFP::VEC3(0); typename PFP::VEC3 Cavg = typename PFP::VEC3(0);
unsigned int degree = m_map.template degree<typename PFP::MAP, VERTEX, VOLUME>(d); unsigned int degree = 0;
Cavg = Algo::Surface::Geometry::volumeCentroid<PFP>(m_map, d, m_position); Traversor3VW<typename PFP::MAP> travVW(m_map, d);
for(Dart dit = travVW.begin() ; dit != travVW.end() ; dit = travVW.next())
{
Cavg += Algo::Surface::Geometry::volumeCentroid<PFP>(m_map, dit, m_position);
++degree;
}
Cavg /= degree; Cavg /= degree;
//P /= degree;
typename PFP::VEC3 Aavg = typename PFP::VEC3(0); typename PFP::VEC3 Aavg = typename PFP::VEC3(0);
degree = m_map.template degree<typename PFP::MAP, VERTEX, FACE>(d); degree = 0;
Aavg = Algo::Surface::Geometry::faceCentroid<PFP>(m_map, d, m_position); Traversor3VF<typename PFP::MAP> travVF(m_map, d);
for(Dart dit = travVF.begin() ; dit != travVF.end() ; dit = travVF.next())
{
Aavg += Algo::Surface::Geometry::faceCentroid<PFP>(m_map, dit, m_position);
++degree;
}
Aavg /= degree; Aavg /= degree;
typename PFP::VEC3 Mavg = typename PFP::VEC3(0); typename PFP::VEC3 Mavg = typename PFP::VEC3(0);
degree = m_map.template degree<typename PFP::MAP, VERTEX, EDGE>(d); degree = 0;
Dart d2 = m_map.phi2(d); Traversor3VE<typename PFP::MAP> travVE(m_map, d);
Mavg = (m_position[d] + m_position[d2]) * typename PFP::REAL(0.5); for(Dart dit = travVE.begin() ; dit != travVE.end() ; dit = travVE.next())
{
Dart d2 = m_map.phi2(dit);
Mavg += (m_position[dit] + m_position[d2]) * typename PFP::REAL(0.5);
++degree;
}
Mavg /= degree; Mavg /= degree;
typename PFP::VEC3 vp = Cavg + Aavg * 3 + Mavg * 3 + P; typename PFP::VEC3 vp = Cavg + Aavg * 3 + Mavg * 3 + P;
...@@ -132,120 +244,13 @@ public: ...@@ -132,120 +244,13 @@ public:
m_map.incCurrentLevel() ; m_map.incCurrentLevel() ;
m_position[d] = P; // += vp; //P; m_position[d] = vp;
} }
return false; return false;
} }
}; };
///* MJ96 basic functions : polyhedral meshes
// *********************************************************************************/
//template <typename PFP>
//class MJ96VertexVertexFunctor : public FunctorType
//{
//protected:
// typename PFP::MAP& m_map ;
// VertexAttribute<typename PFP::VEC3>& m_position ;
//
//public:
// MJ96VertexVertexFunctor(typename PFP::MAP& m, VertexAttribute<typename PFP::VEC3>& p) : m_map(m), m_position(p)
// {}
//
// bool operator() (Dart d)
// {
// if(m_map.isBoundaryVertex(d))
// {
// m_map.decCurrentLevel() ;
//
// Dart db = m_map.findBoundaryFaceOfVertex(d);
//
// typename PFP::VEC3 np1(0) ;
// typename PFP::VEC3 np2(0) ;
// unsigned int degree1 = 0 ;
// unsigned int degree2 = 0 ;
// Dart it = db ;
// do
// {
// ++degree1 ;
// Dart dd = m_map.phi1(it) ;
// np1 += m_position[dd] ;
// Dart end = m_map.phi_1(it) ;
// dd = m_map.phi1(dd) ;
// do
// {
// ++degree2 ;
// np2 += m_position[dd] ;
// dd = m_map.phi1(dd) ;
// } while(dd != end) ;
// it = m_map.phi2(m_map.phi_1(it)) ;
// } while(it != db) ;
//
// float beta = 3.0 / (2.0 * degree1) ;
// float gamma = 1.0 / (4.0 * degree2) ;
// np1 *= beta / degree1 ;
// np2 *= gamma / degree2 ;
//
// typename PFP::VEC3 vp = m_position[db] ;
// vp *= 1.0 - beta - gamma ;
//
// m_map.incCurrentLevel() ;
//
// m_position[d] = np1 + np2 + vp ;
//
// }
// else
// {
// m_map.decCurrentLevel() ;
//
// typename PFP::VEC3 P = m_position[d];
//
// //vertex points
// typename PFP::VEC3 Cavg = typename PFP::VEC3(0);
// unsigned int degree = 0;
// Traversor3VW<typename PFP::MAP> travVW(m_map, d);
// for(Dart dit = travVW.begin() ; dit != travVW.end() ; dit = travVW.next())
// {
// Cavg += Algo::Surface::Geometry::volumeCentroid<PFP>(m_map, dit, m_position);
// ++degree;
// }
// Cavg /= degree;
//
// typename PFP::VEC3 Aavg = typename PFP::VEC3(0);
// degree = 0;
// Traversor3VF<typename PFP::MAP> travVF(m_map, d);
// for(Dart dit = travVF.begin() ; dit != travVF.end() ; dit = travVF.next())
// {
// Aavg += Algo::Surface::Geometry::faceCentroid<PFP>(m_map, dit, m_position);
// ++degree;
// }
// Aavg /= degree;
//
// typename PFP::VEC3 Mavg = typename PFP::VEC3(0);
// degree = 0;
// Traversor3VE<typename PFP::MAP> travVE(m_map, d);
// for(Dart dit = travVE.begin() ; dit != travVE.end() ; dit = travVE.next())
// {
// Dart d2 = m_map.phi2(dit);
// Mavg += (m_position[dit] + m_position[d2]) * typename PFP::REAL(0.5);
// ++degree;
// }
// Mavg /= degree;
//
// typename PFP::VEC3 vp = Cavg + Aavg * 3 + Mavg * 3 + P;
// vp /= 8;
//
// m_map.incCurrentLevel() ;
//
// m_position[d] = vp;
// }
//
// return false;
// }
//};
template <typename PFP> template <typename PFP>
class MJ96EdgeVertexFunctor : public FunctorType class MJ96EdgeVertexFunctor : public FunctorType
{ {
...@@ -368,22 +373,30 @@ public: ...@@ -368,22 +373,30 @@ public:
} }
else else
{ {
Dart df = m_map.phi_1(m_map.phi_1(d)) ; Dart df = m_map.phi1(m_map.phi1(d)) ;
m_map.decCurrentLevel() ; m_map.decCurrentLevel() ;
typename PFP::VEC3 p = Algo::Surface::Geometry::faceCentroid<PFP>(m_map, df, m_position);
//face points
typename PFP::VEC3 C0 = Algo::Surface::Geometry::volumeCentroid<PFP>(m_map, df, m_position);
typename PFP::VEC3 C1 = Algo::Surface::Geometry::volumeCentroid<PFP>(m_map, m_map.phi3(df), m_position);
typename PFP::VEC3 A = Algo::Surface::Geometry::faceCentroid<PFP>(m_map, m_map.phi3(df), m_position);
typename PFP::VEC3 fp = C0 + A * 2 + C1;
fp /= 4;
m_map.incCurrentLevel() ; m_map.incCurrentLevel() ;
m_position[d] = fp; m_position[d] = p ;
// Dart df = m_map.phi_1(m_map.phi_1(d)) ;
//
// m_map.decCurrentLevel() ;
//
// //face points
// typename PFP::VEC3 C0 = Algo::Surface::Geometry::volumeCentroid<PFP>(m_map, df, m_position);
// typename PFP::VEC3 C1 = Algo::Surface::Geometry::volumeCentroid<PFP>(m_map, m_map.phi3(df), m_position);
//
// typename PFP::VEC3 A = Algo::Surface::Geometry::faceCentroid<PFP>(m_map, m_map.phi3(df), m_position);
//
// typename PFP::VEC3 fp = C0 + A * 2 + C1;
// fp /= 4;
//
// m_map.incCurrentLevel() ;
//
// m_position[d] = fp;
} }
return false; return false;
......
...@@ -167,6 +167,8 @@ public: ...@@ -167,6 +167,8 @@ public:
*/ */
unsigned int subdivideVolume(Dart d, bool triQuad = true, bool OneLevelDifference = true); unsigned int subdivideVolume(Dart d, bool triQuad = true, bool OneLevelDifference = true);
unsigned int subdivideHexa(Dart d, bool OneLevelDifference = true);
//! Subdivide the volume of d to hexahedral cells //! Subdivide the volume of d to hexahedral cells
/*! @param d Dart from the volume /*! @param d Dart from the volume
*/ */
......
...@@ -843,98 +843,106 @@ unsigned int Map3MR<PFP>::subdivideVolume(Dart d, bool triQuad, bool OneLevelDif ...@@ -843,98 +843,106 @@ unsigned int Map3MR<PFP>::subdivideVolume(Dart d, bool triQuad, bool OneLevelDif
return vLevel; return vLevel;
} }
template <typename PFP>
unsigned int Map3MR<PFP>::subdivideHexa(Dart d, bool OneLevelDifference)
{
assert(m_map.getDartLevel(d) <= m_map.getCurrentLevel() || !"subdivideVolume : called with a dart inserted after current level") ;
assert(!volumeIsSubdivided(d) || !"Trying to subdivide an already subdivided face") ;
unsigned int vLevel = volumeLevel(d);
Dart old = volumeOldestDart(d);
m_map.pushLevel() ;
m_map.setCurrentLevel(vLevel) ; // go to the level of the volume to subdivide its faces
if(m_map.getCurrentLevel() == m_map.getMaxLevel())
m_map.addLevelBack() ;
//
// Subdivide Faces and Edges
//
Traversor3WF<typename PFP::MAP> traF(m_map, old);
for(Dart dit = traF.begin(); dit != traF.end(); dit = traF.next())
{
//if needed subdivide face
if(!faceIsSubdivided(dit))
subdivideFace(dit,false);
}
std::vector<std::pair<Dart, Dart> > subdividedFaces;
subdividedFaces.reserve(128);
Dart centralDart = NIL;
Traversor3WV<typename PFP::MAP> traWV(m_map, d);
for(Dart ditWV = traWV.begin(); ditWV != traWV.end(); ditWV = traWV.next())
{
m_map.incCurrentLevel() ;
(*vertexVertexFunctor)(ditWV) ;
Dart e = ditWV;
std::vector<Dart> v ;
do
{
v.push_back(m_map.phi1(e));
v.push_back(m_map.phi1(m_map.phi1(e)));
if(!m_map.PFP::MAP::ParentMap::isBoundaryEdge(m_map.phi1(e)))
subdividedFaces.push_back(std::pair<Dart,Dart>(m_map.phi1(e),m_map.phi2(m_map.phi1(e))));
if(m_map.phi1(m_map.phi1(m_map.phi1(e))) != e)
if(!m_map.PFP::MAP::ParentMap::isBoundaryEdge(m_map.phi1(m_map.phi1(e))))
subdividedFaces.push_back(std::pair<Dart,Dart>(m_map.phi1(m_map.phi1(e)),m_map.phi2(m_map.phi1(m_map.phi1(e)))));
e = m_map.phi2(m_map.phi_1(e));
}
while(e != ditWV);
m_map.splitVolume(v);
Dart dd = m_map.phi2(m_map.phi1(ditWV));;
Dart next = m_map.phi1(m_map.phi1(dd)) ;
m_map.splitFace(dd, next) ; // insert a first edge
Dart ne = m_map.phi2(m_map.phi_1(dd)) ;
m_map.cutEdge(ne) ; // cut the new edge to insert the central vertex
centralDart = m_map.phi1(ne);
dd = m_map.phi1(m_map.phi1(next)) ;
while(dd != ne) // turn around the face and insert new edges
{ // linked to the central vertex
Dart tmp = m_map.phi1(ne) ;
m_map.splitFace(tmp, dd) ;
dd = m_map.phi1(m_map.phi1(dd)) ;
}
m_map.decCurrentLevel() ;
}
m_map.incCurrentLevel();
m_map.deleteVolume(m_map.phi3(m_map.phi2(m_map.phi1(d))));
for (std::vector<std::pair<Dart,Dart> >::iterator it = subdividedFaces.begin(); it != subdividedFaces.end(); ++it)
{
Dart f1 = m_map.phi2((*it).first);
Dart f2 = m_map.phi2((*it).second);
if(m_map.isBoundaryFace(f1) && m_map.isBoundaryFace(f2))
{
m_map.sewVolumes(f1, f2);//, false);
}
}
(*volumeVertexFunctor)(centralDart) ;
m_map.popLevel() ;
return vLevel;
}
// Traversor3WV<typename PFP::MAP> traWV(m_map, d);
// for(Dart ditWV = traWV.begin(); ditWV != traWV.end(); ditWV = traWV.next())
// {
// m_map.setCurrentLevel(m_map.getMaxLevel()) ;
//
// Dart e = ditWV;
// std::vector<Dart> v ;
//
// do
// {
// v.push_back(m_map.phi1(e));
// v.push_back(m_map.phi1(m_map.phi1(e)));
//
// if(!m_map.PFP::MAP::ParentMap::isBoundaryEdge(m_map.phi1(e)))
// subdividedFaces.push_back(std::pair<Dart,Dart>(m_map.phi1(e),m_map.phi2(m_map.phi1(e))));
//
// if(!m_map.PFP::MAP::ParentMap::isBoundaryEdge(m_map.phi1(m_map.phi1(e))))
// subdividedFaces.push_back(std::pair<Dart,Dart>(m_map.phi1(m_map.phi1(e)),m_map.phi2(m_map.phi1(m_map.phi1(e)))));
//
// e = m_map.phi2(m_map.phi_1(e));
// }
// while(e != ditWV);
//
// splitSurfaceInVolume(v);
//
// Dart dd = m_map.phi2(m_map.phi1(ditWV));
// Dart next = m_map.phi1(m_map.phi1(dd)) ;
// m_map.PFP::MAP::ParentMap::splitFace(dd, next) ;
//
// Dart ne = m_map.phi2(m_map.phi_1(dd));
// m_map.PFP::MAP::ParentMap::cutEdge(ne) ;
// centralDart = m_map.phi1(ne);
//
// dd = m_map.phi1(m_map.phi1(next)) ;
// while(dd != ne)