Skip to content

Commit

Permalink
Fix some code style issues
Browse files Browse the repository at this point in the history
  • Loading branch information
akuukka committed Jul 4, 2024
1 parent cdc561d commit 4ef66c6
Show file tree
Hide file tree
Showing 5 changed files with 85 additions and 75 deletions.
55 changes: 33 additions & 22 deletions QuickHull.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -44,17 +44,13 @@ namespace quickhull {
template<typename FloatType>
HalfEdgeMesh<FloatType, size_t> QuickHull<FloatType>::getConvexHullAsMesh(const FloatType* vertexData, size_t vertexCount, bool CCW, FloatType epsilon) {
VertexDataSource<FloatType> vertexDataSource((const vec3*)vertexData,vertexCount);
buildMesh(vertexDataSource, CCW, false, epsilon);
buildMesh(vertexDataSource, epsilon);
return HalfEdgeMesh<FloatType, size_t>(m_mesh, m_vertexData);
}

template<typename T>
void QuickHull<T>::buildMesh(const VertexDataSource<T>& pointCloud, bool CCW, bool useOriginalIndices, T epsilon) {
// CCW is unused for now
(void)CCW;
// useOriginalIndices is unused for now
(void)useOriginalIndices;

void QuickHull<T>::buildMesh(const VertexDataSource<T>& pointCloud, T epsilon)
{
if (pointCloud.size()==0) {
m_mesh = MeshBuilder<T>();
return;
Expand All @@ -71,11 +67,14 @@ namespace quickhull {

// Reset diagnostics
m_diagnostics = DiagnosticsData();

// The planar case happens when all the points appear to lie on a two dimensional
// subspace of R^3.
m_planar = false;

m_planar = false; // The planar case happens when all the points appear to lie on a two dimensional subspace of R^3.
createConvexHalfEdgeMesh();
if (m_planar) {
const size_t extraPointIndex = m_planarPointCloudTemp.size()-1;
const size_t extraPointIndex = m_planarPointCloudTemp.size() - 1;
for (auto& he : m_mesh.m_halfEdges) {
if (he.m_endVertex == extraPointIndex) {
he.m_endVertex = 0;
Expand All @@ -88,7 +87,7 @@ namespace quickhull {

template<typename T>
ConvexHull<T> QuickHull<T>::getConvexHull(const VertexDataSource<T>& pointCloud, bool CCW, bool useOriginalIndices, T epsilon) {
buildMesh(pointCloud,CCW,useOriginalIndices,epsilon);
buildMesh(pointCloud, epsilon);
return ConvexHull<T>(m_mesh,m_vertexData, CCW, useOriginalIndices);
}

Expand All @@ -100,7 +99,7 @@ namespace quickhull {

// Compute base tetrahedron
setupInitialTetrahedron();
assert(m_mesh.m_faces.size()==4);
assert(m_mesh.m_faces.size() == 4);

// Init face stack with those faces that have points assigned to them
m_faceList.clear();
Expand All @@ -117,8 +116,10 @@ namespace quickhull {
while (!m_faceList.empty()) {
iter++;
if (iter == std::numeric_limits<size_t>::max()) {
// Visible face traversal marks visited faces with iteration counter (to mark that the face has been visited on this iteration) and the max value represents unvisited faces. At this point we have to reset iteration counter. This shouldn't be an
// issue on 64 bit machines.
// Visible face traversal marks visited faces with iteration counter
// (to mark that the face has been visited on this iteration) and the max
// value represents unvisited faces. At this point we have to reset iteration
// counter. This shouldn't be an issue on 64 bit machines.
iter = 0;
}

Expand All @@ -137,7 +138,9 @@ namespace quickhull {
const vec3& activePoint = m_vertexData[tf.m_mostDistantPoint];
const size_t activePointIndex = tf.m_mostDistantPoint;

// Find out the faces that have our active point on their positive side (these are the "visible faces"). The face on top of the stack of course is one of them. At the same time, we create a list of horizon edges.
// Find out the faces that have our active point on their positive side
// (these are the "visible faces"). The face on top of the stack of course is
// one of them. At the same time, we create a list of horizon edges.
m_horizonEdges.clear();
m_possiblyVisibleFaces.clear();
m_visibleFaces.clear();
Expand Down Expand Up @@ -171,21 +174,28 @@ namespace quickhull {
assert(faceData.m_faceIndex != topFaceIndex);
}

// The face is not visible. Therefore, the halfedge we came from is part of the horizon edge.
// The face is not visible. Therefore, the halfedge we entered from
// is part of the horizon edge.
pvf.m_isVisibleFaceOnCurrentIteration = 0;
m_horizonEdges.push_back(faceData.m_enteredFromHalfEdge);
// Store which half edge is the horizon edge. The other half edges of the face will not be part of the final mesh so their data slots can by recycled.

// Store which half edge is the horizon edge. The other half edges of the face
// will not be part of the final mesh so their data slots can by recycled.
const auto halfEdges = m_mesh.getHalfEdgeIndicesOfFace(m_mesh.m_faces[m_mesh.m_halfEdges[faceData.m_enteredFromHalfEdge].m_face]);
const std::int8_t ind = (halfEdges[0]==faceData.m_enteredFromHalfEdge) ? 0 : (halfEdges[1]==faceData.m_enteredFromHalfEdge ? 1 : 2);
m_mesh.m_faces[m_mesh.m_halfEdges[faceData.m_enteredFromHalfEdge].m_face].m_horizonEdgesOnCurrentIteration |= (1<<ind);
}
const size_t horizonEdgeCount = m_horizonEdges.size();

// Order horizon edges so that they form a loop. This may fail due to numerical instability in which case we give up trying to solve horizon edge for this point and accept a minor degeneration in the convex hull.
// Order horizon edges so that they form a loop. This may fail due to numerical
// inaccuracy in which case we give up trying to solve horizon edge for this point
// and accept a minor degeneration in the convex hull.
if (!reorderHorizonEdges(m_horizonEdges)) {
m_diagnostics.m_failedHorizonEdges++;
std::cerr << "Failed to solve horizon edge." << std::endl;
auto it = std::find(tf.m_pointsOnPositiveSide->begin(),tf.m_pointsOnPositiveSide->end(),activePointIndex);
auto it = std::find(tf.m_pointsOnPositiveSide->begin(),
tf.m_pointsOnPositiveSide->end(),
activePointIndex);
tf.m_pointsOnPositiveSide->erase(it);
if (tf.m_pointsOnPositiveSide->size()==0) {
reclaimToIndexVectorPool(tf.m_pointsOnPositiveSide);
Expand Down Expand Up @@ -307,11 +317,12 @@ namespace quickhull {
*/

template <typename T>
std::array<size_t,6> QuickHull<T>::getExtremeValues() {
std::array<size_t,6> outIndices{0,0,0,0,0,0};
T extremeVals[6] = {m_vertexData[0].x,m_vertexData[0].x,m_vertexData[0].y,m_vertexData[0].y,m_vertexData[0].z,m_vertexData[0].z};
std::array<size_t, 6> QuickHull<T>::getExtremeValues() {
std::array<size_t,6> outIndices{0, 0, 0, 0, 0, 0};
T extremeVals[6] = { m_vertexData[0].x, m_vertexData[0].x, m_vertexData[0].y,
m_vertexData[0].y, m_vertexData[0].z, m_vertexData[0].z };
const size_t vCount = m_vertexData.size();
for (size_t i=1;i<vCount;i++) {
for (size_t i = 1; i < vCount; i++) {
const Vector3<T>& pos = m_vertexData[i];
if (pos.x>extremeVals[0]) {
extremeVals[0]=pos.x;
Expand Down
85 changes: 38 additions & 47 deletions QuickHull.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -43,15 +43,13 @@
* - Each point that was assigned to visible faces is now assigned to at most one of the newly created faces.
* - Those new faces that have points assigned to them are added to the top of Face Stack.
* - M is now the convex hull.
*
* TO DO:
* - Implement a proper 2D QuickHull and use that to solve the degenerate 2D case (when all the points lie on the same plane in 3D space).
* */

namespace quickhull {

struct DiagnosticsData {
size_t m_failedHorizonEdges; // How many times QuickHull failed to solve the horizon edge. Failures lead to degenerated convex hulls.
size_t m_failedHorizonEdges; // How many times QuickHull failed to solve the horizon edge.
// Failures lead to degenerate convex hulls.

DiagnosticsData() : m_failedHorizonEdges(0) { }
};
Expand Down Expand Up @@ -79,51 +77,64 @@ namespace quickhull {
std::vector<size_t> m_horizonEdges;
struct FaceData {
size_t m_faceIndex;
size_t m_enteredFromHalfEdge; // If the face turns out not to be visible, this half edge will be marked as horizon edge
size_t m_enteredFromHalfEdge; // If the face turns out not to be visible,
// this half edge will be marked as horizon edge
FaceData(size_t fi, size_t he) : m_faceIndex(fi),m_enteredFromHalfEdge(he) {}
};
std::vector<FaceData> m_possiblyVisibleFaces;
std::deque<size_t> m_faceList;

// Create a half edge mesh representing the base tetrahedron from which the QuickHull iteration proceeds. m_extremeValues must be properly set up when this is called.
// Create a half edge mesh representing the base tetrahedron from which the QuickHull
// iteration proceeds. m_extremeValues must be properly set up when this is called.
void setupInitialTetrahedron();

// Given a list of half edges, try to rearrange them so that they form a loop. Return true on success.
// Given a list of half edges, try to rearrange them so that they form a loop.
// Return true on success.
bool reorderHorizonEdges(std::vector<size_t>& horizonEdges);

// Find indices of extreme values (max x, min x, max y, min y, max z, min z) for the given point cloud
// Find indices of extreme values (max x, min x, max y, min y, max z, min z) for the
// given point cloud
std::array<size_t,6> getExtremeValues();

// Compute scale of the vertex data.
FloatType getScale(const std::array<size_t,6>& extremeValues);

// Each face contains a unique pointer to a vector of indices. However, many - often most - faces do not have any points on the positive
// side of them especially at the the end of the iteration. When a face is removed from the mesh, its associated point vector, if such
// exists, is moved to the index vector pool, and when we need to add new faces with points on the positive side to the mesh,
// we reuse these vectors. This reduces the amount of std::vectors we have to deal with, and impact on performance is remarkable.
// Each face contains a unique pointer to a vector of indices.
// However, many - often most - faces do not have any points on the positive
// side of them especially at the the end of the iteration. When a face is removed
// from the mesh, its associated point vector, if such exists, is moved to the index
// vector pool, and when we need to add new faces with points on the positive side to the
// mesh, we reuse these vectors. This reduces the amount of std::vectors we have to deal
// with, and impact on performance is remarkable.
Pool<std::vector<size_t>> m_indexVectorPool;
inline std::unique_ptr<std::vector<size_t>> getIndexVectorFromPool();
inline void reclaimToIndexVectorPool(std::unique_ptr<std::vector<size_t>>& ptr);

// Associates a point with a face if the point resides on the positive side of the plane. Returns true if the points was on the positive side.
// Associates a point with a face if the point resides on the positive side of the plane.
// Returns true if the points was on the positive side.
inline bool addPointToFace(typename MeshBuilder<FloatType>::Face& f, size_t pointIndex);

// This will update m_mesh from which we create the ConvexHull object that getConvexHull function returns
// This will update m_mesh from which we create the ConvexHull object that getConvexHull
// function returns
void createConvexHalfEdgeMesh();

// Constructs the convex hull into a MeshBuilder object which can be converted to a ConvexHull or Mesh object
void buildMesh(const VertexDataSource<FloatType>& pointCloud, bool CCW, bool useOriginalIndices, FloatType eps);
// Constructs the convex hull into a MeshBuilder object which can be converted to a
// ConvexHull or Mesh object
void buildMesh(const VertexDataSource<FloatType>& pointCloud, FloatType eps);

// The public getConvexHull functions will setup a VertexDataSource object and call this
ConvexHull<FloatType> getConvexHull(const VertexDataSource<FloatType>& pointCloud, bool CCW, bool useOriginalIndices, FloatType eps);
ConvexHull<FloatType> getConvexHull(const VertexDataSource<FloatType>& pointCloud,
bool CCW, bool useOriginalIndices, FloatType eps);
public:
// Computes convex hull for a given point cloud.
// Params:
// pointCloud: a vector of of 3D points
// CCW: whether the output mesh triangles should have CCW orientation
// useOriginalIndices: should the output mesh use same vertex indices as the original point cloud. If this is false,
// then we generate a new vertex buffer which contains only the vertices that are part of the convex hull.
// eps: minimum distance to a plane to consider a point being on positive of it (for a point cloud with scale 1)
// useOriginalIndices: should the output mesh use same vertex indices as the original point
// cloud. If this is false, then we generate a new vertex buffer which contains only
// the vertices that are part of the convex hull.
// eps: minimum distance to a plane to consider a point being on positive of it
// (for a point cloud with scale 1)
ConvexHull<FloatType> getConvexHull(const std::vector<Vector3<FloatType>>& pointCloud,
bool CCW,
bool useOriginalIndices,
Expand All @@ -133,40 +144,22 @@ namespace quickhull {
// Params:
// vertexData: pointer to the first 3D point of the point cloud
// vertexCount: number of vertices in the point cloud
// CCW: whether the output mesh triangles should have CCW orientation
// useOriginalIndices: should the output mesh use same vertex indices as the original point cloud. If this is false,
// then we generate a new vertex buffer which contains only the vertices that are part of the convex hull.
// eps: minimum distance to a plane to consider a point being on positive side of it (for a point cloud with scale 1)
ConvexHull<FloatType> getConvexHull(const Vector3<FloatType>* vertexData,
size_t vertexCount,
bool CCW,
bool useOriginalIndices,
FloatType eps = defaultEps<FloatType>());

// Computes convex hull for a given point cloud. This function assumes that the vertex data resides in memory
// in the following format: x_0,y_0,z_0,x_1,y_1,z_1,...
// Params:
// vertexData: pointer to the X component of the first point of the point cloud.
// vertexCount: number of vertices in the point cloud
// CCW: whether the output mesh triangles should have CCW orientation
// useOriginalIndices: should the output mesh use same vertex indices as the original point cloud. If this is false,
// then we generate a new vertex buffer which contains only the vertices that are part of the convex hull.
// eps: minimum distance to a plane to consider a point being on positive side of it (for a point cloud with scale 1)
// Computes convex hull for a given point cloud.
// This function assumes that the vertex data resides in memory in the following format:
// x_0,y_0,z_0,x_1,y_1,z_1,...
ConvexHull<FloatType> getConvexHull(const FloatType* vertexData,
size_t vertexCount,
bool CCW,
bool useOriginalIndices,
FloatType eps = defaultEps<FloatType>());

// Computes convex hull for a given point cloud. This function assumes that the vertex data resides in memory
// in the following format: x_0,y_0,z_0,x_1,y_1,z_1,...
// Params:
// vertexData: pointer to the X component of the first point of the point cloud.
// vertexCount: number of vertices in the point cloud
// CCW: whether the output mesh triangles should have CCW orientation
// eps: minimum distance to a plane to consider a point being on positive side of it (for a point cloud with scale 1)
// Returns:
// Convex hull of the point cloud as a mesh object with half edge structure.
// Convex hull of the point cloud as a mesh object with half edge structure.
HalfEdgeMesh<FloatType, size_t> getConvexHullAsMesh(const FloatType* vertexData,
size_t vertexCount,
bool CCW,
Expand All @@ -178,10 +171,6 @@ namespace quickhull {
}
};

/*
* Inline function definitions
*/

template<typename T>
std::unique_ptr<std::vector<size_t>> QuickHull<T>::getIndexVectorFromPool() {
auto r = m_indexVectorPool.get();
Expand All @@ -193,7 +182,9 @@ namespace quickhull {
void QuickHull<T>::reclaimToIndexVectorPool(std::unique_ptr<std::vector<size_t>>& ptr) {
const size_t oldSize = ptr->size();
if ((oldSize+1)*128 < ptr->capacity()) {
// Reduce memory usage! Huge vectors are needed at the beginning of iteration when faces have many points on their positive side. Later on, smaller vectors will suffice.
// Reduce memory usage! Huge vectors are needed at the beginning of iteration when
// faces have many points on their positive side. Later on, smaller vectors will
// suffice.
ptr.reset(nullptr);
return;
}
Expand Down
4 changes: 3 additions & 1 deletion Structs/Plane.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -25,7 +25,9 @@ namespace quickhull {
Plane() = default;

// Construct a plane using normal N and any point P on the plane
Plane(const Vector3<T>& N, const Vector3<T>& P) : m_N(N), m_D(-N.dotProduct(P)), m_sqrNLength(m_N.x*m_N.x+m_N.y*m_N.y+m_N.z*m_N.z) {
Plane(const Vector3<T>& N, const Vector3<T>& P) :
m_N(N), m_D(-N.dotProduct(P)), m_sqrNLength(m_N.x*m_N.x+m_N.y*m_N.y+m_N.z*m_N.z)
{

}
};
Expand Down
5 changes: 4 additions & 1 deletion Structs/Vector3.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -119,7 +119,10 @@ namespace quickhull {
const T dz = z-other.z;
return dx*dx+dy*dy+dz*dz;
}


bool operator==(const Vector3& other) const {
return x == other.x && y == other.y && z == other.z;
}
};

// Overload also << operator for easy printing of debug data
Expand Down
11 changes: 7 additions & 4 deletions Tests/QuickHullTests.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -141,10 +141,13 @@ namespace quickhull {
}

static void testPlanes() {
Vector3<FloatType> N(1,0,0);
Vector3<FloatType> p(2,0,0);
Plane<FloatType> P(N,p);
auto dist = mathutils::getSignedDistanceToPlane(Vector3<FloatType>(3,0,0), P);
vec3 N(1, 0, 0);
vec3 p(2, 0, 0);
Plane<FloatType> P(N, p);
VERIFY(P.isPointOnPositiveSide(p));
VERIFY(P.isPointOnPositiveSide(vec3(3, 0, 0)));
VERIFY(!P.isPointOnPositiveSide(vec3(1, 0, 0)));
FloatType dist = mathutils::getSignedDistanceToPlane(vec3(3,0,0), P);
CHECK_CLOSE(dist,1);
dist = mathutils::getSignedDistanceToPlane(Vector3<FloatType>(1,0,0), P);
CHECK_CLOSE(dist,-1);
Expand Down

0 comments on commit 4ef66c6

Please sign in to comment.