diff --git a/pytorch3d/csrc/utils/geometry_utils.cuh b/pytorch3d/csrc/utils/geometry_utils.cuh index 167463962..2ba60d277 100644 --- a/pytorch3d/csrc/utils/geometry_utils.cuh +++ b/pytorch3d/csrc/utils/geometry_utils.cuh @@ -461,6 +461,27 @@ PointTriangleDistanceBackward( // vec3 utils // // ************************************************************* // +// Computes the area of a triangle (v0, v1, v2). +// +// Args: +// v0, v1, v2: vec3 coordinates of the triangle vertices +// +// Returns +// area: float: The area of the triangle +// +__device__ inline float +AreaOfTriangle(const float3& v0, const float3& v1, const float3& v2) { + float3 p0 = v1 - v0; + float3 p1 = v2 - v0; + + // compute the hypotenus of the scross product (p0 x p1) + float dd = hypot( + p0.y * p1.z - p0.z * p1.y, + hypot(p0.z * p1.x - p0.x * p1.z, p0.x * p1.y - p0.y * p1.x)); + + return dd / 2.0; +} + // Computes the barycentric coordinates of a point p relative // to a triangle (v0, v1, v2), i.e. p = w0 * v0 + w1 * v1 + w2 * v2 // s.t. w0 + w1 + w2 = 1.0 @@ -503,6 +524,7 @@ __device__ inline float3 BarycentricCoords3Forward( // Checks whether the point p is inside the triangle (v0, v1, v2). // A point is inside the triangle, if all barycentric coordinates // wrt the triangle are >= 0 & <= 1. +// If the triangle is degenerate, aka line or point, then return False. // // NOTE that this function assumes that p lives on the space spanned // by (v0, v1, v2). @@ -521,11 +543,16 @@ __device__ inline bool IsInsideTriangle( const float3& v0, const float3& v1, const float3& v2) { - float3 bary = BarycentricCoords3Forward(p, v0, v1, v2); - bool x_in = 0.0f <= bary.x && bary.x <= 1.0f; - bool y_in = 0.0f <= bary.y && bary.y <= 1.0f; - bool z_in = 0.0f <= bary.z && bary.z <= 1.0f; - bool inside = x_in && y_in && z_in; + bool inside; + if (AreaOfTriangle(v0, v1, v2) < 1e-5) { + inside = 0; + } else { + float3 bary = BarycentricCoords3Forward(p, v0, v1, v2); + bool x_in = 0.0f <= bary.x && bary.x <= 1.0f; + bool y_in = 0.0f <= bary.y && bary.y <= 1.0f; + bool z_in = 0.0f <= bary.z && bary.z <= 1.0f; + inside = x_in && y_in && z_in; + } return inside; } diff --git a/pytorch3d/csrc/utils/geometry_utils.h b/pytorch3d/csrc/utils/geometry_utils.h index fd95a815f..ef1f437d3 100644 --- a/pytorch3d/csrc/utils/geometry_utils.h +++ b/pytorch3d/csrc/utils/geometry_utils.h @@ -560,6 +560,26 @@ PointTriangleDistanceBackward( return std::make_tuple(grad_p, grad_v0, grad_v1, grad_v2); } +// Computes the area of a triangle (v0, v1, v2). +// Args: +// v0, v1, v2: vec3 coordinates of the triangle vertices +// +// Returns: +// area: float: the area of the triangle +// +template +T AreaOfTriangle(const vec3& v0, const vec3& v1, const vec3& v2) { + vec3 p0 = v1 - v0; + vec3 p1 = v2 - v0; + + // compute the hypotenus of the scross product (p0 x p1) + float dd = std::hypot( + p0.y * p1.z - p0.z * p1.y, + std::hypot(p0.z * p1.x - p0.x * p1.z, p0.x * p1.y - p0.y * p1.x)); + + return dd / 2.0; +} + // Computes the squared distance of a point p relative to a triangle (v0, v1, // v2). If the point's projection p0 on the plane spanned by (v0, v1, v2) is // inside the triangle with vertices (v0, v1, v2), then the returned value is @@ -604,6 +624,7 @@ vec3 BarycentricCoords3Forward( // Checks whether the point p is inside the triangle (v0, v1, v2). // A point is inside the triangle, if all barycentric coordinates // wrt the triangle are >= 0 & <= 1. +// If the triangle is degenerate, aka line or point, then return False. // // NOTE that this function assumes that p lives on the space spanned // by (v0, v1, v2). @@ -623,11 +644,16 @@ static bool IsInsideTriangle( const vec3& v0, const vec3& v1, const vec3& v2) { - vec3 bary = BarycentricCoords3Forward(p, v0, v1, v2); - bool x_in = 0.0f <= bary.x && bary.x <= 1.0f; - bool y_in = 0.0f <= bary.y && bary.y <= 1.0f; - bool z_in = 0.0f <= bary.z && bary.z <= 1.0f; - bool inside = x_in && y_in && z_in; + bool inside; + if (AreaOfTriangle(v0, v1, v2) < 1e-5) { + inside = 0; + } else { + vec3 bary = BarycentricCoords3Forward(p, v0, v1, v2); + bool x_in = 0.0f <= bary.x && bary.x <= 1.0f; + bool y_in = 0.0f <= bary.y && bary.y <= 1.0f; + bool z_in = 0.0f <= bary.z && bary.z <= 1.0f; + inside = x_in && y_in && z_in; + } return inside; } diff --git a/tests/test_point_mesh_distance.py b/tests/test_point_mesh_distance.py index bef5d9af7..6f35af0fb 100644 --- a/tests/test_point_mesh_distance.py +++ b/tests/test_point_mesh_distance.py @@ -96,7 +96,7 @@ def _point_to_bary(point: torch.Tensor, tri: torch.Tensor) -> torch.Tensor: d20 = v2.dot(v0) d21 = v2.dot(v1) - denom = d00 * d11 - d01 * d01 + denom = d00 * d11 - d01 * d01 + TestPointMeshDistance.eps() s2 = (d11 * d20 - d01 * d21) / denom s3 = (d00 * d21 - d01 * d20) / denom s1 = 1.0 - s2 - s3 @@ -117,6 +117,13 @@ def _is_inside_triangle(point: torch.Tensor, tri: torch.Tensor) -> torch.Tensor: Returns: inside: BoolTensor of shape (1) """ + v0 = tri[1] - tri[0] + v1 = tri[2] - tri[0] + area = torch.cross(v0, v1).norm() / 2.0 + + # check if triangle is a line or a point. In that case, return False + if area < 1e-5: + return False bary = TestPointMeshDistance._point_to_bary(point, tri) inside = ((bary >= 0.0) * (bary <= 1.0)).all() return inside @@ -836,6 +843,28 @@ def test_point_mesh_face_distance(self): ) self.assertClose(pcls.points_list()[i].grad, pcls_op.points_list()[i].grad) + def test_small_faces_case(self): + for device in [torch.device("cpu"), torch.device("cuda:0")]: + mesh_vertices = torch.tensor( + [ + [-0.0021, -0.3769, 0.7146], + [-0.0161, -0.3771, 0.7146], + [-0.0021, -0.3771, 0.7147], + ], + dtype=torch.float32, + device=device, + ) + mesh1_faces = torch.tensor([[0, 2, 1]], device=device) + mesh2_faces = torch.tensor([[2, 0, 1]], device=device) + pcd_points = torch.tensor([[-0.3623, -0.5340, 0.7727]], device=device) + mesh1 = Meshes(verts=[mesh_vertices], faces=[mesh1_faces]) + mesh2 = Meshes(verts=[mesh_vertices], faces=[mesh2_faces]) + pcd = Pointclouds(points=[pcd_points]) + + loss1 = point_mesh_face_distance(mesh1, pcd) + loss2 = point_mesh_face_distance(mesh2, pcd) + self.assertClose(loss1, loss2) + @staticmethod def point_mesh_edge(N: int, V: int, F: int, P: int, device: str): device = torch.device(device)