/******************************************************************************* * CGoGN: Combinatorial and Geometric modeling with Generic N-dimensional Maps * * version 0.1 * * Copyright (C) 2009-2012, 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.unistra.fr/ * * Contact information: cgogn@unistra.fr * * * *******************************************************************************/ #include "Geometry/basic.h" #include "Algo/Geometry/centroid.h" #include "Topology/generic/traversorCell.h" #include "Topology/generic/traversor2.h" namespace CGoGN { namespace Algo { namespace Surface { namespace Geometry { template typename PFP::REAL triangleArea(typename PFP::MAP& map, Dart d, const VertexAttribute& position) { typename PFP::VEC3 p1 = position[d] ; typename PFP::VEC3 p2 = position[map.phi1(d)] ; typename PFP::VEC3 p3 = position[map.phi_1(d)] ; return Geom::triangleArea(p1, p2, p3) ; } template typename PFP::REAL convexFaceArea(typename PFP::MAP& map, Dart d, const VertexAttribute& position) { typedef typename PFP::VEC3 VEC3 ; if(map.faceDegree(d) == 3) return triangleArea(map, d, position) ; else { float area = 0.0f ; VEC3 centroid = faceCentroid(map, d, position) ; Traversor2FE t(map, d) ; for(Dart it = t.begin(); it != t.end(); it = t.next()) { VEC3 p1 = position[it] ; VEC3 p2 = position[map.phi1(it)] ; area += Geom::triangleArea(p1, p2, centroid) ; } return area ; } } template typename PFP::REAL totalArea(typename PFP::MAP& map, const VertexAttribute& position, unsigned int thread) { typename PFP::REAL area(0) ; TraversorF t(map) ; for(Dart d = t.begin(); d != t.end(); d = t.next()) area += convexFaceArea(map, d, position) ; return area ; } template typename PFP::REAL vertexOneRingArea(typename PFP::MAP& map, Dart d, const VertexAttribute& position) { typename PFP::REAL area(0) ; Traversor2VF t(map, d) ; for(Dart it = t.begin(); it != t.end(); it = t.next()) area += convexFaceArea(map, it, position) ; return area ; } template typename PFP::REAL vertexBarycentricArea(typename PFP::MAP& map, Dart d, const VertexAttribute& position) { typename PFP::REAL area(0) ; Traversor2VF t(map, d) ; for(Dart it = t.begin(); it != t.end(); it = t.next()) area += convexFaceArea(map, it, position) / 3 ; return area ; } template typename PFP::REAL vertexVoronoiArea(typename PFP::MAP& map, Dart d, const VertexAttribute& position) { typename PFP::REAL area(0) ; Traversor2VF t(map, d) ; for(Dart it = t.begin(); it != t.end(); it = t.next()) { const typename PFP::VEC3& p1 = position[it] ; const typename PFP::VEC3& p2 = position[map.phi1(it)] ; const typename PFP::VEC3& p3 = position[map.phi_1(it)] ; if(!Geom::isTriangleObtuse(p1, p2, p3)) { typename PFP::REAL a = Geom::angle(p3 - p2, p1 - p2) ; typename PFP::REAL b = Geom::angle(p1 - p3, p2 - p3) ; area += ( (p2 - p1).norm2() / tan(b) + (p3 - p1).norm2() / tan(a) ) / 8 ; } else { typename PFP::REAL tArea = Geom::triangleArea(p1, p2, p3) ; if(Geom::angle(p2 - p1, p3 - p1) > M_PI / 2) area += tArea / 2 ; else area += tArea / 4 ; } } return area ; } template void computeAreaFaces(typename PFP::MAP& map, const VertexAttribute& position, FaceAttribute& face_area) { TraversorF t(map) ; for(Dart d = t.begin(); d != t.end(); d = t.next()) face_area[d] = convexFaceArea(map, d, position) ; } template void computeOneRingAreaVertices(typename PFP::MAP& map, const VertexAttribute& position, VertexAttribute& vertex_area) { TraversorV t(map) ; for(Dart d = t.begin(); d != t.end(); d = t.next()) vertex_area[d] = vertexOneRingArea(map, d, position) ; } template void computeBarycentricAreaVertices(typename PFP::MAP& map, const VertexAttribute& position, VertexAttribute& vertex_area) { TraversorV t(map) ; for(Dart d = t.begin(); d != t.end(); d = t.next()) vertex_area[d] = vertexBarycentricArea(map, d, position) ; } template void computeVoronoiAreaVertices(typename PFP::MAP& map, const VertexAttribute& position, VertexAttribute& vertex_area) { TraversorV t(map) ; for(Dart d = t.begin(); d != t.end(); d = t.next()) vertex_area[d] = vertexVoronoiArea(map, d, position) ; } namespace Parallel { template class FunctorConvexFaceArea: public FunctorMapThreaded { const VertexAttribute& m_position; FaceAttribute& m_area; public: FunctorConvexFaceArea( typename PFP::MAP& map, const VertexAttribute& position, FaceAttribute& area): FunctorMapThreaded(map), m_position(position), m_area(area) { } void run(Dart d, unsigned int threadID) { m_area[d] = convexFaceArea(this->m_map, d, m_position) ; } }; template void computeAreaFaces(typename PFP::MAP& map, const VertexAttribute& position, FaceAttribute& area, unsigned int nbth) { FunctorConvexFaceArea funct(map,position,area); Algo::Parallel::foreach_cell(map, funct, nbth, false); } template class FunctorVertexOneRingArea: public FunctorMapThreaded { const VertexAttribute& m_position; VertexAttribute& m_area; public: FunctorVertexOneRingArea( typename PFP::MAP& map, const VertexAttribute& position, VertexAttribute& area): FunctorMapThreaded(map), m_position(position), m_area(area) { } void run(Dart d, unsigned int threadID) { m_area[d] = vertexOneRingArea(this->m_map, d, m_position) ; } }; template void computeOneRingAreaVertices(typename PFP::MAP& map, const VertexAttribute& position, VertexAttribute& area, unsigned int nbth) { FunctorConvexFaceArea funct(map,position,area); Algo::Parallel::foreach_cell(map, funct, nbth, false); } template class FunctorVertexVoronoiArea: public FunctorMapThreaded { const VertexAttribute& m_position; VertexAttribute& m_area; public: FunctorVertexVoronoiArea( typename PFP::MAP& map, const VertexAttribute& position, VertexAttribute& area): FunctorMapThreaded(map), m_position(position), m_area(area) { } void run(Dart d, unsigned int threadID) { m_area[d] = vertexVoronoiArea(this->m_map, d, m_position) ; } }; template void computeVoronoiAreaVertices(typename PFP::MAP& map, const VertexAttribute& position, VertexAttribute& area, unsigned int nbth) { FunctorConvexFaceArea funct(map,position,area); Algo::Parallel::foreach_cell(map, funct, nbth, false); } } } // namespace Geometry } } // namespace Algo } // namespace CGoGN