From 04e695760ebcf3b0376b625cbbef5c732dbf01e6 Mon Sep 17 00:00:00 2001
From: Karsten Rink <karsten.rink@ufz.de>
Date: Tue, 23 Sep 2014 15:06:33 +0200
Subject: [PATCH] added method for testing if point is inside element for all
 element types

---
 GeoLib/AnalyticalGeometry.h             |  3 +--
 MeshLib/Elements/Element.h              |  6 ++++++
 MeshLib/Elements/TemplateHex-impl.h     | 11 +++++++++++
 MeshLib/Elements/TemplateHex.h          |  3 +++
 MeshLib/Elements/TemplateLine-impl.h    |  9 +++++++++
 MeshLib/Elements/TemplateLine.h         |  3 +++
 MeshLib/Elements/TemplatePrism-impl.h   |  8 ++++++++
 MeshLib/Elements/TemplatePrism.h        |  3 +++
 MeshLib/Elements/TemplatePyramid-impl.h |  7 +++++++
 MeshLib/Elements/TemplatePyramid.h      |  3 +++
 MeshLib/Elements/TemplateQuad.h         |  3 +++
 MeshLib/Elements/TemplateTet-impl.h     |  6 ++++++
 MeshLib/Elements/TemplateTet.h          |  3 +++
 MeshLib/Elements/TemplateTri.h          |  3 +++
 14 files changed, 69 insertions(+), 2 deletions(-)

diff --git a/GeoLib/AnalyticalGeometry.h b/GeoLib/AnalyticalGeometry.h
index c6c403b8289..3c11356b5e9 100644
--- a/GeoLib/AnalyticalGeometry.h
+++ b/GeoLib/AnalyticalGeometry.h
@@ -212,8 +212,7 @@ double orient3d(GeoLib::Point const& p,
  /// Checks if the four given points are located on a plane.
  bool isCoplanar(const GeoLib::Point& a, const GeoLib::Point& b, 
 	 const GeoLib::Point& c, const GeoLib::Point& d);
-
-
+ 
 /**
  * Method first computes the intersection points of line segements of GeoLib::Polyline objects
  * (@see computeIntersectionPoints()) and pushes each intersection point in the GeoLib::PointVec
diff --git a/MeshLib/Elements/Element.h b/MeshLib/Elements/Element.h
index 7e10f07e404..5e402eff727 100644
--- a/MeshLib/Elements/Element.h
+++ b/MeshLib/Elements/Element.h
@@ -19,10 +19,13 @@
 #include <limits>
 #include <boost/optional.hpp>
 
+#include "Point.h"
+
 #include "MeshLib/MeshEnums.h"
 #include "MeshLib/Mesh.h"
 #include "MeshLib/MeshQuality/ElementErrorCode.h"
 
+
 namespace MeshLib {
 
 class Node;
@@ -152,6 +155,9 @@ public:
 	/// Returns true if these two indeces form an edge and false otherwise
 	virtual bool isEdge(unsigned i, unsigned j) const = 0;
 
+    /// Returns true if the given point is not located outside of the element
+    virtual bool isPntInElement(GeoLib::Point const& pnt) const = 0;
+
 	/**
 	 * Tests if the element is geometrically valid.
 	 */
diff --git a/MeshLib/Elements/TemplateHex-impl.h b/MeshLib/Elements/TemplateHex-impl.h
index 3469c18b60a..1b863420a75 100644
--- a/MeshLib/Elements/TemplateHex-impl.h
+++ b/MeshLib/Elements/TemplateHex-impl.h
@@ -139,6 +139,17 @@ bool TemplateHex<NNODES,CELLHEXTYPE>::isEdge(unsigned idx1, unsigned idx2) const
 	return false;
 }
 
