Round-off errors in vtkTriangle::TriangleArea()
This issue was created automatically from an original Mantis Issue. Further discussion may take place here.
vtkTriangle::TrianlgeArea() is implemented by using the Heron's formula, which may cause severe round-off errors at corner cases. The error originates from the fact that sqrt(Machine_Epsilon) = 1e-8, which is far larger than the Machine_Epsilon. The proposed replacement based on the cross-product requires slightly less computation than the original.
----- Sample code in Python ----- import numpy as np import vtk
def TriArea(p1, p2, p3) : n = np.ndarray(3) vtk.vtkTriangle.ComputeNormalDirection(p1, p2, p3, n) return vtk.vtkMath.Norm(n) / 2
Very sharp triangle.
We expect an area to be the order of the machine epsilon.
p1 = np.array([ 0.0, 0.0, 0.0 ]) p2 = np.array([ 1.0e-13, 1.0, 0.0 ]) p3 = np.array([ -1.0e-13, 1.0, 0.0 ])
print "orig: {0:22.16e}".format( vtk.vtkTriangle.TriangleArea(p1, p2, p3) ) print " new: {0:22.16e}".format( TriArea(p1, p2, p3) ) print
Rotate the triangle
tr = vtk.vtkTransform() tr.RotateZ(30) tr.TransformNormal(p1, p1) tr.TransformNormal(p2, p2) tr.TransformNormal(p3, p3)
print "orig: {0:22.16e}".format( vtk.vtkTriangle.TriangleArea(p1, p2, p3) ) print " new: {0:22.16e}".format( TriArea(p1, p2, p3) )
----- Result ----- orig: 0.0000000000000000e+00 new: 1.0000000000000000e-13
orig: 7.4505805969238281e-09 new: 9.9999599856604815e-14
Note: This is a follow-up of the report by Masato, http://public.kitware.com/pipermail/vtkusers/2014-August/084947.html