What is the use of epsilon in the Möller–Trumbore intersection algorithm?

The Möller–Trumbore intersection algorithm compares the determinant with an epsilon. What is this epsilon and how is its value chosen?

The epsilon is there to counter floating-point accuracy problems. Note that the code checks that the algorithm is between $-\epsilon$ and $\epsilon$; that is, close to zero. If the determinant is zero, then the ray lies in the plane of the triangle and can be safely said to not intersect it (as the alternative, an edge-on intersection, should be caught by neighbouring triangles in the mesh).
The value of the epsilon depends on the floating-point precision you're using. It should be small enough that only very grazing intersections are culled, but large enough that it avoids the speckling problem. It can be determined by experiment, or by just eyeballing it. The implementation you linked to uses $\epsilon = 10^{-7}$, which is a reasonable choice for single-precision arithmetic; you may want to use $10^{-9}$ or $10^{-10}$ for double-precision computations.