Skip to content

Commit

Permalink
Fix multiple bugs in SliceTriangleGeometry
Browse files Browse the repository at this point in the history 
Signed-off-by: Michael Jackson <mike.jackson@bluequartz.net>
  • Loading branch information
@imikejackson
imikejackson committed Jan 3, 2025
1 parent 5d109c5 commit 850e97e
Show file tree
Hide file tree
Showing 4 changed files with 138 additions and 80 deletions.
79 changes: 41 additions & 38 deletions src/Plugins/SimplnxCore/src/SimplnxCore/Filters/Algorithms/SliceTriangleGeometry.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,39 +4,9 @@
#include "simplnx/DataStructure/Geometry/EdgeGeom.hpp"
#include "simplnx/DataStructure/Geometry/TriangleGeom.hpp"
#include "simplnx/Utilities/GeometryUtilities.hpp"
#include "simplnx/Utilities/IntersectionUtilities.hpp"

using namespace nx::core;

namespace
{
// -----------------------------------------------------------------------------
char RayIntersectsPlane(const float32 d, const nx::core::Point3Df& q, const nx::core::Point3Df& r, nx::core::Point3Df& p)
{
const float64 rqDelZ = r[2] - q[2];
const float64 dqDelZ = d - q[2];
const float64 t = dqDelZ / rqDelZ;
for(int i = 0; i < 3; i++)
{
p[i] = q[i] + (t * (r[i] - q[i]));
}
if(t > 0.0 && t < 1.0)
{
return '1';
}
if(t == 0.0)
{
return 'q';
}
if(t == 1.0)
{
return 'r';
}

return '0';
}
} // namespace

// -----------------------------------------------------------------------------
SliceTriangleGeometry::SliceTriangleGeometry(DataStructure& dataStructure, const IFilter::MessageHandler& mesgHandler, const std::atomic_bool& shouldCancel,
SliceTriangleGeometryInputValues* inputValues)
Expand Down Expand Up @@ -94,13 +64,13 @@ Result<> SliceTriangleGeometry::operator()()
return MakeErrorResult(-62102, fmt::format("Number of sectioned vertices and edges do not make sense. Number of Vertices: {} and Number of Edges: {}", numVerts, numEdges));
}

auto& edge = m_DataStructure.getDataRefAs<EdgeGeom>(m_InputValues->SliceDataContainerName);
edge.resizeVertexList(numVerts);
edge.resizeEdgeList(numEdges);
INodeGeometry0D::SharedVertexList& verts = edge.getVerticesRef();
INodeGeometry1D::SharedEdgeList& edges = edge.getEdgesRef();
edge.getVertexAttributeMatrix()->resizeTuples({numVerts});
edge.getEdgeAttributeMatrix()->resizeTuples({numEdges});
auto& edgeGeom = m_DataStructure.getDataRefAs<EdgeGeom>(m_InputValues->SliceDataContainerName);
edgeGeom.resizeVertexList(numVerts);
edgeGeom.resizeEdgeList(numEdges);
INodeGeometry0D::SharedVertexList& verts = edgeGeom.getVerticesRef();
INodeGeometry1D::SharedEdgeList& edges = edgeGeom.getEdgesRef();
edgeGeom.getVertexAttributeMatrix()->resizeTuples({numVerts});
edgeGeom.getEdgeAttributeMatrix()->resizeTuples({numEdges});
auto& sliceAM = m_DataStructure.getDataRefAs<AttributeMatrix>(m_InputValues->SliceDataContainerName.createChildPath(m_InputValues->SliceAttributeMatrixName));
sliceAM.resizeTuples({sliceTriangleResult.NumberOfSlices});

Expand Down Expand Up @@ -131,5 +101,38 @@ Result<> SliceTriangleGeometry::operator()()
}
}

return GeometryUtilities::EliminateDuplicateNodes<EdgeGeom>(edge);
Result<> result = GeometryUtilities::EliminateDuplicateNodes<EdgeGeom>(edgeGeom);
if(result.invalid())
{
return result;
}