+template <unsigned NNODES, CellType CELLHEXTYPE>
+bool TemplateHex<NNODES,CELLHEXTYPE>::isPntInElement(GeoLib::Point const& pnt) const
+{
+    return (GeoLib::isPointInTetrahedron(pnt, *_nodes[4], *_nodes[7], *_nodes[5], *_nodes[0]) ||
+            GeoLib::isPointInTetrahedron(pnt, *_nodes[5], *_nodes[3], *_nodes[1], *_nodes[0]) ||
+            GeoLib::isPointInTetrahedron(pnt, *_nodes[5], *_nodes[7], *_nodes[3], *_nodes[0]) ||
+            GeoLib::isPointInTetrahedron(pnt, *_nodes[5], *_nodes[7], *_nodes[6], *_nodes[2]) ||
+            GeoLib::isPointInTetrahedron(pnt, *_nodes[1], *_nodes[3], *_nodes[5], *_nodes[2]) ||
+            GeoLib::isPointInTetrahedron(pnt, *_nodes[3], *_nodes[7], *_nodes[5], *_nodes[2]));
+}
+
 template <unsigned NNODES, CellType CELLHEXTYPE>
 Element* TemplateHex<NNODES,CELLHEXTYPE>::clone() const
 {
diff --git a/MeshLib/Elements/TemplateHex.h b/MeshLib/Elements/TemplateHex.h
index 0fcb9bf951c..33b7f6805d3 100644
--- a/MeshLib/Elements/TemplateHex.h
+++ b/MeshLib/Elements/TemplateHex.h
@@ -104,6 +104,9 @@ public:
 	/// Returns true if these two indices form an edge and false otherwise
 	bool isEdge(unsigned i, unsigned j) const;
 
+    /// Returns true if the given point is not located outside of the hexahedron
+    bool isPntInElement(GeoLib::Point const& pnt) const;
+
 	/**
 	 * Tests if the element is geometrically valid.
 	 * @param check_zero_volume indicates if volume == 0 should be checked
diff --git a/MeshLib/Elements/TemplateLine-impl.h b/MeshLib/Elements/TemplateLine-impl.h
index 875cb8ecf29..b954d3e238d 100644
--- a/MeshLib/Elements/TemplateLine-impl.h
+++ b/MeshLib/Elements/TemplateLine-impl.h
@@ -55,6 +55,15 @@ template <unsigned NNODES, CellType CELLLINETYPE>
 TemplateLine<NNODES,CELLLINETYPE>::~TemplateLine()
 {}
 
+template <unsigned NNODES, CellType CELLLINETYPE>
+bool TemplateLine<NNODES,CELLLINETYPE>::isPntInElement(GeoLib::Point const& pnt) const
+{
+    double tmp;
+    double dist(0);
+    MathLib::calcProjPntToLineAndDists(pnt.getCoords(), _nodes[0]->getCoords(), _nodes[1]->getCoords(), tmp, dist);
+    return (dist < std::numeric_limits<double>::epsilon());
+}
+
 template <unsigned NNODES, CellType CELLLINETYPE>
 ElementErrorCode TemplateLine<NNODES,CELLLINETYPE>::validate() const
 { 
diff --git a/MeshLib/Elements/TemplateLine.h b/MeshLib/Elements/TemplateLine.h
index b3c5bb29e24..858ef66944d 100644
--- a/MeshLib/Elements/TemplateLine.h
+++ b/MeshLib/Elements/TemplateLine.h
@@ -115,6 +115,9 @@ public:
 		return false;
 	}
 
+    /// Returns true if pnt is located on the line segment and false otherwise
+    bool isPntInElement(GeoLib::Point const& pnt) const;
+
 	/**
 		* Tests if the element is geometrically valid.
 		* @param check_zero_volume indicates if area == 0 should be checked
diff --git a/MeshLib/Elements/TemplatePrism-impl.h b/MeshLib/Elements/TemplatePrism-impl.h
index b6c00bf195b..52cf5bb251a 100644
--- a/MeshLib/Elements/TemplatePrism-impl.h
+++ b/MeshLib/Elements/TemplatePrism-impl.h
@@ -148,6 +148,14 @@ bool TemplatePrism<NNODES,CELLPRISMTYPE>::isEdge(unsigned idx1, unsigned idx2) c
 	return false;
 }
 
+template <unsigned NNODES, CellType CELLPRISMTYPE>
+bool TemplatePrism<NNODES,CELLPRISMTYPE>::isPntInElement(GeoLib::Point const& pnt) const
+{
+    return (GeoLib::isPointInTetrahedron(pnt, *_nodes[0], *_nodes[1], *_nodes[2], *_nodes[3]) ||
+            GeoLib::isPointInTetrahedron(pnt, *_nodes[1], *_nodes[4], *_nodes[2], *_nodes[3]) ||
+            GeoLib::isPointInTetrahedron(pnt, *_nodes[2], *_nodes[4], *_nodes[5], *_nodes[3]));
+}
+
 template <unsigned NNODES, CellType CELLPRISMTYPE>
 Element* TemplatePrism<NNODES,CELLPRISMTYPE>::clone() const
 {
diff --git a/MeshLib/Elements/TemplatePrism.h b/MeshLib/Elements/TemplatePrism.h
index d91f2954024..c390e99e9cf 100644
--- a/MeshLib/Elements/TemplatePrism.h
+++ b/MeshLib/Elements/TemplatePrism.h
@@ -102,6 +102,9 @@ public:
 	/// Returns true if these two indeces form an edge and false otherwise
 	bool isEdge(unsigned i, unsigned j) const;
 
+    /// Returns true if the given point is not located outside of the prism
+    bool isPntInElement(GeoLib::Point const& pnt) const;
+
 	/**
 	 * Tests if the element is geometrically valid.
 	 * @param check_zero_volume indicates if volume == 0 should be checked
diff --git a/MeshLib/Elements/TemplatePyramid-impl.h b/MeshLib/Elements/TemplatePyramid-impl.h
index d21b6b0031a..823562a2995 100644
--- a/MeshLib/Elements/TemplatePyramid-impl.h
+++ b/MeshLib/Elements/TemplatePyramid-impl.h
@@ -150,6 +150,13 @@ bool TemplatePyramid<NNODES,CELLPYRAMIDTYPE>::isEdge(unsigned idx1, unsigned idx
 	return false;
 }
 
+template <unsigned NNODES, CellType CELLPYRAMIDTYPE>
+bool TemplatePyramid<NNODES,CELLPYRAMIDTYPE>::isPntInElement(GeoLib::Point const& pnt) const
+{
+    return (GeoLib::isPointInTetrahedron(pnt, *_nodes[0], *_nodes[1], *_nodes[2], *_nodes[4]) ||
+            GeoLib::isPointInTetrahedron(pnt, *_nodes[0], *_nodes[2], *_nodes[3], *_nodes[4]));
+}
+
 template <unsigned NNODES, CellType CELLPYRAMIDTYPE>
 Element* TemplatePyramid<NNODES,CELLPYRAMIDTYPE>::clone() const
 {
diff --git a/MeshLib/Elements/TemplatePyramid.h b/MeshLib/Elements/TemplatePyramid.h
index da8b34f776e..c12f169cc08 100644
--- a/MeshLib/Elements/TemplatePyramid.h
+++ b/MeshLib/Elements/TemplatePyramid.h
@@ -100,6 +100,9 @@ public:
 	/// Returns true if these two indeces form an edge and false otherwise
 	bool isEdge(unsigned i, unsigned j) const;
 
+    /// Returns true if the given point is not located outside of the pyramid
+    bool isPntInElement(GeoLib::Point const& pnt) const;
+
 	/**
 	 * Tests if the element is geometrically valid.
 	 * @param check_zero_volume indicates if volume == 0 should be checked
diff --git a/MeshLib/Elements/TemplateQuad.h b/MeshLib/Elements/TemplateQuad.h
index 21c4fe4ed3e..f7fd8e0889c 100644
--- a/MeshLib/Elements/TemplateQuad.h
+++ b/MeshLib/Elements/TemplateQuad.h
@@ -89,6 +89,9 @@ public:
 	/// Returns true if these two indeces form an edge and false otherwise
 	bool isEdge(unsigned i, unsigned j) const;
 
+    /// Returns true if the given point is not located outside of the quad
+    bool isPntInElement(GeoLib::Point const& pnt) const { return isPntInside(pnt); }
+
 	/**
 	 * Check if the 3d GeoLib::Point is inside of the quad element.
 	 * @param pnt the 3d GeoLib::Point object
diff --git a/MeshLib/Elements/TemplateTet-impl.h b/MeshLib/Elements/TemplateTet-impl.h
index 601f5199bcf..c3c990a4409 100644
--- a/MeshLib/Elements/TemplateTet-impl.h
+++ b/MeshLib/Elements/TemplateTet-impl.h
@@ -128,6 +128,12 @@ bool TemplateTet<NNODES,CELLTETTYPE>::isEdge(unsigned idx1, unsigned idx2) const
 	return false;
 }
 
+template <unsigned NNODES, CellType CELLTETTYPE>
+bool TemplateTet<NNODES,CELLTETTYPE>::isPntInElement(GeoLib::Point const& pnt) const
+{
+    return GeoLib::isPointInTetrahedron(pnt, *_nodes[0], *_nodes[1], *_nodes[2], *_nodes[3]);
+}
+
 template <unsigned NNODES, CellType CELLTETTYPE>
 Element* TemplateTet<NNODES,CELLTETTYPE>::clone() const
 {
diff --git a/MeshLib/Elements/TemplateTet.h b/MeshLib/Elements/TemplateTet.h
index b98c1e79580..530cbf09714 100644
--- a/MeshLib/Elements/TemplateTet.h
+++ b/MeshLib/Elements/TemplateTet.h
@@ -99,6 +99,9 @@ public:
 	/// Returns true if these two indeces form an edge and false otherwise
 	bool isEdge(unsigned i, unsigned j) const;
 
+    /// Returns true if the given point is not located outside of the tetrahedron
+    bool isPntInElement(GeoLib::Point const& pnt) const;
+
 	/**
 	 * Tests if the element is geometrically valid.
 	 * @param check_zero_volume indicates if volume == 0 should be checked
diff --git a/MeshLib/Elements/TemplateTri.h b/MeshLib/Elements/TemplateTri.h
index 465a8daed05..325efe914b3 100644
--- a/MeshLib/Elements/TemplateTri.h
+++ b/MeshLib/Elements/TemplateTri.h
@@ -93,6 +93,9 @@ public:
 	/// Returns true if these two indices form an edge and false otherwise
 	bool isEdge(unsigned idx1, unsigned idx2) const;
 
+    /// Returns true if the given point is not located outside of the triangle
+    bool isPntInElement(GeoLib::Point const& pnt) const { return isPntInside(pnt); }
+
 	/**
 	 * Check if the 3d GeoLib::Point is inside of the element.
 	 * @param pnt the 3d GeoLib::Point object
-- 
GitLab