Commit 80f7c7ec authored by Pierre Kraemer's avatar Pierre Kraemer

Merge cgogn:~cgogn/CGoGN

parents 32117a01 3b87573c
......@@ -49,6 +49,8 @@
#include "Algo/Render/GL2/topo3Render.h"
#include "Utils/Shaders/shaderSimpleColor.h"
#include "Utils/cgognStream.h"
#include "Utils/drawer.h"
......
/*******************************************************************************
* CGoGN: Combinatorial and Geometric modeling with Generic N-dimensional Maps *
* version 0.1 *
* Copyright (C) 2009-2011, IGG Team, LSIIT, University of Strasbourg *
* *
* This library is free software; you can redistribute it and/or modify it *
* under the terms of the GNU Lesser General Public License as published by the *
* Free Software Foundation; either version 2.1 of the License, or (at your *
* option) any later version. *
* *
* This library is distributed in the hope that it will be useful, but WITHOUT *
* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License *
* for more details. *
* *
* You should have received a copy of the GNU Lesser General Public License *
* along with this library; if not, write to the Free Software Foundation, *
* Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. *
* *
* Web site: http://cgogn.u-strasbg.fr/ *
* Contact information: cgogn@unistra.fr *
* *
*******************************************************************************/
#ifndef __ALGO_BOOLEANOPERATOR_VERTICES_H__
#define __ALGO_BOOLEANOPERATOR_VERTICES_H__
#include "Geometry/basic.h"
#include "Geometry/inclusion.h"
#include "Geometry/orientation.h"
namespace CGoGN
{
namespace Algo
{
namespace BooleanOperator
{
template <typename PFP>
void mergeVertex(typename PFP::MAP& map, const typename PFP::TVEC3& positions, Dart d, Dart e);
template <typename PFP>
void mergeVertices(typename PFP::MAP& map, const typename PFP::TVEC3& positions);
}
}
}
#include "mergeVertices.hpp"
#endif
/*******************************************************************************
* CGoGN: Combinatorial and Geometric modeling with Generic N-dimensional Maps *
* version 0.1 *
* Copyright (C) 2009-2011, IGG Team, LSIIT, University of Strasbourg *
* *
* This library is free software; you can redistribute it and/or modify it *
* under the terms of the GNU Lesser General Public License as published by the *
* Free Software Foundation; either version 2.1 of the License, or (at your *
* option) any later version. *
* *
* This library is distributed in the hope that it will be useful, but WITHOUT *
* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License *
* for more details. *
* *
* You should have received a copy of the GNU Lesser General Public License *
* along with this library; if not, write to the Free Software Foundation, *
* Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. *
* *
* Web site: http://cgogn.u-strasbg.fr/ *
* Contact information: cgogn@unistra.fr *
* *
*******************************************************************************/
namespace CGoGN
{
namespace Algo
{
namespace BooleanOperator
{
template <typename PFP>
void mergeVertex(typename PFP::MAP& map, const typename PFP::TVEC3& positions, Dart d, Dart e)
{
assert(Geom::arePointsEquals(positions[d],positions[e]) && !map.sameVertex(d,e));
Dart dd;
do
{
dd = map.alpha1(d);
map.removeEdgeFromVertex(dd);
Dart ee = e;
do
{
if(Geom::testOrientation2D(positions[map.phi1(dd)],positions[ee],positions[map.phi1(ee)])!=Geom::RIGHT
&& Geom::testOrientation2D(positions[map.phi1(dd)],positions[ee],positions[map.phi1(map.alpha1(ee))])==Geom::RIGHT)
{
break;
}
ee = map.alpha1(ee);
} while(ee != e);
map.insertEdgeInVertex(ee,dd);
} while(dd!=d);
}
template <typename PFP>
void mergeVertices(typename PFP::MAP& map, const typename PFP::TVEC3& positions)
{
for(Dart d = map.begin() ; d != map.end() ; map.next(d))
{
CellMarker vM(map,VERTEX);
vM.mark(d);
for(Dart dd = map.begin() ; dd != map.end() ; map.next(dd))
{
if(!vM.isMarked(dd))
{
vM.mark(dd);
if(Geom::arePointsEquals(positions[d],positions[dd]))
{
mergeVertex<PFP>(map,positions,d,dd);
}
}
}
}
}
}
}
}
......@@ -95,7 +95,7 @@ void exportMeshWire(std::ofstream& out, typename PFP::MAP& map, typename PFP::TV
if(position[dd][0]!=position[map.phi1(dd)][0] || position[dd][1]!=position[map.phi1(dd)][1] || position[dd][2]!=position[map.phi1(dd)][2]) {
out << "cylinder{ " << std::endl;
out << "<" << position[dd][0] << "," << position[dd][2] << "," << position[dd][1] << ">," << std::endl;
out << "<" << position[map.phi1(dd)][0] << "," << position[map.phi1(dd)][2] << "," << position[map.phi1(dd)][1] << ">, 0.5" << std::endl;
out << "<" << position[map.phi1(dd)][0] << "," << position[map.phi1(dd)][2] << "," << position[map.phi1(dd)][1] << ">, 1.5" << std::endl;
out << "}" << std::endl;
}
dd = map.phi1(dd);
......
......@@ -24,6 +24,7 @@
#include "Algo/Geometry/normal.h"
#include "Algo/Geometry/centroid.h"
#include "intersection.h"
#include <limits>
......@@ -37,17 +38,18 @@ namespace Geometry
{
template <typename PFP>
bool intersectionLineConvexFace(typename PFP::MAP& map, Dart d, const typename PFP::TVEC3& position, const typename PFP::VEC3& P, const typename PFP::VEC3& Dir, typename PFP::VEC3& Inter)
bool intersectionLineConvexFace(typename PFP::MAP& map, Dart d, const typename PFP::TVEC3& positions, const typename PFP::VEC3& P, const typename PFP::VEC3& Dir, typename PFP::VEC3& Inter)
{
typedef typename PFP::VEC3 VEC3 ;
const float SMALL_NUM = std::numeric_limits<typename PFP::REAL>::min() * 5.0f;
VEC3 p1 = position[d];
VEC3 n = faceNormal<PFP>(map, d, position);
VEC3 p1 = positions[d];
VEC3 n = faceNormal<PFP>(map,d,positions);
VEC3 w0 = P - p1;
float a = -(n * w0);
float b = n * Dir;
float a = -(n*w0);
float b = n*Dir;
if (fabs(b) < SMALL_NUM)
return false;
......@@ -56,7 +58,7 @@ bool intersectionLineConvexFace(typename PFP::MAP& map, Dart d, const typename P
Inter = P + r * Dir; // intersect point of ray and plane
// is I inside the face?
VEC3 p2 = position[map.phi1(d)];
VEC3 p2 = positions[map.phi1(d)];
VEC3 v = p2 - p1 ;
VEC3 vInter = Inter - p1;
float dirV = v * vInter;
......@@ -67,7 +69,7 @@ bool intersectionLineConvexFace(typename PFP::MAP& map, Dart d, const typename P
while(it != d)
{
p1 = p2;
p2 = position[map.phi1(it)];
p2 = positions[map.phi1(it)];
v = p2 - p1;
vInter = Inter - p1;
float dirD = v * vInter;
......@@ -83,22 +85,22 @@ bool intersectionLineConvexFace(typename PFP::MAP& map, Dart d, const typename P
}
template <typename PFP>
bool intersectionSegmentConvexFace(typename PFP::MAP& map, Dart d, const typename PFP::TVEC3& position, const typename PFP::VEC3& PA, const typename PFP::VEC3& PB, typename PFP::VEC3& Inter)
bool intersectionSegmentConvexFace(typename PFP::MAP& map, Dart d, const typename PFP::TVEC3& positions, const typename PFP::VEC3& PA, const typename PFP::VEC3& PB, typename PFP::VEC3& Inter)
{
typename PFP::VEC3 dir = PB - PA;
if (intersectionLineConvexFace(map, d, position, PA, dir, Inter))
if (intersectionLineConvexFace(map,d,positions,PA,dir,Inter))
{
typename PFP::VEC3 dirA = PA - Inter;
typename PFP::VEC3 dirB = PB - Inter;
typename PFP::VEC3 dirB = PB -Inter;
if (dirA * dirB < 0)
if ( (dirA*dirB) < 0 )
return true;
}
return false;
}
template <typename PFP>
bool areTrianglesInIntersection(typename PFP::MAP& map, Dart tri1, Dart tri2, const typename PFP::TVEC3& position)
bool areTrianglesInIntersection(typename PFP::MAP& map, Dart tri1, Dart tri2, const typename PFP::TVEC3& positions)
{
typedef typename PFP::VEC3 VEC3 ;
......@@ -107,8 +109,8 @@ bool areTrianglesInIntersection(typename PFP::MAP& map, Dart tri1, Dart tri2, co
VEC3 tris2[3];
for (unsigned int i = 0; i < 3; ++i)
{
tris1[i] = position[tri1];
tris2[i] = position[tri2];
tris1[i] = positions[tri1];
tris2[i] = positions[tri2];
tri1 = map.phi1(tri1);
tri2 = map.phi1(tri2);
}
......@@ -165,8 +167,8 @@ bool areTrianglesInIntersection(typename PFP::MAP& map, Dart tri1, Dart tri2, co
// return isSegmentInTriangle2D(inter1,inter2,tris2[0],tris2[1],triS2[2],nTri2);
//compute face normal
VEC3 normale1 = faceNormal<PFP>(map, tri1, position);
VEC3 bary1 = faceCentroid<PFP>(map, tri1, position);
VEC3 normale1 = faceNormal<PFP>(map,tri1,positions);
VEC3 bary1 = faceCentroid<PFP>(map,tri1,positions);
int pos = 0;
int neg = 0;
......@@ -186,8 +188,8 @@ bool areTrianglesInIntersection(typename PFP::MAP& map, Dart tri1, Dart tri2, co
return false;
//same for the second triangle
VEC3 normale2 = faceNormal<PFP>(map, tri2, position);
VEC3 bary2 = faceCentroid<PFP>(map, tri2, position);
VEC3 normale2 = faceNormal<PFP>(map,tri2,positions);
VEC3 bary2 = faceCentroid<PFP>(map,tri2,positions);
pos = 0;
neg = 0;
for (unsigned int i = 0; i < 3 ; ++i)
......@@ -208,7 +210,7 @@ bool areTrianglesInIntersection(typename PFP::MAP& map, Dart tri1, Dart tri2, co
for (unsigned int i = 0; i < 3 && !intersection; ++i)
{
VEC3 inter;
intersection = Geom::intersectionSegmentTriangle(tris1[i], tris1[(i+1)%3], tris2[0], tris2[1], tris2[2], inter);
intersection = Geom::intersectionSegmentTriangle(tris1[i],tris1[(i+1)%3],tris2[0],tris2[1],tris2[2],inter);
}
if (intersection)
......@@ -217,34 +219,14 @@ bool areTrianglesInIntersection(typename PFP::MAP& map, Dart tri1, Dart tri2, co
for (unsigned int i = 0; i < 3 && !intersection; ++i)
{
VEC3 inter;
intersection = Geom::intersectionSegmentTriangle(tris2[i], tris2[(i+1)%3], tris1[0], tris1[1], tris1[2], inter);
intersection = Geom::intersectionSegmentTriangle(tris2[i],tris2[(i+1)%3],tris1[0],tris1[1],tris1[2],inter);
}
return intersection;
}
template <typename PFP>
bool intersectionSphereEdge(typename PFP::MAP& map, typename PFP::VEC3& center, typename PFP::REAL radius, Dart d, const typename PFP::TVEC3& position, typename PFP::REAL& alpha)
{
typedef typename PFP::VEC3 VEC3 ;
typedef typename PFP::REAL REAL ;
const typename PFP::VEC3& p1 = position[d];
const typename PFP::VEC3& p2 = position[map.phi1(d)];
if(Geom::isPointInSphere(p1, center, radius) && !Geom::isPointInSphere(p2, center, radius))
{
VEC3 p = p1 - center;
VEC3 qminusp = p2 - center - p;
REAL s = p * qminusp;
REAL n2 = qminusp.norm2();
alpha = (- s + sqrt(s*s + n2 * (radius*radius - p.norm2()))) / n2;
return true ;
}
return false ;
}
} // namespace Geometry
} // namespace Algo
}
} // namespace CGoGN
}
......@@ -240,34 +240,37 @@ inline void ImplicitHierarchicalMap3::next(Dart& d)
inline bool ImplicitHierarchicalMap3::foreach_dart_of_vertex(Dart d, FunctorType& f, unsigned int thread)
{
DartMarkerStore mv(*this); // Lock a marker
DartMarkerStore mv(*this,thread); // Lock a marker
bool found = false; // Last functor return value
std::list<Dart> darts_list; //Darts that are traversed
darts_list.push_back(d); //Start with the dart d
std::list<Dart>::iterator darts;
std::vector<Dart> darts_list; //Darts that are traversed
darts_list.reserve(50);
darts_list.push_back(d); //Start with the dart d
mv.mark(d);
for(darts = darts_list.begin(); !found && darts != darts_list.end() ; ++darts)
for(std::vector<Dart>::iterator darts = darts_list.begin(); !found && darts != darts_list.end() ; ++darts)
{
Dart dc = *darts;
//add phi21 and phi23 successor if they are not marked yet
Dart d2 = phi2(dc);
Dart d21 = phi1(d2); // turn in volume
Dart d23 = phi3(d2); // change volume
if(!mv.isMarked(d21))
if(d2 != dc)
{
darts_list.push_back(d21);
mv.mark(d21);
}
Dart d21 = phi1(d2); // turn in volume
Dart d23 = phi3(d2); // change volume
if((d23!=d2) && !mv.isMarked(d23))
{
darts_list.push_back(d23);
mv.mark(d23);
if(!mv.isMarked(d21))
{
darts_list.push_back(d21);
mv.mark(d21);
}
if((d23!=d2) && !mv.isMarked(d23))
{
darts_list.push_back(d23);
mv.mark(d23);
}
}
found = f(dc);
......
......@@ -62,6 +62,15 @@ void coarsenFace(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& position);
template <typename PFP>
void coarsenVolume(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& position);
/***********************************************************************************
* Raffinement
***********************************************************************************/
/*
* Un brin de la face oppose aux faces a spliter
*/
template <typename PFP>
void splitVolume(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& position);
/*******************************************************
*
*/
......
......@@ -85,7 +85,9 @@ void subdivideFace(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& position
++degree;
p += position[it] ;
if(!map.edgeIsSubdivided(it)) // first cut the edges (if they are not already)
{
Algo::IHM::subdivideEdge<PFP>(map, it, position) ; // and compute the degree of the face
}
it = map.phi1(it) ;
} while(it != old) ;
p /= typename PFP::REAL(degree) ;
......@@ -94,38 +96,38 @@ void subdivideFace(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& position
Dart res;
if(degree == 3 && sType == S_TRI) //subdiviser une face triangulaire
{
Dart dd = map.phi1(old) ;
Dart e = map.phi1(map.phi1(dd)) ;
map.splitFace(dd, e) ; // insert a new edge
unsigned int id = map.getNewEdgeId() ;
map.setEdgeId(map.phi_1(dd), id, EDGE) ; // set the edge id of the inserted edge to the next available id
unsigned int idface = map.getFaceId(old);
map.setFaceId(dd, idface, FACE) ;
map.setFaceId(e, idface, FACE) ;
dd = e ;
e = map.phi1(map.phi1(dd)) ;
map.splitFace(dd, e) ;
id = map.getNewEdgeId() ;
map.setEdgeId(map.phi_1(dd), id, EDGE) ;
map.setFaceId(dd, idface, FACE) ;
map.setFaceId(e, idface, FACE) ;
dd = e ;
e = map.phi1(map.phi1(dd)) ;
map.splitFace(dd, e) ;
id = map.getNewEdgeId() ;
map.setEdgeId(map.phi_1(dd), id, EDGE) ;
map.setFaceId(dd, idface, FACE) ;
map.setFaceId(e, idface, FACE) ;
}
else
{
// if(degree == 3 && sType == S_TRI) //subdiviser une face triangulaire
// {
// Dart dd = map.phi1(old) ;
// Dart e = map.phi1(map.phi1(dd)) ;
// map.splitFace(dd, e) ; // insert a new edge
// unsigned int id = map.getNewEdgeId() ;
// map.setEdgeId(map.phi_1(dd), id, EDGE) ; // set the edge id of the inserted edge to the next available id
//
// unsigned int idface = map.getFaceId(old);
// map.setFaceId(dd, idface, FACE) ;
// map.setFaceId(e, idface, FACE) ;
//
// dd = e ;
// e = map.phi1(map.phi1(dd)) ;
// map.splitFace(dd, e) ;
// id = map.getNewEdgeId() ;
// map.setEdgeId(map.phi_1(dd), id, EDGE) ;
//
// map.setFaceId(dd, idface, FACE) ;
// map.setFaceId(e, idface, FACE) ;
//
// dd = e ;
// e = map.phi1(map.phi1(dd)) ;
// map.splitFace(dd, e) ;
// id = map.getNewEdgeId() ;
// map.setEdgeId(map.phi_1(dd), id, EDGE) ;
//
// map.setFaceId(dd, idface, FACE) ;
// map.setFaceId(e, idface, FACE) ;
// }
// else
// {
Dart dd = map.phi1(old) ;
map.splitFace(dd, map.phi1(map.phi1(dd))) ;
......@@ -162,7 +164,7 @@ void subdivideFace(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& position
dd = map.phi2(map.phi1(dd));
}
while(dd != ne);
}
// }
map.setCurrentLevel(cur) ;
}
......@@ -173,12 +175,14 @@ Dart subdivideVolume(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& positi
assert(map.getDartLevel(d) <= map.getCurrentLevel() || !"Access to a dart introduced after current level") ;
assert(!map.volumeIsSubdivided(d) || !"Trying to subdivide an already subdivided volume") ;
unsigned int vLevel = map.volumeLevel(d) ;
Dart old = map.volumeOldestDart(d) ;
unsigned int vLevel = map.volumeLevel(d);
Dart old = map.volumeOldestDart(d);
unsigned int cur = map.getCurrentLevel();
map.setCurrentLevel(vLevel);
unsigned int cur = map.getCurrentLevel() ;
map.setCurrentLevel(vLevel) ; // go to the level of the face to subdivide its edges
std::cout << "vLevel = " << vLevel << std::endl;
/*
* au niveau du volume courant i
......@@ -187,7 +191,8 @@ Dart subdivideVolume(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& positi
*/
DartMarkerStore mf(map); // Lock a face marker to save one dart per face
DartMarkerStore mv(map); // Lock a vertex marker to compute volume center
//DartMarkerStore mv(map); // Lock a vertex marker to compute volume center
CellMarker mv(map, VERTEX);
typename PFP::VEC3 volCenter;
unsigned count = 0 ;
......@@ -212,9 +217,10 @@ Dart subdivideVolume(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& positi
//compute volume centroid
if(!mv.isMarked(e))
{
//mv.markOrbit(VERTEX, e);
mv.mark(e);
volCenter += position[e];
++count;
mv.markOrbit(VERTEX, e);
oldEdges.push_back(e);
}
......@@ -240,6 +246,7 @@ Dart subdivideVolume(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& positi
std::vector<std::pair<Dart,Dart> > subdividedfaces;
subdividedfaces.reserve(25);
int i = 0;
//First step : subdivide edges and faces
//creates a i+1 edge level and i+1 face level
for (std::vector<Dart>::iterator face = visitedFaces.begin(); face != visitedFaces.end(); ++face)
......@@ -248,12 +255,24 @@ Dart subdivideVolume(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& positi
//if needed subdivide face
if(!map.faceIsSubdivided(d))
{
++i;
Algo::IHM::subdivideFace<PFP>(map, d, position);
}
std::cout << "CHECK FACES " << std::endl;
map.check();
//save a dart from the subdivided face
unsigned int cur = map.getCurrentLevel() ;
unsigned int fLevel = map.faceLevel(d) + 1; //puisque dans tous les cas, la face est subdivisee
map.setCurrentLevel(map.getMaxLevel());
Dart old = map.faceOldestDart(d);
std::cout << " Oldest dart Face = " << old << std::endl;
unsigned int fLevel = map.faceLevel(old); //puisque dans tous les cas, la face est subdivisee
std::cout << "FACE LEVEL = " << fLevel << std::endl;
map.setCurrentLevel(fLevel) ;
......@@ -269,6 +288,8 @@ Dart subdivideVolume(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& positi
map.setCurrentLevel(cur);
}
map.setCurrentLevel(vLevel + 1) ; // go to the next level to perform volume subdivision
std::vector<Dart> newEdges; //save darts from inner edges
......@@ -280,21 +301,33 @@ Dart subdivideVolume(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& positi
Dart e = *edge;
Dart f1 = map.phi1(*edge);
int i = 0;
do
{
map.unsewFaces(map.phi1(map.phi1(e)));
std::cout << "corner #" << i << " degree = " << map.faceDegree(e) << std::endl;
if(map.phi2(map.phi1(map.phi1(e))) != map.phi1(map.phi1(e)))
map.unsewFaces(map.phi1(map.phi1(e)));
//TODO utile ?
//if(map.phi2(map.phi1(e)) != map.phi1(e))
if(map.phi2(map.phi1(e)) != map.phi1(e))
map.unsewFaces(map.phi1(e));
++i;
e = map.phi2(map.phi_1(e));
}
while(e != *edge);
map.closeHole(f1);
std::cout << "CHECK CLOSE HOLE " << std::endl;
map.check();
std::cout << "corner degree = " << i << std::endl;
std::cout << "face degree = " << map.faceDegree(map.phi2(f1)) << std::endl;
// if(map.faceDegree(map.phi2(f1)) > 6)
// return Dart::nil();
Dart old = map.phi2(map.phi1(e));
Dart dd = map.phi1(map.phi1(old)) ;
map.splitFace(old,dd) ;
......@@ -362,6 +395,8 @@ Dart subdivideVolume(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& positi
map.setCurrentLevel(cur) ;
std::cout << std::endl;
return subdividedfaces.begin()->first;
}
......@@ -740,10 +775,81 @@ void coarsenVolume(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& position
}
/***********************************************************************************
* Raffinement
***********************************************************************************/
template <typename PFP>
void splitVolume(typename PFP::MAP& map, Dart d, typename PFP::TVEC3& position)
{