// REMOVE DUPLICATE EDGES FROM THE GENERATED EDGE GEOMETRY
// Remember to also fix up the sliceIds and regionIds arrays
using UniqueEdges = std::set<std::pair<uint64, uint64>>;
UniqueEdges uniqueEdges;
usize currentEdgeListSize = 0;
for(usize edgeIdx = 0; edgeIdx < numEdges; edgeIdx++)
{
uniqueEdges.insert(std::minmax(edges[edgeIdx * 2], edges[edgeIdx * 2 + 1]));
// Did something get inserted
if(currentEdgeListSize < uniqueEdges.size())
{
edges[currentEdgeListSize * 2] = edges[edgeIdx * 2];
edges[currentEdgeListSize * 2 + 1] = edges[edgeIdx * 2 + 1];
sliceId[currentEdgeListSize] = sliceId[edgeIdx];
if(m_InputValues->HaveRegionIds)
{
(*triRegionIds)[currentEdgeListSize] = (*triRegionIds)[edgeIdx];
}
currentEdgeListSize++;
}
}
if(numEdges != uniqueEdges.size())
{
edgeGeom.resizeEdgeList(uniqueEdges.size());
edgeGeom.getEdgeAttributeMatrix()->resizeTuples({uniqueEdges.size()});
}

return result;
}
2 changes: 0 additions & 2 deletions src/Plugins/SimplnxCore/src/SimplnxCore/Filters/Algorithms/SliceTriangleGeometry.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,6 @@ constexpr ChoicesParameter::ValueType k_UserDefinedRange = 1;

struct SIMPLNXCORE_EXPORT SliceTriangleGeometryInputValues
{
// VectorFloat32Parameter::ValueType SliceDirection;
ChoicesParameter::ValueType SliceRange;
float32 Zstart;
float32 Zend;
Expand Down Expand Up @@ -59,7 +58,6 @@ class SIMPLNXCORE_EXPORT SliceTriangleGeometry
protected:
using TriStore = AbstractDataStore<INodeGeometry2D::SharedFaceList::value_type>;
using VertsStore = AbstractDataStore<INodeGeometry0D::SharedVertexList::value_type>;
// usize determineBoundsAndNumSlices(float32& minDim, float32& maxDim, usize numTris, TriStore& tris, VertsStore& triVerts);

private:
DataStructure& m_DataStructure;
Expand Down
135 changes: 96 additions & 39 deletions src/simplnx/Utilities/GeometryUtilities.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,8 @@
#include "simplnx/Common/Result.hpp"
#include "simplnx/Utilities/Math/MatrixMath.hpp"

#include <cstdio>

using namespace nx::core;

namespace
Expand All @@ -12,6 +14,7 @@ constexpr uint64 k_FullRange = 0;
constexpr uint64 k_UserDefinedRange = 1;
constexpr float32 k_PartitionEdgePadding = 0.000001;
const Point3Df k_Padding(k_PartitionEdgePadding, k_PartitionEdgePadding, k_PartitionEdgePadding);

} // namespace

GeometryUtilities::FindUniqueIdsImpl::FindUniqueIdsImpl(VertexStore& vertexStore, const std::vector<std::vector<size_t>>& nodesInBin, nx::core::Int64DataStore& uniqueIds)
Expand Down Expand Up @@ -264,6 +267,9 @@ struct Edge
Point3Df end;
bool valid;
int32 regionId;
uint8 positiveCount = 0;
uint8 negativeCount = 0;
uint8 zeroCount = 0;

// Constructors
Edge()
Expand Down Expand Up @@ -313,13 +319,15 @@ struct Plane
}
};

