/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2018 OpenFOAM Foundation
Copyright (C) 2019 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM 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 General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::triSurfaceTools
Description
A collection of tools for triSurface.
Note
The curvature calculation is an implementation of the algorithm from:
\verbatim
"Estimating Curvatures and their Derivatives on Triangle Meshes"
by S. Rusinkiewicz
3DPVT'04 Proceedings of the 3D Data Processing,
Visualization, and Transmission, 2nd International Symposium
Pages 486-493
http://gfx.cs.princeton.edu/pubs/_2004_ECA/curvpaper.pdf
\endverbatim
SourceFiles
triSurfaceTools.C
triSurfaceCloseness.C
triSurfaceCurvature.C
\*---------------------------------------------------------------------------*/
#ifndef triSurfaceTools_H
#define triSurfaceTools_H
#include "boolList.H"
#include "pointField.H"
#include "vectorField.H"
#include "triadFieldFwd.H"
#include "DynamicList.H"
#include "HashSet.H"
#include "Map.H"
#include "FixedList.H"
#include "Pair.H"
#include "vector2D.H"
#include "triPointRef.H"
#include "surfaceLocation.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// Forward declaration of classes
class boundBox;
class edge;
class labelledTri;
class polyBoundaryMesh;
class plane;
class triSurface;
class face;
class Time;
template<class Face> class MeshedSurface;
/*---------------------------------------------------------------------------*\
Class triSurfaceTools Declaration
\*---------------------------------------------------------------------------*/
class triSurfaceTools
{
// Private Member Functions
// Refinement
enum refineType
{
NONE,
RED,
GREEN
};
static void calcRefineStatus
(
const triSurface& surf,
const label facei,
List<refineType>& refine
);
static void greenRefine
(
const triSurface& surf,
const label facei,
const label edgeI,
const label newPointi,
DynamicList<labelledTri>& newFaces
);
static triSurface doRefine
(
const triSurface& surf,
const List<refineType>& refineStatus
);
// Coarsening
static scalar faceCosAngle
(
const point& pStart,
const point& pEnd,
const point& pLeft,
const point& pRight
);
static void protectNeighbours
(
const triSurface& surf,
const label vertI,
labelList& faceStatus
);
//- Faces to collapse because of edge collapse
static labelHashSet getCollapsedFaces
(
const triSurface& surf,
label edgeI
);
// Return value of faceUsed for faces using vertI (local numbering).
// Used internally.
static label vertexUsesFace
(
const triSurface& surf,
const labelHashSet& faceUsed,
const label vertI
);
// Get new connections between faces (because of edge collapse)
// in form of tables:
// - given edge get other edge
// - given edge get other face
// A face using point v1 on edge will get connected to a face using
// point v2 if they share a common vertex
// (but not a common edge since then the triangles collapse to
// nothing)
static void getMergedEdges
(
const triSurface& surf,
const label edgeI,
const labelHashSet& collapsedFaces,
Map<label>& edgeToEdge,
Map<label>& edgeToFace
);
//- Calculates (cos of) angle across edgeI of facei,
// taking into account updated addressing (resulting from edge
// collapse)
static scalar edgeCosAngle
(
const triSurface& surf,
const label v1,
const point& pt,
const labelHashSet& collapsedFaces,
const Map<label>& edgeToEdge,
const Map<label>& edgeToFace,
const label facei,
const label edgeI
);
//- Calculate minimum (cos of) edge angle using addressing from
// collapsing
// edge to v1 at pt. Returns 1 if v1 is on edge without neighbours
// (and hence no edge angle can be defined)
static scalar collapseMinCosAngle
(
const triSurface& surf,
const label v1,
const point& pt,
const labelHashSet& collapsedFaces,
const Map<label>& edgeToEdge,
const Map<label>& edgeToFace
);
//- Like collapseMinCosAngle but return true for value < minCos
bool collapseCreatesFold
(
const triSurface& surf,
const label v1,
const point& pt,
const labelHashSet& collapsedFaces,
const Map<label>& edgeToEdge,
const Map<label>& edgeToFace,
const scalar minCos
);
////- Checks if edge collapse creates triangles on top of each
//// other
//static bool collapseCreatesDuplicates
//(
// const triSurface& surf,
// const label edgeI,
// const Map<bool>& collapsedFaces
//);
// Tracking
//- Finds the triangle edge/point cut by the plane between
// a point inside/on edge of a triangle and a point outside.
// Returns
// - location on edge/point and hit()
// - or miss() if no intersection found
static surfaceLocation cutEdge
(
const triSurface& s,
const label triI,
const label excludeEdgeI,
const label excludePointi,
const point& triPoint,
const plane& cutPlane,
const point& toPoint
);
//- Checks if current is on the same triangle as the endpoint
// and shifts it there. If so updates current and sets a hit.
static void snapToEnd
(
const triSurface& s,
const surfaceLocation& endInfo,
surfaceLocation& current
);
//- Visits faces eFaces around start. Does not visit triangle
// start.triangle() nor edge excludeEdgeI.
// Returns edge, triangle (if more than one choice) which gets
// us nearer endpoint.
// Returns
// - hit() if triangle contains endpoint
// - triangle()=-1 if no triangle found
// - nearest triangle/edge otherwise
static surfaceLocation visitFaces
(
const triSurface& s,
const labelList& eFaces,
const surfaceLocation& start,
const label excludeEdgeI,
const label excludePointi,
const surfaceLocation& end,
const plane& cutPlane
);
public:
// OBJ writing
//- Write pointField to OBJ format file
static void writeOBJ
(
const fileName& fName,
const pointField& pts
);
//- Write vertex subset to OBJ format file
static void writeOBJ
(
const triSurface& surf,
const fileName& fName,
const boolList& markedVerts
);
// Additional addressing
//- Get all triangles using edge endpoint
static void getVertexTriangles
(
const triSurface& surf,
const label edgeI,
labelList& edgeTris
);
//- Get all vertices (local numbering) connected to vertices of edge
static labelList getVertexVertices
(
const triSurface& surf,
const edge& e
);
//- Get face connected to edge not facei
static label otherFace
(
const triSurface& surf,
const label facei,
const label edgeI
);
//- Get the two edges on facei counterclockwise after edgeI
static void otherEdges
(
const triSurface& surf,
const label facei,
const label edgeI,
label& e1,
label& e2
);
//- Get the two vertices (local numbering) on facei counterclockwise
// vertI
static void otherVertices
(
const triSurface& surf,
const label facei,
const label vertI,
label& vert1I,
label& vert2I
);
//- Get edge opposite vertex (local numbering)
static label oppositeEdge
(
const triSurface& surf,
const label facei,
const label vertI
);
//- Get vertex (local numbering) opposite edge
static label oppositeVertex
(
const triSurface& surf,
const label facei,
const label edgeI
);
//- Returns edge label connecting v1, v2 (local numbering)
static label getEdge
(
const triSurface& surf,
const label vert1I,
const label vert2I
);
//- Return index of triangle (or -1) using all three edges
static label getTriangle
(
const triSurface& surf,
const label e0I,
const label e1I,
const label e2I
);
// Coarsening
//- Create new triSurface by collapsing edges to edge mids.
static triSurface collapseEdges
(
const triSurface& surf,
const labelList& collapsableEdges
);
//- Face collapse status.
// anyEdge: any edge can be collapsed
// noEdge: no edge can be collapsed
// collapsed: already collapsed
// >0: edge label that can be collapsed
static const label ANYEDGE;
static const label NOEDGE;
static const label COLLAPSED;
//- Create new triSurface by collapsing edges to specified
// positions. faceStatus allows
// explicit control over which faces need to be protected (see above).
// faceStatus gets updated to protect collapsing already collapsed
// faces.
static triSurface collapseEdges
(
const triSurface& surf,
const labelList& collapsableEdges,
const pointField& edgeMids,
labelList& faceStatus
);
// Refinement
//- Refine edges by splitting to opposite vertex
static triSurface greenRefine
(
const triSurface& surf,
const labelList& refineEdges
);
//- Refine face by splitting all edges. Neighbouring face is
// greenRefine'd.
static triSurface redGreenRefine
(
const triSurface& surf,
const labelList& refineFaces
);
// Geometric
//- Returns element in edgeIndices with minimum length
static label minEdge
(
const triSurface& surf,
const labelList& edgeIndices
);
//- Returns element in edgeIndices with minimum length
static label maxEdge
(
const triSurface& surf,
const labelList& edgeIndices
);
//- Merge points within distance
static triSurface mergePoints
(
const triSurface& surf,
const scalar mergeTol
);
//- Triangle (unit) normal. If nearest point to triangle on edge use
// edge normal (calculated on the fly); if on vertex use vertex normal.
// Uses planarTol.
static vector surfaceNormal
(
const triSurface& surf,
const label nearestFacei,
const point& nearestPt
);
//- On which side of surface
enum sideType
{
UNKNOWN, // cannot be determined (e.g. non-manifold)
INSIDE, // inside
OUTSIDE // outside
};
//- If nearest point is on edgeI, determine on which side of surface
// sample is.
static sideType edgeSide
(
const triSurface& surf,
const point& sample,
const point& nearestPoint,
const label edgeI
);
//- Given nearest point (to sample) on surface determines which side
// sample is. Uses either face normal, edge normal or point normal
// (non-trivial). Uses triangle::classify.
static sideType surfaceSide
(
const triSurface& surf,
const point& sample,
const label nearestFacei
);
// Triangulation of faces
//- Simple triangulation of (selected patches of) boundaryMesh. Needs
// polyMesh (or polyBoundaryMesh) since only at this level are the
// triangles on neighbouring patches connected.
// Return faceMap from triI to faceI
static triSurface triangulate
(
const polyBoundaryMesh& mBesh,
const labelHashSet& includePatches,
labelList& faceMap,
const bool verbose = false
);
static triSurface triangulate
(
const polyBoundaryMesh& bMesh,
const labelHashSet& includePatches,
const boundBox& bBox,
const bool verbose = false
);
//- Face-centre triangulation of (selected patches of) boundaryMesh.
// Needs
// polyMesh (or polyBoundaryMesh) since only at this level are the
// triangles on neighbouring patches connected.
triSurface triangulateFaceCentre
(
const polyBoundaryMesh& mBesh,
const labelHashSet& includePatches,
const bool verbose = false
);
// Triangulation and interpolation
//- Do unconstrained Delaunay of points. Returns triSurface with 3D
// points with z=0. All triangles in region 0.
static triSurface delaunay2D(const List<vector2D>&);
//- Calculate linear interpolation weights for point (guaranteed to be
// inside triangle)
static void calcInterpolationWeights
(
const triPointRef& tri,
const point& p,
FixedList<scalar, 3>& weights
);
// Calculate weighting factors from samplePts to triangle it is in.
// Uses linear search to find triangle.
// Vertices are:
// (a b c) : vertices of the triangle abc the point is in
// or if the point is outside all triangles:
// (a b -1) : the edge ab the point is nearest to.
// (a -1 -1) : the vertex a the point is nearest to
static void calcInterpolationWeights
(
const triSurface& s,
const pointField& samplePts,
List<FixedList<label, 3>>& allVerts,
List<FixedList<scalar, 3>>& allWeights
);
// Curvature
//- Weighting for normals of faces attached to vertex
static scalar vertexNormalWeight
(
const triFace& f,
const label pI,
const vector& fN,
const UList<point>& points
);
//- Weighted average of normals of attached faces
static tmp<vectorField> vertexNormals(const triSurface& surf);
//- Local coordinate-system for each point normal
static tmp<triadField> vertexTriads
(
const triSurface& surf,
const vectorField& pointNormals
);
//- Surface curvatures at the vertex points
static tmp<scalarField> curvatures
(
const triSurface& surf,
const vectorField& pointNormals,
const triadField& pointTriads
);
//- Surface curvatures at the vertex points
static tmp<scalarField> curvatures
(
const triSurface& surf
);
//- Write surface curvature at the vertex points and return the field
static tmp<scalarField> writeCurvature
(
const Time& runTime,
const word& basename,
const triSurface& surf
);
// Closeness
//- Check and write internal/external closeness fields
static Pair<tmp<scalarField>> writeCloseness
(
const Time& runTime,
const word& basename,
const triSurface& surf,
const scalar internalAngleTolerance = 45,
const scalar externalAngleTolerance = 10
);
// Feature Proximity
//- Calculate feature proximity
scalarField featureProximity
(
const triSurface& surf,
const scalar searchDistance
);
//- Check and write internal/external closeness fields
static void writeFeatureProximity
(
const Time& runTime,
const word& basename,
const triSurface& surf,
const bool writeVTK,
const scalar searchDistance
);
// Surface checking functionality
//- Check single triangle for (topological) validity
static bool validTri
(
const triSurface&,
const label facei,
const bool verbose = true
);
//- Check single triangle for (topological) validity
static bool validTri
(
const MeshedSurface<face>&,
const label facei,
const bool verbose = true
);
// Tracking
//- Test point on plane of triangle to see if on edge or point or inside
static surfaceLocation classify
(
const triSurface&,
const label triI,
const point& trianglePoint
);
//- Track on surface to get closer to point.
// Possible situations:
// - 1. reached endpoint
// - 2. reached edge (normal situation)
// - 3. reached end of surface (edge on single face)
// Input:
// - starting position+triangle/edge/point (so has to be on surface!)
// - (optional) previous position+triangle/edge/point to prevent
// going back. Set index (of triangle/edge/point) to -1 if not
// used.
// - end position+triangle/edge/point (so has to be on surface!)
// - plane to follow. Has to go through end point!
// Returns:
// - true if end point reached (situation 1)
// - new position+triangle/edge/point
// Caller has to check for situation 3 by checking that triangle()
// is not set.
static surfaceLocation trackToEdge
(
const triSurface&,
const surfaceLocation& start,
const surfaceLocation& end,
const plane& cutPlane
);
//- Track from edge to edge across surface. Uses trackToEdge.
// Not really useful by itself, more example of how to use trackToEdge.
// endInfo should be location on surface.
// hitInfo should be initialised to starting location (on surface as
// well). Upon return is set to end location.
static void track
(
const triSurface&,
const surfaceLocation& endInfo,
const plane& cutPlane,
surfaceLocation& hitInfo
);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //