Commit 0acabc38 authored by Pierre Kraemer's avatar Pierre Kraemer

remove normal cycle computation from collectors

parent 49adae74
......@@ -354,6 +354,7 @@ private:
VertexAttribute<VEC3, MAP> normal ;
EdgeAttribute<EdgeInfo, MAP> edgeInfo ;
EdgeAttribute<REAL, MAP> edgeangle ;
EdgeAttribute<REAL, MAP> edgearea ;
VertexAttribute<REAL, MAP> kmax ;
VertexAttribute<REAL, MAP> kmin ;
VertexAttribute<VEC3, MAP> Kmax ;
......@@ -390,6 +391,13 @@ public:
Algo::Surface::Geometry::computeAnglesBetweenNormalsOnEdges<PFP>(m, pos, edgeangle) ;
}
edgearea = m.template getAttribute<REAL, EDGE, MAP>("edgearea") ;
if(!edgearea.isValid())
{
edgearea = m.template addAttribute<REAL, EDGE, MAP>("edgearea") ;
Algo::Surface::Geometry::computeAreaEdges<PFP>(m, pos, edgearea) ;
}
kmax = m.template getAttribute<REAL, VERTEX, MAP>("kmax") ;
kmin = m.template getAttribute<REAL, VERTEX, MAP>("kmin") ;
Kmax = m.template getAttribute<VEC3, VERTEX, MAP>("Kmax") ;
......@@ -404,7 +412,7 @@ public:
Kmax = m.template addAttribute<VEC3, VERTEX, MAP>("Kmax") ;
Kmin = m.template addAttribute<VEC3, VERTEX, MAP>("Kmin") ;
Knormal = m.template addAttribute<VEC3, VERTEX, MAP>("Knormal") ;
Algo::Surface::Geometry::computeCurvatureVertices_NormalCycles<PFP>(m, radius, pos, normal, edgeangle, kmax, kmin, Kmax, Kmin, Knormal) ;
Algo::Surface::Geometry::computeCurvatureVertices_NormalCycles<PFP>(m, radius, pos, normal, edgeangle, edgearea, kmax, kmin, Kmax, Kmin, Knormal) ;
}
edgeInfo = m.template addAttribute<EdgeInfo, EDGE, MAP>("edgeInfo") ;
......@@ -412,6 +420,7 @@ public:
~EdgeSelector_Curvature()
{
// this->m_map.removeAttribute(edgeangle) ;
// this->m_map.removeAttribute(edgearea) ;
// this->m_map.removeAttribute(kmax) ;
// this->m_map.removeAttribute(kmin) ;
// this->m_map.removeAttribute(Kmax) ;
......@@ -452,6 +461,7 @@ private:
EdgeAttribute<EdgeInfo, MAP> edgeInfo ;
EdgeAttribute<REAL, MAP> edgeangle ;
EdgeAttribute<REAL, MAP> edgearea ;
std::multimap<float,Dart> edges ;
typename std::multimap<float,Dart>::iterator cur ;
......@@ -473,11 +483,19 @@ public:
Algo::Surface::Geometry::computeAnglesBetweenNormalsOnEdges<PFP>(m, pos, edgeangle) ;
}
edgearea = m.template getAttribute<REAL, EDGE, MAP>("edgearea") ;
if(!edgearea.isValid())
{
edgearea = m.template addAttribute<REAL, EDGE, MAP>("edgearea") ;
Algo::Surface::Geometry::computeAreaEdges<PFP>(m, pos, edgearea) ;
}
edgeInfo = m.template addAttribute<EdgeInfo, EDGE, MAP>("edgeInfo") ;
}
~EdgeSelector_CurvatureTensor()
{
// this->m_map.removeAttribute(edgeangle) ;
// this->m_map.removeAttribute(edgearea) ;
this->m_map.removeAttribute(edgeInfo) ;
}
SelectorType getType() { return S_CurvatureTensor ; }
......
......@@ -1022,14 +1022,14 @@ void EdgeSelector_Curvature<PFP>::updateAfterCollapse(Dart d2, Dart dd2)
MAP& m = this->m_map ;
normal[d2] = Algo::Surface::Geometry::vertexNormal<PFP>(m, d2, m_position) ;
Algo::Surface::Geometry::computeCurvatureVertex_NormalCycles<PFP>(m, d2, radius, m_position, normal, edgeangle, kmax, kmin, Kmax, Kmin, Knormal) ;
Algo::Surface::Geometry::computeCurvatureVertex_NormalCycles<PFP>(m, d2, radius, m_position, normal, edgeangle, edgearea, kmax, kmin, Kmax, Kmin, Knormal) ;
Dart vit = d2 ;
do
{
Dart nVert = m.phi1(vit) ;
normal[nVert] = Algo::Surface::Geometry::vertexNormal<PFP>(m, nVert, m_position) ;
Algo::Surface::Geometry::computeCurvatureVertex_NormalCycles<PFP>(m, nVert, radius, m_position, normal, edgeangle, kmax, kmin, Kmax, Kmin, Knormal) ;
Algo::Surface::Geometry::computeCurvatureVertex_NormalCycles<PFP>(m, nVert, radius, m_position, normal, edgeangle, edgearea, kmax, kmin, Kmax, Kmin, Knormal) ;
updateEdgeInfo(m.phi1(vit), false) ; // must recompute some edge infos in the
if(vit == d2 || vit == dd2) // neighborhood of the collapsed edge
......@@ -1131,7 +1131,7 @@ void EdgeSelector_Curvature<PFP>::computeEdgeInfo(Dart d, EdgeInfo& einfo)
// compute things on the coarse version of the mesh
normal[newV] = Algo::Surface::Geometry::vertexNormal<PFP>(m, d2, m_position) ;
Algo::Surface::Geometry::computeCurvatureVertex_NormalCycles<PFP>(m, d2, radius, m_position, normal, edgeangle, kmax, kmin, Kmax, Kmin, Knormal) ;
Algo::Surface::Geometry::computeCurvatureVertex_NormalCycles<PFP>(m, d2, radius, m_position, normal, edgeangle, edgearea, kmax, kmin, Kmax, Kmin, Knormal) ;
// VEC3 norm = normal[newV] ;
REAL mCurv = (kmax[newV] + kmin[newV]) / REAL(2) ;
......@@ -1330,7 +1330,7 @@ void EdgeSelector_CurvatureTensor<PFP>::computeEdgeInfo(Dart d, EdgeInfo& einfo)
MATRIX tens1;
Algo::Surface::Selection::Collector_OneRing_AroundEdge<PFP> col1(m);
col1.collectAll(d);
col1.computeNormalCyclesTensor(m_position, edgeangle, tens1); // edgeangle is up to date here
Algo::Surface::Geometry::normalCycles_computeTensor(col1, m_position, edgeangle, tens1); // edgeangle is up to date here
tens1 *= col1.computeArea(m_position); // mean tensor * area = integral of the tensor
Algo::Surface::Geometry::normalCycles_SortTensor<PFP>(tens1);
......@@ -1345,7 +1345,7 @@ void EdgeSelector_CurvatureTensor<PFP>::computeEdgeInfo(Dart d, EdgeInfo& einfo)
MATRIX tens2;
Algo::Surface::Selection::Collector_OneRing<PFP> col2(m);
col2.collectAll(d);
col2.computeNormalCyclesTensor(m_position, tens2); // edgeangle is not up to date here
Algo::Surface::Geometry::normalCycles_computeTensor(col2, m_position, tens2); // edgeangle is not up to date here
tens2 *= col2.computeArea(m_position); // mean tensor * area = integral of the tensor
Algo::Surface::Geometry::normalCycles_SortTensor<PFP>(tens2);
......
......@@ -58,7 +58,8 @@ void computeCurvatureVertices_QuadraticFitting(
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin) ;
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin
);
template <typename PFP>
void computeCurvatureVertex_QuadraticFitting(
......@@ -69,7 +70,8 @@ void computeCurvatureVertex_QuadraticFitting(
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin) ;
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin
);
template <typename PFP>
void vertexQuadraticFitting(
......@@ -78,20 +80,23 @@ void vertexQuadraticFitting(
typename PFP::MATRIX33& localFrame,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
typename PFP::REAL& a, typename PFP::REAL& b, typename PFP::REAL& c, typename PFP::REAL& d, typename PFP::REAL& e) ;
typename PFP::REAL& a, typename PFP::REAL& b, typename PFP::REAL& c, typename PFP::REAL& d, typename PFP::REAL& e
);
template <typename PFP>
void quadraticFittingAddVertexPos(
typename PFP::VEC3& v,
typename PFP::VEC3& p,
typename PFP::MATRIX33& localFrame) ;
typename PFP::MATRIX33& localFrame
);
template <typename PFP>
void quadraticFittingAddVertexNormal(
typename PFP::VEC3& v,
typename PFP::VEC3& n,
typename PFP::VEC3& p,
typename PFP::MATRIX33& localFrame) ;
typename PFP::MATRIX33& localFrame
);
/*
template <typename PFP>
......@@ -118,7 +123,8 @@ void computeCurvatureVertices_NormalCycles(
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal) ;
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal
);
template <typename PFP>
void computeCurvatureVertex_NormalCycles(
......@@ -133,27 +139,8 @@ void computeCurvatureVertex_NormalCycles(
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal) ;
template <typename PFP>
void normalCycles_SortAndSetEigenComponents(
const typename PFP::VEC3& e_val,
const Geom::Matrix<3,3,typename PFP::REAL> & e_vec,
typename PFP::REAL& kmax,
typename PFP::REAL& kmin,
typename PFP::VEC3& Kmax,
typename PFP::VEC3& Kmin,
typename PFP::VEC3& Knormal,
const typename PFP::VEC3& normal) ;
template <typename PFP>
void normalCycles_SortTensor(
Geom::Matrix<3,3,typename PFP::REAL>& tensor) ;
template <typename PFP>
void normalCycles_ProjectTensor(
Geom::Matrix<3,3,typename PFP::REAL>& tensor,
const typename PFP::VEC3& normal_vector) ;
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal
);
template <typename PFP>
void computeCurvatureVertices_NormalCycles_Projected(
......@@ -167,7 +154,8 @@ void computeCurvatureVertices_NormalCycles_Projected(
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal) ;
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal
);
template <typename PFP>
void computeCurvatureVertex_NormalCycles_Projected(
......@@ -182,70 +170,128 @@ void computeCurvatureVertex_NormalCycles_Projected(
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal) ;
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal
);
/* normal cycles with collector as a parameter : not usable in parallel */
template <typename PFP>
void computeCurvatureVertices_NormalCycles(
typename PFP::MAP& map,
Algo::Surface::Selection::Collector<PFP> & neigh,
void normalCycles_computeTensor(
Algo::Surface::Selection::Collector<PFP>& col,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal) ;
Geom::Matrix<3,3,typename PFP::REAL>& tensor
);
template <typename PFP>
void computeCurvatureVertex_NormalCycles(
Vertex v,
Algo::Surface::Selection::Collector<PFP> & neigh,
void normalCycles_computeTensor(
Algo::Surface::Selection::Collector<PFP>& col,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal) ;
Geom::Matrix<3,3,typename PFP::REAL>& tensor
);
template <typename PFP>
void computeCurvatureVertices_NormalCycles_Projected(
typename PFP::MAP& map,
Algo::Surface::Selection::Collector<PFP> & neigh,
void normalCycles_computeTensor(
Algo::Surface::Selection::Collector<PFP>& col,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal) ;
Geom::Matrix<3,3,typename PFP::REAL>& tensor
);
template <typename PFP>
void computeCurvatureVertex_NormalCycles_Projected(
Vertex v,
Algo::Surface::Selection::Collector<PFP> & neigh,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal) ;
void normalCycles_SortAndSetEigenComponents(
const typename PFP::VEC3& e_val,
const Geom::Matrix<3,3,typename PFP::REAL> & e_vec,
typename PFP::REAL& kmax,
typename PFP::REAL& kmin,
typename PFP::VEC3& Kmax,
typename PFP::VEC3& Kmin,
typename PFP::VEC3& Knormal,
const typename PFP::VEC3& normal
);
template <typename PFP>
void normalCycles_SortTensor(Geom::Matrix<3,3,typename PFP::REAL>& tensor);
template <typename PFP>
void normalCycles_ProjectTensor(
Geom::Matrix<3,3,typename PFP::REAL>& tensor,
const typename PFP::VEC3& normal_vector
);
/* normal cycles with collector as a parameter : not usable in parallel */
//template <typename PFP>
//void computeCurvatureVertices_NormalCycles(
// typename PFP::MAP& map,
// Algo::Surface::Selection::Collector<PFP> & neigh,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal) ;
//template <typename PFP>
//void computeCurvatureVertex_NormalCycles(
// Vertex v,
// Algo::Surface::Selection::Collector<PFP> & neigh,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal) ;
//template <typename PFP>
//void computeCurvatureVertices_NormalCycles_Projected(
// typename PFP::MAP& map,
// Algo::Surface::Selection::Collector<PFP> & neigh,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal) ;
//template <typename PFP>
//void computeCurvatureVertex_NormalCycles_Projected(
// Vertex v,
// Algo::Surface::Selection::Collector<PFP> & neigh,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal) ;
namespace Parallel
{
//template <typename PFP>
//void computeCurvatureVertices_QuadraticFitting(
// typename PFP::MAP& map,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin);
template <typename PFP>
void computeCurvatureVertices_NormalCycles(
typename PFP::MAP& map,
......@@ -274,16 +320,6 @@ void computeCurvatureVertices_NormalCycles_Projected(
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal) ;
template <typename PFP>
void computeCurvatureVertices_QuadraticFitting(
typename PFP::MAP& map,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin);
} // namespace Parallel
......
......@@ -49,22 +49,19 @@ void computeCurvatureVertices_QuadraticFitting(
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin)
{
// TODO: nl not thread safe
if (CGoGN::Parallel::NumberOfThreads > 1)
{
Parallel::computeCurvatureVertices_QuadraticFitting<PFP>(map, position, normal, kmax, kmin, Kmax, Kmin);
return;
}
// if (CGoGN::Parallel::NumberOfThreads > 1)
// {
// Parallel::computeCurvatureVertices_QuadraticFitting<PFP>(map, position, normal, kmax, kmin, Kmax, Kmin);
// return;
// }
foreach_cell<VERTEX>(map, [&] (Vertex v)
{
computeCurvatureVertex_QuadraticFitting<PFP>(map, v, position, normal, kmax, kmin, Kmax, Kmin) ;
}
, FORCE_CELL_MARKING);
// TraversorV<typename PFP::MAP> t(map) ;
// for(Vertex v = t.begin(); v != t.end(); v = t.next())
// computeCurvatureVertex_QuadraticFitting<PFP>(map, v, position, normal, kmax, kmin, Kmax, Kmin) ;
}
template <typename PFP>
......@@ -352,7 +349,7 @@ void computeCurvatureVertex_NormalCycles(
neigh.collectAll(v) ;
MATRIX tensor(0) ;
neigh.computeNormalCyclesTensor(position, edgeangle, edgearea, tensor);
normalCycles_computeTensor(neigh, position, edgeangle, edgearea, tensor);
// solve eigen problem
Eigen::SelfAdjointEigenSolver<E_MATRIX> solver(Utils::convertRef<E_MATRIX>(tensor));
......@@ -376,11 +373,11 @@ void computeCurvatureVertices_NormalCycles_Projected(
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal)
{
// if (CGoGN::Parallel::NumberOfThreads > 1)
// {
// Parallel::computeCurvatureVertices_NormalCycles_Projected<PFP>(map, radius, position, normal, edgeangle, kmax, kmin, Kmax, Kmin, Knormal);
// return;
// }
if (CGoGN::Parallel::NumberOfThreads > 1)
{
Parallel::computeCurvatureVertices_NormalCycles_Projected<PFP>(map, radius, position, normal, edgeangle, edgearea, kmax, kmin, Kmax, Kmin, Knormal);
return;
}
foreach_cell<VERTEX>(map, [&] (Vertex v)
{
......@@ -415,7 +412,7 @@ void computeCurvatureVertex_NormalCycles_Projected(
neigh.collectAll(v) ;
MATRIX tensor(0) ;
neigh.computeNormalCyclesTensor(position, edgeangle, edgearea, tensor);
normalCycles_computeTensor(neigh, position, edgeangle, edgearea, tensor);
// project the tensor
normalCycles_ProjectTensor<PFP>(tensor, normal[v]);
......@@ -428,113 +425,196 @@ void computeCurvatureVertex_NormalCycles_Projected(
normalCycles_SortAndSetEigenComponents<PFP>(ev,evec,kmax[v],kmin[v],Kmax[v],Kmin[v],Knormal[v],normal[v]);
}
template <typename PFP>
void computeCurvatureVertices_NormalCycles(
typename PFP::MAP& map,
Algo::Surface::Selection::Collector<PFP>& neigh,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal)
{
TraversorV<typename PFP::MAP> t(map) ;
for(Vertex v = t.begin(); v != t.end(); v = t.next())
computeCurvatureVertex_NormalCycles<PFP>(map, v, neigh, position, normal, edgeangle, edgearea, kmax, kmin, Kmax, Kmin, Knormal) ;
}
//template <typename PFP>
//void computeCurvatureVertices_NormalCycles(
// typename PFP::MAP& map,
// Algo::Surface::Selection::Collector<PFP>& neigh,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal)
//{
// TraversorV<typename PFP::MAP> t(map) ;
// for(Vertex v = t.begin(); v != t.end(); v = t.next())
// computeCurvatureVertex_NormalCycles<PFP>(map, v, neigh, position, normal, edgeangle, edgearea, kmax, kmin, Kmax, Kmin, Knormal) ;
//}
//template <typename PFP>
//void computeCurvatureVertex_NormalCycles(
// Vertex v,
// Algo::Surface::Selection::Collector<PFP>& neigh,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal)
//{
// typedef typename PFP::REAL REAL ;
// typedef typename PFP::VEC3 VEC3 ;
// typedef Geom::Matrix<3,3,REAL> MATRIX;
// typedef Eigen::Matrix<REAL,3,1> E_VEC3;
// typedef Eigen::Matrix<REAL,3,3,Eigen::RowMajor> E_MATRIX;
// // collect the normal cycle tensor
// neigh.collectAll(v) ;
// MATRIX tensor(0) ;
// neigh.computeNormalCyclesTensor(position, edgeangle, edgearea, tensor);
// // solve eigen problem
// Eigen::SelfAdjointEigenSolver<E_MATRIX> solver(Utils::convertRef<E_MATRIX>(tensor));
// const VEC3& ev = Utils::convertRef<VEC3>(solver.eigenvalues());
// const MATRIX& evec = Utils::convertRef<MATRIX>(solver.eigenvectors());
// normalCycles_SortAndSetEigenComponents<PFP>(ev,evec,kmax[v],kmin[v],Kmax[v],Kmin[v],Knormal[v],normal[v]);
//}
//template <typename PFP>
//void computeCurvatureVertices_NormalCycles_Projected(
// typename PFP::MAP& map,
// Algo::Surface::Selection::Collector<PFP>& neigh,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal)
//{
// TraversorV<typename PFP::MAP> t(map) ;
// for(Vertex v = t.begin(); v != t.end(); v = t.next())
// computeCurvatureVertex_NormalCycles_Projected<PFP>(map, v, neigh, position, normal, edgeangle, edgearea, kmax, kmin, Kmax, Kmin, Knormal) ;
//}
//template <typename PFP>
//void computeCurvatureVertex_NormalCycles_Projected(
// Vertex v,
// Algo::Surface::Selection::Collector<PFP>& neigh,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
// const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
// const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
// VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
// VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal)
//{
// typedef typename PFP::REAL REAL ;
// typedef typename PFP::VEC3 VEC3 ;
// typedef Geom::Matrix<3,3,REAL> MATRIX;
// typedef Eigen::Matrix<REAL,3,1> E_VEC3;
// typedef Eigen::Matrix<REAL,3,3,Eigen::RowMajor> E_MATRIX;
// // collect the normal cycle tensor
// neigh.collectAll(v) ;
// MATRIX tensor(0) ;
// neigh.computeNormalCyclesTensor(position, edgeangle, edgearea, tensor);
// // project the tensor
// normalCycles_ProjectTensor<PFP>(tensor, normal[v]);
// // solve eigen problem
// Eigen::SelfAdjointEigenSolver<E_MATRIX> solver(Utils::convertRef<E_MATRIX>(tensor));
// const VEC3& ev = Utils::convertRef<VEC3>(solver.eigenvalues());
// const MATRIX& evec = Utils::convertRef<MATRIX>(solver.eigenvectors());
// normalCycles_SortAndSetEigenComponents<PFP>(ev,evec,kmax[v],kmin[v],Kmax[v],Kmin[v],Knormal[v],normal[v]);
//}
template <typename PFP>
void computeCurvatureVertex_NormalCycles(
Vertex v,
Algo::Surface::Selection::Collector<PFP>& neigh,
void normalCycles_computeTensor(
Algo::Surface::Selection::Collector<PFP>& col,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& position,
const VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& normal,
const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgeangle,
const EdgeAttribute<typename PFP::REAL, typename PFP::MAP>& edgearea,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmax,
VertexAttribute<typename PFP::REAL, typename PFP::MAP>& kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmax,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Kmin,
VertexAttribute<typename PFP::VEC3, typename PFP::MAP>& Knormal)
Geom::Matrix<3,3,typename PFP::REAL>& tensor)
{
typedef typename PFP::REAL REAL ;
typedef typename PFP::VEC3 VEC3 ;
typedef Geom::Matrix<3,3,REAL> MATRIX;
typedef Eigen::Matrix<REAL,3,1> E_VEC3;
typedef Eigen::Matrix<REAL,3,3,Eigen::RowMajor> E_MATRIX;
// collect the normal cycle tensor
neigh.collectAll(v) ;
tensor.zero();
MATRIX tensor(0) ;
neigh.computeNormalCyclesTensor(position, edgeangle, edgearea, tensor);
// collect edges inside the neighborhood
for (Edge e : col.getInsideEdges())
{
typename PFP::REAL edgeangle = Algo::Surface::Geometry::computeAngleBetweenNormalsOnEdge<PFP>(col.getMap(), e, position);
typename PFP::VEC3 ev = Algo::Surface::Geometry::vectorOutOfDart<PFP>(col.getMap(), e, position);
tensor += Geom::transposed_vectors_mult(ev,ev) * edgeangle * (1 / ev.norm());
}
// solve eigen problem
Eigen::SelfAdjointEigenSolver<E_MATRIX> solver(Utils::convertRef<E_MATRIX>(tensor));
const VEC3& ev = Utils::convertRef<VEC3>(solver.eigenvalues());
const MATRIX& evec = Utils::convertRef<MATRIX>(solver.eigenvectors());
// collect edges on the border
for (Dart d : col.getBorder())
{
typename PFP::REAL edgeangle = Algo::Surface::Geometry::computeAngleBetweenNormalsOnEdge<PFP>(col.getMap(), d, position);
typename PFP::REAL alpha = col.borderEdgeRatio(d, position);
typename PFP::VEC3 ev = Algo::Surface::Geometry::vectorOutOfDart<PFP>(col.getMap(), d, position);
tensor += Geom::transposed_vectors_mult(ev,ev) * edgeangle * (1 / ev.norm()) * alpha;
}
normalCycles_SortAndSetEigenComponents<PFP>(ev,evec,kmax[v],kmin[v],Kmax[v],Kmin[v],Knormal[v],normal[v]);
tensor /= col.computeArea(position);
}
template <typename PFP>
void computeCurvatureVertices_NormalCycles_Projected(
typename PFP::MAP& map,
Algo::Surface::Selection::Collector<PFP>& neigh,