diff --git a/src/pe/raytracing/Raytracer.cpp b/src/pe/raytracing/Raytracer.cpp
index 5f230265024392603f907ecaf0ff33feb5fe6476..6996b7483f056f6c92c3910b821d52c1f3f08968 100644
--- a/src/pe/raytracing/Raytracer.cpp
+++ b/src/pe/raytracing/Raytracer.cpp
@@ -30,6 +30,29 @@ real_t deg2rad(real_t deg) {
 namespace walberla {
 namespace pe {
 namespace raytracing {
+   
+void BodyIntersectionInfo_Comparator_MPI_OP( BodyIntersectionInfo *in, BodyIntersectionInfo *inout, int *len, MPI_Datatype *dptr) {
+   for (int i = 0; i < *len; ++i) {
+      if (in->bodySystemID != 0 && inout->bodySystemID != 0) {
+         WALBERLA_ASSERT(in->imageX == inout->imageX && in->imageY == inout->imageY, "coordinates of infos do not match: " << in->imageX << "/" << in->imageY << " and " << inout->imageX << "/" << inout->imageY);
+      }
+      
+      if ((in->t < inout->t && in->bodySystemID != 0) || (inout->bodySystemID == 0 && in->bodySystemID != 0)) {
+         // info in "in" is closer than the one in "inout" -> update inout to values of in
+         inout->imageX = in->imageX;
+         inout->imageY = in->imageY;
+         inout->bodySystemID = in->bodySystemID;
+         inout->t = in->t;
+         inout->r = in->r;
+         inout->g = in->g;
+         inout->b = in->b;
+      }
+      
+      in++;
+      inout++;
+   }
+}
+   
 /*!\brief Instantiation constructor for the Raytracer class.
  *
  * \param forest BlockForest the raytracer operates on.
@@ -71,6 +94,7 @@ Raytracer::Raytracer(const shared_ptr<BlockStorage> forest, const BlockDataID st
    
    setupView_();
    setupFilenameRankWidth_();
+   setupMPI_();
 }
 
 /*!\brief Instantiation constructor for the Raytracer class using a config object for view setup.
@@ -105,6 +129,8 @@ Raytracer::Raytracer(const shared_ptr<BlockStorage> forest, const BlockDataID st
       setTBufferOutputEnabled(true);
       setTBufferOutputDirectory(config.getParameter<std::string>("tbuffer_output_directory"));
       WALBERLA_LOG_INFO_ON_ROOT("t buffers will be written to " << getTBufferOutputDirectory() << ".");
+   } else {
+      setTBufferOutputEnabled(false);
    }
    
    setLocalImageOutputEnabled(config.getParameter<bool>("local_image_output_enabled", false));
@@ -126,9 +152,17 @@ Raytracer::Raytracer(const shared_ptr<BlockStorage> forest, const BlockDataID st
    backgroundColor_ = config.getParameter<Color>("backgroundColor", Vec3(0.1, 0.1, 0.1));
 
    blockAABBIntersectionPadding_ = config.getParameter<real_t>("blockAABBIntersectionPadding", real_t(0.0));
+
+   std::string reductionMethod = config.getParameter<std::string>("reductionMethod", "MPI_REDUCE");
+   if (reductionMethod == "MPI_REDUCE") {
+      reductionMethod_ = MPI_REDUCE;
+   } else if (reductionMethod == "MPI_GATHER") {
+      reductionMethod_ = MPI_GATHER;
+   }
       
    setupView_();
    setupFilenameRankWidth_();
+   setupMPI_();
 }
 
 /*!\brief Utility function for setting up the view plane and calculating required variables.
@@ -157,6 +191,39 @@ void Raytracer::setupFilenameRankWidth_() {
    filenameRankWidth_ = uint8_c(log10(numProcesses))+1;
 }
 
+void Raytracer::setupMPI_() {
+   MPI_Op_create((MPI_User_function *)BodyIntersectionInfo_Comparator_MPI_OP, true, &bodyIntersectionInfo_reduction_op);
+   
+   const int nblocks = 7;
+   const int blocklengths[nblocks] = {1,1,1,1,1,1,1};
+   MPI_Datatype types[nblocks] = {
+      MPI_UNSIGNED, // for coordinate
+      MPI_UNSIGNED, // for coordinate
+      MPI_UNSIGNED_LONG_LONG, // for id
+      MPI_DOUBLE, // for distance
+      MPI_DOUBLE, // for color
+      MPI_DOUBLE, // for color
+      MPI_DOUBLE // for color
+   };
+   MPI_Aint displacements[nblocks];
+   displacements[0] = offsetof(BodyIntersectionInfo, imageX);
+   displacements[1] = offsetof(BodyIntersectionInfo, imageY);
+   displacements[2] = offsetof(BodyIntersectionInfo, bodySystemID);
+   displacements[3] = offsetof(BodyIntersectionInfo, t);
+   displacements[4] = offsetof(BodyIntersectionInfo, r);
+   displacements[5] = offsetof(BodyIntersectionInfo, g);
+   displacements[6] = offsetof(BodyIntersectionInfo, b);
+   
+   MPI_Datatype tmp_type;
+   MPI_Type_create_struct(nblocks, blocklengths, displacements, types, &tmp_type);
+   
+   MPI_Aint lb, extent;
+   MPI_Type_get_extent( tmp_type, &lb, &extent );
+   MPI_Type_create_resized( tmp_type, lb, extent, &bodyIntersectionInfo_mpi_type );
+   
+   MPI_Type_commit(&bodyIntersectionInfo_mpi_type);
+}
+   
 /*!\brief Generates the filename for output files.
  * \param base String that precedes the timestap and rank info.
  * \param timestep Timestep this image is from.
@@ -286,6 +353,64 @@ void Raytracer::writeImageBufferToFile(const std::vector<Color>& imageBuffer, co
       WALBERLA_LOG_WARNING("lodePNG error " << error << " when trying to save image file to " << fullPath.string() << ": " << lodepng_error_text(error));
    }
 }
+   
+void Raytracer::syncImageUsingMPIReduce(std::vector<BodyIntersectionInfo>& intersectionsBuffer, WcTimingTree* tt) {
+   WALBERLA_MPI_BARRIER();
+   if (tt != NULL) tt->start("Reduction");
+   int rank = mpi::MPIManager::instance()->rank();
+
+   const int recvRank = 0;
+   if( rank == recvRank ) {
+      MPI_Reduce(MPI_IN_PLACE,
+                 &intersectionsBuffer[0], int_c(intersectionsBuffer.size()),
+                 bodyIntersectionInfo_mpi_type, bodyIntersectionInfo_reduction_op,
+                 recvRank, MPI_COMM_WORLD);
+   } else {
+      MPI_Reduce(&intersectionsBuffer[0], 0, int_c(intersectionsBuffer.size()),
+                 bodyIntersectionInfo_mpi_type, bodyIntersectionInfo_reduction_op,
+                 recvRank, MPI_COMM_WORLD);
+   }
+   
+   WALBERLA_MPI_BARRIER();
+   if (tt != NULL) tt->stop("Reduction");
+}
+   
+void Raytracer::syncImageUsingMPIGather(std::vector<BodyIntersectionInfo>& intersections, std::vector<BodyIntersectionInfo>& intersectionsBuffer, WcTimingTree* tt) {
+   WALBERLA_MPI_BARRIER();
+   if (tt != NULL) tt->start("Reduction");
+   
+   mpi::SendBuffer sendBuffer;
+   for (auto& info: intersections) {
+      sendBuffer << info.imageX << info.imageY
+      << info.bodySystemID << info.t
+      << info.r << info.g << info.b;
+   }
+
+   mpi::RecvBuffer recvBuffer;
+   mpi::gathervBuffer(sendBuffer, recvBuffer, 0);
+   
+   WALBERLA_ROOT_SECTION() {
+      BodyIntersectionInfo info;
+      while (!recvBuffer.isEmpty()) {
+         recvBuffer >> info.imageX;
+         recvBuffer >> info.imageY;
+         recvBuffer >> info.bodySystemID;
+         recvBuffer >> info.t;
+         recvBuffer >> info.r;
+         recvBuffer >> info.g;
+         recvBuffer >> info.b;
+         
+         size_t i = coordinateToArrayIndex(info.imageX, info.imageY);
+         
+         if (intersectionsBuffer[i].bodySystemID == 0 || info.t < intersectionsBuffer[i].t) {
+            intersectionsBuffer[i] = info;
+         }
+      }
+   }
+   
+   WALBERLA_MPI_BARRIER();
+   if (tt != NULL) tt->stop("Reduction");
+}
 }
 }
 }
diff --git a/src/pe/raytracing/Raytracer.h b/src/pe/raytracing/Raytracer.h
index 2e9727287edbe99520acbbc046bfa7b2c149583e..2901fca22659526a341292bc6ec8ba057a5236bf 100644
--- a/src/pe/raytracing/Raytracer.h
+++ b/src/pe/raytracing/Raytracer.h
@@ -35,6 +35,11 @@
 #include "pe/ccd/ICCD.h"
 #include <pe/ccd/HashGrids.h>
 
+#include "core/mpi/MPIManager.h"
+#include "core/mpi/MPIWrapper.h"
+#include "core/mpi/all.h"
+#include <stddef.h>
+
 using namespace walberla;
 using namespace walberla::pe;
 using namespace walberla::timing;
@@ -46,14 +51,18 @@ namespace raytracing {
 /*!\brief Contains information about a ray-body intersection.
  */
 struct BodyIntersectionInfo {
-   size_t imageX;                //!< viewing plane x pixel coordinate the ray belongs to.
-   size_t imageY;                //!< viewing plane y pixel coordinate the ray belongs to.
-   walberla::id_t bodySystemID;  //!< system ID of body which was hit.
-   real_t t;                     //!< distance from camera to intersection point on body.
-   Vec3 color;                   //!< color computed for the pixel.
+   unsigned int imageX;          //!< Viewing plane x pixel coordinate the ray belongs to.   -> MPI_UNSIGNED
+   unsigned int imageY;          //!< Viewing plane y pixel coordinate the ray belongs to.   -> MPI_UNSIGNED
+   walberla::id_t bodySystemID;  //!< System ID of body which was hit.                       -> MPI_UNSIGNED_LONG_LONG
+   double t;                     //!< Distance from camera to intersection point on body.    -> MPI_DOUBLE
+   double r;                     //!< Red value for the pixel.                               -> MPI_DOUBLE
+   double g;                     //!< Green value for the pixel.                             -> MPI_DOUBLE
+   double b;                     //!< Blue value for the pixel.                              -> MPI_DOUBLE
 };
-
+   
 class Raytracer {
+   enum ReductionMethod { MPI_REDUCE, MPI_GATHER };
+   
 public:
    enum Algorithm { RAYTRACE_HASHGRIDS, RAYTRACE_NAIVE, RAYTRACE_COMPARE_BOTH };
    
@@ -127,6 +136,9 @@ private:
    real_t pixelHeight_;       //!< The height of a pixel of the generated image in the viewing plane.
    //@}
    
+   MPI_Op bodyIntersectionInfo_reduction_op;
+   MPI_Datatype bodyIntersectionInfo_mpi_type;
+   
 public:
    /*!\name Get functions */
    //@{
@@ -164,6 +176,7 @@ public:
    
    void setupView_();
    void setupFilenameRankWidth_();
+   void setupMPI_();
    
 private:
    std::string getOutputFilename(const std::string& base, size_t timestep, bool isGlobalImage) const;
@@ -172,6 +185,9 @@ private:
    void writeImageBufferToFile(const std::vector<Color>& imageBuffer, size_t timestep, bool isGlobalImage = false) const;
    void writeImageBufferToFile(const std::vector<Color>& imageBuffer, const std::string& fileName) const;
    
+   void syncImageUsingMPIReduce(std::vector<BodyIntersectionInfo>& intersectionsBuffer, WcTimingTree* tt = NULL);
+   void syncImageUsingMPIGather(std::vector<BodyIntersectionInfo>& intersections, std::vector<BodyIntersectionInfo>& intersectionsBuffer, WcTimingTree* tt);
+   
    inline bool isPlaneVisible(const PlaneID plane, const Ray& ray) const;
    inline size_t coordinateToArrayIndex(size_t x, size_t y) const;
    
@@ -557,12 +573,16 @@ void Raytracer::generateImage(const size_t timestep, WcTimingTree* tt) {
             
             imageBuffer[coordinateToArrayIndex(x, y)] = color;
             BodyIntersectionInfo intersectionInfo = {
-               x,
-               y,
+               uint32_t(x),
+               uint32_t(y),
                body_closest->getSystemID(),
                t_closest,
-               color
+               color[0],
+               color[1],
+               color[2]
             };
+            
+            intersectionsBuffer[coordinateToArrayIndex(x, y)] = intersectionInfo;
             intersections.push_back(intersectionInfo);
          }
          
@@ -634,10 +654,6 @@ void Raytracer::generateImage(const size_t timestep, WcTimingTree* tt) {
          gatheredIntersectionCount++;
       }
    }
-   if (tt != NULL) tt->stop("Reduction");
-
-   //WALBERLA_LOG_INFO("#particles visible: " << visibleBodyIDs.size());
-   WALBERLA_LOG_INFO_ON_ROOT("#gathered intersections: " << gatheredIntersectionCount);
    
    if (tt != NULL) tt->start("Output");
    if (getImageOutputEnabled()) {
@@ -645,6 +661,12 @@ void Raytracer::generateImage(const size_t timestep, WcTimingTree* tt) {
          writeImageBufferToFile(imageBuffer, timestep);
       }
       WALBERLA_ROOT_SECTION() {
+         std::vector<Color> fullImageBuffer(pixelsHorizontal_ * pixelsVertical_, backgroundColor_);
+
+         for (auto& info: intersectionsBuffer) {
+            fullImageBuffer[coordinateToArrayIndex(info.imageX, info.imageY)] = Color(info.r, info.g, info.b);
+         }
+         
          writeImageBufferToFile(fullImageBuffer, timestep, true);
       }
    }
@@ -652,9 +674,16 @@ void Raytracer::generateImage(const size_t timestep, WcTimingTree* tt) {
    if (getTBufferOutputEnabled()) {
       writeTBufferToFile(tBuffer, timestep);
       WALBERLA_ROOT_SECTION() {
+         std::vector<real_t> fullTBuffer(pixelsHorizontal_ * pixelsVertical_, inf);
+
+         for (auto& info: intersectionsBuffer) {
+            fullTBuffer[coordinateToArrayIndex(info.imageX, info.imageY)] = info.t;
+         }
+         
          writeTBufferToFile(fullTBuffer, timestep, true);
       }
    }
+   
    if (tt != NULL) tt->stop("Output");
    if (tt != NULL) tt->stop("Raytracing");
 }