/**
* @brief
*/
struct PointInfo
{

float SignedDistance;
uint8 location;

PointInfo(float signedDistance)
explicit PointInfo(float signedDistance)
: SignedDistance(signedDistance)
, location(3) // Default the point is on the plane
{
Expand Down Expand Up @@ -355,95 +363,141 @@ struct PointInfo
// Function to compute the intersection between a triangle and a plane
Edge IntersectTriangleWithPlane(const Point3Df& v0, const Point3Df& v1, const Point3Df& v2, const Plane& plane)
{
PointInfo p0 = {plane.signedDistance(v0)};
PointInfo p1 = {plane.signedDistance(v1)};
PointInfo p2 = {plane.signedDistance(v2)};
PointInfo p0{plane.signedDistance(v0)};
PointInfo p1{plane.signedDistance(v1)};
PointInfo p2{plane.signedDistance(v2)};

// Count the number of vertices on each side of the plane
int positiveCount = p0.positive() + p1.positive() + p2.positive();
int negativeCount = p0.negative() + p1.negative() + p2.negative();
int zeroCount = p0.onPlane() + p1.onPlane() + p2.onPlane();

// No intersection if all vertices are on one side of the plane
if(positiveCount == 3 || negativeCount == 3)
// Handle case where the triangle lies entirely on the plane
if(positiveCount == 3 || negativeCount == 3 || zeroCount == 3)
{
return {}; // Invalid edge because triangle is completely above or below the plane
Edge e;
e.positiveCount = positiveCount;
e.negativeCount = negativeCount;
e.zeroCount = zeroCount;
return std::move(e); // Invalid edge because triangle is completely above or below the plane
}

// Edge to store the intersection line segment
Edge intersectionEdge;
// int intersectionCount = 0;
intersectionEdge.positiveCount = positiveCount;
intersectionEdge.negativeCount = negativeCount;
intersectionEdge.zeroCount = zeroCount;

// Helper lambda to compute intersection point
auto computeIntersection = [](const Edge& e, float _dist1, float _dist2) -> Point3Df {
float t = _dist1 / (_dist1 - _dist2); // Compute interpolation parameter
return e.start + (e.end - e.start) * t; // Return the interpolated point
};

// Handle case where the triangle lies entirely on the plane
if(zeroCount == 3)
// Handle cases where only one intersection point is found (vertex lies on plane)
// and the other two vertices are either both above or below the plane
// Find the vertex that lies on the plane
if(zeroCount == 1 && (positiveCount == 2 || negativeCount == 2))
{
// Return any edge of the triangle
// return {v0, v1};
return {}; // Return invalid
Edge e;
e.positiveCount = positiveCount;
e.negativeCount = negativeCount;
e.zeroCount = zeroCount;
return std::move(e);
}

// Handle cases where only one intersection point is found (vertex lies on plane)
// Find the vertex that lies on the plane
if(zeroCount == 1)
if(positiveCount == 1 && negativeCount == 1 && zeroCount == 1)
{
return {};
Edge e;
if(p0.onPlane())
{
e.start = v0;
e.end = computeIntersection({v1, v2}, p1.SignedDistance, p2.SignedDistance);
}
else if(p1.onPlane())
{
e.start = v1;
e.end = computeIntersection({v0, v2}, p0.SignedDistance, p2.SignedDistance);
}
else if(p2.onPlane())
{
e.start = v2;
e.end = computeIntersection({v0, v1}, p0.SignedDistance, p1.SignedDistance);
}
e.positiveCount = positiveCount;
e.negativeCount = negativeCount;
e.zeroCount = zeroCount;
e.valid = true;
return std::move(e);
}
// if(p0.onPlane() && zeroCount == 1)
// {
// return {v0, v0};
// }
// else if(p1.onPlane() && zeroCount == 1)
// {
// return {v1, v1};
// }
// else if(p2.onPlane() && zeroCount == 1)
// {
// return {v2, v2};
// }

// Check edges for coincidence with plane
if(p0.onPlane() && p1.onPlane() && zeroCount == 2)
{
return {v0, v1};
Edge e{v0, v1};
e.positiveCount = positiveCount;
e.negativeCount = negativeCount;
e.zeroCount = zeroCount;
return std::move(e);
}
if(p1.onPlane() && p2.onPlane() && zeroCount == 2)
{
return {v1, v2};
Edge e{v1, v2};
e.positiveCount = positiveCount;
e.negativeCount = negativeCount;
e.zeroCount = zeroCount;
return std::move(e);
}
if(p0.onPlane() && p2.onPlane() && zeroCount == 2)
{
return {v0, v2};
Edge e{v0, v2};
e.positiveCount = positiveCount;
e.negativeCount = negativeCount;
e.zeroCount = zeroCount;
return std::move(e);
}

if(p0.planeSplitsEdge(p1) && p0.planeSplitsEdge(p2))
{
auto intersectionPoint0 = computeIntersection({v0, v1}, p0.SignedDistance, p1.SignedDistance);
auto intersectionPoint1 = computeIntersection({v0, v2}, p0.SignedDistance, p2.SignedDistance);
return {intersectionPoint0, intersectionPoint1};

Edge e{intersectionPoint0, intersectionPoint1};
e.positiveCount = positiveCount;
e.negativeCount = negativeCount;
e.zeroCount = zeroCount;
return std::move(e);
}

if(p0.planeSplitsEdge(p1) && p1.planeSplitsEdge(p2))
{
auto intersectionPoint0 = computeIntersection({v0, v1}, p0.SignedDistance, p1.SignedDistance);
auto intersectionPoint1 = computeIntersection({v1, v2}, p1.SignedDistance, p2.SignedDistance);
return {intersectionPoint0, intersectionPoint1};
Edge e{intersectionPoint0, intersectionPoint1};
e.positiveCount = positiveCount;
e.negativeCount = negativeCount;
e.zeroCount = zeroCount;
return std::move(e);
}

if(p1.planeSplitsEdge(p2) && p2.planeSplitsEdge(p0))
{
auto intersectionPoint0 = computeIntersection({v1, v2}, p1.SignedDistance, p2.SignedDistance);
auto intersectionPoint1 = computeIntersection({v2, v0}, p2.SignedDistance, p0.SignedDistance);
return {intersectionPoint0, intersectionPoint1};
Edge e{intersectionPoint0, intersectionPoint1};
e.positiveCount = positiveCount;
e.negativeCount = negativeCount;
e.zeroCount = zeroCount;
return std::move(e);
}

// No valid intersection found
return {}; // Invalid edge
Edge e;
e.positiveCount = positiveCount;
e.negativeCount = negativeCount;
e.zeroCount = zeroCount;
return std::move(e); // Invalid edge
}
} // namespace slice_helper

Expand Down Expand Up @@ -513,14 +567,15 @@ GeometryUtilities::SliceTriangleReturnType GeometryUtilities::SliceTriangleGeome
std::vector<int32> regionIds;

int32 edgeCounter = 0;
int32 sliceIndex = 0;
int32 sliceIndex = -1;
// Loop over each slice plane
for(float zValue = zStart; zValue <= zEnd; zValue = zValue + sliceSpacing)
{
if(shouldCancel)
{
break;
}
sliceIndex++;
d = zValue;

// Define a plane with a normal vector and a point on the plane
Expand All @@ -529,6 +584,7 @@ GeometryUtilities::SliceTriangleReturnType GeometryUtilities::SliceTriangleGeome

// Create the plane
slice_helper::Plane plane(planeNormal, pointOnPlane);

// Loop over each Triangle and get edges/vertices of any intersection
for(usize triIdx = 0; triIdx < numTris; triIdx++)
{
Expand Down Expand Up @@ -559,9 +615,10 @@ GeometryUtilities::SliceTriangleReturnType GeometryUtilities::SliceTriangleGeome
}
edgeCounter++;
}
}
} // END TRIANGLE LOOP

sliceIndex++;
}
} // END SLICE LOOP

return {std::move(slicedVerts), std::move(sliceIds), std::move(regionIds), numberOfSlices};
}
2 changes: 1 addition & 1 deletion src/simplnx/Utilities/ImageRotationUtilities.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -278,7 +278,7 @@ class FilterProgressCallback
auto now = std::chrono::steady_clock::now();
if(std::chrono::duration_cast<std::chrono::milliseconds>(now - m_InitialTime).count() > 1000)
{
m_MessageHandler(IFilter::ProgressMessage{IFilter::Message::Type::Progress, fmt::format("Nodes Completed: {}", m_Progcounter), m_Progcounter});
m_MessageHandler(IFilter::Message{IFilter::Message::Type::Info, fmt::format("Nodes Completed: {}", m_Progcounter)});
m_InitialTime = std::chrono::steady_clock::now();
}
}
Expand Down

0 comments on commit 850e97e

Please sign in to comment.