diff --git a/src/pe/raytracing/Intersects.h b/src/pe/raytracing/Intersects.h
index f1dd5d0017407e96ff72ba10f9bf99428f93ae4e..8e5e38cf5802994b55a405f97be64bf7b890c33b 100644
--- a/src/pe/raytracing/Intersects.h
+++ b/src/pe/raytracing/Intersects.h
@@ -28,15 +28,12 @@
 #include "pe/rigidbody/Sphere.h"
 #include "pe/rigidbody/Union.h"
 #include "pe/utility/BodyCast.h"
+#include <core/math/Utility.h>
 #include <boost/math/special_functions/sign.hpp>
 #include <boost/tuple/tuple.hpp>
 
 #include <pe/raytracing/Ray.h>
 
-#define EPSILON real_t(1e-4)
-
-using namespace boost::math;
-
 namespace walberla {
 namespace pe {
 namespace raytracing {
@@ -49,7 +46,7 @@ inline bool intersects(const BodyID body, const Ray& ray, real_t& t, Vec3& n);
    
 inline bool intersects(const AABB& aabb, const Ray& ray, real_t& t, real_t padding = real_t(0.0), Vec3* n = NULL);
 inline bool intersectsSphere(const Vec3& gpos, real_t radius, const Ray& ray, real_t& t0, real_t& t1);
-   
+
 struct IntersectsFunctor
 {
    const Ray& ray_;
@@ -64,6 +61,8 @@ struct IntersectsFunctor
    }
 };
 
+const real_t discriminantEps = real_t(1e-4);
+   
 inline bool intersects(const SphereID sphere, const Ray& ray, real_t& t, Vec3& n) {
    real_t inf = std::numeric_limits<real_t>::max();
    
@@ -90,12 +89,12 @@ inline bool intersects(const PlaneID plane, const Ray& ray, real_t& t, Vec3& n)
    const Vec3& origin = ray.getOrigin();
    const Vec3& planeNormal = plane->getNormal();
    real_t denominator = planeNormal * direction;
-   if (std::abs(denominator) > EPSILON) {
+   if (std::abs(denominator) > discriminantEps) {
       real_t t_;
       t_ = ((plane->getPosition() - origin) * planeNormal) / denominator;
-      if (t_ > EPSILON) {
+      if (t_ > 0) {
          t = t_;
-         n = planeNormal * sign(-denominator);
+         n = planeNormal * walberla::math::sign(-denominator);
          return true;
       } else {
          t = inf;
@@ -196,7 +195,7 @@ inline bool intersects(const CapsuleID capsule, const Ray& ray, real_t& t, Vec3&
    real_t b = real_t(2)*origin[2]*direction[2] + real_t(2)*origin[1]*direction[1];
    real_t c = origin[2]*origin[2] + origin[1]*origin[1] - capsule->getRadius()*capsule->getRadius();
    real_t discriminant = b*b - real_t(4.)*a*c;
-   if (discriminant >= EPSILON) {
+   if (std::abs(discriminant) >= discriminantEps) {
       // With discriminant smaller than 0, cylinder is not hit by ray (no solution for quadratic equation).
       // Thus only enter this section if the equation is actually solvable.
       
diff --git a/src/pe/raytracing/Raytracer.cpp b/src/pe/raytracing/Raytracer.cpp
index 3ae2e9454a036bf8a864dca4d41d2c2d1fe38415..1808bc242ba86daef11c08fdfabe2d3fb4f29993 100644
--- a/src/pe/raytracing/Raytracer.cpp
+++ b/src/pe/raytracing/Raytracer.cpp
@@ -21,12 +21,6 @@
 #include "Raytracer.h"
 #include "geometry/structured/extern/lodepng.h"
 
-using namespace walberla;
-
-real_t deg2rad(real_t deg) {
-   return deg * math::M_PI / real_t(180.0);
-}
-
 namespace walberla {
 namespace pe {
 namespace raytracing {
@@ -194,7 +188,8 @@ void Raytracer::setupView_() {
    // viewing plane setup
    d_ = (cameraPosition_ - lookAtPoint_).length();
    aspectRatio_ = real_t(pixelsHorizontal_) / real_t(pixelsVertical_);
-   viewingPlaneHeight_ = real_c(tan(deg2rad(fov_vertical_)/real_t(2.))) * real_t(2.) * d_;
+   real_t fov_vertical_rad = fov_vertical_ * math::M_PI / real_t(180.0);
+   viewingPlaneHeight_ = real_c(tan(fov_vertical_rad/real_t(2.))) * real_t(2.) * d_;
    viewingPlaneWidth_ = viewingPlaneHeight_ * aspectRatio_;
    viewingPlaneOrigin_ = lookAtPoint_ - u_*viewingPlaneWidth_/real_t(2.) - v_*viewingPlaneHeight_/real_t(2.);